Differential evolution particle swarm optimization algorithm based on good point set for computing Nash equilibrium of finite noncooperative game

1.   Introduction

Game theory is a mathematical theory that studies interactions between decision makers. It mainly studies how players use the information they have to make decisions. Game theory is widely used in many fields, such as economics, political science, biology and other fields ^{[1,2,3,4,5,6]}. Noncooperative game is the main component of game theory, and the Nash equilibrium is the core concept of noncooperative game. However, achieving equilibrium requires players to make predictions according to certain steps in the game. Therefore, it can be reduced to a calculation problem to a certain extent.

At present, many numerical methods have been proposed to solve the Nash equilibrium, such as the Lemke-Howson algorithm ^[7], global Newton algorithm ^[8], projection-like method ^[9], trust region algorithm ^[10] and other traditional calculation methods. In recent years, with the increased complexity of game problems, traditional numerical methods are almost impossible to use in solving due to the increasing difficulty of finding the solution and the computation time. Increasing numbers of scholars have focused on intelligent algorithms for biological simulation. Because of the rationality of imitating biological behaviors, such algorithms have the implicit rational characteristics of games. By such means, such algorithms can be considered as some sort of paths to achieve the Nash equilibrium.

Pavlidis ^[11] verified the effectiveness of three computational intelligence algorithms, namely, covariance matrix adaptation evolution strategies, particle swarm optimization (PSO) and differential evolution (DE), to compute the Nash equilibrium of finite strategic games. Boryczka ^[12] compared DE with two well-known algorithms: simplicial subdivision and the Lemke-Howson algorithm. It has also been proven that the DE method can obtain an approximate Nash equilibrium solution. Chen ^[13] utilized the genetic algorithm (GA) to acquire the Nash equilibrium of an $N$-person noncooperative game. With the examples of double matrix games, the implementation of the genetic algorithm was discussed and proven to be effective. Qiu ^[14] proposed an immune algorithm (IA) for solving game equilibrium. The advantages of this algorithm in solving game problems and its stable convergence were verified through examples. Franken ^[15] investigated the application of coevolutionary training techniques based on PSO to evolve the iterated prisoner's dilemma (IPD). Jia ^[16] proposed an immune particle swarm optimization (IPSO) to compute finite noncooperative games among $N$ people. The results show that this algorithm is superior to the immune algorithm and original swarm algorithm in solving game problems. Yang ^[17] proposed the fireworks algorithm (FWA) to compute the Nash equilibrium of $N$-person noncooperative game. Computer simulation results demonstrate that the proposed algorithm is effective and superior to IPSO.

The above studies all demonstrate the advantages of intelligent algorithms in solving Nash equilibrium problems. However, several issues remain, such as the higher complexity of the algorithm, its slow convergence speed and its low accuracy; in particular, the convergence of the algorithm is not proven theoretically. The goal of this paper is to propose an efficient hybrid of GPDEPSO to solve the game problem and then prove its convergence. The paper is organized as follows: In Section 2, we summarize the concept of a noncooperative $N$-person game. In Section 3, we propose a hybrid algorithm, GPDEPSO, to compute the Nash equilibrium of noncooperative $N$-person games. First, we initialize the population with a good point set to ensure the initial particle distribution is global, which will help the algorithm avoid local convergence. Then, the position updating formula of PSO is simplified in a new form that does not have a velocity term, and it is used as a variation operator to perform variation crossover selection with DE. The convergence of GPDEPSO is proved by stochastic functional analysis in Section 4. Section 5 is devoted to computational experiments, and by comparing the algorithm proposed in this paper with other algorithms, its superiority is proven.

2.   Game and Nash equilibrium Definitions ^[18]

Definition 1. We consider an $N$-person finite strategic game defined by $\Gamma = ((N, S_{i}, U_{i}, X_{i}, f_{i}), i = 1, \ldots, n)$,

where

(1). $N = \{1, \ldots, n\}$ is the set of players and $n$ is the number of players;

(2). $S_{i} = \{s_{i1}, \ldots, s_{im_{i}}\} \ \forall i\in N$ is the pure strategy set of player $i$, $m_{i}$ represents the number of strategies available to player $i$, $S = \prod_{i = 1}^nS_{i}$ is the Cartesian product of pure strategy sets of all players, and each pure strategy profile meets $(S_{1}, S_{2}, \ldots, S_{n})\in S$;

(3). $U_{i}:S\rightarrow \mathbb{R}\ \forall i\in N$ represents the payoff function;

(4). $X_{i} = \{x_{i} = (x_{i1}, \ldots, x_{im_{i}}):x_{ik}\geq0, k = 1, \ldots, m_{i}, \sum_{k_{i} = 1}^{m_{i}}x_{ik_{i}} = 1\}\ \forall i\in N$ is the set of mixed strategies, where $x_{ij}$ is the probability that player $i$ adopts $s_{ij}$ for $j = 1, \cdots, m_i$. $X = \prod_{i = 1}^nX_{i}$ is the Cartesian product of mixed strategy sets of all players, and each mixed strategy profile meets $(x_{1}, x_{2}, \ldots, x_{n})\in X$;

(5). $f_{i}:X\rightarrow \mathbb{R} \; \forall i\in N$ represents the expected payoff function.

$f_{i}(x_{1}, \ldots, x_{n}) = \sum_{k_{1} = 1}^{m_{1}}\ldots\sum_{k_{n} = 1}^{m_{n}}U_{i}(s_{1k_{1}}, \ldots, s_{nk_{n}})\prod_{i = 1}^{n}x_{ik_{i}}$ represents player $i$ getting the expected payoff value when he chooses a mixed strategy $x_{i} = (x_{i1}, \ldots, x_{im_{i}})\in X_{i}$. where $U_{i}(s_{ik_{i}}, \ldots, s_{nk_{n}})$ represents player $i$ getting the payoff value when each player $i$ chooses the pure strategy $s_{ik_{i}}\in S_{i}, i = 1, \ldots, n$.

Definition 2. If there is $x^{\ast} = (x_{1}^{\ast}, \ldots, x_{n}^{\ast})$, such that $f_{i}(x_{i}^{\ast}, x_{i^{\wedge}}^{\ast}) = \max\limits_{u_i\in X_i}f_{i}(u_{i}, x_{i^{\wedge}}^{\ast}), \; \forall i\in N$, then $x^{\ast}$ is the Nash equilibrium point of an $N$-person finite noncooperative game, where $i^{\wedge} = N \backslash \{i\}, \ \forall i\in N$.

Conclusion 1. Mixed strategy $x^{\ast}$ is the Nash equilibrium point if and only if every pure strategy $s_{ik_{i}}(1\leq k_{i}\leq m_{i})$ of each player $i$ has $f_{i}(x^{\ast})\geq f_{i}(x^{\ast} \shortparallel s_{ik_{i}})$, where $(x^{\ast} \shortparallel s_{ik_{i}})$ is only the player $i$ replacing their own strategy with $s_{ik_{i}}$, and the other players do not change their own strategy under the condition of equilibrium solution $x^{\ast}$.

In particular, the $\Gamma$ is a bimatrix game when $N = 2$, $(x^{\ast}, y^{\ast})$ is the Nash equilibrium solution if and only if $\left\{\begin{aligned} &x^{\ast} Ay^\ast\geq xAy^{\ast T}, \ \ \forall x, \\ &x^{\ast}By^{\ast}\geq x^{\ast}By^{ T}, \ \ \forall y, \end{aligned}\right.$ where $A$ and $B$ are payoff matrices for each player.

Theorem 1. A mixed strategy $x^{\ast}\in X$ is the Nash equilibrium point of a game $\Gamma$ if and only if $x^{\ast}$ is an optimal solution to the following optimization problem, and the optimal value is 0:

Proof Necessity: Suppose that $x^{\ast}$ is the Nash equilibrium point. According to Conclusion 1, we have $f_{i}(x^{\ast}) \geq f_{i}(x^{\ast} \shortparallel s_{ik_{i}}), \forall i\in N, \forall s_{ik_{i}}\in S_i$, then

Because of $x^{\ast}\in X = \prod_{i = 1}^nX_{i}$, the constraint condition of (2.1) hold. In formula (2.1), $f(x)$ is nonnegative, so $f(x)$ gets the minimum value 0 at $x^{\ast}$.

Sufficiency: Assume that $x^{\ast}$ is the solution of problem (2.1). According to $x^{\ast}$ satisfies conditions of (2.1), we can know that $x^{\ast}\in X = \prod_{i = 1}^nX_{i}$. And because of $f(x^{\ast}) = 0$, so $f_{i}(x^{\ast} \shortparallel s_{ik_{i}}) -f_i(x^{\ast})\leq 0$, that is $f_i(x^{\ast})\geq f_{i}(x^{\ast} \shortparallel s_{ik_{i}}), \forall i\in N, \forall s_{ik_{i}}\in S_i$. So $x^{\ast}$ is the Nash equilibrium point of the game $\Gamma$.

For the two-person double matrix game, the above optimization problem can be simplified as:

where $A_{i}$ is the $i$th row of matrix $A$ and $B_{j}$ is the $j$th column of matrix $B$. Then, $(x^{\ast}, y^{\ast})$ is a Nash equilibrium solution of a two-person noncooperative game, which is also $f(x^{\ast}, y^{\ast}) = 0$.

3.   Differential evolution particle swarm optimization algorithm based on good point set (GPDEPSO)

3.1.   Good point set

Good point set was originally proposed by Hua Luogeng et al. ^[19] and defined as follows:

(1). Let $G_{s}$ be a unit cube in $s$-dimensional Euclidean space. Let $x\in G_{s}$

(2). Let $G_{s}$ have a set of points $P_n(i)$ with $n$ points:

where $\{\cdot\}$ represents the decimal part of the value.

(3). For any given point $r = (r_1, r_2, \cdots, r_s)\in G_s$, let $N_n(r) = N_n(r_1, r_2, \cdots, r_s)$ represents the number of points in $P_n(i)$ that meet the following inequality:

$\varphi(n) = sup|(N_n(r)/n)-|r||$, where $|r| = r_1r_2\cdot \cdots \cdot r_s$, is called a deviation of the point set $P_n(i)$. If $\forall n, \varphi(n) = O(1)$, then $P_n(i)$ is said to be uniformly distributed on $G_s$ and the deviation is $\varphi(n)$.

(4). Let $r\in G_s$, and $P_n(i) = \{(r_1\ast i, r_2\ast i, \cdots, r_s\ast i, ), i = 1, \cdots, n\}$, the deviation $\varphi(n)$ meets $\varphi(n) = C(r, \varepsilon)\cdot n^{-1+\varepsilon}$, where $C(r, \varepsilon)$ is a constant related only to $r$, $\varepsilon$($\varepsilon$ is an arbitrarily small positive number), then, the $P_n(i)$ is called the good point set, $r$ is called the good point.

In this paper, we take

In the following, we generate two distribution maps (Figures 1 and 2) of 500 populations with the random point method and the exponential good point set method, respectively.

As shown, the good point sequence is more uniform and global than the random point method. In addition, the good point method is independent with spatial dimensions, thus it can be well adapted to the high-dimensional problems. It is also stable that the distribution results are the same every time when the number of points is the same.

3.2.   Differential evolution algorithm (DE)

DE is a new simple and robust evolutionary algorithm that was first introduced by Storm and Price ^[20]. There are four operations of DE: initialization, variation, crossover, and selection operation.

(1). Initialization

Let each individual in a population be

where $N$ and $D$ represent, respectively, the population size and space dimension. In the study of DE, it is generally assumed that the initial population conforms to uniform probability distribution, and its form is as follows:

where $rand[0,1]$ represents random values in the range [0,1], $x_{j}^{U}$ and $x_{j}^{L}$ represent, respectively, the upper and lower bounds of parameter variables.

(2). Variation

The variation operation is mainly executed to distinguish DE from other evolutionary algorithms. The variation of individual $V = (v_{i1}, \ldots, v_{iD})$ is generated by the following equation:

where $r1$, $r2$, and $r3$ are different integers between 1 and $N$, and they are also different from $i$; $F$ is a constriction factor to control the size of difference of two individuals, and $t$ is the current iterate point.

(3). Crossover

We use the crossover between the parent and offspring with the given probability for generating new individual $U = (u_{i1}, \ldots, u_{iD})$:

where $rand(j)$ is random value in the range [0,1], $CR$ is crossover operator in the range [0,1], and $rnbr(i)\in\{1, \ldots, D\}$ is a randomly selected sequence, which ensures that a new individual gets at least one component value from the variation vector.

(4). Selection

The offspring $X_{i}^{t+1}$ is generated by selecting the individual and parent according to the following equation:

where $f(\cdot)$ is the fitness function value.

3.3.   Particle swarm optimization algorithm (PSO)

PSO was proposed by Eberhart and Kennedy ^[21], and is a random search algorithm, which is inspired by the activities of flocking birds. PSO uses a population of individuals called particles, and there are two main operations in PSO, speed updating and position updating:

The $\omega$ is the inertia weight, $c_{1}$ and $c_{2}$ are acceleration constants, $r_{1}$ and $r_{2}$ are random values in the interval (0, 1), $v_{ij}^{t}$ is the $i$th particle's velocity in generation $t$, $x_{ij}^{t}$ is the $i$th particle's position in generation $t$, $x_{pbest}^{t}$ is the personal best position of particle $i$ before generation $t$, and $x_{gbest}^{t}$ is the global best position in the searching history^[22].

Speed updating can be further explained as follows: $\omega v_{ij}^{t}$ is called the current state of the particle with the ability to balance global and local search. $c_{1}r_{1}(x_{pbest}^{t}-x_{ij}^{t})$ is the cognitive modal of the particle, which represents the ability of learning from itself and endows particles with strong local search capabilities. $c_{2}r_{2}(x_{gbest}^{t}-x_{ij}^{t})$ represents the social cognition modal of the particle and information sharing among particles, that is, the ability to learn from the entire population and endow particles with strong global search ability. Then, the position updating makes the particles reach the new position.

A simplified position transformation formula without velocity term is expressed as follows ^[23]:

3.4.   GPDEPSO experimental steps and its implementation

The steps of GPDEPSO is described as follows :

Step 1: Set the parameters of GPDEPSO, such as $N$, $D$, $CR$, $F_0$, $\omega_{min}$, $\omega_{max}$, $x^L$, $x^U$, $c_1$, $c_2$, and set the maximum number of iteratios $T$, accuracy $\varepsilon$, where

Step 2: Randomly generate $N$ initial populations $P(0)$ by using a good point set, and $\sum_{k_{i} = 1}^{m_{i}}x_{ik_{i}} = 1, x_{ik_{i}}\geq0, x_{ik_{i}}\in X_{i}, \ \ i = 1, \ldots, N;k_{i} = 1, \ldots, m_{i}$.

Step 3: Calculate the fitness function value $f(x)$ of each individual in population $P(t)$ and determine the $x_{pbest}^{t}$ and $x_{gbest}^{t}$.

Step 4: The next generation population $P_{1}(t)$ is generated by variation of formula (3.10), and population $P_{2}(t)$ is generated by variation of formula (3.6).

Step 5: The population $P^{'}{(t)}$ is generated by crossover of formula (3.7).

Step 6: According to formula (3.8), populations $P(t)$ and $P^{'}{(t)}$ are selected to generate offspring population $P(t+1)$ and the fitness function value of population $P(t+1)$ is calculated.

Step 7: Determine whether to end according to the accuracy and the maximum number of iterations, and output the optimal value; otherwise, turn to step 3.

The pseudo code of GPDEPSO is as follows:

4.   Convergence analysis of GPDEPSO

The main operations of GPDEPSO are a variation operation based on PSO, a crossover operation, and a selection operation of DE. Considering variation crossover operations as a difference operator based on PSO (DOPSO), which is essentially a variation operation, the individual $X_i$ generates an intermediate individual $V_i$ with at least $1-(1-CR-1/D)^D$ probability. The selection operation is regarded as a selection operator (SO), which is strictly based on the strategy of survival of the fittest. It produces a better new generation of individuals by eliminating the inferior individuals in $X_i$ and $V_i$.

For the minimal optimization problem, GPDEPSO evaluation function $\{f(X_i^t) :1 < t < T\}$ is a monotonous nonincremental sequence. The convergence of GPDEPSO is illustrated by He's stochastic functional analysis method in the literature ^[24]. An iteration process of GPDEPSO is abstracted as a composite random mapping of a DOPSO and a SO.

Definition 3. The process of generating test vector $U$ by DOPSO is to recombine and transform everyone dimensional component of the target vector according to probability $\theta = CR+1/D$. In addition, the process can be described as a random mapping of the solution space $\Psi_1:\Omega\times S\rightarrow S^2$, which is defined as:

where $(\Omega, A, \mu)$ is the complete probability measure space and $\Omega$ is the nonempty abstract set, and its element $\omega$ is a basic event. $A$ is the $\sigma-algebra$ composed of some subsets of $\Omega$, $\mu$ is the probability measure on $A$, and $S$ is the solution space: $S = \{X|X = (x_1, x_2, \cdots, x_D), x_j^L\leq x_j\leq x_j^U, j = 1, 2, \cdots, D\}$.

Definition 4. The selection operator SO is the process of selecting the optimal individual from the test vector $U_i$ and the target vector $X_i$ according to the greedy selection method. It is a mapping on the solution space $\Psi_2:S^2\rightarrow S$:

Combining the above two definitions, it can be seen that one iteration of GPDEPSO is equivalent to mapping $\Psi = (\Psi_2\circ \Psi_1):\Omega \times S\rightarrow S(\Omega\times S\rightarrow S^2\rightarrow S)$ to the current population $P(t)$, where $\Psi$ is a reverse-order synthesis mapping corresponding to DOPSO and SO mapping. Then, for the new population $P(t+1)$ it can be expressed as $P(t+1) = \Psi(\omega, P(t)) = \Psi_2(\Psi_1(\omega, P(t))), 0\leq t\leq T-1$.

Let $f(X_{best}^t)$ be the fitness function value of the best individual $X_{best}^t$ in $P(t)$. Under the effect of $\Psi$, the new generation population generated by GPDEPSO will be superior to the previous one. Therefore, the sequence $\{f(X_{best}^t)\}_{1\leq t\leq T}$ is necessarily a monotonous nonincremental sequence (assuming that $f(X)$ is the minimum optimization function in this paper). In addition, the evolutionary process of GPDEPSO can also be characterized by the optimal individual, and the mapping $\Psi$ can be redefined as a mapping corresponding to the process of generating the optimal individual, that is, $X_{best}^{t+1} = \Psi(\omega, X_{best}^t) = \Psi_2(\Psi_1(\omega, X_{best}^t))$.

Lemma 1. Let $\lambda:S\times S\rightarrow R$ be the distance defined on $S$ and satisfy $\lambda(X_i, X_j) = |f(X_i)-f(X_j)|, \forall X_i, X_j\in S$; then, $(S, \lambda)$ is a complete separable metric space.

Similar to the method presented in Ref. ^[25], Lemma 1 is established.

Theorem 2. The random mapping $\Psi = (\Psi_2\circ \Psi_1):\Omega \times S\rightarrow S$ formed by one iteration of GPDEPSO is a random contraction operator.

Proof. According to the definition of DOPSO and SO operation, it can be seen that the new population generated by GPDEPSO in each iteration is better than the previous one. Therefore, for random mapping $\Psi = (\Psi_2\circ \Psi_1):\Omega \times S\rightarrow S$, there exists a random variable with nonnegative real value $0\leq K(\omega) < 1$, which makes the following formula hold:

letting

then

Therefore, the mapping formed by GPDEPSO, $\Psi:\Omega \times S \rightarrow S$ is a random contraction operator.

Lemma 2. (Random Contraction Mapping Theorem) ^[26]

Letting $\Psi:\Omega \times S \rightarrow S$ be a random operator, satisfying that for almost all $\omega\in \Omega$ and $\Psi(\omega)$ is contractive operator, then for $\Psi(\omega)$ there exists a unique random fixed point $\xi(\omega)$, and $\Psi(\omega, \xi(\omega)) = \xi(\omega)$.

According to Theorem 2, the random mapping $\Psi$ is a random contraction operator. By using the conclusion of Lemma 2, for $\Psi(\omega)$ there must exist a unique random fixed point. Then, the convergence criterion of GPDEPSO can be derived, which means that GPDEPSO is asymptotically convergent.

5.   Experimental design and results

5.1.   Numerical examples

We take the classical games given in ^[16] and ^[17] as our numerical examples, and solve them separately by using the GPDEPSO given in this paper and comparing with other algorithms. (For Examples 1-3, the parameters of the algorithm are set as $N = 50$, $T = 100$, $\varepsilon = 10^{-4}$, etc. For examples 4 and 5 of high dimension, the parameters setting are, respectively, $N = 50$, $T = 100$, $\varepsilon = 10^{-5}$, etc., and $N = 100$, $T = 300$, $\varepsilon = 10^{-5}$, etc.).

Example 1 Prisoner's Dilemma Game $\Gamma(X_{1}, Y_{1}, A_{1}, B_{1})$:

Example 2 Guessing Game $\Gamma(X_{2}, Y_{2}, A_{2}, B_{2})$:

Example 3 Supervisory Game $\Gamma(X_{3}, Y_{3}, A_{3}, B_{3})$:

Example 4 Consider a Three-Dimensional Payoff Matrix Game $\Gamma(X_{4}, Y_{4}, A_{4}, B_{4})$:

Example 5 Consider a 10-Dimensional Payoff Matrix Game $\Gamma(X_{5}, Y_{5}, A_{5}, B_{5})$:

5.2.   Results and discussion

The following Tables 1-4 present computational results of the above examples by GPDEPSO.

From above four tables, we can see that the Nash equilibrium can be obtained under the given parameters by using GPDEPSO proposed in this paper. In addition, the accuracy of fitness function values is higher approximately $10^4$ times and $10^2$ times than that of Refs. ^[16] and ^[17]. It can be seen that the speed of this algorithm does not be affected under the condition of increasing the population size (N = 50). On the contrary, according to the results of Examples 1, 3, and 4, it is show that the number of iterations is reduced significantly, and the Example 4 is the most obvious. Compared with Refs.^[16], the number of iteration of Example 4 is reduced by approximately 16 times, which is approximately 8 times compared with Ref.^[17]. Through the above analysis, GPDEPSO is superior to the algorithms presented in the existing literature in terms of iteration times and accuracy of results. Furthermore, although the number of populations is lager than that in the previous literature, it does not affect the speed and accuracy of the algorithm, but plays a powerful role in finding the global optimal solution.

The following are two comparison figures for solving Example 5, a high-dimensional payoff matrix game. Figure 3 is a comparison between GPDEPSO and hybrid differential evolution particle swarm optimization algorithm (DEPSO). Figure 4 is a comparison of GPDEPSO, differential evolution algorithm based on good point set (GPDE), and particle swarm optimization algorithm based on good point set (GPPSO). The Nash equilibrium solution calculated by all methods is $x^{\ast} = y^{\ast} = (0.1, \ldots, 0.1)_{1\times10}$.

As shown in Figures 3 and 4, the high-dimensional payoff matrix game can be solved better by GPDEPSO. Compared with Ref. ^[17], it can be seen that not only the calculation speed is greatly improved, but the accuracy of the fitness function value is also better. It can be seen from Figure 3 that GPDEPSO converges faster and more smoothly than DEPSO, which indicates that the good point set can avoid falling into a local optimum during the calculation process. From Figure 4, comparing GPDEPSO, GPPSO, with GPDE, and find that GPDEPSO combines the ability of PSO to find the global optimum quickly with the fast convergence ability of DE. In Example 5, GPDEPSO can quickly find the global optimal solution within the first 50 iterations, and the approximate exact solution can be obtained after 150 iterations. Therefore, GPDEPSO has a good advantage in solving Nash equilibrium problems with a high-dimensional payoff matrix.

6.   Conclusions

We propose GPDEPSO from the point of view of DE, considering the different characteristics of DE and PSO, and we use a good point set to make the initial data more uniform. This algorithm combines the advantages of DE and PSO, which not only ensures the simple operation, easy implementation, and fast convergence of DE, but also enhances its global optimization ability. By solving the Nash equilibrium of noncooperative games, we find that the proposed algorithm is superior to the comparative algorithms in terms of the calculation accuracy and convergence. In particular, for a high-dimensional payoff matrix game, the efficiency of the algorithm is remarkable. In the future, since a Nash equilibrium problem is a complex, NP-hard problem, it would be interesting to consider the influence of different variation operations and selection strategies on the solution of the Nash equilibrium, as well as to consider the Nash equilibrium problems for more complex multiobjective games and multiple games.

Acknowledgments

This study was supported by the National Natural Science Foundation of China (71961003), the Qian Jiaohe YJSCXJH ((2019)029), the Qian Kehe Platform Talents ((2017)5788), and the Qian Kehe LH Zi ((2017)7223).

Conflict of interest

The authors declare no conflict of interest.