egttools.analytical.PairwiseComparison¶
- class PairwiseComparison(self: egttools.numerical.numerical.PairwiseComparison, population_size: int, game: egttools.numerical.numerical.games.AbstractGame)¶
Bases:
pybind11_object
A class containing methods to study analytically the evolutionary dynamics using the Pairwise comparison rule.
This class defines methods to compute fixation probabilities, transition matrices in the Small Mutation Limit (SML), gradients of selection, and the full transition matrices of the system when considering mutation > 0.
- Parameters
population_size (int) – Size of the population.
game (egttools.games.AbstractGame) – A game object which must implement the abstract class
egttools.games.AbstractGame
. This game will contain the expected payoffs for each strategy in the game, or at least a method to compute it, and a method to calculate the fitness of each strategy for a given population state.
See also
egttools.numerical.PairwiseComparisonNumerical
,egttools.analytical.StochDynamics
,egttools.games.AbstractGame
Note
Analytical computations should be avoided for problems with very large state spaces. This means very big populations with many strategies. The bigger the state space, the more memory and time these methods will require!
Also, for now it is not possible to update the game without having to instantiate PairwiseComparison again. Hopefully, this will be fixed in the future.
Methods
Calculates the fixation probability of an invading strategy in a population o resident strategy.
Calculates the gradient of selection without mutation for the given state.
Calculates the transition matrix of the reduced Markov Chain that emerges when assuming SML.
Computes the transition matrix of the Markov Chain which defines the population dynamics.
- __init__(self: egttools.numerical.numerical.PairwiseComparison, population_size: int, game: egttools.numerical.numerical.games.AbstractGame) None ¶
A class containing methods to study analytically the evolutionary dynamics using the Pairwise comparison rule.
This class defines methods to compute fixation probabilities, transition matrices in the Small Mutation Limit (SML), gradients of selection, and the full transition matrices of the system when considering mutation > 0.
- Parameters
population_size (int) – Size of the population.
game (egttools.games.AbstractGame) – A game object which must implement the abstract class
egttools.games.AbstractGame
. This game will contain the expected payoffs for each strategy in the game, or at least a method to compute it, and a method to calculate the fitness of each strategy for a given population state.
See also
egttools.numerical.PairwiseComparisonNumerical
,egttools.analytical.StochDynamics
,egttools.games.AbstractGame
Note
Analytical computations should be avoided for problems with very large state spaces. This means very big populations with many strategies. The bigger the state space, the more memory and time these methods will require!
Also, for now it is not possible to update the game without having to instantiate PairwiseComparison again. Hopefully, this will be fixed in the future.
- __new__(**kwargs)¶
- calculate_fixation_probability(self: egttools.numerical.numerical.PairwiseComparison, invading_strategy_index: int, resident_strategy_index: int, beta: float) float ¶
Calculates the fixation probability of an invading strategy in a population o resident strategy.
This method calculates the fixation probability of one mutant of the invading strategy in a population where all other individuals adopt the resident strategy.
- Parameters
- Returns
The fixation probability of one mutant of the invading strategy in a population where all other members adopt the resident strategy.
- Return type
See also
egttools.analytical.StochDynamics
,egttools.analytical.StochDynamics.fixation_probability
,egttools.analytical.PairwiseComparison.calculate_transition_matrix
,egttools.analytical.PairwiseComparison.calculate_transition_and_fixation_matrix_sml
,egttools.analytical.PairwiseComparison.calculate_gradient_of_selection
,egttools.numerical.PairwiseComparisonNumerical
,egttools.numerical.PairwiseComparisonNumerical.estimate_stationary_distribution
,egttools.numerical.PairwiseComparisonNumerical.estimate_stationary_distribution_sparse
,egttools.numerical.PairwiseComparisonNumerical.estimate_fixation_probability
- calculate_gradient_of_selection(self: egttools.numerical.numerical.PairwiseComparison, beta: float, state: numpy.ndarray[numpy.uint64[m, 1]]) numpy.ndarray[numpy.float64[m, 1]] ¶
Calculates the gradient of selection without mutation for the given state.
This method calculates the gradient of selection (without mutation), which is, the most likely direction of evolution of the system.
- Parameters
beta (float) – Intensity of selection
state (numpy.ndarray) – Vector containing the counts of each strategy in the population.
- Returns
Vector of shape (nb_strategies,) containing the gradient of selection, i.e., The most likely path of evolution of the stochastic system.
- Return type
See also
egttools.analytical.StochDynamics
,egttools.analytical.StochDynamics.full_gradient_selection
,egttools.analytical.PairwiseComparison.calculate_transition_matrix
,egttools.analytical.PairwiseComparison.calculate_transition_and_fixation_matrix_sml
,egttools.numerical.PairwiseComparisonNumerical
,egttools.numerical.PairwiseComparisonNumerical.estimate_stationary_distribution
,egttools.numerical.PairwiseComparisonNumerical.estimate_stationary_distribution_sparse
- calculate_transition_and_fixation_matrix_sml(self: egttools.numerical.numerical.PairwiseComparison, beta: float) Tuple[numpy.ndarray[numpy.float64[m, n]], numpy.ndarray[numpy.float64[m, n]]] ¶
Calculates the transition matrix of the reduced Markov Chain that emerges when assuming SML.
By assuming the limit of small mutations (SML), we can reduce the number of states of the dynamical system to those which are monomorphic, i.e., the whole population adopts the same strategy.
Thus, the dimensions of the transition matrix in the SML is (nb_strategies, nb_strategies), and the transitions are given by the normalized fixation probabilities. This means that a transition where i neq j, T[i, j] = fixation(i, j) / (nb_strategies - 1) and T[i, i] = 1 - sum{T[i, j]}.
This method will also return the matrix of fixation probabilities, where fixation_probabilities[i, j] gives the probability that one mutant j fixates in a population of i.
- Parameters
beta (float) – Intensity of selection
- Returns
A tuple including the transition matrix and a matrix with the fixation probabilities. Both matrices have shape (nb_strategies, nb_strategies).
- Return type
Tuple[numpy.ndarray, numpy.ndarray]
See also
egttools.analytical.StochDynamics
,egttools.analytical.StochDynamics.fixation_probability
,egttools.analytical.StochDynamics.transition_and_fixation_matrix
,egttools.analytical.PairwiseComparison.calculate_fixation_probability
,egttools.analytical.PairwiseComparison.calculate_transition_matrix
,egttools.analytical.PairwiseComparison.calculate_transition_and_fixation_matrix_sml
,egttools.analytical.PairwiseComparison.calculate_gradient_of_selection
,egttools.numerical.PairwiseComparisonNumerical
,egttools.numerical.PairwiseComparisonNumerical.estimate_stationary_distribution
,egttools.numerical.PairwiseComparisonNumerical.estimate_stationary_distribution_sparse
,egttools.numerical.PairwiseComparisonNumerical.estimate_fixation_probability
- calculate_transition_matrix(self: egttools.numerical.numerical.PairwiseComparison, beta: float, mu: float) scipy.sparse.csr_matrix[numpy.float64] ¶
Computes the transition matrix of the Markov Chain which defines the population dynamics.
It is not advisable to use this method for very large state spaces since the memory required to store the matrix might explode. In these cases you should resort to dimensional reduction techniques, such as the Small Mutation Limit (SML).
- Parameters
- Returns
Sparse vector containing the transition probabilities from any population state to another. This matrix will be of shape nb_states x nb_states.
- Return type
scipy.sparse.csr_matrix
See also
egttools.analytical.StochDynamics
,egttools.analytical.StochDynamics.calculate_full_transition_matrix
,egttools.analytical.PairwiseComparison.calculate_transition_and_fixation_matrix_sml
,egttools.numerical.PairwiseComparisonNumerical
,egttools.numerical.PairwiseComparisonNumerical.estimate_stationary_distribution
,egttools.numerical.PairwiseComparisonNumerical.estimate_stationary_distribution_sparse
- game(self: egttools.numerical.numerical.PairwiseComparison) egttools.numerical.numerical.games.AbstractGame ¶
- nb_states(self: egttools.numerical.numerical.PairwiseComparison) int ¶
- nb_strategies(self: egttools.numerical.numerical.PairwiseComparison) int ¶
- population_size(self: egttools.numerical.numerical.PairwiseComparison) int ¶
- update_population_size(self: egttools.numerical.numerical.PairwiseComparison, arg0: int) None ¶