Mean Field Multi-Agent Reinforcement Learning

Yaodong Yang, Rui Luo, Minne Li, Ming Zhou, Weinan Zhang, Jun Wang

Introduction

Multi-agent reinforcement learning (MARL) is concerned with a set of autonomous agents that share a common environment (Busoniu et al., 2008). Learning in MARL is fundamentally difficult since agents not only interact with the environment but also with each other. Independent QQ-learning (Tan, 1993) that considers other agents as a part of the environment often fails as the multi-agent setting breaks the theoretical convergence guarantee and makes the learning unstable: changes in the policy of one agent will affect that of the others, and vice versa (Matignon et al., 2012).

Instead, accounting for the extra information from conjecturing the policies of other agents is beneficial to each single learner (Foerster et al., 2017; Lowe et al., 2017a). Studies show that an agent who learns the effect of joint actions has better performance than those who do not in many scenarios, including cooperative games (Panait & Luke, 2005), zero-sum stochastic games (Littman, 1994), and general-sum stochastic games (Littman, 2001; Hu & Wellman, 2003).

The existing equilibrium-solving approaches, although principled, are only capable of solving a handful of agents (Hu & Wellman, 2003; Bowling & Veloso, 2002). The computational complexity of directly solving (Nash) equilibrium would prevent them from applying to the situations with a large group or even a population of agents. Yet, in practice, many cases do require strategic interactions among a large number of agents, such as the gaming bots in Massively Multiplayer Online Role-Playing Game (Jeong et al., 2015), the trading agents in stock markets (Troy, 1997), or the online advertising bidding agents (Wang et al., 2017).

In this paper, we tackle MARL when a large number of agents co-exist. We consider a setting where each agent is directly interacting with a finite set of other agents; through a chain of direct interactions, any pair of agents is interconnected globally (Blume, 1993). The scalability is solved by employing Mean Field Theory (Stanley, 1971) – the interactions within the population of agents are approximated by that of a single agent played with the average effect from the overall (local) population. The learning is mutually reinforced between two entities rather than many entities: the learning of the individual agent’s optimal policy is based on the dynamics of the agent population, meanwhile, the dynamics of the population is updated according to the individual policies. Based on such formulation, we develop practical mean field QQ-learning and mean field Actor-Critic algorithms, and discuss the convergence of our solution under certain assumptions. Our experiment on a simple multi-agent resource allocation shows that our mean field MARL is capable of learning over many-agent interactions when others fail. We also demonstrate that with temporal-difference learning, mean field MARL manages to learn and solve the Ising model without even explicitly knowing the energy function. At last, in a mixed cooperative-competitive battle game, we show that the mean field MARL achieves high winning rates against other baselines previously reported for many agent systems.

Preliminary

MARL intersects between reinforcement learning and game theory. The marriage of the two gives rise to the general framework of stochastic game (Shapley, 1953).

The agents choose actions according to their policies, also known as strategies. For agent jj, the corresponding policy is defined as πj:S→Ω(Aj)\pi^{j}:\mathcal{S}\to\Omega(\mathcal{A}^{j}), where Ω(Aj)\Omega(\mathcal{A}^{j}) is the collection of probability distributions over agent jj’s action space Aj\mathcal{A}^{j}. Let \vectorsymπ≜[π1,…,πN]\vectorsym{\pi}\triangleq[\pi^{1},\dots,\pi^{N}] denote the joint policy of all agents; we assume, as one usually does, \vectorsymπ\vectorsym{\pi} to be time-independent, which is referred to be stationary. Provided an initial state ss, the value function of agent jj under the joint policy \vectorsymπ\vectorsym{\pi} is written as the expected cumulative discounted future reward:

where s′s^{\prime} is the state at the next time step. The value function v\vectorsymπjv^{j}_{\vectorsym{\pi}} can be expressed in terms of the QQ-function in Eq. (2) as

The QQ-function for NN-agent game in Eq. (2) extends the formulation for single-agent game by considering the joint action taken by all agents \vectorsyma≜[a1,…,aN]\vectorsym{a}\triangleq[a^{1},\dots,a^{N}], and by taking the expectation over the joint action in Eq. (3).

We formulate MARL by the stochastic game with a discrete-time non-cooperative setting, i.e. no explicit coalitions are considered. The game is assumed to be incomplete but to have perfect information (Littman, 1994), i.e. each agent knows neither the game dynamics nor the reward functions of others, but it is able to observe and react to the previous actions and the resulting immediate rewards of other agents.

2 Nash Q𝑄Q-learning

Here we adopt a compact notation for the joint policy of all agents except jj as \vectorsymπ∗−j≜[π∗1,…,π∗j−1,π∗j+1,…,π∗N]\vectorsym{\pi}^{-j}_{*}\triangleq[\pi^{1}_{*},\dots,\pi^{j-1}_{*},\pi^{j+1}_{*},\dots,\pi^{N}_{*}].

In a Nash equilibrium, each agent acts with the best response π∗j\pi^{j}_{*} to others, provided that all other agents follow the policy \vectorsymπ∗−j\vectorsym{\pi}^{-j}_{*}. It has been shown that, for a NN-agent stochastic game, there is at least one Nash equilibrium with stationary policies (Fink et al., 1964). Given a Nash policy \vectorsymπ∗\vectorsym{\pi}_{*}, the Nash value function \vectorsymv ⁣Nash(s)≜[v\vectorsymπ∗1(s),…,v\vectorsymπ∗N(s)]\vectorsym{v}^{\mathop{}\!\mathtt{Nash}}(s)\triangleq[v^{1}_{\vectorsym{\pi}_{*}}(s),\dots,v^{N}_{\vectorsym{\pi}_{*}}(s)] is calculated with all agents following \vectorsymπ∗\vectorsym{\pi}_{*} from the initial state ss onward.

Nash QQ-learning (Hu & Wellman, 2003) defines an iterative procedure with two alternating steps for computing the Nash policy: 1) solving the Nash equilibrium of the current stage game defined by {\vectorsymQt}\{\vectorsym{Q}_{t}\} using the Lemke-Howson algorithm (Lemke & Howson, 1964), 2) improving the estimation of the QQ-function with the new Nash equilibrium value. It can be proved that under certain assumptions, the Nash operator H ⁣Nash\mathcal{H}^{\mathop{}\!\mathtt{Nash}} defined by the following expression

forms a contraction mapping, where \vectorsymQ≜[Q1,…,QN]\vectorsym{Q}\triangleq[Q^{1},\dots,Q^{N}], and \vectorsymr(s,\vectorsyma)≜[r1(s,\vectorsyma),…,rN(s,\vectorsyma)]\vectorsym{r}(s,\vectorsym{a})\triangleq[r^{1}(s,\vectorsym{a}),\dots,r^{N}(s,\vectorsym{a})]. The QQ-function will eventually converge to the value received in a Nash equilibrium of the game, referred to as the Nash QQ-value.

Mean Field MARL

The dimension of joint action \vectorsyma\vectorsym{a} grows proportionally w.r.t. the number of agents NN. As all agents act strategically and evaluate simultaneously their value functions based on the joint actions, it becomes infeasible to learn the standard QQ-function Qj(s,\vectorsyma)Q^{j}(s,\vectorsym{a}). To address this issue, we factorize the QQ-function using only the pairwise local interactions:

where N(j)\mathcal{N}(j) is the index set of the neighboring agents of agent jj with the size Nj=∣N(j)∣N^{j}=|\mathcal{N}(j)| determined by the settings of different applications. It is worth noting that the pairwise approximation of the agent and its neighbors, while significantly reducing the complexity of the interactions among agents, still preserves global interactions between any pair of agents implicitly (Blume, 1993). Similar approaches can be found in factorization machine (Rendle, 2012) and learning to rank (Cao et al., 2007).

The pairwise interaction Qj(s,aj,ak)Q^{j}(s,a^{j},a^{k}) as in Eq. (5) can be approximated using the mean field theory (Stanley, 1971). Here we consider discrete action spaces, where the action aja^{j} of agent jj is a discrete categorical variable represented as the one-hot encoding with each component indicating one of the DD possible actions: aj≜[a1j,…,aDj]a^{j}\triangleq[a^{j}_{1},\dots,a^{j}_{D}]. We calculate the mean action aˉj\bar{a}^{j} based on the neighborhood N(j)\mathcal{N}(j) of agent jj, and express the one-hot action aka^{k} of each neighbor kk in terms of the sum of aˉj\bar{a}^{j} and a small fluctuation δaj,k\delta{a^{j,k}} as

where aˉj≜[aˉ1j,…,aˉDj]\bar{a}^{j}\triangleq[\bar{a}^{j}_{1},\dots,\bar{a}^{j}_{D}] can be interpreted as the empirical distribution of the actions taken by agent jj’s neighbors. By Taylor’s theorem, the pairwise QQ-function Qj(s,aj,ak)Q^{j}(s,a^{j},a^{k}), if twice-differentiable w.r.t. the action aka^{k} taken by neighbor kk, can be expended and expressed as

As illustrated in Fig. 1, with the mean field approximation, the pairwise interactions Qj(s,aj,ak)Q^{j}(s,a^{j},a^{k}) between agent jj and each neighboring agent kk are simplified as that between jj, the central agent, and the virtual mean agent, that is abstracted by the mean effect of all neighbors within jj’s neighborhood. The interaction is thus simplified and expressed by the mean field QQ-function Qj(s,aj,aˉj)Q^{j}(s,a^{j},\bar{a}^{j}) in Eq. (8). During the learning phase, given an experience e=\big{(}s,\{a^{k}\},\{r^{j}\},s^{\prime}\big{)}, the mean field QQ-function is updated in a recurrent manner as

where αt\alpha_{t} denotes the learning rate, and aˉj\bar{a}^{j} is the mean action of all neighbors of agent jj as defined in Eq. (6). The mean field value function vtj(s′)v^{j}_{t}(s^{\prime}) for agent jj at time tt in Eq. (9) is

As shown in Eqs. (9) and (10), with the mean field approximation, the MARL problem is converted into that of solving for the central agent jj’s best response πtj\pi^{j}_{t} w.r.t. the mean action aˉj\bar{a}^{j} of all jj’s neighbors, which represents the action distribution of all neighboring agents of the central agent jj.

We introduce an iterative procedure in computing the best response πtj\pi^{j}_{t} of each agent jj. In the stage game {\vectorsymQt}\{\vectorsym{Q}_{t}\}, the mean action aˉj\bar{a}^{j} of all jj’s neighbors is first calculated by averaging the actions aka^{k} taken by jj’s NjN^{j} neighbors from the policies πtk\pi^{k}_{t} parametrized by their previous mean actions aˉ−k\bar{a}^{k}_{-}

With each aˉj\bar{a}^{j} calculated as in Eq. (11), the policy πtj\pi^{j}_{t} changes consequently due to the dependence on the current aˉj\bar{a}^{j}. The new Boltzmann policy is then determined for each jj that

By iterating Eqs. (11) and (12), the mean actions aˉj\bar{a}^{j} and the corresponding policies πtj\pi^{j}_{t} for all agents improves alternatively. In spite of lacking an intuitive impression of being stationary, in the following subsections, we will show that the mean action aˉj\bar{a}^{j} will be equilibrated at an unique point after several iterations, and hence the policy πtj\pi^{j}_{t} converges.

To distinguish from the Nash value function \vectorsymv ⁣Nash(s)\vectorsym{v}^{\mathop{}\!\mathtt{Nash}}(s) in Eq. (4), we denote the mean field value function in Eq. (10) as \vectorsymv ⁣MF(s)≜[v1(s),…,vN(s)]\vectorsym{v}^{\mathop{}\!\mathtt{MF}}(s)\triangleq[v^{1}(s),\dots,v^{N}(s)]. With \vectorsymv ⁣MF\vectorsym{v}^{\mathop{}\!\mathtt{MF}} assembled, we now define the mean field operator H ⁣MF\mathcal{H}^{\mathop{}\!\mathtt{MF}} in the form of

In fact, we can prove that H ⁣MF\mathcal{H}^{\mathop{}\!\mathtt{MF}} forms a contraction mapping; that is, one updates \vectorsymQ\vectorsym{Q} by iteratively applying the mean field operator H ⁣MF\mathcal{H}^{\mathop{}\!\mathtt{MF}}, the mean field QQ-function will eventually converge to the Nash QQ-value under certain assumptions.

2 Implementation

We can implement the mean field QQ-function in Eq. (8) by universal function approximators such as neural networks, where the QQ-function is parameterized with the weights ϕ\phi. The update rule in Eq. (9) can be reformulated as weights adjustment. For off-policy learning, we exploit either standard QQ-learning (Watkins & Dayan, 1992) for discrete action spaces or DPG (Silver et al., 2014) for continuous action spaces. Here we focus on the former, which we call MF-QQ.

In MF-QQ, agent jj is trained by minimizing the loss function

where yj=rj+γ vϕ−j ⁣MF(s′)y^{j}=r^{j}+\gamma\,v^{\mathop{}\!\mathtt{MF}}_{\phi^{j}_{-}}(s^{\prime}) is the target mean field value calculated with the weights ϕ−j\phi^{j}_{-}. Differentiating L(ϕj)\mathcal{L}(\phi^{j}) gives

which enables the gradient-based optimizers for training.

Instead of setting up Boltzmann policy using the QQ-function as in MF-QQ, we can explicitly model the policy by neural networks with the weights θ\theta, which leads to the on-policy actor-critic method (Konda & Tsitsiklis, 2000) that we call MF-AC. The policy network πθj\pi_{\theta^{j}}, i.e. the actor, of MF-AC is trained by the sampled policy gradient:

The critic of MF-AC follows the same setting for MF-QQ with Eq. (14). During the training of MF-AC, one needs to alternatively update ϕ\phi and θ\theta until convergence. We illustrate the MF-QQ iterations in Fig. 2, and present the pesudocode for both MF-QQ and MF-AC in Appendix A.

3 Proof of Convergence

We now prove the convergence of \vectorsymQt≜[Qt1,…,QtN]\vectorsym{Q}_{t}\triangleq[Q^{1}_{t},\dots,Q^{N}_{t}] to the Nash QQ-value \vectorsymQ∗=[Q∗1,…,Q∗N]\vectorsym{Q}_{*}=[Q^{1}_{*},\dots,Q^{N}_{*}] as the iterations of MF-QQ is applied. The proof is presented by showing that the mean field operator H ⁣MF\mathcal{H}^{\mathop{}\!\mathtt{MF}} in Eq. (13) forms a contraction mapping with the fixed point at \vectorsymQ∗\vectorsym{Q}_{*} under the main assumptions. We start from introducing the assumptions:

Each action-value pair is visited infinitely often, and the reward is bounded by some constant KK.

Agent’s policy is Greedy in the Limit with Infinite Exploration (GLIE). In the case with the Boltzmann policy, the policy becomes greedy w.r.t. the QQ-function in the limit as the temperature decays asymptotically to zero.

For each stage game [Qt1(s),...,QtN(s)][Q_{t}^{1}(s),...,Q_{t}^{N}(s)] at time tt and in state ss in training, for all tt, ss, j∈{1,…,N}j\in\{1,\dots,N\}, the Nash equilibrium \vectorsymπ∗=[π∗1,…,π∗N]\vectorsym{\pi}_{*}=[\pi^{1}_{*},\dots,\pi^{N}_{*}] is recognized either as 1) the global optimum or 2) a saddle point expressed as:

Note that Assumption 3 imposes a strong constraint on every single stage game encountered in training. In practice, however, we find this constraint appears not to be a necessary condition for the learning algorithm to converge. This is in line with the empirical findings in Hu & Wellman (2003).

Our proof is also built upon the two lemmas as follows:

Under Assumption 3, the Nash operator H ⁣Nash\mathcal{H}^{\mathop{}\!\mathtt{Nash}} in Eq. (4) forms a contraction mapping on the complete metric space from Q\mathcal{Q} to Q\mathcal{Q} with the fixed point being the Nash QQ-value of the entire game, i.e. Ht ⁣Nash\vectorsymQ∗=\vectorsymQ∗\mathcal{H}^{\mathop{}\!\mathtt{Nash}}_{t}\vectorsym{Q}_{*}=\vectorsym{Q}_{*}.

converges to zero with probability 11 (w.p.11) when

0≤αt(x)≤10\leq\alpha_{t}(x)\leq 1, ∑tαt(x)=∞\sum_{t}\alpha_{t}(x)=\infty, ∑tαt2(x)<∞\sum_{t}\alpha_{t}^{2}(x)<\infty;

x∈Xx\in\mathcal{X}, the set of possible states, and ∣X∣<∞|\mathcal{X}|<\infty;

var[Ft(x)∣Ft]≤K(1+∥Δt∥W2)\mathbf{var}[F_{t}(x)|\mathcal{F}_{t}]\leq K(1+\|\Delta_{t}\|_{W}^{2}) with constant K>0K>0.

Here Ft\mathcal{F}_{t} denotes the filtration of an increasing sequence of σ\sigma-fields including the history of processes; αt,Δt,Ft∈Ft\alpha_{t},\Delta_{t},F_{t}\in\mathcal{F}_{t} and ∥⋅∥W\|\cdot\|_{W} is a weighted maximum norm (Bertsekas, 2012).

See Theorem 1 in Jaakkola et al. (1994) and Corollary 5 Szepesvári & Littman (1999) for detailed derivation. We include it here to stay self-contained. ∎

By subtracting \vectorsymQ∗(s,\vectorsyma)\vectorsym{Q}_{*}(s,\vectorsym{a}) on both sides of Eq. (9), we present the relation from the comparison with Eq. (15) such that

where x≜(st,\vectorsymat)x\triangleq(s_{t},\vectorsym{a}_{t}) denotes the visited state-action pair at time tt. In Eq. (15), α(t)\alpha(t) is interpreted as the learning rate with αt(s′,\vectorsyma′)=0\alpha_{t}(s^{\prime},\vectorsym{a}^{\prime})=0 for any (s′,\vectorsyma′)≠(st,\vectorsymat)(s^{\prime},\vectorsym{a}^{\prime})\neq(s_{t},\vectorsym{a}_{t}); this is because that each agent only updates the QQ-function with the state sts_{t} and actions \vectorsymat\vectorsym{a}_{t} visited at time tt. Lemma 2 suggests Δt(x)\Delta_{t}(x)’s convergence to zero, which means, if it holds, the sequence of QQ’s will asymptotically tend to the Nash QQ-value \vectorsymQ∗\vectorsym{Q}_{*}.

One last piece to establish the main theorem is the below:

See details in Appendix D due to the space limit. ∎

In a finite-state stochastic game, the \vectorsymQt\vectorsym{Q}_{t} values computed by the update rule of MF-QQ in Eq. (9) converges to the Nash QQ-value \vectorsymQ∗=[Q∗1,…,Q∗N]\vectorsym{Q}_{*}=[Q_{*}^{1},\dots,Q_{*}^{N}], if Assumptions 1, 2 & 3, and Lemma 2’s first and second conditions are met.

Let Ft\mathcal{F}_{t} denote the σ\sigma-field generated by all random variables in the history of the stochastic game up to time tt: (st,αt,\vectorsymat,rt−1,...,s1,α1,\vectorsyma1,\vectorsymQ0)(s_{t},\alpha_{t},\vectorsym{a}_{t},r_{t-1},...,s_{1},\alpha_{1},\vectorsym{a}_{1},\vectorsym{Q}_{0}). Note that \vectorsymQt\vectorsym{Q}_{t} is a random variable derived from the historical trajectory up to time tt. Given the fact that all \vectorsymQτ\vectorsym{Q}_{\tau} with τ<t\tau<t are Ft\mathcal{F}_{t}-measurable, both \vectorsymΔt\vectorsym{\Delta}_{t} and \vectorsymFt−1\vectorsym{F}_{t-1} are therefore also Ft\mathcal{F}_{t}-measurable, which satisfies the measurability condition of Lemma 2.

To apply Lemma 2, we need to show that the mean field operator H ⁣MF\mathcal{H}^{\mathop{}\!\mathtt{MF}} meets Lemma 2’s third and fourth conditions. For Lemma 2’s third condition, we begin with Eq. (16) that

Note the fact that \vectorsymFt ⁣Nash\vectorsym{F}_{t}^{\mathop{}\!\mathtt{Nash}} in Eq. (17) is essentially the \vectorsymFt\vectorsym{F}_{t} in Lemma 2 in proving the convergence of the Nash QQ-learning algorithm. From Lemma 1, it is straightforward to show that \vectorsymFt ⁣Nash\vectorsym{F}_{t}^{\mathop{}\!\mathtt{Nash}} forms a contraction mapping with the norm ∥⋅∥∞\|\cdot\|_{\infty} being the maximum norm on \vectorsyma\vectorsym{a}. We thus have for all tt that

In meeting the third condition, we obtain from Eq. (17) that

We are left to prove that ct=∥\vectorsymCt(st,\vectorsymat)∣Ft∥c_{t}=\|\vectorsym{C}_{t}(s_{t},\vectorsym{a}_{t})|\mathcal{F}_{t}\| converges to zero w.p.11. With Assumption 3, for each stage game, all the globally optimal equilibrium(s) share the same Nash value, so does the saddle-point equilibrium(s). Each of the two following results is essentially associated with one of the two mutually exclusive scenarios in Assumption 3:

For globally optimal equilibriums, all players obtain the joint maximum values that are unique and identical for all equilibriums according to the definition;

Suppose that the stage game {\vectorsymQt}\{\vectorsym{Q}_{t}\} has two saddle-point equilibriums, \vectorsymπ\vectorsym{\pi} and \vectorsymρ\vectorsym{\rho}. It holds for agent jj that

By combing the above inequalities, we obtain

Given Proposition 1 that the policy based on the mean field QQ-function forms a contraction mapping, and that all optimal/saddle points share the same Nash value in each stage game, with the homogeneity of agents, \vectorsymv ⁣MF\vectorsym{v}^{\mathop{}\!\mathtt{MF}} will asymptotically converges to \vectorsymv ⁣Nash\vectorsym{v}^{\mathop{}\!\mathtt{Nash}}, the third condition is thus satisfied.

For the fourth condition, we exploit the conclusion that is proved above that H ⁣MF\mathcal{H}^{\mathop{}\!\mathtt{MF}} forms a contraction mapping, i.e. H ⁣MF\vectorsymQ∗=\vectorsymQ∗\mathcal{H}^{\mathop{}\!\mathtt{MF}}\vectorsym{Q}_{*}=\vectorsym{Q}_{*}, and it follows that

In the last step of Eq. (19), we employ Assumption 1 that the reward \vectorsymrt\vectorsym{r}_{t} is always bounded by some constant. Finally, with all conditions met, it follows Lemma 2 that \vectorsymΔt\vectorsym{\Delta}_{t} converges to zero w.p.11, i.e. \vectorsymQt\vectorsym{Q}_{t} converges to \vectorsymQ∗\vectorsym{Q}_{*} w.p.11. ∎

Apart from being convergent to the Nash QQ-value, MF-QQ is also Rational (Bowling & Veloso, 2001, 2002). We leave the corresponding discussion in Appendix D for details.

Related Work

We continue our discussion on related work from Introduction and make comparisons with existing techniques in a greater scope. Our work follows the same direction as Littman (1994); Hu & Wellman (2003); Bowling & Veloso (2002) on adapting a Stochastic Game (van der Wal et al., 1981) into the MARL formulation. Specifically, Littman (1994) addressed two-player zero-sum stochastic games by introducing a “minimax” operator in QQ-learning, whereas Hu & Wellman (2003) extended it to the general-sum case by learning a Nash equilibrium in each stage game and considering a mixed strategy. Nash-Q learning is guaranteed to converge to Nash strategies under the (strong) assumption that there exists an equilibrium for every stage game. In the situation where agents can be identified as either "friends" or "foes" (Littman, 2001), one can simply solve it by alternating between fully cooperative and zero-sum learning. Considering the convergence speed, Littman & Stone (2005) and de Cote & Littman (2008) draw on the folk theorem and acquired a polynomial-time Nash equilibrium algorithm for repeated stochastic games, while Bowling & Veloso (2002) tried varying the learning rate to improve the convergence.

The recent treatment of MARL was using deep neural networks as the function approximator. In addressing the non-stationary issue in MARL, various solutions have been proposed including neural-based opponent modeling (He & Boyd-Graber, 2016), policy parameters sharing (Gupta et al., 2017), etc. Researchers have also adopted the paradigm of centralized training with decentralized execution for multi-agent policy-gradient learning: BICNET (Peng et al., 2017), COMA (Foerster et al., 2018) and MADDPG (Lowe et al., 2017a), which allows the centralized critic QQ-function to be trained with the actions of other agents, while the actor needs only local observation to optimize agent’s policy.

The above MARL approaches limit their studies mostly to tens of agents. As the number of agents grows larger, not only the input space of QQ grows exponentially, but most critically, the accumulated noises by the exploratory actions of other agents make the QQ-function learning no longer feasible. Our work addresses the issue by employing the mean field approximation (Stanley, 1971) over the joint action space. The parameters of the QQ-function is independent of the number of agents as it transforms multiple agents interactions into two entities interactions (single agent v.s. the distribution of the neighboring agents). This would effectively alleviate the problem of the exploratory noise (Colby et al., 2015) caused by many other agents, and allow each agent to determine which actions are beneficial to itself.

Our work is also closely related to the recent development of mean field games (MFG) (Lasry & Lions, 2007; Huang et al., 2006; Weintraub et al., 2006). MFG studies population behaviors resulting from the aggregations of decisions taken from individuals. Mathematically, the dynamics are governed by a set of two stochastic differential equations that model the backward dynamics of individual’s value function, and the forward dynamics of the aggregate distribution of agent population. Despite that the backward equation equivalently describes what Bellman equation indicates in the MDP, the primarily goal for MFG is rather for a model-based planning and to infer the movements of the individual density through time. The mean field approximation (Stanley, 1971) in also employed in physics, but our work is different in that we focus on a model-free solution of learning optimal actions when the dynamics of the system and the reward function are unknown. Very recently, Yang et al. (2017) built a connection between MFG and reinforcement learning. Their focus is, however, on the inverse RL in order to learn both the reward function and the forward dynamics of the MFG from the policy data, whereas our goal is to form a computable QQ-learning algorithm under the framework of temporal difference learning.

Experiments

We analyze and evaluate our algorithms in three different scenarios, including two stage games: the Gaussian Squeeze and the Ising Model, and the mixed cooperative-competitive battle game.

Environment. In the Gaussian Squeeze (GS) task (HolmesParker et al., 2014), NN homogeneous agents determine their individual action aja^{j} to jointly optimize the most appropriate summation x=∑j=1Najx=\sum_{j=1}^{N}{a^{j}}. Each agent has 10 action choices – integers to 99. The system objective is defined as G(x)=xe−(x−μ)2σ2G(x)=xe^{{-(x-\mu)^{2}\over\sigma^{2}}}, where μ\mu and σ\sigma are the pre-defined mean and variance of the system. In the scenario of traffic congestion, each agent is one traffic controller trying to send aja^{j} vehicles into the main road. Controllers are expected to coordinate with each other to make the full use of the main route while avoiding congestions. The goal of each agent is to learn to allocate system resources efficiently, avoiding either over-use or under-use. The GS problem here sits ideally as an ablation study on the impact of multi-agent exploratory noises toward the learning (Colby et al., 2015).

Model Settings. We implement MF-QQ and MF-AC following the framework of centralized training (shared critic) with decentralized execution (independent actor). We compare against 4 baseline models: (1) Independent Learner (IL), a traditional QQ-Learning algorithm that does not consider the actions performed by other agents; (2) Frequency Maximum QQ-value (FMQ) (Kapetanakis & Kudenko, 2002), a modified IL which increases the QQ-values of actions that frequently produced good rewards in the past; (3) Recursive Frequency Maximum QQ-value (Rec-FMQ) (Matignon et al., 2012), an improved version of FMQ that recursively computes the occurrence frequency to evaluate and then choose actions; (4) Multi-agent Actor-Critic (MAAC), a variant of MADDPG architecture for the discrete action space (see Eq. (4) in Lowe et al. (2017b)). All models use the multilayer perception as the function approximator. The detailed settings of the implementation are in the Appendix C.1.

Results. Figure. 3 illustrates the results for the GS environment of μ=400\mu=400 and σ=200\sigma=200 with three different numbers of agents (N=100,500,1000N=100,500,1000) that stand for 33 levels of congestions. In the smallest GS setting of Fig. 3(a), all models show excellent performance. As the agent number increases, Figs. 3(b) and 3(c) show MF-QQ and MF-AC’s capabilities of learning the optimal allocation effectively after a few iterations, whereas all four baselines fail to learn at all. We believe this advantage is due to the awareness of other agents’ actions under the mean field framework; such mechanism keeps the interactions among agents manageable while reducing the noisy effect of the exploratory behaviors from the other agents. Between MF-QQ and MF-AC, MF-QQ converges faster. Both FMQ and Rec-FMQ fail to reach pleasant performance, it might be because agents are essentially unable to distinguish the rewards received for the same actions, and are thus unable to update their own QQ-values w.r.t. the actual contributions. It is worth noting that MAAC is surprisingly inefficient in learning when the number of agents becomes large; it simply fails to handle the non-aggregated noises due to the agents’ explorations.

2 Model-free MARL for Ising Model

Environment. In statistical mechanics, the Ising model is a mathematical framework to describe ferromagnetism (Ising, 1925). It also has wide applications in sociophysics (Galam & Walliser, 2010). With the energy function explicitly defined, mean field approximation (Stanley, 1971) is a typical way to solve the Ising model for every spin jj, i.e. ⟨aj⟩=∑aajP(a)\langle a^{j}\rangle=\sum_{a}{a^{j}}{P(a)}. See the Appendix C.2 for more details.

In addition to the reward, the order parameter (OP) (Stanley, 1971) is a traditional measure of purity for the Ising model. OP is defined as ξ=∣N↑−N↓∣N\xi={|N_{\uparrow}-N_{\downarrow}|\over N}, where N↑N_{\uparrow} represents the number of up spins, and N↓N_{\downarrow} for the down spins. The closer the OP is to 11, the more orderly the system is.

Model Settings. To validate the correctness of the MF-QQ learning, we implement MCMC methods (Binder et al., 1993) to simulate the same Ising model and provide the ground truth for comparison. The full settings of MCMC and MF-QQ for Ising model are provided in the Appendix C.2. One of the learning goals is to obtain the accurate approximation of ⟨aj⟩\langle a^{j}\rangle. Notice that agents here do not know exactly the energy function, but rather use the temporal difference learning to approximate ⟨aj⟩\langle a^{j}\rangle during the learning procedure. Once this is accurately approximated, the Ising model as a whole should be able to converge to the same simulation result suggested by MCMC.

Correctness of MF-QQ. Figure. 5 illustrates the relationship between the order parameter at equilibrium under different system temperatures. MF-QQ converges nearly to the exact same plot as MCMC, this justifies the correctness of our algorithms. Critically, MF-QQ finds a similar Curie temperature (the phase change point) as MCMC that is τ=1.2\tau=1.2. As far as we know, this is the first work that manages to solve the Ising model via model-free reinforcement learning methods. Figure. 5 illustrates the mean squared error between the learned QQ-value and the reward target. MF-QQ is shown in Fig. 4(a) to be able to learn the target well under low temperature settings. When it comes to the Curie temperature, the environment enters into the phase change when the stochasticity dominates, resulting in a lower OP and higher MSE observed in Fig. 4(b). We visualize the equilibrium in Fig. 6. The equilibrium points from MF-QQ in fact match MCMC’s results under three types of temperatures. The spins tend to stay aligned under a low temperature (τ=0.9\tau=0.9). As the temperature rises (τ=1.2\tau=1.2), some spins become volatile and patches start to form as spontaneous magnetization. This phenomenon is mostly observed around the Curie temperature. After passing the Curie temperature, the system becomes unstable and disordered due to the large thermal fluctuations, resulting in random spinning patterns.

3 Mixed Cooperative-Competitive Battle Game

Environment. The Battle game in the Open-source MAgent system (Zheng et al., 2018) is a Mixed Cooperative-Competitive scenario with two armies fighting against each other in a grid world, each empowered by a different RL algorithm. In the setting of Fig. 7(a), each army consists of 6464 homogeneous agents. The goal of each army is to get more rewards by collaborating with teammates to destroy all the opponents. Agent can takes actions to either move to or attack nearby grids. Ideally, the agents army should learn skills such as chasing to hunt after training. We adopt the default reward setting: −0.005-0.005 for every move, 0.20.2 for attacking an enemy, 55 for killing an enemy, −0.1-0.1 for attacking an empty grid, and −0.1-0.1 for being attacked or killed.

Model Settings. Our MF-QQ and MF-AC are compared against the baselines that are proved successful on the MAgent platform. We focus on the battles between mean field methods (MF-QQ, MF-AC) and their non-mean field counterparts, independent QQ-learning (IL) and advantageous actor critic (AC). We exclude MADDPG/MAAC as baselines, as the framework of centralized critic cannot deal with the varying number of agents for the battle (simply because agents could die in the battle). Also, as we demonstrated in the previous experiment of Fig. 3, MAAC tends to scale poorly and fail when the agent number is in hundreds.

Results and Discussion. We train all four models by 2000 rounds self-plays, and then use them for comparative battles. During the training, agents can quickly pick up the skills of chasing and cooperation to kill in Fig. 7(a). The Fig. 8 shows the result of winning rate and the total reward over 2000 rounds cross-comparative experiments. It is evident that on all the metrics mean field methods, MF-QQ largely outperforms the corresponding baselines, i.e. IL and AC respectively, which shows the effectiveness of the mean field MARL algorithms. Interestingly, IL performs far better than AC and MF-AC (2nd block from the left in Fig. 8(a)), although it is worse than the mean field counterpart MF-QQ. This might imply the effectiveness of off-policy learning with shuffled buffer replay in many-agent RL towards a more stable learning process. Also, the QQ-learning family tends to introduce a positive bias (Hasselt, 2010) by using the maximum action value as an approximation for the maximum expected action value, and such overestimation can be beneficial for each single agent to find the best response to others even though the environment itself is still changing. On the other hand, On-policy methods need to comply with the GLIE assumption (Assumption 2 in Sec 3.3) so as to converge properly to the optimal value (Singh et al., 2000), which is in the end a greedy policy as off-policy methods. Figure. 7(b) further shows the self-play learning curves of MF-AC and MF-QQ. MF-QQ presents a faster convergence speed than MF-AC, which is consistent with the findings in the Gaussian Squeeze task (see Fig. 3(b) & 3(c)). Apart from 64, we further test the scenarios when the agent size is 8, 144, 256, the comparative results keep the same relativity as Fig. 8; we omit the presentations for clarity.

Conclusions

In this paper, we developed mean field reinforcement learning methods to model the dynamics of interactions in the multi-agent systems. MF-QQ iteratively learns each agent’s best response to the mean effect from its neighbors; this effectively transform the many-body problem into a two-body problem. Theoretical analysis on the convergence of the MF-QQ algorithm to Nash QQ-value was provided. Three types of tasks have justified the effectiveness of our approaches. In particular, we report the first result to solve the Ising model using model-free reinforcement learning methods.

Acknowledgement

We sincerely thank Ms. Yi Qu for her generous help on the graphic design.

References

Appendix A Detailed mean field reinforcement learning algorithms

We published the code at https://github.com/mlii/mfrl.

Appendix B Proof of the bound for the remainder term in Eq. 7

Recall Eq. (8) that we approximate the action aka^{k} taken by the neighboring agent kk with the mean action aˉ\bar{a} calculated from the neighborhood N(j)\mathcal{N}(j). The state ss and the action aja^{j} of the central agent jj can be considered as fixed parameters; the indices j,kj,k of agents are essentially irrelevant to the derivation. With those omitted for simplicity, We rewrite the expression of the pairwise QQ-function as Q(a)≜Qj(s,aj,ak)Q(a)\triangleq Q^{j}(s,a^{j},a^{k}).

Suppose that QQ is MM-smooth, where its gradient ∇Q\nabla Q is Lipschitz-continuous with constant MM such that for all a,aˉa,\bar{a}

With the Lagrange’s mean value theorem, we have

Define δa≜a−aˉ\delta{a}\triangleq a-\bar{a} and the normalized vector δa^≜\nicefraca−aˉ∥a−aˉ∥2\delta{\hat{a}}\triangleq\nicefrac{{a-\bar{a}}}{{\|a-\bar{a}\|_{2}}} with ∥δa^∥2=1\|\delta{\hat{a}}\|_{2}=1, it follows from the above inequality

By arbitrary choice of (the unnormalized vector) δa\delta{a} such that the magnitude ∥δa∥2→0\|\delta{a}\|_{2}\to 0, it follows from above that

By aligning (the normalized vector) δa^\delta{\hat{a}} in the direction of the eigenvectors of the Hessian matrix ∇2Q\nabla^{2}Q, we can obtain for any eigenvalue λ\lambda of ∇2Q\nabla^{2}Q that

which indicates that all eigenvalues of ∇2Q\nabla^{2}Q can be bounded in the symmetric interval [−M,M][-M,M].

As the Hessian matrix ∇2Q\nabla^{2}Q is real symmetric and hence diagonalizable, there exist an orthogonal matrix UU such that U⊤[∇2Q]U=Λ≜diag⁡[λ1,…,λD]U^{\top}[\nabla^{2}Q]U=\Lambda\triangleq\operatorname{diag}[\lambda_{1},\dots,\lambda_{D}]. It then follows that

Recall the definition δa=a−aˉ\delta{a}=a-\bar{a} in Eq. (6), where aa is the one-hot encoding for DD actions, and aˉ\bar{a} is a DD-dimensional multinomial distribution. It can be shown that

where ii represents the specific action aa has represented such that ai′=0a_{i^{\prime}}=0 for i′≠ii^{\prime}\neq i.

With all elements assembled, we have proved that each single remainder term Rs,ajj(ak)R^{j}_{s,a^{j}}(a^{k}) in Eq. (8) is bounded in [−2M,2M][-2M,2M].

Appendix C Experiment details

IL, FMQ, Rec-FMQ and MF-QQ all use a three-layer MLP to approximate QQ-value. All agents share the same QQ-network for each experiment. The shared QQ-network takes an agent embedding as input and computes QQ-value for each candidate action. For MF-QQ, we also feed in the action approximation aˉ\bar{a}. We use the Adam optimizer with a learning rate of 0.00001 and ϵ\epsilon-greedy exploration unless otherwise specified. For FMQ, we set the exponential decay rate s=0.006s=0.006, start temperature max_temp=1000 and FMQ heuristic c=5c=5. For Rec-FMQ, we set the frequency learning rate αf=0.01\alpha_{f}=0.01.

MAAC and MF-AC use the Adam optimizer with a learning rate of 0.001 and 0.0001 for Critics and Actors respectively, and τ=0.01\tau=0.01 for updating the target networks. We share the Critic among all agents in each experiment and feed in an agent embedding as extra input. Actors are kept separate. The discounted factor γ\gamma is set to be 0.95 and the mini-batch size is set to be 200. The size of the replay buffer is 10610^{6} and we update the network parameters after every 500 samples added to the replay buffer.

For all models, we use the performance of the joint-policy learned up to that point if learning and exploration were turned off (i.e., take the greedy action w.r.t. the learned policy) to compare our method with the above baseline models.

C.2 Ising Model

An Ising model is defined as a stateless system with NN homogeneous sites on a finite square lattice. Each site determines their individual spin aja^{j} to interact with each other and aims to minimize the system energy for a more stable environment. The system energy is defined as

where τ\tau is the system temperature. When the temperature rises beyond a certain point (the Curie temperature), the system can no longer keep a stable form and a phase transition happens. As the ground-truth is known, we would be able to evaluate the correctness of the QQ-function learning when there is a large body of agents interacted.

The mean field theory provides an approximate solution to ⟨aj⟩=∑aajP(a)\langle a^{j}\rangle=\sum_{a}{a^{j}}{P(a)} through a set of self-consistent mean field equations

where tt represents the number of iterations.

To learn an optimal joint policy \vectorsymπ∗\vectorsym{\pi}^{*} for Ising model, we use the stateless QQ-learning with mean field approximation (MF-QQ), defined as

where the mean aˉj\bar{a}^{j} is given as the mean ⟨aj⟩\langle a^{j}\rangle from the last time step, and the individual reward is

To balance the trade-off between exploration and exploitation under low temperature settings, we use a policy with Boltzmann exploration and a decayed exploring temperature. The temperature for Boltzmann exploration of MF-QQ is multiplied by a decay factor exponentially through out the training process.

Without lost of generality, we assume λ>0\lambda>0, thus neighboring sites with the same action result in lower energy (observe higher reward) and are more stable. Each site should also align with the sign of external field hjh^{j} to reduce the system energy. For simplification, we eliminate the effect of external fields and assume the model to be discrete, i.e., ∀j∈N,hj=0,aj∈{−1,1}\forall j\in N,h^{j}=0,a^{j}\in\{-1,1\}.

We simulate the Ising model using Metropolis Monte Carlo methods (MCMC). After initialization, we randomly change a site’s spin state and calculate the energy change, select a random number between 0 and 1, and accept the state change only if the number is less than e(Ej−E−j)τe^{{(E^{j}-E^{j}_{-})\over\tau}}. This is called the Metropolis technique, which saves computation time by selecting the more probable spin states.

C.3 Battle Game

IL and MF-QQ have almost the same hyper-parameters settings. The learning rate is α=10−4\alpha=10^{-4}, and with a dynamic exploration rate linearly decays from γ=1.0\gamma=1.0 to γ=0.05\gamma=0.05 during the 2000 rounds training. The discounted factor γ\gamma is set to be 0.95 and the mini-batch size is 128. The size of replay buffer is 5×1055\times 10^{5}.

AC and MF-AC also have almost the same hyper-parameters settings. The learning rate is α=10−4\alpha=10^{-4}, the temperature of soft-max layer in actoractor is τ=0.1\tau=0.1. And the coefficient of entropy in the total loss is 0.08, the coefficient of value in the total loss is 0.1.

Appendix D Further details towards the theoretical guarantee of MF-Q𝑄Q

Following the contraction mapping theorem (Kreyszig, 1978), in order to be a contraction, the operator has to satisfy:

where 0≤α<10\leq\alpha<1 and B(\vectorsyma)≜[B(a1),…,B(aN)]\mathcal{B}(\vectorsym{a})\triangleq[\mathcal{B}(a^{1}),\dots,\mathcal{B}(a^{N})].

Here we start from binomial case and then adapt to the multinomial case in general. We first rewrite B(aj)\mathcal{B}(a^{j}) as

where ΔQ(s,aj,aˉ)=Q(s,a¬j,aˉ)−Q(s,aj,aˉ)\Delta Q(s,a^{j},\bar{a})=Q(s,{a^{\neg}}^{j},\bar{a})-Q(s,a^{j},\bar{a}).

In the second equation, we apply the mean value theorem in calculus: ∃x0∈[x1,x2]\exists x_{0}\in[x_{1},x_{2}], s.t., f(x1)−f(x2)=f′(x0)(x1−x2)f(x_{1})-f(x_{2})=f^{\prime}(x_{0})(x_{1}-x_{2}). In the third equation we use the maximum value for e−βΔQ0/(1+e−βΔQ0)2=1/4e^{-\beta\Delta Q_{0}}/({1+e^{-\beta\Delta Q_{0}}})^{2}=1/4 when Q0=0Q_{0}=0. In the last equation we apply the Lipschitz constraint in the assumption where constant K≥0K\geq 0. Finally, we have:

In order for the contraction to hold, T>K2T>{K\over 2}. In other words, when the action space is binary for each agent, and the temperature is sufficiently large, the mean field procedure converges.

This proposition can be easily extended to multinomial case by replacing binary variable aja^{j} by a multi-dimensional binary indicator vector \vectorsymaj\vectorsym{a}^{j}, on each dimension, the rest of the derivations would remain essentially the same. ∎

In line with (Bowling & Veloso, 2001, 2002), we argue that to better evaluate a multi-agent learning algorithm, on top of the convergence guarantee, discussion on property of Rationality is also needed.

(also see (Bowling & Veloso, 2001, 2002)) In an NN-agent stochastic game defined in this paper, given all agents converge to stationary policies, if the learning algorithm converges to a policy that is a best response to the other agents’ policies, then the algorithm is Rationale.

Our mean field QQ-learning is rational in that Eq. (5) converts many agents interactions into two-body interactions between a single agent and the distribution of other agents actions. When all agents follow stationary policies, their policy distribution would be stationary too. As such the two-body stochastic game becomes an MDP, and the agent would choose a policy (based on Assumption 2) which is the best response to the distribution of other stationary policies. As agents are symmetric in our case, they all show the best response to the distributions, and are therefore rational.

Appendix E Proof of Mean Field Reinforcement Learning with Function Approximation

Previous convergence results in Theorem.1 has shown that the Mean Field Q-learning algorithm will converge when the Q function is in a tabular case. We now move onto the proof that the MF-Q algorithm will converge when the Q function is represented by some functional approximations.

An NN-agent (or, NN-player) stochastic game Γ\Gamma is formalized by the tuple \Gamma\triangleq\big{(}\mathcal{S},\mathcal{A}^{1},\dots,\mathcal{A}^{N},r^{1},\dots,r^{N},p,\gamma\big{)}. The state-space S\mathcal{S} is finite. Let (S,p\vectorsymπ)(\mathcal{S},p_{\vectorsym}{\pi}) be the Markov chain induced by the joint policy \vectorsymπ\vectorsym{\pi}, and we assume it to be uniformly ergodic.

In the functional approximation setting, we can apply the update rules:

In the above, Δt\Delta_{t} is the temporal difference at time tt.

And the goal is to derive the parameter vector \vectorsymϕ={ϕj}\vectorsym{\phi}=\{\phi^{j}\} such that \vectorsymω⊤\vectorsymϕ{\vectorsym{\omega}}^{\top}\vectorsym{\phi} approximates the (local) Nash Q-values. At each time step, the learning policy \vectorsymπϕt\vectorsym{\pi}_{\phi_{t}} is the Botlzmann policy with respect to \vectorsymω⊤\vectorsymϕ{\vectorsym{\omega}}^{\top}\vectorsym{\phi}. Give the Proposition 1, we know that the policy \vectorsymπϕt\vectorsym{\pi}_{\phi_{t}} is K2T{K\over 2T} Lipschitz continuous with respect to ϕt\phi_{t}.

Similar to the framework used in the convergence proof of Q-learning with function approximation (Melo et al., 2008), we establish convergence of Eq. (38) by adopting an ordinary differentiable equation (ODE) with a globally asymptotically stable equilibrium point where the trajectories closely follow.

Given the MDP Γ\Gamma, \vectorsymπϕt\vectorsym{\pi}_{\phi_{t}}, {\vectorsymωp,p=1,...,P}\{\vectorsym{\omega}_{p},p=1,...,P\}, and the learning policy \vectorsymπϕt\vectorsym{\pi}_{\phi_{t}} that is K2T{K\over 2T} Lipschitz continuous with respect to ϕt\phi_{t}, if the Assumptions 1, 2 & 3, and Lemma 2’s first and second conditions are met, then there exists C0C_{0} such that the algorithm in Eq. (38) converges w.p.1 if K2T<C0{K\over 2T}<C_{0}.

We first re-write the Eq. (38) as on ODE:

Notice that we use a vector for considering the updating rule for the Q function of each agent. We can easily know that necessity condition of the equilibrium is that it must follow \vectorsymϕ∗=\vectorsymAϕ∗−1\vectorsymbϕ∗\vectorsym{\phi}^{*}=\vectorsym{A}_{\phi^{*}}^{-1}\vectorsym{b}_{\phi^{*}}. The existence of the such equilibrium has been restricted in the scenario that meets Assumption 33. In the proof of Theorem 1 we have already pointed out that under the Assumption 33, the existing equilibrium, either in the form of global equilibrium or in the form of saddle-point equilibrium, is unique.

As we know that the policy \vectorsymπϕt\vectorsym{\pi}_{\phi_{t}} is Lipschitz w.r.t ϕt\phi_{t}, this implies that \vectorsymAϕ\vectorsym{A}_{\phi} and \vectorsymvϕ\vectorsym{v}_{\phi} are also Lipschitz continuous w.r.t to ϕ\phi. In other words, if K2T≤C0{K\over 2T}\leq C_{0} is sufficiently small and close to zero, then the norm term of (sup⁡\vectorsymϕ∣∣\vectorsymAϕ−\vectorsymAϕ∗∣∣2+sup⁡\vectorsymϕ∣∣\vectorsymbϕ−\vectorsymbϕ∗∣∣2∣∣\vectorsymϕ−\vectorsymϕ∗∣∣2)\left(\sup_{\vectorsym{\phi}}||\vectorsym{A}_{\phi}-\vectorsym{A}_{\phi^{*}}||_{2}+\sup_{\vectorsym{\phi}}\dfrac{||\vectorsym{b}_{\phi}-\vectorsym{b}_{\phi}^{*}||_{2}}{||\vectorsym{\phi}-\vectorsym{\phi}^{*}||_{2}}\right) goes to zero. Considering near the equilibrium point ϕ∗\phi^{*}, \vectorsymAϕ∗\vectorsym{A}_{\phi^{*}} is a negative definite matrix, the Eq. (42) tends to be negative definite as well, so the ODE in Eq.(41) is globally asymptotically stable and the conclusion of the theorem follows. ∎