### Multilocus model and computer simulations

Figure 1 illustrates the 'mating pool' mode of reproduction on which our model is based. A haploid finite population of size *N* = *m* × *n* is randomly subdivided into *m* groups with *n* individuals per group. Offspring genotypes at each generation are sampled from a common population and randomly subdivided into temporary groups where individuals representing a finite sample from the pooled distribution reproduce proportional to their fitness. Productivity enhancing (offspring number) genes when in genetically diverse groups coded by *g* = 1, ⋯, *G* loci were stored as lists on a single chromosome. The haploid population was initially polymorphic for *i* = 3 alleles at equal frequencies. With random grouping the distributions of the different compositions of groups for each of the *i* alleles will be binomial with parameters (*p*_{
i
}, *n*), where *n* is the group size. Assuming that allele 1 is the fittest one, within-group selection at the gth locus was modelled as

\begin{array}{l}\text{Homogeneousgroupsat}\\ \text{locus}g\{\begin{array}{l}{W}_{1}^{g}=1\\ {W}_{j\left(j\ne 1\right)}^{g}={f}_{s}\text{;}{f}_{s}\le 1\text{;}\end{array}\\ \text{Heterogeneousgroupsat}\\ \text{locus}g\{\begin{array}{l}{W}_{1\cup j}^{g}=1+y\left({a}^{z}-1\right)\\ {W}_{j\cup k}^{g}={f}_{d}\left[1+y\left({a}^{z}-1\right)\right]\text{;}\forall k\ne j,{f}_{d}\le 1\end{array}\end{array}

(3)

where *y* is a constant for the proportional increase of group's fitness as a function of allele diversity *a* (i.e. number of alleles in a group); and *z* models the benefits of allele diversity upon fitness. Therefore, without genetic diversity increasing group's fitness directional selection drives allele 1 to fixation whenever *f*_{
s
}or *f*_{
d
}<1. We set *y* = 0.5, *z* = 1, and *f*_{
s
}= *f*_{
d
}= 0.95. Allelic diversity can easily persist in this model under a wide range of parameter values [16, 21], and there are up to 3^{G}possible segregating haplotypes.

After an initial period of 500 generations to allow the population to reach a quasi-equilibrium state under selection and free recombination between adjacent loci, the focal strongly altruistic locus was introduced in the centre of the set of those whole-group beneficial loci with alleles sampled from a binomial distribution with *p* = 0.05 and *q* = 0.95 for the *A*- and *S*-types, respectively. Selection at the focal locus was modelled as in Hamilton's [1] original model

\begin{array}{l}{W}_{A}=1-k+\left(v-1\right)\frac{K}{n-1}\text{;}\\ {W}_{S}=1+v\frac{K}{n-1},\end{array}

(4)

where *v* is the number of altruistic members in a group of size *n*, and *k* is the units of fitness given by an *A*-type in order to add *K* units to the joint fitness of his (*n* -1) companions. The degree of positive assortment (*F*) required for altruism to evolve is [1]

\frac{K}{k}=\frac{1}{F}.

(5)

In other words, the critical *c/b* ratio for positive selection of altruism is *k*/*K*.

With *N* = 240 and n = 4 as assumed in most simulations (but see below) only 60 groups are formed, and the hypergeometric probability (i.e., sampling is without replacement) of randomly forming homogeneous groups of altruists with initial *p* = 0.05 is 3.67 × 10^{-6}. Random group formation events will result in relatedness (*r*^{o}) values which fluctuate around – 1/(*N-* 1) for the altruistic locus. With *k* = 0.1 and *K* = 0.2 as assumed here strong altruism is expected to be lost a few generations after its introduction into the population since homogeneous groups of altruists will be surely absent.

With multiple loci the probability of allele fixation increases in the whole-group beneficial loci [21], and mutation was introduced whenever a locus *g* was fixed for one allele. This is, however, a relatively minor problem with less than *G* = 40 loci or so [21], and the only point to introducing mutation was to keep constant the 'effective' number of segregating *G* loci in the population.

We simulated repeated introductions of the *A* allele into the population whenever the *S*-type reached fixation. The lists of loci could undergo a recombination process following a stochastic multilocus method [22] with recombination frequency between adjacent loci denoted by *R* (no interference was assumed). In each generation, the sequence of events was mating pool formation, mutation, selection, recombination, and subdivision. We assumed multiplicativity of selective effects. Roulette selection operator, in which the chance of a chromosome getting selected is proportional to its fitness, was used. Selected chromosomes were randomly paired for recombination and randomly assigned to a group. Up to seven computer-simulation trials with different random number seeds were run for each combination of parameter values, and each trial was run for up to 50,000 generations. By following the fate of the *A*-type variant introduced repeatedly into the population, until it is fixed or lost, we can estimate the mean times to fixation or loss, respectively. The simulations were implemented in MATLAB version 7.2 [23]. An analytical treatment of the model in a simplified situation with group size *n* = 2 and no recombination is given after the numerical results.

### Numerical results

In order to numerically demonstrate how strong altruism can evolve in single-generation, randomly formed groups, consider a focal locus with two alleles (*A*: altruistic; *S*: selfish types) embedded in a multilocus genome containing *g* = 1, ⋯, *G* selected loci stored as lists on a single chromosome with whole-group beneficial traits in quasi-equilibrium state under selection (see above). Group members interact for one generation and because of behaviours individuals perform during such single-generation associations, group members affect each other fitness and may also enjoy certain individual-level benefits expressed by groups composed of genetically different individuals before the population is pooled back again and randomly forms new groups in the next generation. We have chosen this model primarily because it explicitly incorporates the key characteristics of the Hill-Robertson effect – random sampling and selection at multiple (independent) loci. Here we explore the extent to which recombination influences the evolution of an altruistic allele surrounded by other loci subject to selection.

The population size was fixed to *N* = 240 individuals. This introduces high stochasticity (but see below for the effect of both population size and group size), but parameters in the *c/b* ratio were chosen as to make it very unlikely that the *A*-type could be fixed by genetic drift in the single locus case (ultimate loss of the altruistic allele was constantly observed after > 2.5 × 10^{5} introductions into the population with *n* = 4). An initial survey of a range of number of loci and values for the proportional increase of group's fitness as a function of allele diversity at locus *g* was made, and here we present summary results with group size *n* = 4; *G* = 16, 32, 48, 64 loci. As expected, the numerical outcomes show that the mean persistence time of the altruistic *A*-type is influenced by recombination between adjacent linked loci (Figure 2a). When embedded in a chromosome, or chromosome fragment, with tight linkage *R* ≤ 0.001 (notice that an *R* value of 10^{-3} between adjacent loci would correspond to a region approximately (*G* + 1)× 100 × 10^{-3} centimorgans in length and roughly equivalent to a chromosome with 1,000 genes and total map distance 100 centimorgans), the altruistic allele could easily persist in the population for quite a long time because the genome eventually crystallized in a few segregating haplotypes. This effect was already clear with only *G* = 8 loci, a situation that allows tracking the haplotype dynamics and the eventual trapping of the *A*-type by complementary segregating haplotypes (Figure 2b).

The results can be understood as follows. Deterministic haplotype dynamics assuming linkage equilibrium at generation *t*_{0} guarantees that any *A*-type mutant introduced into the population will be damned to eventual loss no matter how many *G* loci are segregating in the population (see analytical treatment below). Variance in linkage disequilibrium is generated directly by random drift, and the *A*-type can be kept in the population when statistically associated (hitchhike) with haplotypes that remain segregating at quasi-equilibrium state. Stochasticity in this model is obviously a function of haplotype space that a given population size can explore. Since there are up to 3^{G}possible segregating haplotypes for the *G* loci, no real population will eventually behave according to the deterministic dynamics and hitchhiking effects related to the action of social heterosis may prevent the altruistic allele from going quickly extinct. With low recombination neighbouring genes tend to be inherited together, and those three segregating haplotypes in the sampled simulation (Figure 2b) can be thought of as alleles of a whole-group effect supergene, with all group members experiencing a net benefit from genetic diversity. Relatedness in this case is [13] *r*^{w}= 1/*n* (superscript ^{w}stands for whole-group relatedness) and is obviously positive. This does not, however, invalidate our claim that strong altruism can invade without kin selection because relatedness at the focal locus still remains *r*^{o}= -1/(*N*-1), so there is no kin selection effect involved. The role of group size *n* can also be easily visualized. A large group size will tend to have low across-group variance in genetic diversity, which eventually drives haplotype [1211*S*2112] in Figure 2b to fixation since we have assumed that at each *g*th locus allele 1 is the fittest one (eq. 3).

Recombination breaks down the statistical associations between alleles at linked sites (i.e. detaches *r*^{o}from *r*^{w}) but, as we might expect, the hitchhiking effect is most marked when there are many segregating loci. With intermediate linkage (*R* = 0.005) fixation of strong altruism was a sporadic outcome for *G* = 32 and *G* = 48 loci, and happened with a relatively high probability for *G* = 64 loci (Figure 3). In these cases it was computationally unfeasible to keep track of all segregating haplotypes in the population, but the results clearly suggest that the eventually successful *A*-allele was locked around a haplotype region that went to fixation or, alternatively, to various regions that happened to reach a high enough frequency (see below). Finally, swift loss of the *A*-type was always the final outcome with *R* ≥ 0.10.

### Analytical dynamics of whole-group genes

In order to understand the genetic hitchhiking effect responsible for the spread of strong altruism in our group selection model, it is important to first understand the dynamics of whole-group genes since linkage disequilibria can also be generated by deterministic forces. To clarify matters, we deal here with the simplest situation where the *g* = 1, ⋯, *G* loci with whole-group beneficial traits have two alleles each (e.g., 1^{1} and 2^{1}; 1^{2} and 2^{2} for a 2-locus genome). For each pair, allele 1 is the fittest one as assumed in Methods. This creates up to 2^{G}possible haplotypes. We set group size to *n* = 2 and assume no recombination (*R* = 0). It is then possible to easily follow the dynamics of all segregating haplotypes in the population through time. Assuming that the benefits of allele diversity upon fitness are linear (*z* = 1 in eq. 3), the payoff matrix for the *g*th locus is

{P}^{g}\text{:}\left|\begin{array}{cc}{1}^{g}& {2}^{g}\\ {1}^{g}& {w}_{11}=1& {w}_{12}=\left(1+y\right)\\ {2}^{g}& {w}_{21}={f}_{d}\left(1+y\right)& {w}_{22}={f}_{s}\end{array}\right|

(6)

where element *w*_{
ij
}is the fitness of individual *i* (1^{g}, 2^{g) }in a group with individual *j*. With *G* loci and multiplicativity of selective effects the payoff matrix is simply

**Γ** = **P**^{1} ⊗ **P**^{2} ⊗ ⋯ ⊗**P**^{G}

where ⊗ is the Kronecker tensor product. For *G* = 2 this gives the matrix

\begin{array}{cc}\Gamma \text{:}& \left|\begin{array}{cccc}{1}^{1}{1}^{2}& {1}^{1}{2}^{2}& {2}^{1}{1}^{2}& {2}^{1}{2}^{2}\\ {1}^{1}{1}^{2}& 1& \left(1+y\right)& \left(1+y\right)& {\left(1+y\right)}^{2}\\ {1}^{1}{2}^{2}& {f}_{d}\left(1+y\right)& {f}_{s}& {f}_{d}{\left(1+y\right)}^{2}& {f}_{s}\left(1+y\right)\\ {2}^{1}{1}^{2}& {f}_{d}\left(1+y\right)& {f}_{d}{\left(1+y\right)}^{2}& {f}_{s}& {f}_{s}\left(1+y\right)\\ {2}^{1}{2}^{2}& {f}_{d}^{2}{\left(1+y\right)}^{2}& {f}_{s}{f}_{d}\left(1+y\right)& {f}_{s}{f}_{d}\left(1+y\right)& {f}_{s}^{2}\end{array}\right|\end{array}

(8)

It is clear that group mean is

\Psi =\frac{1}{2}\left(\Gamma +{\Gamma}^{T}\right)

(9)

where superscript ^{T}signifies matrix transposition. Now we have to generate the appropriate matrix of haplotypes, which can be done as follows. At generation *t*_{0} segregating alleles have frequencies {p}_{0}^{g} and 1 - {p}_{0}^{g} for locus *g*, and linkage disequilibrium is absent (*D* = 0). Assume for simplicity *G* = 2 loci; then the 4 × 4 haplotype matrix is obtained in two steps. First, we multiply the allele frequency vectors to obtain the haplotype frequencies

h=\left[{p}_{0}^{1}\left(1-{p}_{0}^{1}\right)\right]\otimes \left[{p}_{0}^{2}\left(1-{p}_{0}^{2}\right)\right],

(10)

and second, we multiply the resulting **h** vector

**v** = **h** ⊗ **h**

This row vector can now be rearranged by noting that its first 2^{G}elements are the corresponding random group frequencies of haplotype 1^{1}1^{2} in a group with haplotype 1^{1}1^{2}, ⋯, 2^{1}2^{2}; the second 2^{G}elements the corresponding frequencies of haplotype1^{1}2^{2} in a group with haplotype 1^{1}1^{2}, ⋯, 2^{1}2^{2}; etc. After rearranging this vector we obtain the suitable matrix of haplotypes (**H**). Since *R* = 0, the recurrence relations for the haplotype frequencies are

\overline{w}{h}_{i}={\displaystyle \sum _{i}{H}_{i}\odot {\Gamma}_{i}\odot {\Psi}_{i}},

(12)

where

\overline{w}={\displaystyle \sum _{i,j}{H}_{i,j}\odot {\Gamma}_{i,j}\odot {\Psi}_{i,j}}

is the average fitness, *h*_{
i
}is the corresponding row for the *i*th haplotype (individual), and ⊙ denotes element-by-element multiplication of the corresponding *i*, *j* elements in the three matrices (also known as a Hadamard product). Letting *f*_{
s
}= *f*_{
d
}= 0.95 in payoff matrix (6), Figure 4 shows the haplotype dynamics for different values of the proportional increase of group's fitness *y* as a function of allele diversity.

It is clear that with the assumed starting conditions the deterministic short-term evolution drives the population to end up with only coupling haplotypes at equilibrium state. We can also define ESS (Evolutionary Stable Strategy) conditions for the payoff matrix (8) by noting that the game theoretical model is isomorphic to a trait group model with group size *n* = 2. Inspection of (8), by application of standard game theory [24] and dynamics [25], reveals that there are only two ESSs, corresponding to complete coupling {1^{1}1^{2}, 2^{1}2^{2}} and repulsion {1^{1}2^{2}, 2^{1}1^{2}} haplotypes. All other equilibria are unstable. This holds if *f*_{
s
}, *f*_{
d
}, *y* < 1; *y* > (1- *f*_{
d
})/*f*_{
d
}. The domains of attraction of the two ESSs are not equal, hence from initial linkage equilibrium the system converges to the coupling ESS (Figure 4). It is straightforward to show that by increasing the number of loci by one, the number of mixed ESSs doubles, i.e. the latter scales with 2^{G-1 }if the number of alleles per locus is 2. These considerations clearly suggest that there is a balance between recombination (always decreases linkage disequilibrium) and selection (convergence to coupling or repulsion haplotypes) in the system.

The preceding analysis can explain some previous numerical results [26]. Computer trials using bi-allelic loci with all haplotypes initially at equal frequencies and finite population size *N* = 120, with group size set to *n* = 2 and parameter values *f*_{
s
}= *f*_{
d
}= 0.95 and *y* = 0.2, had found that coupling haplotypes were the most common ones left in the population assuming no recombination, but in 42.8% of the trials with *G* = 2 loci, and 56.5% with *G* = 3 loci, the fittest haplotype was lost from the population. There was, however, strong selection for allelic diversity at each locus to be maintained. This was most evident when repulsion haplotypes were the only ones left, as they remained segregating in the population at equilibrium state. Finally, it should be stated here that by letting *f*_{
d
}= 1 in payoff matrix (6) the system converges to the standard heterotic case [17, 18, 27].

### Fate of the altruistic allele in the analytical dynamics

It is straightforward to incorporate the altruistic locus in the former analytical treatment for the whole-group genes. The payoff matrix is now (eq. 4)

{P}^{a}\text{:}\left|\begin{array}{cc}S& A\\ S& 1& \left(1+K\right)\\ A& \left(1-k\right)& \left(1-k+K\right)\end{array}\right|

(13)

where *S* and *A* stand for the selfish and (strongly) altruistic types, respectively; and *k/K* is the critical ratio for positive selection of altruism. With multiplicativity of selective effects and *G* whole-group loci the multilocus payoff matrix is

**Γ** = **P**^{a}⊗ **P**^{1} ⊗ ⋯ ⊗ **P**^{G}

which, with *G* = 1 becomes

\begin{array}{cc}\Gamma \text{:}& |\begin{array}{cc}S{1}^{1}& S{2}^{1}\\ S{1}^{1}& 1& \left(1+y\right)\\ S{2}^{1}& {f}_{d}\left(1+y\right)& {f}_{s}\\ A{1}^{1}& \left(1-k\right)& \left(1-k\right)\left(1+y\right)\\ A{2}^{1}& \left(1-k\right){f}_{d}\left(1+y\right)& {f}_{s}\left(1-k\right)\end{array}\cdots \\ \begin{array}{cc}A{1}^{1}& A{2}^{1}\\ S{1}^{1}& \left(1+K\right)& \left(1+K\right)\left(1+y\right)\\ S{2}^{1}& \left(1+K\right){f}_{d}\left(1+y\right)& {f}_{s}\left(1+K\right)\\ A{1}^{1}& \left(1-k+K\right)& \left(1-k+K\right)\left(1+y\right)\\ A{2}^{1}& \left(1-k+K\right){f}_{d}\left(1+y\right)& {f}_{s}\left(1-k+K\right)\end{array}|\end{array}

(15)

It is easy to see that {S1^{1}, S2^{1}} is the only (mixed) ESS, thus selfishness prevails. With *G* ≥ 2 the problem is also straightforward. If the altruistic allele is fixed in the population, but otherwise all haplotypes are present, then the metastable equilibria exactly mirror those of the system without the altruistic locus, i.e. for *G* = 2 loci the equilibria are the same as for matrix (8). Suppose now that the system is in the equilibrium {*A*1^{1}1^{2}, *A*2^{1}2^{2}}. Which haplotypes can invade the population? Inspection of the **Γ** matrix (15) reveals that only haplotypes *S*1^{1}1^{2} and *S*2^{1}2^{2} can invade, but *S*1^{1}2^{2} and S2^{1}1^{2} cannot, *even though the latter also carry the selfish allele*, provided *f*_{
d
}>(1 + *K*)/[1- *k* + *K*)(1 + *y*)]. This means that the pair {*A*1^{1}1^{2}, A2^{1}2^{2}} is unstable *in the direction of the same heterotic haplotypes only*, carrying the selfish allele. This generalizes to more loci in an important way. As the number of loci increases, it remains true that the resident, altruistic pair of haplotypes can be invaded by two appropriate haplotypes only; all the other selfish haplotypes are repelled.

Deterministic haplotype dynamics assuming linkage equilibrium at generation *t*_{0} also shows that any *A*-type mutant introduced into the population is damned to eventual loss no matter how many *G* loci are segregating in the population. The only possibility for the *A*-type to be kept in the population is to statistically associate itself (hitchhike) with haplotypes that remain segregating at quasi-equilibrium state.

In finite populations allelic diversity can be lost when the fittest haplotype is present since its eventual fixation can occur [26], a result which strongly depends on the number of *G* loci initially introduced in the population and the strength of selection [21]. This can be understood by comparing Figures 4b and 4c, where the equilibrium frequency of the fittest haplotype is higher as the number of loci increase. Because in the model for whole-group genes there is also strong selection to maintain allelic diversity it is, therefore, not surprising that a mutant *A*-type could eventually hitchhike a haplotype that is complementary to those already present in the population. The interaction between linkage, selection, and sampling in structured-deme models is precisely the phenomenon that can be described by the 'Hill-Robertson' effect [19] or genetic hitchhiking [20].

We can now envisage the effects of group size *n* and population size *N* on the eventual fate of allele *A* under a given set of conditions. Increasing group size decreases the across-group variance in allelic diversity at the *G* loci, which in turns increases the within-group selection and reduces the likelihood of polymorphism since the fittest allele *1*^{g}at the *g*th locus goes to fixation (eqs. 3 and 6). Hence, hitchhiking may not be possible and the eventual dynamics of the *A*-type converges to the single locus case. In other words, the initial spread of altruism in this model is conditional on a large enough among-group genetic variance and will not happen in genetically homogeneous settings. Simulations with *N* = 240, *G* ≥ 16, and parameter values as above indicate that group size *n* ≥ 16 is large enough as to dramatically reduce or even prevent the eventual trapping of the *A*-type by complementary segregating haplotypes when *R* = 0.

The effect of population size *N* is relatively more complex since stochasticity in this model is obviously a function of the strength of selection at each locus and the haplotype space that a given population can explore. For instance, in simulations with *G* = 8 (haplotype space 2 × 3^{8} = 13, 122; including the focal altruistic locus), *n* = 4, *R* = 0 and strong selection quasi-stable polymorphism of the *A*-type was not detected when *N* ≥ 500 since haplotype evolution tends to approach the deterministic situation; but was a frequent outcome with *G* ≥ 32 loci even when *N* = 2, 000. On the other hand, reducing the strength of selection by 60% (i.e., *y* = 0.2, *k* = 0.04 and *K* = 0.08) allowed to increase population size to *N* = 500 as to qualitatively obtain similar results for *G* = 8 as those obtained with strong selection and *N* = 240. In summary, the hitchhiking effect is robust for large population sizes as long as group size *n* is small enough, so that among-group allelic diversity at the *G* loci can be selectively kept.

In the simulations with *G* = 64 loci and *R* = 0.005 we detected a relatively high probability of fixation of strong altruism (5.78 × 10^{-4}). Figure 5 also shows a sample simulation where a snapshot of segregating haplotypes right after fixation of the *A*-type was obtained. As expected from the analytical treatment above, it is clear that the successful *A*-allele was locked around a complementary region that reached a high enough frequency because it allowed group members to experience a net benefit from genetic diversity at the *G* loci.