Fixation dynamics on hypergraphs

Hypergraphs have been a useful tool for analyzing population dynamics such as opinion formation and the public goods game occurring in overlapping groups of individuals. In the present study, we propose and analyze evolutionary dynamics on hypergraphs, in which each node takes one of the two types of different but constant fitness values. For the corresponding dynamics on conventional networks, under the birth-death process and uniform initial conditions, most networks are known to be amplifiers of natural selection; amplifiers by definition enhance the difference in the strength of the two competing types in terms of the probability that the mutant type fixates in the population. In contrast, we provide strong computational evidence that a majority of hypergraphs are suppressors of selection under the same conditions by combining theoretical and numerical analyses. We also show that this suppressing effect is not explained by one-mode projection, which is a standard method for expressing hypergraph data as a conventional network. Our results suggest that the modeling framework for structured populations in addition to the specific network structure is an important determinant of evolutionary dynamics, paving a way to studying fixation dynamics on higher-order networks including hypergraphs.


Introduction
Populations of individuals, cells, and habitats, on which evolutionary processes take place often have structure that may be described by networks.Evolutionary graph theory enables us to mathematically and computationally investigate population dynamics in which multiple types of different fitness compete on networks under a selection pressure [1][2][3][4].A minimal evolutionary dynamics model on graphs, or networks, is to assume that there are two types of different fitness values, 1 and r, which are constant over time (i.e., constant selection), and that each node is occupied by either type, which can change over time.A question then is which type eventually occupies the entire population in the absence of mutation during evolutionary dynamics, i.e., which type fixates [5].In physics and mathematics literature, the same model is often called the biased voter model, and the fixation is often called the consensus [6][7][8][9].Crucially, the fixation probability of each type, as well as other properties of the evolutionary dynamics such as the average time to fixation, depends on the network structure in addition to the value of r.Application of evolutionary graph theory to social dilemma games has also been successful in giving analytical insights into various conditions under which cooperation is favored over defection [1-3, 10, 11].
Evolutionary graph theory with conventional networks assumes that interaction between nodes occurs pairwise because each edge in a conventional network represents direct connectivity between two nodes.However, in reality, more than two individuals may simultaneously interact and compete in evolutionary dynamics.For example, an evolutionary dynamics under constant selection can be regarded as a model of opinion dynamics in a population of human or animal individuals, in which they may meet in groups with more than two individuals for opinion formation.As another example, the public goods game is naturally defined for group interaction; each individual in the group decides whether or not to contribute to a collective good, and all individuals will receive a share of an augmented amount of the collective good regardless of whether or not they contributed.Hypergraphs are a natural extension of conventional networks to the case of group interactions [12][13][14].In a hypergraph, the concept of edge is extended to the hyperedge, which represents a unit of interaction and may contain more than two nodes.Population dynamics on hypergraphs such as evolutionary dynamics of social dilemma games [15,16] and opinion formation [17,18] have been investigated (see [13,19] for reviews including different types of dynamics).
A key question in the study of constant-selection dynamics on networks is whether the given network is amplifier or suppressor of natural selection [1,5].Suppose that the resident and mutant type have constant fitness values 1 and r, respectively.We anticipate that the mutant type is more likely to fixate than the resident type under the same condition if r > 1 and vice versa.In fact, how much the fixation probability of the mutant type increases as one increases r depends on the network structure.Some networks are amplifiers of selection; on these networks, a single mutant has a higher probability of fixation than in the well-mixed population, corresponding to the Moran process, at any r > 1 and a smaller fixation probability than in the Moran process at any r < 1.Other networks are suppressors of selection; on these networks, a single mutant has a lower fixation probability than in the Moran process at any r > 1 and a higher fixation probability than in the Moran process at any r < 1.Under the socalled birth-death processes, which we focus on, most networks are known to be amplifiers of selection, at least when the initial mutant is located on a node selected uniformly at random [20][21][22].Research has discovered various classes of amplifiers of selection [5,[23][24][25][26][27] while few for suppressors of selection [28].
In the present study, we study constant-selection evolutionary dynamics on hypergraphs without mutation.We propose two such models and formulate the fixation probability of mutants using theory based on Markov chains, which extends the same approach for conventional networks to the case of hypergraphs.Then, we examine whether hypergraphs are amplifiers of selection, suppressors of selection, or neither.We show to our numerical efforts that, contrary to the known results for networks, all the hypergraphs that we have investigated are suppressors of selection even under the birth-death process and uniform initial condition.We reach this conclusion by semi-analytical investigations for hypergraphs with high symmetry and numerical simulations on empirical hypergraphs.

Models
We introduce two models of evolutionary dynamics on undirected hypergraphs.Let H be an undirected hypergraph with node set V = {1, . .., N}, where N is the number of nodes, and hyperedge set E. Each e 2 E, where e � V and e 6 ¼ ;, is a hyperedge, intuitively representing group interaction among the nodes belonging to e.If each e 2 E is a set containing exactly two nodes, H is a conventional undirected network.We assume that H is connected in the sense that there is a hyperedge intersecting both W and V − W for every non-empty proper subset W of V.
We define model 1 by the following discrete-time evolutionary dynamics on hypergraph H, which extends the birth-death process for conventional networks (precisely speaking, birthdeath process with selection on the birth, or the Bd rule [29][30][31]) and the Moran process in well-mixed populations.We assume that there are two types of individuals, referred to as A and B, and that A and B have fitness r and 1, respectively.We refer to A as the mutant type and B as the resident type.In each time step, we select one node for reproduction with the probability proportional to its fitness.We call this node the parent.Then, the parent node selects one of the hyperedges to which it belongs, denoted by e, with the equal probability, regardless of the size of e (i.e., the number of nodes contained in e).Finally, the parent converts the type (i.e., A or B) of all other nodes belonging to e into the parent's type.We repeat this process until all nodes in V have the same type.Once this unanimity configuration is reached, which is the fixation, no node will change its type even if one further runs the evolutionary dynamics.We examine the probability that every node eventually has type A, which we call the fixation probability of type A.
Model 2 is the same as model 1 in that there are two types of individuals with fitness r and 1, that we select a parent node in each time step with the probability proportional to its fitness, and that the parent selects one of its hyperedges, e. Differently from the model 1, the parent node converts all the other nodes belonging to e into the parent's type if and only if the parent's type is the majority in e, i.e., if more than half of the nodes in e including the parent have the same type as the parent's type.Nothing occurs in a time step if the parent's type is not the majority in e, including in the case of a tie when the size of e is even.By assumption, a type has to be a majority in an hyperedge for it to spread even when it is stronger than the other type (i.e., type A stronger than type B when r > 1, and vice versa when r < 1).Model 2 is an evolutionary dynamics variant of opinion formation models under the majority rule [9,32,33].
The assumption that all nodes belonging to a hyperedge e copy the parent's type, made for both models, may seem strong, especially when the size of e is large.In fact, if we allow only one node belonging to e to copy the parent's type in a time step, then it is straightforward to show that model 1 is equivalent to the birth-death process on the conventional network called the one-mode projection (see sections 3.1.4and 4 for results for the one-mode projection).Therefore, by assuming that all nodes in e update their type in a time step, we examine effects of an extreme case of allowing simultaneous updating of nodes constrained by the hypergraph structure.
When the network is the complete graph, which is a conventional undirected network and therefore a hypergraph, model 1 is called the Moran process.For the Moran process, the fixation probability of type A when there are initially i individuals of type A, denoted by x i , is given by [1] 3 Results for synthetic hypergraphs In this section, for three model hypergraphs that are mathematically convenient, we calculate the fixation probability for mutant type A when there are initially i nodes of type A and N − i nodes of type B, for both models 1 and 2. For the cyclic 3-uniform hypergraph, which we introduce in section 3.1.2,we impose i = 1 and i = 2 for models 1 and 2, respectively, for a technical reason to be explained.The fixation probability depends on the initial condition.We select the i mutant nodes from the N nodes uniformly at random, which is called the uniform initialization [34].
For any hypergraph of size N, at any time step, either type A or B inhabits each node.Therefore, there are 2 N states in total.The fixation probability for type A depends on each state, not just on i.To know the fixation probability for type A for the initial states with i mutants, we need to solve a linear system of 2 N − 2 unknowns, where each unknown is the fixation probability for an initial state.Note that we safely excluded the two unknowns corresponding to the two initial states in which all nodes unanimously have type A or B. We did so because the fixation probability for type A is trivially equal to 1 and 0 for these two states.To solve a linear system with 2 N − 2 unknowns is daunting except when N is small.Therefore, we analyze three types of symmetric hypergraphs in which all or most nodes are structurally equivalent to each other.For these hypergraphs, we only need to track the number of nodes with type A among the structurally equivalent nodes.By doing so, we can drastically reduce the dimension of the linear system to be analytically or numerically solved.In this section, we denote the number of nodes of type A in the entire hypergraph by i 2 {0, . .., N}.The fixation of type A and B corresponds to i = N and i = 0, respectively.

Model 1
3.1.1Complete 3-uniform hypergraph.We have mentioned that our model 1 on the complete graph is equivalent to the Moran process.To investigate whether model 1 on counterparts of the complete graph for hypergraphs is equivalent to the Moran process, we consider the complete 3-uniform hypergraph [35].A complete 3-uniform hypergraph on node set V is defined by hyperedge set E, which is the set of all subsets of V containing just three nodes.In other words, In this section, we refer to i, i.e., the number of nodes of type A, as the state of the evolutionary dynamics.Note that knowing the dynamics of i is enough for completely understanding the evolutionary dynamics on the complete 3-uniform hypergraph owing to its symmetry.Under model 1, state i either remains unchanged or moves to i − 2, i − 1, i + 1, or i + 2 in a single time step of the evolutionary dynamics.This is because every hyperedge of the complete 3-uniform hypergraph has three nodes, and therefore there are at most two nodes that flip their type in a time step.
We denote the (N + 1) × (N + 1) transition probability matrix by P = [p i,j ], where p i,j is the probability that the state moves from i to j in a time step.At state i, the probability that a node of type A and B is selected as parent is equal to ri/(ri + N − i) and (N − i)/(ri + N − i), respectively.If the selected parent node is of type A, a hyperedge containing the parent and two nodes of type B is used for reproduction with probability NÀ i , where n k À � represents the binomial coefficient "n choose k".In this case, i increases by 2. Alternatively, a hyperedge containing the parent, a different node of type A, and a node of type B is used for reproduction with probability iÀ À � .In this case, i increases by 1. Otherwise, a hyperedge containing the parent and two other nodes of type A is used for reproduction.In this case, i does not change.If the parent is of type B, a hyperedge containing the parent and two nodes of type A is used for reproduction with probability i 2 À � = NÀ 1 2 À � .In this case, i decreases by 2. Alternatively, a hyperedge containing the parent, a node of type A, and a different node of type B is used for reproduction with probability i In this case, i decreases by 1. Otherwise, a hyperedge containing the parent and two other nodes of type B is selected.In this case, i does not change.Therefore, the transition probabilities are given by All the other entries of P are equal to zero.Therefore, P is a pentadiagonal matrix.States i = 0 and i = N are absorbing states.
Denote by x i the probability of ending up in state N, i.e., fixation probability of the mutant type A, when the initial state is i.We obtain In vector notation, we can concisely write Eqs ( 11) to (15) as where x = (x 0 , x 1 , x 2 , . .., x N ) > , and > represents the transposition.Eq ( 16) is equivalent to where I is the identity matrix.Using Eqs ( 11), ( 15) and ( 17), we obtain

PLOS COMPUTATIONAL BIOLOGY
where b = (0, 0, . .., 0, 1) > , and Like P, the (N + 1) × (N + 1) matrix M is a pentadiagonal matrix.PTRANS-I and PTRANS-II are numerical algorithms for efficiently solving pentadiagonal linear systems [36].Note that we are calculating x although we only need x 1 .These two algorithms run in O(log n) time, where n is the number of unknowns, and, for n = 4000, they are about ten times faster than the SciPy algorithm for banded matrices, scipy.linalg.solve_banded[37].We use the PTRANS-II built in a Python package pentapy [38] to calculate the fixation probability as a function of r unless N is small.For a given value of r, the calculation requires O (log N) time, which is much faster than solving a full linear system of 2 N − 2 unknowns.
For small complete 3-uniform hypergraphs, one can analytically calculate x 1 by directly solving Eq (18).We obtain for N = 4 and for N = 5.We show x i for i 2 {2, . .., N − 1}, N 2 {4, 5} in Text A in S1 Text.We compare Eqs (20) and ( 21) with x 1 for the Moran process (i.e., Eq (1) with i = 1) with N = 4 and N = 5 as a function of r in Fig 2A and 2B, respectively.Although the results are shown by blue lines, they are hidden behind the orange lines.The figure indicates that these complete 3-uniform hypergraphs are suppressors of selection because their x 1 is smaller than that for the Moran process for r > 1 and larger for r < 1.We also obtained x 1 by numerically solving Eq (18) for N = 20 and N = 200.The results shown in Fig 2C and 2D for N = 20 and N = 200, respectively, indicate that these larger complete 3-uniform hypergraphs are also suppressors of selection.Therefore, we conclude that the complete 3-uniform hypergraphs are suppressors of selection under the evolutionary dynamics described by model 1.

Cyclic 3-uniform hypergraph.
The fixation probability as a function of r for undirected cycle graphs is the same as that for the Moran process because the cycle graphs are socalled isothermal graphs under the birth-death process [1,5,31].To examine whether the same equivalence result holds true for hypergraphs, we consider an extension of the cycle graph to the hypergraph, which we call the cyclic 3-uniform hypergraph.A cyclic 3-uniform hypergraph consists of a node set V = {1, . .., N} and a hyperedge set E = {{1, 2, 3}, {2, 3, 4}, . .., {N − 2, N − 1, N}, {N − 1, N, 1}, {N, 1, 2}}, i.e., any three consecutively numbered nodes with a periodic boundary condition form a hyperedge.The cyclic 3-uniform hypergraph with N = 3 nodes is the complete 3-uniform hypergraph.Therefore, we consider cyclic 3-uniform hypergrapys with N � 4. We show the cyclic 3-uniform hypergraph on 8 nodes in Fig 1B.
We assume that there is initially one individual of type A. Then, under model 1, all the nodes of type A are consecutively numbered at any time, without being interrupted by nodes of type B. Therefore, similarly to the analysis of evolutionary dynamics on cycles [39,40], it suffices to track the number of nodes of type A, which we again denote by i, to understand the evolutionary dynamics on this hypergraph.It should be noted that this technique can only be used when i = 1 or when the i (� 2) nodes of type A are consecutively located on the cycle in the initial state.
When i = 1, there are three types of events that can occur next.Without loss of generality, we assume that the ℓth node is the only node of type A (see Fig 3A for a schematic).The state moves from i = 1 to i = 0 in one time step such that B fixates in the following two cases.In the first case, either node v ℓ−2 or v ℓ+2 , which is of type B, is selected as parent.If N � 5, these two nodes are distinct.Therefore, this event occurs with probability 2/(r + N − 1).Then, the hyperedge that contains the parent and v ℓ is used for reproduction, which occurs with probability 1/ 3.For example, if v ℓ−2 is selected as parent and hyperedge {ℓ − 2, ℓ − 1, ℓ} is used for reproduction, then the state moves from 1 to 0. In the second case, either v ℓ−1 or v ℓ+1 is selected as parent, which occurs with probability 2/(r + N − 1).Then, one of the two hyperedges that contain the parent and v ℓ is used for reproduction, which occurs with probability 2/3.For example, if v ℓ−1 is selected as parent and hyperedge {ℓ − 2, ℓ − 1, ℓ} or {ℓ − 1, ℓ, ℓ + 1} is used for reproduction, then the state moves from 1 to 0. By summing up these probabilities, we obtain In fact, Eq (22) also holds true for N = 4 although v ℓ−2 and v ℓ+2 are identical nodes when N = 4; see Text B in S1 Text for the proof.Alternatively, the state moves from i = 1 to i = 3 whenever v ℓ is selected as parent, which occurs with probability r/(r + N − 1).Therefore, we obtain If any other event occurs, then i = 1 remains unchanged.Therefore, we obtain and All the other entries of transition probability matrix P are equal to zero.By solving Eq (18), we obtain for N = 4 and for N = 5.We show x i for i 2 {2, . .., N − 1}, N 2 {4, 5} in Text A in S1 Text.We compare Eqs (35) and ( 36) with x 1 for the Moran process with N = 4 and N = 5 in Fig 2A and 2B, respectively.We find that these cyclic 3-uniform hypergraphs are suppressors of selection.The result for N = 4 given by Eq (35) coincides with that for the complete 3-uniform hypergraph given by Eq (20).For larger N, we use the same numerical method for solving Eq (18) as that for the complete 3-uniform hypergraph because P is a pentadiagonal matrix.We show the thus calculated x 1 for N = 20 and N = 200 in Fig 2C and 2D, respectively.These figures confirm that these larger cyclic 3-uniform hypergraphs are also suppressors of selection.Therefore, we conclude that the cyclic 3-uniform hypergraph is a suppressor of selection.

Star 3-uniform hypergraph.
The conventional star graphs are strong amplifiers of selection [5,26].To examine whether counterparts of the star graph for hypergraphs are also amplifiers of selection, we define the star 3-uniform hypergraph as follows and examine the fixation probability of the mutant on it.A star 3-uniform hypergraph consists of a node set V with a single hub node and N − 1 leaf nodes, and hyperedges each of which is of size three and consists of the hub node and a pair of N − 1 leaf nodes.We use all the NÀ 1 2 À � pairs of leaf nodes to form hyperedges such that the hub node belongs to Owing to the structural equivalence among the N − 1 leaf nodes, the state of the evolutionary dynamics on the star 3-uniform hypergraph is completely determined by the type on the hub (i.e., either A or B) and the number of leaf nodes of type A, which ranges between 0 and N − 1.Therefore, we denote the state by (i 1 , i 2 ), where i 1 = 1 or 0 if the hub node is of type A or B, respectively, and i 2 2 {0, 1, . .., N − 1} represents the number of the leaf nodes of type A. The total number of nodes of type A is given by i = i 1 + i 2 .The fixation of type A and B corresponds to (i 1 , i 2 ) = (1, N − 1) and (0, 0), respectively.
We denote the probability that type A fixates starting with state (i where the 2N × 2N stochastic matrix P is given by and C, D, E, and F are N × N matrices given in Text G in S1 Text.The transition matrix P is a sparse matrix.We use the scipy implementation of the DGESV routine of LAPACK, scipy.linalg.solve, to numerically solve Eq (37) to obtain xð0;iÞ and xð1;iÀ 1Þ .We remind that x i is the probability that type A fixates when there are initially i nodes of type A that are selected uniformly at random.There are NÀ 1 iÀ 1

À �
states in which i − 1 nodes have type A and the hub has type A. There are NÀ 1 i À � states in which i nodes have type A and the hub has type B. Therefore, we obtain We analytically solve Eqs ( 37) and ( 39) to obtain for N = 4 and for N = 5.We show x i for i 2 {2, . .., N − 1}, N 2 {4, 5} in Text A in S1 Text.We compare Eqs (40) and ( 41) with x 1 for the Moran process with N = 4 and N = 5 in Fig 2A and 2B, respectively.We find that these star 3-uniform hypergraphs are suppressors of selection.We then computed x 1 by numerically solving Eq (37)  We also verified that the star 3-uniform hypergraph with N = 1500 nodes is also a suppressor of selection although the suppressing effect is even weaker (see Text H in S1 Text).This result is consistent with the analytically known result that larger conventional star graphs are stronger amplifiers of selection [5,26].Therefore, we conclude based on our numerical results that the star 3-uniform hypergraph is likely to be a suppressor of selection for any N.

Comparison with conventional networks obtained by one-mode projection.
A lossy representation of a hypergraph as a conventional network is the one-mode projection, in which two nodes are adjacent by an edge if and only if they belong to at least one common hyperedge.A weighted network version of the one-mode projection defines the edge weight by the number of hyperedges that the two nodes share.The one-mode projection is a common approach for analyzing dynamics on hypergraphs including evolutionary dynamics [4].The one-mode projections of the complete and cyclic 3-uniform hypergraphs are regular networks (i.e., networks in which all the nodes have the same degree), in the case of both unweighted and weighted variants of the one-mode projection.Then, the isothermal theorem [5] guarantees that the birth-death process on the one-mode projections of the complete and cyclic 3-uniform hypergraphs is equivalent to the Moran process.Therefore, our result that the complete and cyclic 3-uniform hypergraphs are suppressors of selection is not an artifact of onemode projection.
The same argument does not apply for the star 3-uniform hypergraph because the onemode projection of the star 3-uniform hypergraph is not a regular network.Furthermore, the star 3-uniform hypergraph are only analogously similar to the conventional star graph.In fact, leaf nodes are adjacent to each other in the star 3-uniform hypergraph, whereas they are not in the conventional star graph.Then, the direct connection between leaf nodes might be a reason why the star 3-uniform hypergraph is a suppressor of selection.To exclude this possibility, using similar analytical techniques, we investigated the fixation probability for the weighted one-mode projection of the star 3-uniform hypergraph, which is a weighted complete graph.Note that the unweighted one-mode projection of the star 3-uniform hypergraph is the unweighted complete graph, which is trivially equivalent to the Moran process.As we show in Text I in S1 Text, we have found that the obtained weighted complete graph is an amplifier of selection although the amplifying effect is weak.Therefore, our result that the star 3-uniform hypergraph is a suppressor of selection is not expected from the one-mode projection.

Model 2
In this section, we analyze the fixation probability for the birth-death process governed by model 2 on the complete, cyclic, and star 3-uniform hypergraphs.

Complete 3-uniform hypergraph.
On the complete 3-uniform hypergraph, the state i either remains unchanged or moves to i − 1 or i + 1 in a single time step because all the hyperedges are composed of three nodes such that there are at most one node that changes the state under the majority rule given by model 2. At state i, the probability that a node of type A and B is selected as parent is equal to ri/(ri + N − i) and (N − i)/(ri + N − i), respectively.If the parent is of type A, then the following three types of events are possible.First, if a hyperedge containing the parent and two nodes of type B is used for reproduction with probability , then the state does not change.Second, if a hyperedge containing the parent, a different node of type A, and a node of type B is used for reproduction with probability , then the state moves from i to i + 1.Third, if a hyperedge containing the parent and two other nodes of type A is used for reproduction with the remaining probability, then the state does not change.If the parent is of type B, then the following three types of events are possible.First, if a hyperedge containing the parent and two nodes of type A is used for reproduction with probability i 2 À � = NÀ 1 2 À � , then the state does not change.Second, if a hyperedge containing the parent, a node of type A, and a different node of type B is used for reproduction with probability i , then the state moves from i to i − 1.Third, if a hyperedge containing the parent and two other nodes of type B is used for reproduction with the remaining probability, then the state does not change.Therefore, the transition probabilities are given by All the other entries of P are equal to zero.Therefore, P is a tridiagonal matrix.We note that the states i = 0 and i = N are absorbing states.
The fixation probability of type A starting from state i, i.e., x i , satisfies Similar to the analysis of the fixation probability for the Moran process [1], we set . Eq (48) leads to y i+1 = y i γ i with i 2 {2, . .., N − 2}.Therefore, we obtain y 1 = 0, By summing all these expressions and using y N = 0, which Eq (49) implies, we obtain Therefore, we obtain and This expression does not simplify further.
For model 2, it always holds true that x 1 = 0 because the evolutionary dynamics is driven by a majority rule.Therefore, we compare x 2 , instead of x 1 , as a function of r with x 2 for the Moran process, to examine whether a given hypergraph is an amplifier of selection, suppressor of selection, equivalent to the Moran process, or neither.We compare x 2 calculated from Eq (52) with that for the Moran process at four values of N in Fig 4 .The figure indicates that the complete 3-uniform hypergraph with 4 nodes is a suppressor of selection.However, the complete 3-uniform hypergraph with N > 4 is neither amplifier nor suppressor of selection.This is because Eq (52) implies that x 2 = 2 3−N when r = 1, which is different from the result for the Moran process, i.e., x 2 = 2/N, when N > 4.

Cyclic 3-uniform hypergraph.
Consider the evolutionary dynamics under model 2 on the cyclic 3-uniform hypergraph.Assume that there are initially two nodes of type A, which are distributed uniformly at random, and N − 2 nodes of type B. We derived in Text J in S1 Text the fixation probability for type A, x 2 , using techniques similar to those used for the complete 3-uniform hypergraph (section 3.2.1).We obtain We show x 3 for N = 5 in Text A in S1 Text; note that x 2 for N 2 {4, 5} and x 3 for N = 5 are the only nontrivial x i values for N 2 {4, 5} because x 1 = 0 and x N−1 = 1 for model 2.
We compare x 2 given by Eq (54) with x 2 for the Moran process at four values of N in

Star 3-uniform hypergraph.
We calculate the fixation probability for model 2 on the star 3-uniform hypergraph.As in the case of model 1, we exploit the fact that the combination of the type of the hub (i.e., either A or B) and the number of leaf nodes of type A, which ranges between 0 and N − 1, completely specifies the state of the evolutionary dynamics on this hypergraph.We show the derivation of the fixation probability in Text K in S1 Text.
We obtain for N = 4 and for N = 5.We show x 3 for N = 5 in Text A in S1 Text.We compare Eqs ( 55) and ( 56) with x 2 for the Moran process with N = 4 and N = 5 in Fig 4A and 4B, respectively.We find that the star 3-uniform hypergraph with N = 4 is a suppressor of selection.In contrast, the star 3-uniform hypergraph with N = 5 is neither an amplifier nor suppressor of selection because x 2 = 9/35 at r = 1, which is different from the corresponding result for the Moran process, i.e., x 2 = 2/5.We also obtained x 2 for N = 20 and N = 200 by numerically solving a system of linear equations with 2N − 2 unknowns (see Text K in S1 Text).Fig 4C and 4D, which compare the obtained x 2 values with x 2 for the Moran process for N = 20 and N = 200, respectively, indicate that these larger star 3-uniform hypergraphs are neither amplifiers nor suppressors of selection.Therefore, we conclude that star 3-uniform hypergraphs are neither amplifiers nor suppressors of selection for N � 5 under the evolutionary dynamics described by model 2.

Neutral drift
The discussion of amplifier and suppressor of selection requires x 1 = 1/N (or x 2 = 2/N in the case of our model 2) when r = 1.In this section, we prove the following theorem, which justifies such discussion for model 1 and not for model 2.
Theorem 1.Consider the neutral drift, i.e., r = 1.When there are initially i mutants (i.e., type A) that are distributed uniformly at random over the N nodes, the fixation probability for the mutant is equal to i/N.
Proof.Consider model 1', which is defined by the same evolutionary dynamics as that for model 1 with the initial condition in which each node i is occupied by a distinct neutral type Ãi .Therefore, there are initially N types, all of which have the constant fitness equal to 1.The evolutionary dynamics terminates when a single type fixates.We denote by q i the fixation probability for type Ãi under the aforementioned initial condition.It should be noted that P N i¼1 q i ¼ 1.Now we consider the original model 1 with the initial condition in which there is just one mutant located at node i.Then, the fixation probability for the mutant is equal to q i .This is because model 1' is reduced to model 1 with the present initial condition if we regard type Ãi as the mutant type and all the N − 1 types Ãj , where j 6 ¼ i, as the resident type.We recall that x 1 is the average of the fixation probability over the N initial conditions in which there is just one mutant.Therefore, we obtain x 1 ¼ P N i¼1 q i =N ¼ 1=N.Next we consider model 1 with the initial condition in which there are two mutants whose locations are selected uniformly at random.Assume that the two mutants are initially located at nodes 1 and 2. The mutant type fixates with probability q 1 + q 2 because model 1' is reduced to model 1 with the present initial condition if we regard type Ã1 and Ã2 as the mutant type and all the other N − 2 types Ãj , where j 6 ¼ 1, 2, as the resident type.Therefore, we obtain Similarly, we obtain Remark 2. This theorem is known for the birth-death process on conventional networks with the proof being unchanged [41][42][43].
Remark 3. The theorem does not hold true for model 2. This is because, in model 2, whether the parent node u propagates its type to the other nodes in the selected hyperedge e containing u depends on the other nodes belonging to e.As an example, consider the case in which e contains three nodes, u, v 1 , and v 2 .Parent u imposes its type on v 1 and v 2 only when either v 1 or v 2 has the same type as u.Otherwise, i.e., if both v 1 and v 2 are of the opposite type as u's, then no state change occurs after u is selected as parent and hyperedge e is selected.Model 1' cannot handle this situation.In model 1', the parent disseminates its type to all the other nodes in the selected hyperedge e, as in model 1, regardless of the type of the other nodes belonging to e. Therefore, we cannot map model 1' to model 2. In fact, Fig 4 shows that x 2 6 ¼ 2/N in most cases, verifying that Theorem 1 does not hold true in general.

Numerical results for empirical hypergraphs
In this section, we carry out numerical simulations for empirical hypergraphs.The present study is motivated by the result that most networks are amplifiers of selection under the birthdeath process [20][21][22].In section 3.1, we showed that three model hypergraphs are suppressors of selection under model 1.In section 3.2, we argued that one cannot discuss whether the same hypergraphs are amplifiers or suppressors of selection under model 2 except for the complete, cyclic, and star 3-uniform hypergraphs with N = 4 nodes.This is because model 2 does not respect x 2 = 2/N in general (see Remark 3).Therefore, we focus on model 1 in this section and examine whether empirical hypergraphs tend to be suppressors of selection, as is the case for the complete, cyclic, and star 3-uniform hypergraphs.
Empirical, or general, hypergraphs are distinct from the complete, cyclic, and star 3-uniform hypergraphs in two main aspects.First, in general, empirical hypergraphs do not have much symmetry that we can exploit to simplify the probability transition matrix P. Therefore, pursuing analytical solutions, which would involve the solution of a linear system with 2 N − 2 unknowns, is formidable unless N is small.Second, empirical hypergraphs contain hyperedges of different sizes, whereas the 3-uniform hypergraphs only have hyperedges of size 3. Therefore, in this section, we numerically examine the stochastic evolutionary dynamics on empirical hypergraphs.In each time step, we select a parent node with the probability proportional to its fitness from all the N nodes.Then, the parent propagates its type to all the other nodes in a hyperedge e to which the parent belongs.We select e uniformly at random from all the hyperedges to which the parent belongs regardless of the size of the hyperedge.
To calculate the fixation probability of a single mutant of type A for an arbitrary hypergraph and for each value of r, we run the birth-death process until type A or B fixates for each node v initially occupied by type A. Note that all the N − 1 nodes except v are initially of type B. For each v, we run 3 × 10 3 simulations for r � 1 and 4 × 10 4 simulations for r < 1.We use a substantially larger number of simulations for r < 1 than r � 1 because the fixation probability is small when r is small and therefore it can be estimated with a higher accuracy with more simulations.For a given value of r, we estimate the fixation probability of type A as a fraction of the runs in which type A has fixated among the 3 × 10 3 × N or 4 × 10 4 × N simulations; the factor N is due to the N different choices of v.
We examine the fixation probability on four empirical hypergraphs.The corporate club membership hypergraph (corporate hypergraph for short) contains 25 nodes and 15 hyperedges [44,45].Each node represents a corporate executive officer.Each hyperedge represents a social organization such as clubs and boards of which some officers are members.The Enron hypergraph is an email communication network and has 143 nodes and 10,885 hyperedges [46,47].Each node represents an email address at Enron.Each hyperedge is comprised of all recipient addresses of an email.The Senate committee hypergraph (Senate hypergraph for short) has 282 nodes and 315 hyperedges.Each node is a member of the US Senate.Each hyperedge represents committee memberships [48,49].The high-school hypergraph is a social contact network with 327 nodes and 7,818 hyperedges.Each node is a student.Each hyperedge represents to a group of students that were in proximity of one another [49,50].All the four hypergraphs are connected hypergraphs.
In Fig 5, we show the relationships between the fixation probability for a single mutant as a function of r, for the four empirical hypergraphs, one per panel.We find that all the hypergraphs are suppressors of selection.We also simulate the birth-death process on the weighted one-mode projection of each empirical hypergraph.We find that the fixation probability for the obtained weighted network is almost the same as that for the Moran process for the corporate (Fig 5A  that for the Moran process.In contrast, the original Enron hypergraph is a much stronger suppressor of selection.Therefore, our result that the four empirical hypergraphs are suppressors of selection is not expected from the one-mode projection. To further examine the result that the empirical hypergraphs are suppressors of selection, we now simulate the same evolutionary dynamics on the randomized hypergraph.We obtained the randomized hypergraph for each empirical hypergraph by randomly shuffling the hyperedges of the original hypergraph.In the random shuffling, we preserved the degree of each node and the size of each hyperedge [51,52].We show the fixation probability for the randomized hypergraphs by the green circles in Fig 5 .We find that the randomized hypergraph is also a suppressor of selection for all the four hypergraphs.The randomized hypergraph is less suppressing than the original hypergraph for three empirical hypergraphs (see Fig 5A -5C) and vice versa for the other one empirical hypergraph (see Fig 5D).Therefore, how the features of empirical hypergraphs except for the distribution of the node's degree and hyperedge's size affects the fixation probability is non-equivocal.However, these results for the randomized hypergraphs further strengthen our main claim that the hypergraphs that we have investigated are suppressors of selection under model 1.

Discussion
We have proposed two models of evolutionary dynamics on hypergraphs that are extensions of the birth-death process on networks.For both models of evolutionary dynamics, we semianalytically derived the fixation probability for the mutant type under the constant selection for three synthetic hypergraphs with high symmetry, which generalize the complete graph, cycle graph, and star graph.For model 1, which is appropriate for discussing the strength of natural selection, we showed that these synthetic hypergraphs are suppressors of selection.Furthermore, by numerical simulations, we have shown that four empirical hypergraphs of our arbitrary choices are also suppressors of selection.Our results are in stark contrast to the known result that most networks are amplifiers of selection under the birth-death updating rule and the uniform initialization [20][21][22], which we also assumed.It is often the case that interaction among individuals is often of higher order, and hypergraphs rather than conventional networks are more direct representation of many empirical data of social and ecological interactions [13,14].Therefore, amplification of natural selection by birth-death processes on networks may not be as universal as it has been suggested once we expand the class of network models to be considered.
For conventional networks, finding suppressors of selection under the birth-death process is difficult [28].In contrast, under the death-birth process, most conventional networks are suppressors of selection [20], and amplification of selection is bounded and transient [53].Our main result that all the hypergraphs that we have investigated are suppressors of selection under model 1 begs various related research questions.Is there a theoretical bound on amplification of selection in birth-death processes on hypergraphs?Is there a systematically constructed class of amplifying hypergraphs even if they are rare?Are hypergraphs suppressors of selection under the death-birth process?If it is the case, hypergraphs are even more suppressing under the death-birth than birth-death processes?Can we find an optimal amplifier [24,26,27] or suppressor of selection for hypergraphs?We save these questions for the future.
As has been done in other studies, we used the Moran process as the baseline for judging the amplifying and suppressing effects of selection.An alternative is to use a hypergraph equivalent of the complete graph as the baseline.The complete k-uniform hypergraph is a candidate.).Therefore, the cyclic 3-uniform hypergraph with N 2 {5, 20, 200} is not an isothermal graph relative to the complete 3-uniform hypergraph, whereas the conventional circle graph is an isothermal graph relative to the Moran process.The result may change if we use a different baseline hypergraph such as the one that contains all the possible hyperedges whose size is at most k or the complete k-uniform hypergraph with a different k value.We may want to use a baseline satisfying that cyclic uniform hypergraphs are isothermal graphs relative to the baseline.Our choice of the Moran process as the baseline is to avoid these difficulties.How to determine the baseline hypergraph for discussing amplification and suppression of selection is an open question.
Metapopulation models assume networks in which a node is a patch containing individuals.Individuals either travel from a patch to an adjacent one to another according to a mobility rule such as a simple random walk or do not move but interact with the other individuals in the same patch more frequently than those in adjacent patches.Extending the island models in which the patches are connected as the complete graph [54], constant-selection dynamics in metapopulation network models have been investigated [55].Metapopulation models and hypergraphs are apparently similar because the entire network is composed of interacting groups of individuals (i.e., patches or hyperedges) in both models.Recent studies showed that metapopulation models can be suppressors of selection [56,57].However, constant-selection dynamics in metapopulation models and hypergraphs are crucially different in the following two aspects.First, interaction between individuals is pairwise in metapopulation models.In contrast, our hypergraph models assume group interaction, allowing multiple individuals to simultaneously change their type.Second, metapopulation models assume that within-group interaction between individuals is more likely than between-group interaction.This is not the case for hypergraphs in general although one can enforce it by, for example, allowing duplicated or weighted hyperedges representing a relatively high frequency of within-group interaction.Therefore, we consider that the hypergraphs used in the present study, and potentially many others, are suppressors of selection for different reasons from those for suppressing metapopulation networks, such as the possibility of multiple individuals to be simultaneously updated in the case of the hypergraph.
Evolutionary set theory is a mathematical framework with which to analyze coevolution of the strategy (i.e., type) of the node and the membership of the node to sets (i.e., groups) [58][59][60].It assumes that each node belongs to different groups and play games with other nodes belonging to the same group.A node v with a large fitness obtained from playing the game disseminates v's group membership as well as v's type to other nodes with a high probability.Therefore, evolutionary set theory is a dynamic graph theory.Evolutionary set theory is distinct from evolutionary dynamics on hypergraphs studied in the present study in the following aspects.First, in evolutionary set theory, interaction between players is pairwise by default.In contrast, in evolutionary dynamics on hypergraphs, the outcome of an interaction in a group (i.e., hyperedge) may not be decomposed into the summation of pairwise interaction in the group.For example, we assumed that a parent node v simultaneously imposes its type on all the other nodes in the selected hyperedge.Second, the group membership evolves in evolutionary set theory.In contrast, it is fixed in the evolutionary dynamics considered in this study, which is a simplification assumption.Third, reproduction in evolutionary set theory occurs globally, not limited to between nodes in the same group.In contrast, we have assumed that reproduction occurs between nodes in the same hyperedge.These differences create opportunities for future work.For example, extending the present model to the case of dynamic hypergraphs may be interesting.For conventional networks, the time to fixation has been shown to be in a tradeoff relationship with the fixation probability [61][62][63][64][65].In other words, a strongly amplifying network tends to accompany a large mean fixation time.Our mathematical and computational framework is readily adaptable to the examination of fixation times.The fixation time for hypergraphs including the comparison with the case of conventional networks warrants future work.Other topics for further investigations include the effects of different initial conditions [7,34,42,66], mathematical analyses of other symmetric hypergraphs, weak selection expansion of fixation probability [22] in the case of hypergraphs, amplification and suppression of selection in the mutation-selection equilibrium under a small mutation rate [67], and fixation probability of cooperation in evolutionary game dynamics on hypergraphs [15,16].Focusing on fixation of a mutant (or resident) type in general allows us to use Markov chain theory to reach various mathematical insights.We believe that the present study is an important step toward deploying this powerful mathematical and computational machinery to evolutionary dynamics on higher-order networks.
We show the complete 3-uniform hypergraph on four nodes in Fig 1A as an example.

Fig 2 .
Fig 2. Fixation probability for different hypergraph models under model 1.We compare the Moran process, which is the baseline, complete 3-uniform hypergraphs, cyclic 3-uniform hypergraphs, and star 3-uniform hypergraphs.(A) N = 4. (B) N = 5. (C) N = 20.(D) N = 200.The insets in (C) and (D) magnify the results for r values smaller than and close to r = 1.In these insets and main panel (B), the results for the complete 3-uniform hypergraph are close to those for the cyclic hypergraph such that the blue lines are almost hidden behind the orange lines.In (A), the two results are exactly the same such that the blue line is completely hidden behind the orange line.https://doi.org/10.1371/journal.pcbi.1011494.g002

NÀ 1 2 À
� hyperedges.Each leaf node belongs to N − 1 hyperedges and are structurally equivalent to each other.The star 3-uniform hypergraph with N = 3 is the complete 3-uniform hypergraph.Therefore, we consider star 3-uniform hypergraphs with N � 4. We show the 3-uniform star hypergraph on 5 nodes in Fig 1C as an example.
for N = 20 and N = 200.Fig 2C and 2D compare the obtained x 1 values with those for the Moran process with N = 20 and N = 200, respectively.These figures indicate that these larger star 3-uniform hypergraphs are also suppressors of selection although the degree of suppression is smaller for larger star 3-uniform hypergraphs; Fig 2D shows that the difference between the star 3-uniform hypergraph and the Moran process in terms of x 1 is tiny when N = 200.Fig 2 suggests that larger star 3-uniform hypergraphs are weaker suppressors of selection.
Fig 4. The figure indicates that the cyclic 3-uniform hypergraph with N = 4 is a suppressor of selection.However, the cyclic 3-uniform hypergraph with N > 4 is neither amplifier nor suppressor of selection under model 2. This is because x 2 = 14/[5(N − 1)(N − 2)] when r = 1 for cyclic 3-uniform hypergraphs.This x 2 value is different from the value for the Moran process, i.e., x 2 = 2/N.

Fig 4 .
Fig 4. Fixation probability for different hypergraph models under model 2. We compare the Moran process, which is the baseline, complete 3-uniform hypergraphs, cyclic 3-uniform hypergraphs, and star 3-uniform hypergraphs.(A) N = 4. (B) N = 5. (C) N = 20.(D) N = 200.In (A), the result for the complete 3-uniform hypergraph is exactly the same as that for the cyclic 3-uniform hypergraph such that the blue line is completely hidden behind the orange line.In (C) and (D), the results for the complete 3-uniform hypergraph (shown by the blue lines) are not identical but close to those for the star hypergraph (shown by the green lines) such that the former are hidden behind the latter.https://doi.org/10.1371/journal.pcbi.1011494.g004 ), Senate (Fig 5C), and high-school (Fig 5D) hypergraphs.The one-mode projection of the Enron hypergraph is a suppressor of selection (Fig 5B).However, the fixation probability for the one-mode projection of the Enron hypergraph as a function of r is close to

Fig 5 .
Fig 5. Fixation probability for empirical hypergraphs, their one-mode projection, and the randomized hypergraphs.(A) Corporate.(B) Enron.(C) Senate.(D) High-school.The insets of (B), (C), and (D) magnify the results for r values smaller than and close to r = 1.We estimated the fixation probability at each value of r as the fraction of runs in which fixation of type A is reached.We find that empirical and randomized hypergraphs are suppressors of selection for all the four data sets.https://doi.org/10.1371/journal.pcbi.1011494.g005