Abstract
This document gives a few use cases for the EvolutionaryGames package. EvolutionaryGames provides basic concepts of evolutionary game theory, like e.g. finding evolutionary stable strategies and computing and drawing evolutionarily stable sets as well as phase diagrams for various evolutionary dynamics for single-population games with two, three and four different phenotypes.
Keywords: Evolutionary game theory, evolutionary stable strategies, evolutionarily stable sets, evolutionary game dynamics.The field of evolutionary game theory is concerned with interactions among large populations of strategically interacting agents. These agents may adapt their behaviour according to current payoff opportunities. Evolutionary game theory originates from biology (J. M. Smith and Price 1973), but nowadays has become an established tool in other disciplines, in particular in economics (see e.g. (Sandholm 2010)) and computer science (see e.g. (Suri 2007)).
The purpose of our package EvolutionaryGames is to enhance the R ecosystem by tools to plot and present evolutionary game dynamics in a both informative and attractive way. One similar tool in the public domain is William Sandholm’s package dynamo (Franchetti and Sandholm 2013) written in Mathematica. Dynamo (Franchetti and Sandholm 2013) frequently served as an inspiration for our implementation. Still, even though dynamo (Franchetti and Sandholm 2013) itself is free, Mathematica is proprietary software and it is nowhere near as widely used as R. Also, the package structures in R and the widespread availability of packages via CRAN make it much easier for the user to employ parts of our package EvolutionaryGames in his or her own software.
The following very brief presentation of concepts in evolutionary game theory and how EvolutionaryGames addresses them mainly follows the two books by Weibull (Weibull 1997) and Sandholm (Sandholm 2010). Furthermore, our vignette provides detailed references to the respective original articles. Our point is to show how various concepts from evolutionary game theory can be computed via the package EvolutionaryGames and where to find the conceptual, technical and mathematical details. In other words, we would like to stress that this little document is by no means intended to serve as an introduction to the fascinating field of evolutionary game theory itself. Again, for the latter the reader is referred to the the two books by Weibull (Weibull 1997) and Sandholm (Sandholm 2010). For an in-depth mathematical discussion we also recommend the book by Hofbauer and Sigmund (Hofbauer and Sigmund 1998).
The package EvolutionaryGames focusses on single-population games with two, three or four phenotypes. (The first author hopes to address multi-population games in a separate package based on EvolutionaryGames in the future.) For single-population games we regard EvolutionaryGames as a powerful and attractive alternative to William Sandholm’s Mathematica package dynamo (Franchetti and Sandholm 2013).
The concept of evolutionary stable strategies (ESS) was first proposed by Maynard Smith and Price in (J. M. Smith and Price 1973). An incumbent strategy \(x\) is said to be evolutionary stable if
If \(x\) is a strict Nash
equilibrium, then \(x\) must also be an
ESS.
Using our package EvolutionaryGames one can compute ESS for both \(2\) x \(2\) and \(3\) x \(3\) symmetric matrix games. The function
ESS.R receives three input arguments:
A
: A matrix specifying the symmetric gamestrategies
: A vector of strings of length \(2\) or \(3\) specifying the names of all
strategiesfloats
: A logical value that handles number
representation. If floats is set to TRUE
(default), a
floating-point number will be used for the output, otherwise the output
will be specified as fractions.Let us first take a look at a classical Hawk-Dove game.
library(EvolutionaryGames)
ESS(matrix(c(-1, 4, 0, 2), 2, byrow=TRUE), c("Hawk", "Dove"))
## Hawk Dove
## [1,] 0.6666667 0.3333333
ESS(matrix(c(-1, 4, 0, 2), 2, byrow=TRUE), c("Hawk", "Dove"), FALSE)
## Hawk Dove
## [1,] 2/3 1/3
It is well known that there are games which do not possess an ESS. A classical example is the game Rock-Scissors-Paper. In such a case our function ESS.R will return NULL.
library(EvolutionaryGames)
<- matrix(c(1, 2, 0, 0, 1, 2, 2, 0, 1), 3, byrow=T)) (A
## [,1] [,2] [,3]
## [1,] 1 2 0
## [2,] 0 1 2
## [3,] 2 0 1
ESS(A, c("Rock", "Scissors", "Paper"))
## NULL
The concept of evolutionarily stable sets was first discussed by B.
Thomas in (Thomas 1985). A very readable
introduction can also be found in the book by Weibull (Weibull 1997), section 2.4.1.
An evolutionarily stable set is a nonempty closed set of symmetric Nash
equilibrium strategies \(X\) such that
each strategy \(x \in X\) earns at
least the same payoff against any nearby alternative best reply \(y\) as \(y\) earns against itself with equal payoffs
limited to the case \(y \in X\). Note
that any evolutionary stable strategy (ESS) constitutes an
evolutionarily stable set and that the union of evolutionarily stable
sets is again an evolutionarily stable set. See the book by Weibull
(Weibull 1997), section 2.4.1, for further
details.
Evolutionarily stable sets are not easy to compute and to plot. The
package EvolutionaryGames computes evolutionarily stable sets of a game
with two players and three strategies in the case that the game has an
evolutionary stable strategy (ESS). If the two player three strategy
game has no ESS, then the code returns a message stating that our
algorithm cannot calculate evolutionarily stable sets for models that do
not have a proper ESS. The authors are very well aware that there are
games having evolutionarily stable sets but no proper ESS. Still, as our
package is not devoted to finding all symmetric Nash equilibria of a
game and as there currently is no package for this task on CRAN, we
decided only to handle the case of games with two players and three
strategies possessing at least one proper ESS. Within our algorithm, we
need a proper ESS as a starting point for our computations. The authors
feel that the possibility of computing and drawing evolutionarily stable
sets in such cases is precious for teaching (and research) purposes.
However, note that computing and drawing evolutionarily stable sets is
time consuming and is clearly the most elaborate task currently
performed by the package EvolutionaryGames. Finally, here is an example
of an evolutionarily stable set together with a plot including sample
trajectories (– the example itself is taken from the book (Broom and Rychtár 2013), p. 52):
library(EvolutionaryGames)
<- matrix(c(-2, 5, 10/9, 0, 5/2, 10/9, -10/9, 35/9, 10/9), 3, byrow=TRUE)
A <- c("Hawk", "Dove", "Mixed ESS")
strategies ESS(A, strategies)
## Hawk Dove Mixed ESS
## [1,] 0 0 1
ESset(A, strategies)
## [,1] [,2] [,3]
## [1,] 0.0000000 0.0000000 1
## [2,] 0.5555556 0.4444444 0
The basic game-theoretic model of biological natural selection is the replicator dynamic (Taylor and Jonker 1978). Still, various alternative dynamics have been proposed and investigated for different applications. Our package currently focusses on continuous dynamics only. The following itemization states which evolutionary dynamics are currently available in our package EvolutionaryGames:
eta
.eta
.Drawing phase diagrams for single-population games with two, three or four phenotypes with different dynamics is the main feature of the package EvolutionaryGames. We start with a rather simple Hawk-Dove game we studied before when analyzing for its ESS. Using phaseDiagram2S.R we obtain a phase diagram for the population share of hawks invading a population of doves under the replicator dynamics.
library(EvolutionaryGames)
<- matrix(c(-1, 4, 0, 2), 2, byrow=TRUE)
A phaseDiagram2S(A, Replicator, strategies = c("Hawk", "Dove"))
For similar phase diagrams, see e.g. the book by Peters (Peters 2015).
In the phase diagrams for the three strategies the user may specify the following parameters for phaseDiagram3S.R:
A
: A Numeric matrix of size \(3\) x \(3\) representing the number of strategies
of a symmetric matrix game.dynamic
: A function representing an evolutionary
dynamic.params
: A numeric vector with additional parameters for
the evolutionary dynamic, like e.g. in the cases of the Logit or ILogit
dynamics.trajectories
: A numeric matrix of size \(m\) x \(3\). Each row represents the initial values
for the trajectory to be examined.contour
: A logical value that handles contour diagram
presentation. By default FALSE
, a nicely coloured contour
plot will only be shown if the user sets
contour = TRUE
.vectorField
: A logical value that handles vector field
presentation. By default FALSE
, a vector field will only be
shown if the user sets vectorField = TRUE
.strategies
: A vector of strings of length \(3\) specifying the names of all strategies.
By default strategies = c("1","2","3")
.In the following we show three plots of the game Rock-Scissors-Paper under the Replicator dynamics with the same initial state:
library(EvolutionaryGames)
<- matrix(c(1, 2, 0, 0, 1, 2, 2, 0, 1), 3, byrow=T)
A <- matrix(c(0.7, 0.2, 0.1), 1, 3, byrow=TRUE)
state <- c("Rock", "Scissors", "Paper")
RSP phaseDiagram3S(A, Replicator, NULL, state, FALSE, FALSE, strategies = RSP)
library(EvolutionaryGames)
<- matrix(c(1, 2, 0, 0, 1, 2, 2, 0, 1), 3, byrow=T)
A <- matrix(c(0.7, 0.2, 0.1), 1, 3, byrow=TRUE)
state <- c("Rock", "Scissors", "Paper")
RSP phaseDiagram3S(A, Replicator, NULL, state, TRUE, FALSE, strategies = RSP)
library(EvolutionaryGames)
<- matrix(c(1, 2, 0, 0, 1, 2, 2, 0, 1), 3, byrow=T)
A <- matrix(c(0.7, 0.2, 0.1), 1, 3, byrow=TRUE)
state phaseDiagram3S(A, Replicator, NULL, state, TRUE, TRUE, strategies = RSP)
Finally, there is also the function phaseDiagram4S.R
for the case of a symmetric matrix game with four strategies. Its
parameters are rather similar to those of
phaseDiagram3S.R except for the facts that
trajectory
is a numeric vector of size \(4\) (rather than a matrix) and that there
is no possibility to enable contour plots or to plot vector fields.
Instead, there is an additional logical value noRGL
handling diagram rotation. By default TRUE
, the diagram
will only be constructed using rgl and thus become rotatable if the user
specifies noRGL = FALSE
. Finally, we present an example
using the Smith dynamic.
library(EvolutionaryGames)
<- matrix(c(5, -9, 6, 8, 20, 1, 2, -18, -14, 0, 2, 20, 13, 0, 4, -13), 4, 4, byrow=TRUE)
A <- c(0.6, 0.15, 0.1, 0.15)
state phaseDiagram4S(A, Smith, NULL, state, noRGL=TRUE)
EvolutionaryGames offers you to create your own dynamics. In particular, it is easy to write your own continuous dynamics. First of all, a dynamic is nothing other than a function that is passed as a parameter to the corresponding function for creating phase diagrams. The following code fragment shows you the minimum necessary structure of an arbitrary dynamic:
<- function (time, state, parameters) {
MyDynamic
#...
return(list (dX))
}
Any dynamic requires the parameters time
,
state
and parameters
. While the parameter
time
is used internally by deSolve to solve the initial
value problem, the other two parameters state
and
parameters
are used to specify the desired dynamic and are
available as numeric vectors. In this context, state
stands
for the desired initial state under which the model is to be simulated
and parameters
contains further parameters, such as the
symmetric matrix specifying the game and, depending on the dynamic,
noise levels or similar parameters.
As can be seen in the example, the return value of a specified dynamic has to be a numerical list. Each component represents the corresponding rate of change of a phenotype under the respective dynamic.
We now show how to implement a very well known dynamic, the replicator dynamics. Our function definition is as follows:
<- function (time, state, parameters) {
Replicator
#...
return(list (dX))
}
The above game and an initial value is passed as a parameter to a function for the generation of phase diagrams and made retrievable within the dynamic. A small limitation, however, is that our matrix needs to be converted into a vector due to certain constraints in our internal usage of the package deSolve. We recommend that you first transfer this vector back into a matrix and maintain the original game. This can be done as follows:
<- function (time, state, parameters) {
Replicator <- parameters
a <- sqrt(length(a))
states <- matrix(a, states, byrow = TRUE)
A <- t(A) # original symmetric game
A
return(list (dX))
}
We now come to the actual part of the implementation. We first calculate the rate of change of each phenotype depending on the other phenotypes. Immediately afterwards, we calculate the average fitness of each phenotype and then set up the actual replicator dynamics:
<- function(time, state, parameters) {
Replicator <- parameters
a <- sqrt(length(a))
states <- matrix(a, states, byrow = TRUE)
A <- t(A)
A
<- c()
dX
for(i in 1:states) {
<- sum(state * A[i, ])
dX[i]
}
<- sum(dX * state)
avgFitness
for(i in 1:states) {
<- state[i] * (dX[i] - avgFitness)
dX[i]
}
return(list(dX))
}
Our dynamics is now applicable within the predefined functions for generating phase diagrams.
We made use of various packages of the CRAN ecosystem. In particular, it was of paramount importance not to write any differential equation solvers ourselves, but to make use of an established solver, instead.
The R-package deSolve (Soetaert, Petzoldt, and Setzer 2010) has long been established as a powerful solver for differential equations. In our context, deSolve is used to generate the data needed in order to visualize the time evolution of a given game under a certain dynamic, i.e. we use it for obtaining the correct input for the visualization of trajectories in our phase diagrams.
We use the package rgl (Adler, Murdoch, et al. 2014), i.e. the R interface to OpenGL, for the four strategies case, because here the phase diagram represents a three-dimensional simplex. Our users have the possibility to follow the development of a game under a given dynamic for all phenotypes by rotating the three-dimensional simplex generated by rgl.
The R package geometry (Habel et al. 2015) helped us with the conversion from barycentric to cartesian coordinates for drawing trajectories. Solving the initial value problem returns the rates of change of the model, which have to be converted into cartesian coordinates before they can be drawn into the phase diagram in order to finally form a trajectory.