The effects of linkage and gene flow on local adaptation: a two-locus continent-island model.

Population subdivision and migration are generally considered to be important causes of linkage disequilibrium (LD). We explore the combined effects of recombination and gene flow on the amount of LD, the maintenance of polymorphism, and the degree of local adaptation in a subdivided population by analyzing a diploid, deterministic continent-island model with genic selection on two linked loci (i.e., no dominance or epistasis). For this simple model, we characterize explicitly all possible equilibrium configurations. Simple and intuitive approximations for many quantities of interest are obtained in limiting cases, such as weak migration, weak selection, weak or strong recombination. For instance, we derive explicit expressions for the measures D(=p(AB)-p(A)p(B)) and r(2) (the squared correlation in allelic state) of LD. They depend in qualitatively different ways on the migration rate. Remarkably high values of r(2) are maintained between weakly linked loci, especially if gene flow is low. We determine how the maximum amount of gene flow that admits preservation of the locally adapted haplotype, hence of polymorphism at both loci, depends on recombination rate and selection coefficients. We also investigate the evolution of differentiation by examining the invasion of beneficial mutants of small effect that are linked to an already present, locally adapted allele. Mutants of much smaller effect can invade successfully than predicted by naive single-locus theory provided they are at least weakly linked. Finally, the influence of linkage on the degree of local adaptation, the migration load, and the effective migration rate at a neutral locus is explored. We discuss possible consequences for the evolution of genetic architecture, in particular, for the emergence of clusters of tightly linked, slightly beneficial mutations and the evolution of recombination and chromosome inversions.

Introduction

Evolution in geographically structured populations is governed by two potentially conflicting forces: selection for improved local adaptation may be opposed by migration and concomitant maladaptive gene flow. For one-locus models, both general theory as well as the study of numerous particular models have provided considerable insight into the evolutionary consequences of this interaction (reviewed by Karlin, 1982; Lenormand, 2002; Nagylaki and Lou, 2008).

Frequently, however, selection acts on multiple loci which may be linked. For this situation much less is known. The recently developed multilocus theory for the evolution in subdivided populations focuses mainly on limiting or special cases, such as weak or strong migration (Bürger, 2009a,b) or the Levene model (Nagylaki, 2009; Bürger, 2009c, 2010), in which LD is weak or absent. The applications treated there concern primarily conditions for the maintenance of multilocus polymorphism.

The main goal of this paper is to explore and quantify the role of selection on linked loci for local adaptation in the presence of maladaptive gene flow. To this aim, we study a diploid continent–island model with genic selection (i.e., no dominance or epistasis) acting on two recombining loci. Continent–island models seem to be the simplest scenario to investigate this topic, and they have been justified by and applied in studies of local adaptation in natural populations (e.g., King and Lawson, 1995). A key role in our analysis will be played by the LD generated by the interaction of selection and migration.

Population subdivision and admixture are known to be important agents in generating LD (Charlesworth and Charlesworth, 2010), and recent studies revealed high levels of LD in geographically structured species as well as substantial differences among species and subpopulations (e.g., Remington et al., 2001; Conrad et al., 2006; Cutter et al., 2006; Hernandez et al., 2007). However, analytical estimates for the expected amount of LD have been obtained mainly for neutral loci. In particular, the expectation and variance of measures of pairwise LD, such as or (the squared correlation in allelic state), have been derived for various island models (Ohta, 1982; McVean, 2002; Wakeley and Lessard, 2003; De and Durrett, 2007). Such estimates are particularly useful in testing hypotheses about neutral evolution (Nordborg and Tavaré, 2002; Slatkin, 2008).

For large subdivided populations in which selection and migration are the dominating forces, little quantitative information is available on the extent of LD. Li and Nei (1974) and Christiansen and Feldman (1975) investigated a two-deme model with symmetric migration and genic selection on two dialellic loci. They showed that for weak migration an equilibrium exhibiting LD exists. Li and Nei also provided numerical examples showing that considerable LD can be maintained. Slatkin (1975) studied gene flow and selection on two linked loci in a cline. He found that substantial LD is maintained if the recombination rate is of the same order of magnitude as the selection coefficients or smaller. In addition to studying the influence of linkage among multiple loci on the width of a cline, Barton (1983) derived a weak-migration approximation for the allele frequencies and the linkage disequilibria in a two-deme model with heterozygote inferiority at finitely many equivalent loci. Spichtig and Kawecki (2004) observed high levels of LD in a model of antagonistic directional selection on a quantitative trait in two demes. However, in general, their fitness functions induce epistasis and dominance, and the dependence of LD on the degree of epistasis was not reported.

For our model, we can quantify explicitly the role of linkage in generating LD and in maintaining genetic polymorphism and local adaptation in the presence of gene flow. On the island, locally adapted alleles and may occur at two linked loci. Evolution occurs in continuous time and we ignore dominance and epistasis. Thus, the model also describes selection on haploids. All immigrants from the continent carry the (maladapted) haplotype . Mutation and random drift are ignored. Our parameters are the selection coefficients at the two loci, the recombination rate, and the migration rate. Section 2 contains a detailed description of the model.

In Section 3, the complete equilibrium and stability structure is derived. Below a critical migration rate, which can be determined explicitly and depends strongly on the recombination rate and the selection coefficients, there is a unique, stable, fully polymorphic equilibrium. If recombination is approximately as strong as selection, an additional, unstable fully polymorphic equilibrium can exist. In this case, the stable polymorphic equilibrium is simultaneously stable with a boundary equilibrium.

In Section 4, we consider various applications. (i) We derive simple and intuitive approximations for the stable, fully polymorphic equilibrium for several limiting cases. (ii) We obtain explicit conditions for the maintenance of both locally adapted alleles in the presence of gene flow. These reveal the role of linkage, and illustrate how it interacts with migration and selection in maintaining polymorphism. (iii) From an explicit analytical formula for the measure of LD at the stable fully polymorphic equilibrium, we derive simple and informative approximations for and in important limiting cases. Interestingly, and show qualitatively different dependence on the migration rate. In particular, may be high for a large range of parameters. (iv) Assuming that a locally adapted allele () is in migration–selection balance on the island, we study the invasion condition for and the final frequency of a linked beneficial mutation of smaller effect (). In the absence of LD, only mutations can invade whose selective advantage exceeds the immigration rate of the maladapted (continental) allele. We show that mutants having a much smaller effect than the immigration rate can invade if they are sufficiently (but not extremely) tightly linked to the already polymorphic locus. Such an invasion always entails an increase of ; the tighter the linkage, the higher this increase. (v) We derive explicit formulas for the migration load and examine how the load and the degree of local adaptation depend on the parameters of the model. (vi) We derive an approximation for the effective migration rate at a linked neutral marker locus.

In Section 5, we summarize and discuss these results. In particular, we elaborate on their consequences for the evolution of the genetic architecture, such as linkage groups or chromosome inversions, in the presence of gene flow.

The model

We extend the classical continent–island model (Haldane, 1930; Nagylaki, 1992) to two loci. We consider a sexually reproducing population of monoecious diploid individuals in which two diallelic loci are under selection. On the continent, the alleles and are fixed. On the island, the alleles and are selectively favored. We call them the island alleles, and and the continental alleles. We ignore mutation and random genetic drift and employ a deterministic continuous-time model to describe evolution on the island. We assume one-way migration from the continent to the island at rate . Thus, all immigrants on the island are of haplotype . The recombination rate between the two loci is denoted by .

The population can be described in terms of , the frequencies of the four gametes on the island. The state space is the probability simplex . Its vertices correspond to the monomorphic states at which one of the gametes has frequency 1. The edges connecting the states and , or and , correspond to the marginal single-locus systems in which the island allele , or , respectively, is absent.

We assume absence of dominance and of epistasis, hence genic selection, and assign the Malthusian parameters and to the alleles and , and and to and . The resulting fitness matrix for the genotypes reads

Assuming continuous time, standard population genetics modeling (e.g. Bürger, 2000, Chap. II.1) yields the following system of differential equations for the evolution of gamete frequencies on the island: where denotes the classical measure of LD.

For our purposes, it is convenient to describe the population composition by the frequencies and of the alleles and and by the LD measure . The gamete frequencies are calculated from , and by The constraints and transform into and It follows that , and evolve according to We note that by rescaling time, for instance to units of or (provided or is positive), the number of independent parameters can be reduced by one (i.e., to or to , respectively). We refrain from doing so because in our applications it will be illuminating to use either or as an independent parameter.

We have also studied the corresponding discrete-time model which allows for strong evolutionary forces and converges to the continuous-time model in the limit of weak selection, recombination, and migration. It yields more complicated expressions for most quantities of interest and, therefore, is less amenable to an explicit analysis. Qualitatively, however, it appears to give analogous results.

Equilibrium structure and stability

To present all possible equilibrium and stability configurations, it is most instructive to describe the equilibrium and stability properties as functions of the migration rate . Thus, we use as a bifurcation parameter. We treat the boundary equilibria first, then the internal, fully polymorphic equilibria. Section 3.3 contains the main results. In Section 3.4, we prove global stability results for a number of special cases. The results derived below provide the basis for the applications treated in Section 4.

Throughout, we posit without loss of generality that selection on the -locus is at least as strong as on the -locus, i.e.,  Therefore, we call and the minor and major locus, respectively.

Boundary equilibria

Monomorphic equilibria

(i) The state at which both island alleles are fixed, is an equilibrium if and only if . In this case, it is globally asymptotically stable (in the sense that it attracts all solutions with the property that both island alleles are initially present; see Section 3.4.1). If , there also exist the three other monomorphic equilibria. They are unstable and will not be needed.

(ii) The equilibrium , at which both continental alleles are fixed, exists always: A linear stability analysis yields the eigenvalues , and . Therefore, is asymptotically stable if and only if , where Thus, if linkage is sufficiently tight (), the continental haplotype can become fixed if the critical migration rate is exceeded. It will be proved in Section 3.4.3 that is globally asymptotically stable if , independently of .

Single-locus polymorphisms

There may exist two single-locus polymorphisms.

(i) If , there exists the equilibrium at which the -locus is polymorphic and the -locus is fixed for the continental allele . As approaches . Although is globally attracting within its one-locus marginal system whenever it is admissible (Nagylaki, 1975, 1992, p. 128–132), it is always unstable with respect to the full two-locus system. The latter statement follows because the eigenvalues are Since admissibility of requires , we infer immediately that Hence, all eigenvalues are real and a simple calculation invoking shows that .

It follows that is not hyperbolic (i.e., no eigenvalue has a vanishing real part) if and only if , where This observation will be of importance below. If , we have , which is in if and only if i.e., if and only if .

(ii) If , there exists the equilibrium at which the -locus is polymorphic and the -locus is fixed for the continental allele . As approaches . If is admissible, it is globally attracting within its one-locus marginal system. Moreover, the eigenvalues are Because is admissible only if , we infer immediately that Hence, all eigenvalues are real. A simple calculation shows that this implies . Finally, we find that if and only if , where and if . We note that holds always, and is fulfilled if and only if

Thus, we have shown that is asymptotically stable if and only if and is hyperbolic except when . If , then , and if and only if (3.12) holds, i.e., if and only if .

No recombination

If , there is the following boundary equilibrium at which only the two gametes and are present: The equilibrium is admissible (and polymorphic) if and only if . (3.14) shows that if , and if . Because the eigenvalues are , and is asymptotically stable whenever it is admissible. If, initially, only the two gametes and are present, single-locus theory (Nagylaki, 1992, p. 128–132) implies convergence to . We shall prove global convergence from arbitrary initial conditions in Section 3.4.4.

Internal equilibria

We assume . Somewhat unexpectedly, the internal equilibria can be calculated explicitly. This is shown by substituting the value obtained from solving (2.5a) into the right-hand side of (2.5c), which remains linear in . Now, solving for , substituting the solution into the right-hand side of (2.5b) and solving the resulting quadratic equation in yields the following coordinates for fully polymorphic equilibria: where We call these equilibria and . The calculations show that no other internal equilibrium can exist. It is straightforward to show that converges to as and converges to as .

Our first task is to determine when the equilibria and are biologically feasible, i.e., in the state space . Obviously, the equilibria can exist only if the root is real, which is the case if and only if , where We have if and only if . The following proposition states the admissibility conditions for the fully polymorphic equilibria.

(a)  is an admissible internal equilibrium if and only if one of the following three, mutually exclusive conditions is satisfied:Therefore, the equilibrium is admissible whenever is positive and below a critical value.

(b)  is an admissible internal equilibrium if and only iforTherefore, the equilibrium is never admissible if is below or above a critical value. It is also not admissible if , i.e., if the locus effects differ by more than a factor of two.

The tedious proof is given in Appendix A.1. The cases in Proposition 1(a), as well as those in (b), are indeed mutually exclusive and cover the whole parameter range for because of the following simple observation: Therefore,

Now we briefly describe the bifurcation patterns underlying Proposition 1. If (3.17a) applies, then leaves the state space through as increases above . If (3.17c) applies, then leaves the state space through as increases above . The condition (3.17b) can be satisfied only if . In this case, coexists with . Then enters the state space through (if (3.18a) holds) or through (if (3.18b) holds). As increases further and approaches from below, and merge and get extinguished at . In Fig. 1, the possible equilibrium configurations and bifurcation patterns for are displayed as a function of .

Fig. 1

Bifurcation diagrams for the case . Diagrams (a)–(g) represent all possible equilibrium configurations, corresponding to the cases (a)–(g) in Theorem 2. Each diagram displays the possible equilibria as a function of the migration rate . Each line indicates one equilibrium. The lines are drawn such that intersections occur if and only if the corresponding equilibria bifurcate. Solid lines represent asymptotically stable equilibria, dashed lines unstable equilibria. Equilibria are shown if and only if they are admissible. (The curves in this diagrams can be produced by using the function . For optimal visibility, diagrams have been drawn for different parameter choices and on different scales.)

We note that neither nor can leave or enter the state space through . This is easy to understand because is the major locus and exists only if , hence . Therefore, whenever exists, selection on locus is strong enough to increase the frequency of near (as is reflected by the fact that the eigenvalue at is always positive). Hence, no admissible equilibrium can be arbitrarily close to .

It remains to explore the stability properties of the internal equilibria. Apparently, and also supported by numerical results, the equilibrium is asymptotically stable whenever it exists, and is always unstable. We could prove the following: (i) is globally asymptotically stable if is small (Section 3.4.2) or if is small (Section 3.4.5). (ii) is asymptotically stable if is large (Appendix A.2.2) or if is close to (and below) a bifurcation value of (Appendix A.2.1). (iii) is unstable if is close to (and above) a bifurcation value of (Appendix A.2.3). As shown by diagrams (b), (d)–(f) of Fig. 1 (which are based on these results), may be simultaneously stable with either or for sufficiently large migration rates, intermediate recombination rates, and similar locus effects. The above results support our conjecture that is globally asymptotically stable whenever it is the only fully polymorphic equilibrium. In Section 3.4.3, we will prove that is globally asymptotically stable if .

Main results

We now formulate our main theorem from which most further conclusions will be derived. It is an immediate consequence of the above results about existence and stability of equilibria.

Let . Then all possible equilibrium configurations, displayed as functions of the migration rate, are given by the schematic bifurcation diagrams (a)–(g) of Fig. 1.

Diagram (a) applies if and only iforDiagram (b) applies if and only ifDiagram (c) applies if and only ifDiagram (d) applies if and only ifDiagram (e) applies if and only ifDiagram (f) applies if and only ifDiagram (g) applies if and only ifor

As a brief guide to these results and Fig. 1, we note that, by (3.19), diagrams (a) or (c) of Fig. 1 apply if (3.17a) holds. This is the case whenever . Similarly, diagram (g) applies if (3.17c) holds, which is satisfied whenever . Thus, diagrams (a) or (g) apply whenever linkage is sufficiently tight or sufficiently weak, respectively. If and are fixed and , then for increasing the diagrams (a), (c), and (g) apply in this order. If , then for increasing diagrams apply in the order (a), (b), (d)–(g). Diagrams (b), (d)–(f) are covered by (3.17b). We distinguished cases (e) and (f) because can be simultaneously stable with or in (e), but only with in (f).

The following arguments offer an explanation why qualitatively different equilibrium configurations can occur for and . At , the order of fitnesses of diploid genotypes changes. If (i.e., the selection coefficients are very different), the three genotypes of highest fitness are those composed by the two gametes and . If (i.e., similar locus effects), the three genotypes of highest fitness are , and . Hence, experiences a higher (relative) fitness advantage than in the other case and, by symmetry, a higher fitness disadvantage. Recombination between the island and the continental haplotype always produces both gametes and . If , moderately strong recombination can maintain the fully polymorphic equilibrium even if one of the boundary equilibria or is already (locally) stable. In addition, if is given, the overall strength of selection increases with increasing , thus facilitating the maintenance of polymorphism in the presence of gene flow.

Since the case of equivalent loci is degenerate, we formulate it separately.

Let . Then all possible equilibrium configurations, displayed as functions of the migration rate, are given by the schematic bifurcation diagrams (a)–(d) of Fig. 2.

Fig. 2

Bifurcation diagrams for the case . Diagrams (a)–(d) represent all possible equilibrium configurations, corresponding to the cases (a)–(d) listed in Theorem 3. Note that in this case we have . Everything else is as in Fig. 1.

Diagram (a) applies if and only ifDiagram (b) applies if and only ifDiagram (c) applies if and only ifDiagram (d) applies if and only if

Because holds in Theorem 3, and bifurcate through simultaneously. Hence, coexistence of and occurs if and only if and . Using (A.13) and (3.18), it is straightforward to show that the interval (of values ) of coexistence of and is maximized if . Then, .

Global convergence results

Here, we prove global convergence results for some limiting cases.

We show that is globally asymptotically stable. For additive fitnesses, mean fitness increases strictly monotonically in time, except at equilibria (Ewens, 1969). In our special case, mean (Malthusian) fitness on the island is given by which is maximized at . Employing the well-known inequality , we obtain where holds only at equilibrium. Because the only equilibria are the four monomorphic states, and is the only stable equilibrium, all trajectories such that initially both alleles are present converge to (Theorem 6.4 and Corollary 6.5 in LaSalle, 1976). This conclusion can also be deduced from Theorem 3 of Karlin and Feldman (1970).

Small

We prove global convergence to if is sufficiently small. The proof follows immediately from Theorem 5.4 in Bürger (2009a), which is a general global perturbation result for weak migration. This applies here because, for , all equilibria (which are exactly the four monomorphic states) are hyperbolic, is globally asymptotically stable (Section 3.4.1), and is the perturbation of if is small. We note that the results by Karlin and McGregor (1972) imply only asymptotic stability.

Large

We demonstrate global convergence to if . Because holds always, we obtain from (2.5a) and if and only if . Analogously, implies with equality if and only if . Thus, we have two global Lyapunov functions, and all trajectories must converge to , i.e., to (LaSalle, 1976).

We prove that is globally asymptotically stable whenever it exists (i.e., if ). Asymptotic stability follows immediately from the eigenvalues, which are given by , and . We prove global convergence to if initially the island haplotype is present, i.e., if . (Note that for every if .) We observe that and Thus, we infer that the -limit of every trajectory with is contained in the invariant edge (LaSalle, 1976). On this edge, single-locus theory guarantees global convergence to .

If, initially, and , then (provided ) and convergence to occurs whenever it is admissible; otherwise, all trajectories converge to . Of course, if holds initially, then all trajectories converge to (if it is admissible) or to . The simple case is left to the reader.

We note that for , the dynamics (2.2) has the global Lyapunov function This follows immediately from Eqs. (2.13)–(2.16) on p. 103 in Bürger (2000) because, if , the system (2.2) is formally equivalent to a one-locus selection-mutation model with four alleles and so-called house-of-cards mutation. To see this, set all mutation rates to types 1, 2, and 3 (gametes ) zero, and assume that each of types 1, 2, and 3 mutates to type 4 () at rate . Therefore (p. 103 in Bürger, 2000), (2.2) is a generalized gradient system.

We prove global convergence of all trajectories to if is sufficiently small and . The crucial point is that the results in Section 3.4.4 imply that, if , the chain-recurrent points of (2.2) are exactly the equilibria (see Lemma 2.2 in Nagylaki et al., 1999). The only possible equilibria are , and , and they are hyperbolic if is small (except when or ). As a consequence, the proof of Theorem 2.3 in Nagylaki et al. (1999) applies unaltered and yields global convergence to for sufficiently small . Asymptotic (local) stability of can also be inferred from Theorem 4.4 of Karlin and McGregor (1972).

Applications

We now treat a number of applications. We begin by deriving simple approximations for the stable fully polymorphic equilibrium in several limiting cases. Then we discuss the potential of maintaining both loci polymorphic on the island. Next, we explore the LD that is maintained by calculating and approximating the measures and . We continue by studying the invasion of a new beneficial mutation of small effect which is linked to a locus that is maintained polymorphic on the island by migration–selection balance. Then we investigate the degree of local adaptation and the migration load. Finally, we compute an approximation for the effective migration rate at a neutral locus that is linked to both selected loci.

Limiting cases

The following limiting cases are highly instructive and will be used in the applications treated below. In each case, we give the coordinates of the stable internal equilibrium . They are readily deduced from Eqs. (3.15) by using a formula manipulation program, such as Mathematica.

Weak migration

We assume that migration is much weaker than selection and recombination, i.e., . Then is given by where the second-order terms will be needed to derive the series expansion of the LD measure . We note that for equivalent loci (), the first-order terms can be obtained from Barton’s (1983) weak-migration approximation (8) for multiple equivalent loci.

Tight linkage

If we assume that recombination is much weaker than selection and migration, i.e., , then is given by It is a straightforward exercise to show that this is admissible if .

Strong recombination

If recombination is strong relative to selection and migration, i.e., , then the internal equilibrium is given by To leading order in , this equilibrium is admissible if In the limit , this condition simplifies to . Obviously, as or . (4.3) may be called the quasi-linkage equilibrium approximation.

Maintenance of polymorphism

The above results admit some worthwhile conclusions about the maintenance of polymorphism at two linked loci. It is well known that at an isolated single locus, i.e., if LD is ignored, a polymorphism is maintained on the island if and only if , where the fitnesses of the three genotypes are , and (Haldane, 1930; Nagylaki, 1975). Naive application of this result to our model would imply that both loci can be maintained polymorphic only if . However, as we show below, the effects of linkage, even if weak, cannot be ignored. Our analysis is simplified by the facts that (i) is the only stable equilibrium at which the island haplotype, hence a polymorphism at both loci, is maintained, and (ii) is (apparently) asymptotically stable whenever it exists.

If (3.27) holds (Fig. 1(g)), which is always the case if recombination is at least twice as strong as selection on the minor locus, the full polymorphism is maintained if As (3.11) shows, we have if (and ), where the upper bound tends to as , and tends to as . Therefore, if recombination is not much stronger than selection on the minor locus, a fully polymorphic equilibrium can be maintained for considerably higher migration rates than predicted by single-locus theory, namely up to values of that are only slightly smaller than .

If linkage is tight, i.e., if , then a two-locus polymorphism is maintained if (Section 4.1.2) In the limit , this yields , the bound obtained from a diallelic one-locus model if the alleles are identified with the non-recombining haplotypes and .

For moderately strong linkage and if the effects of the two loci are similar, polymorphism at both loci can be maintained for migration rates higher than . This occurs if and only if diagrams (b), (d), or (e) of Fig. 1 apply (corresponding to the conditions (3.22), (3.24), or (3.25)), or if diagrams (b) or (c) of Fig. 2 apply (corresponding to conditions (3.29) or (3.30)). Then the polymorphic equilibrium is simultaneously stable with the equilibrium at which the continental haplotype is fixed. Although, in general, this range of parameters is not large, it becomes significant for loci of approximately equal effects. As noted above, this range is maximized if . Then and a full polymorphism is maintained if .

In summary, polymorphism at both loci can be maintained for higher migration rates than expected from simple one-locus considerations. Once , migration is strong enough to fix the continental haplotype independently of linkage and initial conditions (Section 3.4.3).

Linkage disequilibrium

The amount of LD maintained by migration–selection balance on the island can be computed readily because we have the explicit expression (3.15c) for the measure of LD maintained at . We also know that this is the only fully polymorphic equilibrium that can be stable. In Appendix A.3.1 it is proved that the measure is always positive at . Thus, the interaction of nonepistatic selection on the island and immigration from the continent generates positive LD.

For statistical purposes, the following measure of LD is often used (e.g., Nordborg and Tavaré, 2002; Slatkin, 2008): It can be interpreted as the squared correlation in allelic state between the two loci. A general, though complicated, expression for can be obtained directly from (3.15). We investigate and analytically for limiting cases and numerically otherwise.

We start by examining weak migration, i.e., . At the polymorphic equilibrium , we obtain from Eqs. (4.1) and, by series expansion (for which the second-order terms in Eqs. (4.1) are needed), Notably, does not vanish as . In fact, is a decreasing function of both and . It decays rapidly as a function of , and very slowly as a function of if is small. Thus, for tightly linked loci, is close to its maximum value one. Even if selection and recombination are of comparable strength, may assume relatively high values if migration is weak. This is illustrated by Fig. 3(b), in which is assumed, and by the following examples. If , the leading order term in (4.9) becomes , which is between (if ) and (if ). The leading-order term in (4.9) is small if is small, and very small if is small. If , the expression (4.8) is also obtained from Barton’s (1983) weak-migration approximation (8b); cf. Section 4.1.1.

Fig. 3

Linkage disequilibrium at the stable two-locus polymorphism. Respectively, diagrams (a) and (b) display the measures and of LD as a function of . is calculated from (3.15c) () and from and (4.7). We note that and are represented in terms of the compound parameters , and , which determine the model and all derived quantities uniquely. In both figures, is assumed and is as indicated in the legend. Diagram (b) might give the impression that LD decreases with decreasing recombination rate . This, however, is not the case (e.g., (4.9)) because the compound parameters , and are used. Hence, a decrease in implies a decrease in and .

For weak recombination, i.e., if , we obtain from Section 4.1.2 the expansion where the first-order term is independent of . Setting , we observe that . Thus, for small decreases from unity at a linear rate between and .

In the limit of large , we obtain from Section 4.1.3 that decays quadratically in : This approximation requires . Obviously, we have as .

Fig. 3 displays and as functions of for various parameter combinations.

It seems remarkable that and show such different behavior if considered as functions of . Whereas always assumes its maximum at intermediate migration rates, the squared correlation always decreases as increases. The latter fact is clearly indicated by the approximations (4.9) and (4.11).

That tends to a positive value (as opposed to zero) as follows from its definition because, as shown by (4.8), at migration–selection balance is of order and, of course, the frequencies of both continental alleles are of order . That is of order , and not smaller, follows because immigration of the continental haplotype increases the rate of change of by , whereas the contributions of selection and recombination are proportional to , hence of order only if ; see (2.5c). A more intuitive explanation is that for arbitrarily weak migration and tight linkage, the population exists predominantly of and haplotypes because the other haplotypes are produced only at the small rate . The stronger selection on is insufficient to decrease the frequency of to values below those of the recombinant haplotypes. Clearly, this leads to a high correlation in allelic state.

We note without proof that Lewontin’s (1964) measure behaves qualitatively similar to . In summary, migration may be a very efficient agent in inducing quite strong LD, even if selection is genic and weak.

Invasion of a linked beneficial mutation

In accordance with the general model, immigrants from the continent are of gametic type . Let us assume that at locus a locally adapted allele, , is segregating at migration–selection balance on the island. By single-locus theory, this requires . The corresponding equilibrium in our model is . Now suppose that on the island a favorable mutant, , occurs at the linked locus . Because our model is deterministic, the mutant can invade if and only if is unstable. By (3.13) this is the case if and only if Mutants of effect can invade for every recombination rate (which is consistent with the fact ). In addition, mutants of effect can invade for every migration rate because this implies , hence (4.12). This reflects the intuition that sufficiently strong linkage leads to reinforcement between the loci, thereby enabling invasion of .

However, condition (4.12) implies that even mutants of effect can invade. Indeed, straightforward rearrangement shows that holds if and only if Hence, is the minimum selective advantage required for invasion. Fig. 4 displays the dependence of on and .

Fig. 4

Regions of invasion of mutants of small effect. The shaded regions display the scaled fitness advantage admitting invasion of a mutant linked to a beneficial allele that is in migration–selection balance on the island. The curves are obtained from (4.13) and show as a function of for different values of . As explained in the main text, mutants of effect can invade whenever (black region).

If , the condition (4.13) for the invasion of becomes approximately where the right hand side is much smaller than . The invasion condition (4.13) can also be rewritten as where the approximation holds if . Thus, for given effect , (4.15) puts an upper bound on the recombination rate that permits invasion.

After invasion of , the asymptotically stable equilibrium is approached. Because if is unstable (Appendix A.3.2), the frequencies of both island alleles increase. This is demonstrated by Fig. 5, in which the frequencies of the alleles before and after invasion are displayed as functions of .

Fig. 5

Frequencies of island alleles after invasion of a new mutant. The solid lines show the frequencies and of the island alleles reached at the equilibrium , (3.15), after invasion of a mutant () of small effect . The frequency of decreases faster as a function of and vanishes if is such that equality is attained in (4.13). The dashed lines signify the equilibrium frequency of before invasion of . The selection coefficients are and . Qualitatively similar graphs are obtained for other choices of selection coefficients.

Not surprisingly, the largest increase occurs if the new mutant is completely linked (). Then the frequencies of both island alleles increase to because is the (unique) stable equilibrium. The figure shows that the frequencies of the island alleles at the polymorphic equilibrium decrease with increasing recombination rate, and this is proved in Appendix A.3.2. Fig. 5 demonstrates that the invasion of a tightly linked mutant of smaller effect can elevate the frequency of the previously present island allele () considerably. In the displayed example, if is small, the frequency of is pushed above 50%. Thus, becomes more frequent than the continental allele .

Local adaptation and migration load

Local adaptation is impossible if gene flow is so strong that the equilibrium , at which the continental haplotype is fixed, is globally asymptotically stable. By Theorem 2, this is the case if or if Because the mean (Malthusian) fitness on the island is , the mean fitness at is , whereas in the absence of gene flow (at ) it is .

However, local adaptation may already be severely reduced if gene flow is strong enough to eradicate the island haplotype (but not the island allele ). Then the single-locus polymorphism is asymptotically stable. By Theorem 2, this occurs if holds, but also if (3.25) or (3.26) are satisfied and (Fig. 1). The mean fitness at is .

If gene flow is sufficiently weak to maintain both alleles in the population, the degree of local adaptation can be measured by the migration load, by the frequencies of the island alleles, or by the frequency of the island haplotype.

First, we briefly investigate the migration load at the fully polymorphic equilibrium . Because and are decreasing functions of and of (Appendix A.3.2), the same is true for mean fitness. For the mean fitness at , we obtain where the approximation holds for weak migration (Section 4.1.1). Thus, the migration load increases from (approximately) for completely linked loci to for independent loci. Simple approximations for other limiting cases are readily derived from Section 4.1.

For weak gene flow, the frequencies of the island alleles are given by (4.1). From (4.1) and (2.3), the following first-order approximation in is obtained for the frequency of the island haplotype: This, as well as (4.1), is accurate for a fairly large range of parameters. Remarkably, if , the linear approximations for the allele and gamete frequencies give the exact solutions, i.e., the equilibrium values at . Therefore, (4.1) and (4.20) provide highly accurate approximations if and . Fig. 6 displays the exact value and the weak-migration approximation of for various parameter combinations. As this figure shows, (4.20) is also very accurate if and . In general, (4.20) provides a useful approximation if .

Fig. 6

Frequency of the island haplotype as a function of the migration rate for various parameter combinations. The exact value of is displayed by black lines, the corresponding approximation (4.20) by gray lines of the same line type. Dashing indicates different recombination rates, as given in the legends. The selection coefficients are and .

Finally, we explore the case when the island alleles are maintained at frequencies greater than , so that they are more frequent than the alternative continental alleles, and when the island haplotype has a frequency greater than . As is suggested by Figs. 5 and 6, and proved in Appendix A.3.2, the frequencies of both island alleles at the stable fully polymorphic equilibrium increase with increasing linkage and are maximized if . A rough estimate for the maximum value that guarantees that is obtained from (4.20):

For small , the following slightly more accurate approximation for the critical value of is obtained directly from the tight-linkage approximation (4.2): Here, the error is . If , the approximation (4.22) fails, and (4.21) is more accurate.

We treat two examples. If , we compute from (3.15) that both island alleles are more frequent than the respective continental alleles if , and if . The approximation (4.21) produces .

In the extreme case that LD is negligible, i.e., if selection is very weak relative to recombination, the more advantageous allele () has a frequency greater than if . Both island alleles are more frequent than the respective continental alleles if . The island haplotype frequency exceeds if . If , this gives . From (4.21), one obtains and , respectively.

In summary, these results show that the linkage disequilibria generated by at least moderately tight linkage improve local adaptation much more than would occur under linkage equilibrium.

Effective migration rate

The concept of the effective migration rate was introduced as a measure of gene flow at a neutral marker locus because the flow of alleles may be impeded or enhanced by linkage to selected loci (Bengtsson, 1985; Barton and Bengtsson, 1986). Recently, Kobayashi et al. (2008) proposed a precise definition of the effective migration rate for quite general models and derived approximate formulas that facilitate the computation of the effective migration rate under a variety of scenarios.

Here, we calculate the effective migration rate for a neutral locus that is linked to the two selected loci based on a continuous-time variant of the treatments of Bengtsson (1985), Barton and Bengtsson (1986), and Kobayashi et al. (2008). Since we assume one-way migration, the rate of change of allele frequency at an isolated neutral locus is , where denotes the (relative) frequency of the allele on the island and is the (constant) fraction of alleles on the continent. Thus, the fraction of alternative alleles at the neutral locus on the continent, subsumed under , is . In this simple scenario, the migration rate satisfies for every . Thus, determines the rate of approach to equilibrium and can be interpreted as the leading eigenvalue of the (one-dimensional) Jacobian at equilibrium.

If the neutral locus is linked to one or more selected loci, the right side of (4.23) is no longer constant because evolution at the neutral locus will depend on allele frequencies at the selected loci and on the linkage disequilibria. If the selected loci are in equilibrium, the asymptotic rate of convergence of can be defined by the leading eigenvalue of the Jacobian of the system that describes evolution of the frequency of and the corresponding linkage disequilibria. Therefore, we define the effective migration rate as (cf. Bengtsson, 1985; Kobayashi et al., 2008). We explicate this approach below.

It is straightforward to derive the differential equations for the gamete frequencies. We order the gametes as , and denote their frequencies by , respectively. The rates of change caused by selection and immigration are immediately obtained from (2.2) because the locus is neutral. We denote the recombination rate between loci and by and that between and by . Because our model is continuous in time, in the absence of interference we have . The rate of change of due to recombination is and analogous expressions are obtained for the other gametes. The total rate of change is the sum of those caused by selection, migration, and recombination.

As above, for our purposes it is more convenient to work with allele frequencies and linkage disequilibria. Let , and denote the linkage disequilibria between the indicated pairs of loci, and let denote the three-way linkage disequilibrium.

Of course, the equations for the change of , and are given by (2.5). The allele frequency at the neutral locus evolves according to and the equations for the linkage disequilibria involving the neutral locus are

This system has an (asymptotically stable) equilibrium such that the selected loci are at the equilibrium (3.15), and and hold. We denote this equilibrium by . We use the order for the allele frequencies and linkage disequilibria. With the help of Mathematica, it is straightforward to show that the Jacobian at the equilibrium has the block structure Here, is the Jacobian approximating convergence of to , and is that approximating convergence of to : is given in Box I.

The rate of convergence to equilibrium at the neutral locus is determined by the leading eigenvalue of . In the limit of weak migration, the following approximation for this leading eigenvalue is deduced (as can be checked easily with Mathematica): Since we define the effective migration rate as , we obtain the approximation This shows that the factor by which the effective migration rate is reduced relative to the actual migration rate, , is the product of the factors resulting from linkage to a single selected locus (cf. Petry, 1983). Bengtsson (1985) called the gene-flow factor. It is a measure of the penetrability of a genetic barrier.

We note that also approximations for the other eigenvalues can be derived. As a consequence, asymptotic stability of this equilibrium has been shown for weak migration.

A simple calculation yields that is maximized at the recombination rate and is the corresponding maximum value of . In the simple case , we obtain and Therefore, the effective migration rate is massively reduced unless is much larger than . For instance, if , it is reduced to at least of . Fig. 7 displays the gene-flow factor for various parameter combinations.

Fig. 7

Gene-flow factor. The lines in the figure display the gene-flow factor (4.30) as a function of for and . The lines from top to bottom are for , 1, and . The gene-flow factor is maximized at which is given by (4.31). For the given parameters, these values are between 0.415 and 0.484 (from top to bottom).

Discussion

For a diploid two-locus continent–island model with genic selection on two loci, we investigated the combined effects of gene flow and linkage on the maintenance of polymorphism, the amount of LD, the degree of local adaptation, and the fate of newly arising slightly beneficial mutations that would be eliminated by gene flow in the absence of LD. Such a model appears to be an appropriate choice for studying evolution in a derived (island) population that experiences altered environmental conditions and maladaptive gene flow from the ancestral (continental) population.

Our prime analytical achievement (Theorem 2 and Fig. 1) is the complete characterization of all equilibrium configurations. In particular, Eqs. (3.15) provide explicit analytical formulas for the feasible fully polymorphic equilibria, and Eqs. (3.17) and (3.18) give their admissibility conditions. From these results, simple and informative approximations are derived for the measures and of LD and for the migration load, as well as simple criteria for the invasion of slightly beneficial mutants and expressions for their final frequency. In the following we discuss the implications on a number of biological issues.

Maintenance of genetic variation

Theorem 2 and the concomitant evaluation of the potential for the maintenance of polymorphism show that both selectively favored alleles can be maintained in the population for higher levels of gene flow than would be concluded by ignoring LD. Under the assumption of independent loci, both loci can be maintained polymorphic if and only if the migration rate is less than the smaller of the selection coefficients (). The actual upper bounds are given by Eqs. (3.17a)–(3.17c), which show that, essentially, three scenarios have to be distinguished: low, intermediate, and strong recombination. In the first case, the approximate condition (4.6) is obtained, in the last case, the upper bound, (4.5) and (3.11), is between and .

The case of moderately strong linkage is more complicated, although the upper bound is always (3.16). Interestingly, if the selection coefficients of the alleles differ by a factor of less than two and if the recombination rate is about as large as the selection coefficient at the minor locus, two fully polymorphic equilibria ( and ) may exist, one stable, the other unstable (Fig. 1, diagrams (b), (d)–(f); see the discussion above Theorem 3). Then the stable fully polymorphic equilibrium coexists with a stable boundary equilibrium (either or ).

This bistable situation occurs for migration rates slightly below the maximum rate that admits coexistence of island and continental alleles at both loci. Which equilibrium is approached depends on the initial conditions. If, initially the island haplotype () prevails, then the fully polymorphic equilibrium is approached. If, initially, the island haplotype is rare, then recombination breaks up this haplotype sufficiently quickly so that it is lost and, depending on the parameters, the continental haplotype becomes fixed or the single-locus polymorphism , at which is fixed, is approached.

Linkage disequilibrium under gene flow

In large panmictic populations, LD can be maintained only if there is epistasis. This changes in spatially structured populations. Li and Nei (1974), Christiansen and Feldman (1975) and Slatkin (1975) demonstrated that nonepistatic selection can maintain LD within local subpopulations. For sufficiently weak migration between two demes, so that each subpopulation is close to fixation, Barton (1983) derived a simple recursion for the linkage disequilibria of arbitrary order among a finite set of equivalent loci carrying a deleterious mutant. For our continent–island model we were able to derive the explicit analytical formula (3.15c) for the amount of LD maintained under migration–selection balance. It holds for arbitrary relative strength of selection, recombination, and migration, and admits simple approximations in several limiting cases (Eqs. (4.8)–(4.11) and Section 4.1).

Our results show that the measure of LD is maximized at intermediate migration rates, whereas the measure is a monotone decreasing function of the migration rate (Fig. 3). As shown by (4.9) and explained at the end of Section 4.3, decays very slowly as a function of if is small.

Both and decay slowly as functions of if is small (Eqs. (4.8)–(4.10)), but faster for large (4.11). For instance, if the recombination rate is equal to the larger of the two selection coefficients (), then is reduced by less than 50% of its maximum value attained at . The measure is somewhat more sensitive to changes in , but even for loci of equal effects and is reduced only by a factor of relative to the case . Thus, LD is high if the recombination rate is of the same order of magnitude as the larger of the selection coefficients.

Consequences for linked neutral variation

The results above provide information about the LD to be expected among selected loci but not about the consequences for linked neutral variation. Therefore, we performed a preliminary analysis of a model with a third neutral locus that is located between the two selected loci. It shows that, at least for weak migration, linkage equilibrium between the neutral locus and the selected loci will be established. For weak migration, we derived the simple formula (4.30) for the so-called effective migration rate (Bengtsson, 1985; Barton and Bengtsson, 1986; Kobayashi et al., 2008). This is a measure for the strength of gene flow at a neutral locus embedded in a selected background. From this formula and Fig. 7, the extent of impediment of gene flow caused by linked selected loci becomes apparent. The formula shows that the linked loci contribute multiplicatively to the reduction in migration rate and thereby build an efficient barrier to neutral gene flow. A detailed study of the effects on linked neutral variation will have to examine the effects of random genetic drift and is an interesting topic for future research.

Dominance

Our analysis is based on the assumption of genic selection because this is the simplest case and almost all inference methods of selection are based on it. However, we also briefly investigated the influence of dominance (results not shown). In short, small deviations from no dominance lead to the same qualitative model behavior. If both island alleles are partially dominant, a two-locus polymorphism can be maintained for higher migration rates, whereas the opposite is true if they are partially recessive. For weak migration, to leading order in , the equilibrium values , and at the fully polymorphic equilibrium are obtained from Eqs. (4.1) by the substitutions and , where and are the selection coefficients of the single-locus genotypes and , respectively.

If the island alleles are completely dominant (), the maximum migration rate admitting polymorphism is about twice as high as without dominance. If the island alleles are highly or completely recessive, the equilibrium structure of the model is much more complicated because at least four fully polymorphic equilibria can coexist. This occurs for weak migration and strong recombination. In view of the fact that in the one-locus model, two polymorphic equilibria can coexist if the island allele is sufficiently strongly recessive (Nagylaki, 1992, Chap. 6.1), this is not too surprising, but the model behavior could become more complicated at higher migration rates (as it is the case for genic selection).

Local adaptation

The present model enables us to study several aspects of local adaptation in the presence of gene flow. Simple measures of the degree of local adaptation are the migration load and the frequencies of the island alleles or the island haplotype. The effects of linkage on the migration load are relatively straightforward. If gene flow is low, so that the island haplotype is maintained at appreciable frequency, Eq. (4.19) shows that the load is twice as high for very strong recombination as for no recombination. This is easy to understand because recombination predominantly breaks up the well adapted island haplotype if it is frequent. Tighter linkage always increases the frequency of the island alleles and of the island haplotype; see Figs. 5 and 6, Appendix A.3.2, and Eq. (4.20).

Simple instructive approximations are derived for the maximum rate of gene flow such that the frequency of the island haplotype exceeds ; e.g. (4.21). Fig. 6 illustrates the finding that, in the presence of gene flow, already moderately strong linkage enhances the degree of local adaptation substantially. Thus, although gene flow diminishes local adaptation, linkage among advantageous genes may reduce the adverse effects of gene flow considerably in comparison to conclusions from single-locus models (Lenormand, 2002).

For a given rate of gene flow, and disregarding LD, a mutant can invade in our deterministic model if and only if its selective advantage exceeds . The present analysis demonstrates that mutants of much smaller effect can invade successfully provided they occur close to a locus at which an advantageous allele is segregating in migration–selection balance (which requires a selective advantage greater than ). This may be interpreted as a form of hitch-hiking. Eqs. (4.13)–(4.15) provide simple explicit estimates on the relation between linkage and the selective advantage required for invasion. Fig. 4 documents the extent to which linkage decreases the minimum selective advantage that admits invasion. In addition, invasion of a linked mutant of small effect increases the frequency of the already previously present island allele at the major locus (Fig. 5). It would be a challenging, though worthwhile, enterprise to derive the probability of invasion in a population of finite size. Then it will matter whether the allele occurs on a or background.

Evolution of genetic architecture

We showed that weakly advantageous alleles can invade if and only if they are sufficiently tightly linked to an already segregating mutant of large effect. This suggests that long-term evolution under continued gene flow should lead to the emergence of clusters, or linkage groups, of locally adapted alleles. A similar conclusion was drawn recently by Yeaman and Whitlock (2011), who performed a simulation study on the evolution of genetic architecture when two populations connected by migration experience stabilizing selection on a polygenic trait towards different optima. Griswold (2006) investigated a model similar to that of Yeaman and Whitlock, but (mainly) with divergent directional selection. He explored the distribution of allelic effects but presented no information about linkage relations among loci causing phenotypic differentiation.

Our model can also be interpreted in terms of selection on a quantitative trait if the two loci contribute additively to the trait, i.e., without epistasis and dominance, and if on the island linear directional selection acts on the trait (cf. Bürger, 2010). Therefore, clusters, or linkage groups, of adaptive loci can also evolve under directional selection. Our results not only support the main finding of Yeaman and Whitlock (2011), but they provide quantitative relations between the selective advantage and the degree of linkage to the major locus required for invasion and, hence, the emergence of such clusters (Eqs. (4.13)–(4.15), Fig. 4). It would be worthwhile to study the generality of the hypothesis that evolution in spatially structured populations leads to the emergence of clustered genetic architectures. For frequency-dependent disruptive selection, similar hypotheses have been raised (Kopp and Hermisson, 2006; van Doorn and Dieckmann, 2006), but slightly different model assumptions may lead to deviating conclusions (Schneider, 2007). For a discussion of the evidence from and the consequences for QTL analysis, we refer to Yeaman and Whitlock (2011); see also Griswold (2006).

The expression (4.19) shows that the mean fitness at the polymorphic equilibrium increases with decreasing recombination rate between the selected loci. This suggests that reduced recombination rates will be evolutionarily favored in populations subject to new selective challenges but suffering from gene flow through maladapted individuals. As shown by Lenormand and Otto (2000), recombination can also be favored in spatially structured populations, and predicting whether increased or decreased recombination is evolutionary advantageous requires detailed information on the pattern and intensity of migration as well as on the variability of selection across environments.

An efficient mechanism in suppressing recombination among a set of genes is provided by chromosome inversions. Indeed, Kirkpatrick and Barton (2006) suggested that chromosome inversions that capture locally adapted alleles should spread in a population exposed to gene flow by maladapted genotypes. Their arguments are based on the calculation of the invasion rate of such an inversion using a continent–island framework similar to ours. Their derivation assumed weak migration and strong recombination, and either multiplicative effects among a set of loci in a haploid population or additive effects in a diploid population. Below, we show that the true invasion rate may be higher than their approximation suggests.

In our continuous-time model, the invasion rate of the inversion is , which is calculated from (2.2a) by omitting the term and evaluating at the fully polymorphic equilibrium . The result is If migration is weak, the invasion rate is approximated by This is a highly accurate approximation if is small or if is either much smaller or much larger than both selection coefficients. If is approximately as large as the selection coefficients and if the migration rate is not much smaller than the maximum migration rate under which the island haplotype (with recombination) can be maintained, then may underestimate the actual fitness gain substantially. For instance, if and (which, as we showed below Theorem 3, is the maximum migration rate), then (5.2) predicts , whereas the precise value is . Similarly, if , and , then (5.2) predicts , whereas the exact rate is . Again, the approximation (5.2) may underestimate the true value up to a factor of two.

We also note that the gain in mean fitness resulting from the establishment of a chromosomal inversion capturing the locally advantageous alleles is twice the invasion rate . This follows because the (mean) fitness after the inversion has spread is (independently of the recombination rate between the loci since the subpopulation carrying the inversion is fixed for the island alleles); hence, the fitness gain can be approximated by , where is evaluated at , and (4.19a) yields the result.

For two loci of equal effects, is (after correcting for haploidy and a misprint) equivalent to Eq. (3) in Kirkpatrick and Barton (2006), who assumed strong recombination. From the above arguments, we conclude that if recombination rates are of the same order of magnitude as the selective coefficients and migration is moderately strong (such that the island genotype is maintained), the simple approximation may underestimate the true invasion rate by up to a factor of two. Therefore, chromosome inversions may spread even if they capture relatively tightly linked, locally adapted alleles. Thus, this mechanism may be even more general and efficient than suggested by Kirkpatrick and Barton (2006).