Monotone and nonmonotone clines with partial panmixia across a geographical barrier

The number of clines (i.e., nonconstant equilibria) maintained by viability selection, migration, and partial global panmixia in a step-environment with a geographical barrier is investigated. Our results extend the results of T. Nagylaki (2016, Clines with partial panmixia across a geographical barrier, Theor. Popul. Biol. 109) from the no dominance case to arbitrary dominance and to various other selection functions. Unexpectedly, besides the usual monotone clines, we discover nonmonotone clines with both equal and unequal limits at \begin{document}$ \pm\infty $\end{document} .

1. Introduction. The term cline was first coined by the English evolutionary biologist Julian Huxley in 1938, which describes a continuous gradient of gene frequencies, genotypic frequencies, or phenotypic frequencies of a species across its geographical range. Since the formation of clines is closely related to species adaptation, speciation, and the maintenance of genetic variability, the study of clines has been an important research subject in population genetics, ecology, and related fields.
A geographical barrier that divides a habitat into two often occurs in nature such as a mountain chain, a river, or a railway. However, relatively little cline theory has considered this factor; see, e.g., Nagylaki [13] and the references therein. An obvious distinction is that clines will be discontinuous at a barrier. The spread of an advantageous allele may also be delayed by a barrier [25].
In the pioneer work [16], Nagylaki derived a new model that includes all of the above features, namely, a single-locus migration-selection model with partial panmixia in the presence of a geographical barrier. In particular, he deduced an explicit formula for the unique cline, which is monotone, under the assumptions of two alleles without dominance, a step-environment, and homogeneous and isotropic local migration on the entire line.
The purpose of this paper is to generalize Nagylaki's analysis in [16] from the no dominance case to arbitrary dominance, and to various other selection functions. In previously studied models, clines maintained by a step-environment are usually monotone [4,27,16]. Surprisingly, for a certain range of degree of dominance, we discover the existence of nonmonotone clines maintained by a step-environment. Our results shed some light on the interaction of these evolutionary factors on clines, especially the joint effects of degree of dominance, partial panmixia, and geographical barrier.
In Section 2, we formulate the model (1) and introduce some preliminary properties. We present our results on clines of (1) in Section 3 for the case p − < p + ; the main results are Theorems 3.1, 3.4, 3.7, 3.8, 3.11, and 3.12. Section 4 is devoted to the case p − = p + ; the main results are Theorems 4.1 and 4.3. In particular, the existence of nonmonotone clines is established in Theorems 3.11 and 4.3. We summarize our results and discuss open problems in Section 5.
2. The model problem and preliminary properties. Here, we briefly recapitulate the formulation in [16,Sect. 4]. Assume that the gene under consideration is at a single-locus with two alleles A 1 and A 2 . The diploid population occupies the entire line R. Under the joint action of viability selection, local adult migration, and partial panmixia with a geographical barrier at x = 0, Nagylaki derived the following equation for the frequency of A 1 at equilibrium [16]: In (1a), g(x)f (p) designates the effect of selection on allele A 1 . If the fitness of genotype A i A j , denoted by r ij for i, j = 1, 2, is independent of gene frequency, one may assume where the constant h ∈ [−1, 1] is the degree of dominance. When h = 0, h = −1, and h = 1, there is no dominance, A 2 is completely dominant to A 1 , and A 2 is recessive, respectively. Under (2), the selection function f reads The spatial factor g(x) specifies the direction of selection, i.e., at location x, allele A 1 (A 2 ) is selectively favored if g(x) > 0 (< 0).
In this paper, following [16], we focus on a step-environment, i.e., where the parameter α > 0 measures the relative strength of negative selection to positive selection. Note that this step-environment applies also to a plane habitat of two types that is divided by a linear boundary, and x is the directed distance of any point from this boundary. In general, directional selection functions satisfy In [19] and [27], the unimodality of f emerged as a crucial simplifying assumption, i.e., ∃p ∈ (0, 1) s.t. f is strictly increasing & decreasing in (0,p) & (p, 1), resp., (6) which is applied also in this paper. In particular, the cubic f in (3) is unimodal. The term β[p − p(x)] in (1a) describes the effect of long-distance migration. It means that at location x, a portion of population is replaced by the "averaged" population over the habitat due to panmixia. The parameter β > 0 is the scaled panmictic rate; a larger β means a larger portion of the population is pamictic. The term p represents the effect of population local migration, whose rate is scaled out.
The discontinuities of p(x) and p (x) at x = 0 are due to a geographical barrier there. They satisfy the transmission conditions (1d,f). The parameters θ ± > 0 are the scaled transmissivities crossing the barrier from left to right and from right to left, respectively. Smaller θ ± corresponds to a stronger barrier. In particular, θ ± → 0 and θ ± → ∞ show that the barrier becomes impenetrable and the barrier disappears, respectively.
In the absence of both frequency dependence and dominance, i.e., h = 0 in (3), the selection function has the simplest form For every α > 0, define the critical panmictic rate as If α > 0, β ∈ (0, β 0 ), and θ ± > 0, then (1) has a unique cline p(x) ∈ C 2 (R \ {0}), which satisfies (i) The existence and uniqueness in Theorem 2.1 were proved by tedious calculations that verify the sign and number of roots of many involved polynomials, relying crucially on the special form (8). In the next section, we mainly use phase-plane analysis to generalize Theorem 2.1 from (8) to a much wider class of selection functions.
3. Clines with p − < p + . We establish the existence, uniqueness, multiplicity, monotonicity, and non-monotonicity of clines for problem (1) with various selection functions f . According to Proposition 1(iii), the gene frequencies at ±∞ satisfy either 0 < p − < p + < 1, p − = p + = 0, or p − = p + = 1. In this section, we will focus on clines with p − < p + , and treat the latter two cases in Section 4. We state our results for α = 1 and α > 1 in Sections 3.1 and 3.2, respectively. The case α < 1 can be converted to the case α > 1 by considering the equation of 1 − p(−x); see [27,Sect.7.2.3] for more details regarding this transformation.
The following functions play important roles in our proofs.

Figure 1
which is Hamiltonian with energy function H + (p, q) as in (11). From (11) and (14a,b) we have the following observations. (a) H + is an even function of q.
Then one can easily draw in the pq-plane the phase portrait for system (15). Fig.1(A) shows the solutions that lie in the level H + (p + , 0). In light of the boundary conditions (1b,c), the solid arcs tending to (p + , 0) as x → ∞, as indicated by the arrow heads, are of special interest. They are labeled by C + in Fig.1(A).
Then one can easily draw in the pq-plane the phase portrait for system (18). Fig.1(B) shows the solutions that lie in the level H − 1 (p − , 0). In light of the boundary conditions (1b,c), the solid arcs emanating from (p − , 0) as x → −∞, as indicated by the arrow heads, are of special interest. They are labeled by C − in Fig. 1(B).
We draw a typical cline p(x) guaranteed by Theorem 3.1 in Fig. 3(A), in which p ± are the p-coordinate of points B and A in Fig. 2(B), respectively. 3.2. α > 1. In this section, we focus on unimodal f with a maximum value achieved atp, which is either concave down on the left ofp or has one point of inflection there. Define Note that the β 0 in (9) agrees with β * when α > 1 and f (p) = p(1 − p).
The assumptions α > 1, (5), (6), and f (p) ∈ C 2 ([0,p]) with (23) imply that system (7) has a unique solution pair p ± with 0 < p − < p + < 1 for every β ∈ (0, β * ) and only trivial solutions for every β ≥ β * . Condition (23) says that f is concave down on the left of its maximump. Next, we discuss the case that f has one inflection point on the left ofp, i.e., there existš p ∈ (0,p) such that In this case, as explained in Remark 3, the phase portrait of (18) remains the same as in Fig. 1(B). However, we will see that the phase portrait of (15) can be much more complicated than in Fig. 1 Remark 3. For α > 1 and unimodal f with maximum value achieved atp, the proof of Lemma 7.12 (especially, Cases 1 and 3) in [27] shows that p − is the unique zero of [αf (p) + β(p −p)] in (0, 1), which implies that (14c,d) hold. Then the arguments for x < 0 in the proof of Theorem 3.1 apply, and the phase portrait of system (18) is as in Fig. 1(B). (i) m = 1: The unique intersection point has to be p + and it is clear that (14a,b) hold from the geometry of f and L. Then the arguments for x > 0 in the proof of Theorem 3.1 apply, and the corresponding phase portrait of system (15) is as in Fig. 1(A). (ii) m = 2: We order the two intersection points as p 1 and p 2 with p 1 < p 2 ; then p 1 <p < p 2 and L is tangent to f at either p 1 or p 2 . Figs. 4 and 5 show L and f with the corresponding phase portraits of system (15) for the tangent point being p 1 and p 2 , respectively. There are four possibilities.
(a) Fig. 4 with p + = p 1 : The orbit tending to p + in Fig. 4(B) is qualitatively the same as in Fig. 1 Fig. 4 with p + = p 2 : Same as (a). (c) Fig. 5 with p + = p 1 : Same as (a). (d) Fig. 5 with p + = p 2 : The orbit tending to p + in Fig. 5(B) is below q = 0, which has no point of intersection with Φ θ−,θ+ (C − ) (see Fig. 2(B)); thus no cline exists.  In light of Remarks 3 and 4, the same arguments as in Theorem 3.1 immediately yield the following Theorems 3.7 and 3.8.
Proof. The proof is straightforward by (10) and (11): Lemma 3.10. Assume that the assumptions in Lemma 3.9 hold. (i) If I < II, then the phase portrait of (15) is as in Fig. 7(A); there exists a homoclinic orbit connecting the equilibrium (p 1 , 0) with itself.
(ii) If I = II, then the phase portrait of (15) is as in Fig. 7(B); there exist heteroclinic orbits connecting the equilibria (p 1 , 0) and (p 3 , 0). (iii) If I > II, then the phase portrait of (15) is as in Fig. 7(C); there exists a homoclinic orbit connecting the equilibrium (p 3 , 0) with itself.

Proof. Observations (a) and (b) in the proof of Theorem 3.1 still hold, and observation (c) becomes
Then according to the relation between H + (p 1 , 0) and H + (p 3 , 0), or the equivalent relation between I and II, in each case, one can easily draw in the pq-plane the solution portraits that lie in these two levels for system (15). Note that the point (p 0 , 0) in Figs. 7(A,C) is not an equilibrium.
We first consider the case (A) I < II in Fig. 7. In fact, this is the only case that may produce nonmonotone clines with p − < p + . We make the following assumption.
two intersection points D 1 and D 2 for example. For each such point of intersection D i , we denote its unique preimage via Φ θ−,θ+ by E i , then the corresponding cline p(x) of (1) has the path Moreover, p (x) > 0 until p(x) reaches its maximum value p 0 , and p (x) < 0 afterwards. This completes the proof of Part (ii).
A typical cline in Theorem 3.11(i) is qualitatively the same as in Fig. 3(A). We exhibit a typical nonmonotone cline in Theorem 3.11(ii) in Fig. 3(B).
Next, we consider the cases (B) I = II and (C) I > II in Figs. 7.
Proof. From Figs. 7(B,C), we see that the phase portrait of the orbit that tends to p + is qualitatively the same as in Fig. 1(A), whence with Remark 3, the conclusion in Theorem 3.12 follows from the same arguments as in Theorem 3.1.
Remark 6. We observe in either case of Fig. 7, that the level curve of H + (p 2 , 0) contains a single point p 2 . Therefore, if p + = p 2 , there is no cline p(x) connecting p ± as x → ±∞.
Remark 7. In either case of Fig. 7, if p + = p 3 , only monotone clines may exist. (A) I < II: We denote the curve with energy H + (p 3 , 0) by C + and decompose it into Recall Remark 3, for every θ ± > 0, the phase portrait of (18) and the image Φ θ−,θ+ (C − ) are as in Figs.1(B) and 2(A). Note that {C 1,+ } is not monotone now, which may have multiple points of intersection with Φ θ−,θ+ (C − ), each of which will produce a monotone cline connecting p ± . (B) I = II: We denote the heteroclinic orbit emanating from (p 1 , 0) to (p 3 , 0) by C + . Now for each θ − > 0, in light of (13), one may raise or lower the curve Φ θ−,θ+ (C − ) by varying θ + . Let θ * be the critical value such that Φ θ−,θ+ (C − ) does not intersect C + if θ + > θ * and intersects C + at two or more points if θ + < θ * . Each such point of intersection will produce a monotone cline connecting p ± . (C) I > II: We denote the orbit from (p 0 , 0) to (p 3 , 0) by C + . Then the situation is very similar to (B).
Lastly, we apply the above results on clines to the cubic f in (3) with h ∈ [−1, −1/3). In this case, the cubic f satisfies (5), (6), and (24) witȟ We shall find out for what values of h, the assumption (H1) may hold so that nonmonotone clines may exist. Lemma 3.13. Assume that f is as in (3) with h ∈ [−1, −1/3) and that [f (p) − β(p −p)] has three zeros p 1 < p 2 < p 3 in (0, 1) as in Fig. 6. Then I II if and only if p 2 p.
Proof. We first show that if p 2 =p, then I = II. Define the line in which the "=" is due to the fact that Sincě p is the unique critical point of the quadratic polynomial f (p), we see that Let φ(p) = f (p) − L(p), then (31) and (32) show that We conclude from (33) that φ(p) is an odd function about p =p, which with φ(p 1 ) = φ(p 3 ) = 0 reveals that p 1 + p 3 = 2p and I = II. Next, suppose that p 2 <p. We denote the strict line that passes through (p 1 , 0) and (p, f (p)) by L 1 (p), and the two areas enclosed by L 1 and f by I and II .
Then I = II as above and it is clear that I < I = II < II from the geometry of L, L 1 , and f . If p 2 >p, we see that I > II similarly. Thus, Lemma 3.13 is demonstrated.
Proof. The proof of this theorem is similar to the one of Theorem 3.17, except that for h ∈ [−1, −(5 + √ 17)/12), there exist α > 1, β > 0, p ± , andp such that Case (ii) in Remark 8 happens by Remark 9 and Lemma 3.14. Once this happens, the existence of nonmonotone clines in addition to a unique monotone one follows from Theorem 3.11.
We focus on unimodal selection functions f described by (6) that satisfy either (23) or (24) as in Section 3. We first consider nonmonotone clines of (1) with p − = p + = 0.
Proof. We use the phase-plane analysis again. The phase portrait for x < 0 is similar to the one in Fig. 1(B) with p − being the origin now. We observe that under the conditions in either (i) or (ii), the line βp may have at most one point of intersection with the graph of f (p) besides the origin, which determines the phase portrait for x > 0 in a simple way. Consequently, in either case, the trivial solution p(x) ≡ 0 is the only solution with p(±∞) = 0. We omit the details here, since the method is essentially the same as the proof of Theorem 3.1.
Remark 10. If f satisfies (5), (6), and (24), then by the definition of β 1 and β 2 in Theorem 4.1(ii), we see that for every β ∈ (β 1 , β 2 ), the line βp and the graph of f (p) have two points of intersection besides the origin. The situation is exactly the same as in Fig. 6 with p − being the origin now. I and II stand for the two areas bounded by the line βp and the graph of f (p) as in Fig. 6. It is clear that there exists a uniqueβ ∈ (β 1 , β 2 ) such that I ≶ II if and only β ≶β. Similarly to Lemma 3.10, we have the following lemma. (i) If β ∈ (β 1 ,β), then the phase portrait of (15) withp = 0 is as in Fig. 7(A) with p 1 being the origin; there exists a homoclinic orbit joining the equilibrium (p 1 , 0) and itself.
The proof of Theorems 4.3(i) and 4.3(ii) is essentially the same as the proof of Theorems 3.11 and 3.12, respectively, and hence is omitted here.
unimodal f which is concave on the left of its maximum for α > 1 (Theorem 3.4); and to unimodal f which has only one point of inflection on the left of its maximum for α > 1 (theorems summarized in Remark 8). The results on clines for α < 1 can be obtained from the results for α > 1 through the transformation 1 − p(−x).
In particular, we obtain the complete cline structure for the biological important cubic f in (3): The configuration of clines with p − < p + is established in Remark 1 for α = 1, and in Remark 2 and Theorems 3.17 and 3.18 for α > 1, respectively. The configuration of clines for α > 0 with p − = p + = 0 and 1 is established in Theorems 4.4 and 4.5, respectively. Now we posit some unsolved problems. First, the frequencies p(0±) at the barrier are important characteristics of a cline p(x). Therefore, it is desirable to do the approximations at the barrier in various limiting cases of the parameters α, β, and θ ± for the selection functions that we investigated in this paper as in [16,Sect. 4.3] for f (u) = u(1 − u).
Second, it will be much more challenging to study the stability of the clines with respect to the corresponding time-dependent model.
Lastly, if we allow more points of inflection on the left ofp instead of only one as in (24), how complex will the cline structure be? Another generalization is g(x) having two steps or being a continuous monotone function instead of a single step as in (4).