Journal of Mathematical Biology 3,263--269 (1976) 9 by Springer-Verlag 1976
Stochastic Selection in Both Haplophase and Diplophase* D. L. Hartl, West Lafayette, Indiana Received December 9, 1975
Summary A population genetic model incorporating the evolutionary forces of zygotic selection, gametic selection and non-Mendelian segregation has been analyzed for the case in which all selection coefficients and the segregation parameter are assumed to be r a n d o m variables that are uncorrelated from generation to generation. The diffusion approximation of the model is developed, and the subsequent analysis shows that one of four limiting outcomes of the stochastic process may obtain - - an allele m a y be fixed or lost almost surely and irrespective of the initial gene frequency, the gene frequency may converge to a unique stationary distribution, or an allele m a y be fixed or lost with probabilities depending on the initial gene frequency. These outcomes correspond rather closely with the possible outcomes of the deterministic model - - fixation or loss of an allele, convergence to a stable equilibrium, or the existence of an unstable equilibrium.
This paper presents an analysis of a population genetic model that combines the forces of zygotic (diplophase) selection, gametic (haplophase) selection, and nonMendelian segregation, in which all the parameters are assumed to vary randomly from generation to generation. Interest in the model actually centers on two special cases: the case of zygotic selection combined with gametic selection, on the one hand, and the case of zygotic selection combined with non-Mendelian segregation, on the other hand. Scudo (1967) has studied a deterministic model of zygotic selection combined with gametic selection; Hiraizumi, Sandler and Crow (1960) have analyzed a deterministic model of zygotic selection combined with non-Mendelian segregation. Although there is no experimental precedent for combining the two sorts of models, they are conveniently combined for theoretical purposes and the results of the analysis specialized to whichever particular model may be of interest. The whole question of the effects of stochastic variation in selection intensity has come under close scrutiny recently, especially in diploid, one-locus, diallelic models with variation in the intensity of zygotic selection (Hartl and Cook 1973, 1976; Gillespie 1973; Karlin and Lieberman 1974, 1975; Karlin and Levikson 1974; Levikson 1974; Levikson and Karlin 1975; Cook and Hartl 1975; Norman * Work supported by N. I. H. grants GM21732 and GM21623. The author is supported by Research Career Award GM2301.
264
D.L. Hartl:
1975). A principal motivation has been to determine whether the nuisance caused by variation in selection intensity renders the results of classical constant fitness models totally inappropriate in the sense that variation in selection intensity generates qualitatively new phenomena that are not encompassed by the deterministic theory9 This is an important question because selective values in nature are almost certainly not constants, at least not for any appreciable length of time. (How much selective values may vary is, of course, unknown; autocorrelations between the fitnesses in successive generations are also unknown.) The most optimistic hope for the classical deterministic models has been that their results, suitably re-interpreted, would be valid for the limiting behaviour of stochastic models9 So far, in infinite population models, this has been the case, although some care must be exercised in the re-interpretation. It seems likely, however, that the deterministic results will not carry over straightforwardly in the cases of finite population models9 The results of the present heuristic analysis have been briefly summarized elsewhere (Hartl 1975), and the details of the analysis will be presented here. Throughout the analysis I will use the methods employed by Levikson (1974) and Levikson and Karlin (1975), methods which are originally due largely to Feller (1954)9 The Deterministic Model
Consider a locus having two alleles, A and a, in a large, randomly-mating, diploid population. Let the frequency of A among newly-formed zygotes be denoted p, and let q = 1 - p. The frequencies of A A, A a and a a genotypes among newly-formed zygotes will therefore be p2, 2 p q and q2, respectively9 Suppose the egg-to-adult viabilities of AA, A a and a a zygotes are w11, w12 and wz2, where wis>O. Among adults, the frequencies of AA, Aa and aa will be p2w11/~,, 2 p q w 12/# and q2 w22/w, respectively, where # = p2 Wll + 2 p q w 12 + q2 w22. Suppose now that the fraction of functional A-bearing gametes produced by an A a heterozygote is k, where 0 < k < 1. (This is the main assumption in the model of Hiraizumi et al. 1960). Then, at the moment of gamete formation, the frequency of functional A-bearing gametes will be (p2 Wl1 "at-2 p q k Wl2)/W.Suppose that these gametes are subject to further selection, and that A- and a-bearing gametes survive with probabilities vl and v2, respectively9 (This is the main assumption in the model of Scudo 1967.) After the haplophase selection, let the frequency of functional A-bearing gametes be denoted by p'. This frequency will perforce be the frequency of A among newly-formed zygotes of the next generation and
\W12/')2,]
pr \
12 2 /
\ /)2 ,//
\Wl2J
It will be useful in the sequel to have an expression for A p = p' - p :
Stochastic Selection in Both Haplophase and Diplophase
w11/)1
2(1-k)
PLwI2/)2
Pq Ap=
+q v2
PZ(W w12 ~ l /v2 ) l /) + p q [
265
w12
.
(1)
/)2 + 2 ( l - k ) ] +qz(w22]kw12 J
This deterministic model has a fixed point,/3, at 2 k W12 /)1 --W22 /)2
/3= 2 k w12/)1 + 2 (1 -k)
WI.2 I)2 - - W l l /)1 --W22 /)2 '
2 (1 - k) v z , or if (b) w22 > 2 k/)1 and 0 W12
then p--*l, whereas if (d) w2z > W12
2kvl /)2
and
wll
3. Thus the process converges weakly to the diffusion process with drift and diffusion coefficients in (3) and (4).
Analysis of the Stochastic Model A biologically pathological case arises if there exists a value of p, call it ~, such that # (~) = o-z (~) = 0. Inspection of (3) and (4) shows that this can happen only if the random variables are linearly dependent in such a way that $ 4 - $ 2 = ~ ( $ 1 - $ 2 - $ 3 + $ 4 ) , where e is some constant, 0 < c ~ < l . In this case, /3=e. Assume (for the moment) that the random variables are not linearly dependent in this manner. Then the diffusion process has only two boundaries, p = 0 and p = 1, and the limiting behaviour of the process is largely determined by the stochastic stability of these boundaries. Local stochastic stability of a boundary can be investigated using the function S (p)= ~f,os (z) dz, where po=p (0) and s (p)=exp ( - 2 ~o # (x)/0-2 (x) dx). The integral in S (p) converges near a boundary if and only if the boundary is locally stochastically stable. Partial fraction expansion of #(p)/o-2 (p) and
M2-M4-R24+
V4
integration yields s(p)=p-2a(1-p)2Bg(p) where A V2_2R24+V4 , M 1 - M 3 +R13 - V1 B, and where 9 (p) is a function that is continuous in
V1-2R13+V3
[0, 1] and bounded away from zero. Therefore, around p = 0, S (p) converges iff -2A>-1, which is to say iff M 2 - V2/2 - 1, which is to say M1 - V1/2 > M3 - V3/2; this is a necessary and sufficient condition for p-- 1 to be locally stochastically stable. A key theorem in working out the global behavior of the diffusion process is the following: Given p (0) = Po and two constants, ~ and 6, with 0 < e < Po -< 1 - 6 < 1, then
s (po)- s (~)
Prob {p (t) reaches 1 - 6 before reachin9 e I Po, 0 < e_< Po < 1 - 6 < 1} - S (1 - 6 ) - S (~) = L (e, 6) where S (p) is as defined above. We are interested in the behavior of L (e, 6) as ~, 6 4 0 . Consider first the case M 2 - V z / 2 < M 4 - V 4 / 2 and M 1 - V~/2 M 4 - V4/2 and M ~ - V I / 2 > M 3 - 1/3/2, the case in which p = 0 is stochastically unstable and p = 1 is locally stochastically stable. Then S ( 1 - 6) remains bounded as 6--.0 whereas S (e)~oo as e~0. Therefore L (e, 6)--,1 as e, 6--*0, which implies p ( 0 - ' 1 almost surely. Thus, in this case, p = 1 is globally stochastically stable. The third case to consider is M 2 - V2/2 > M 4 - V4/2 and M 1 - V1/2 < M 3 - V3/2, in which case both p = 0 and p = 1 are stochastically unstable. In this case, a necessary and sufficient condition for the gene frequency to converge to a unique stationary distribution is that 1/0-2 (p) s (p) be integrable on [0, 1]. Now, 1/0- 2 (p) s ( p ) = p -2+2A ( 1 - - p ) - Z - 2 ~ g * (p), where A and B are as defined above and g* (P) is continuous on [0, 1] and bounded away from zero. Obviously 1/o-2 (p)S (/9) is integrable on the open interval (0, 1). The function is integrable near 0 iff - 2 + 2 A > - 1, which is to say M 2 - - V 2 / 2 > M 4 - V 4 / 2 . The function is integrable near 1 iff - 2 - 2 B > - 1 , which is to say M 1 - V 1 / 2 < M 3 - V3/2. These are the identical conditions that make p = 0 and p = 1 stochastically unstable. Therefore, in this case, the gene frequency always converges to a unique stationary distribution. The final case to consider is M 2 - 1/2/2 < M 4 - V4/2 and m l - V~/2 > M 3 - V3/2 , in which case both p = 0 and p = 1 are locally stochastically stable. In this case both S ( 1 - 6 ) and S(e) remain finite as 5, 6--.0. Thus L (e, 6) remains finite as e, 6"~0. In this case, as a consequence, some sample paths will converge to p = 0 and the rest will converge to p = 1, with the probability of convergence to a given boundary being determined by the initial gene frequency. In summary, presenting the cases in the same order as the analogous cases in the deterministic model, we have if (a) M 2 - V2/2> M 4 - 1/4/2 and M 1 - VI/2 < M 3 - V3/2, then p converges to a unique stationary distribution; if (b) M 2 - V2/2 < M , - V4/2 and M 1 - I/1/2 > M 3 - V3/2 , then, with some probability dependent on the initial gene frequency, p will converge to 1, and with the remaining probability it will converge to 0; if (c) M 2 - V2/2 > M 4 - V4/2 and MI - V~/2 > M 3 - 1/3/2, then p ~ 1 almost surely; and if (d) M 2 - V2/2 < M 4 - I/-J2 and m 1 - V1/2 < M 3 - V3/2 , then p-~0 almost surely. Returning now to the one case excluded at the beginning of this section, assume that $ 4 - Sz = c~( S I - $ 2 - $3 + $4) where e is a constant and 0 < c~< 1. Since the case is of little interest biologically, local stochastic stability analysis will suffice. Partial fraction expansion of #(p)/0-2(p) by tedious algebra reveals that s ( p ) = p -2A ( 1 - p ) 2 B ( p - ~ ) - 2 c where A and B are as previously defined and C-M1-Mz-M3+M4-R14+R23 R 1 4 - R 1 2 - - R 3 4 q-R23 Journ. Math. Biol. 3 / 3 ~
. Thus the conditions for p = 0 and p = 1 to 19
268
D.L. Hartl:
be locally stochastically stable are as stated above, and the point p--~ is locally stochastically stable if and only if C < 1/2.
Discussion
To the order of the diffusion approximation, E In(1 + St) = M i - V~/2.Thus the four cases of interest may be restated in termsofthe original random variables as follows: if (a) Eln (W2g/W12)Eln [2 (1 - k ) v2/vi], then p-o0 or p ~ l with probabilities depending on the initial gene frequency; if (c) Eln (w22/w12)< Eln (2 k vi/v2) and Eln (wii/Wig)>Eln [2 (1 -k) v2/vl], then p--,1 almost surely; and if (d) Eln (w22/w12)>Eln (2 k vl/v2) and Eln (wli/wl2)
Stochastic Selection in Both Haplophase and Diplophase* D. L. Hartl, West Lafayette, Indiana Received December 9, 1975
Summary A population genetic model incorporating the evolutionary forces of zygotic selection, gametic selection and non-Mendelian segregation has been analyzed for the case in which all selection coefficients and the segregation parameter are assumed to be r a n d o m variables that are uncorrelated from generation to generation. The diffusion approximation of the model is developed, and the subsequent analysis shows that one of four limiting outcomes of the stochastic process may obtain - - an allele m a y be fixed or lost almost surely and irrespective of the initial gene frequency, the gene frequency may converge to a unique stationary distribution, or an allele m a y be fixed or lost with probabilities depending on the initial gene frequency. These outcomes correspond rather closely with the possible outcomes of the deterministic model - - fixation or loss of an allele, convergence to a stable equilibrium, or the existence of an unstable equilibrium.
This paper presents an analysis of a population genetic model that combines the forces of zygotic (diplophase) selection, gametic (haplophase) selection, and nonMendelian segregation, in which all the parameters are assumed to vary randomly from generation to generation. Interest in the model actually centers on two special cases: the case of zygotic selection combined with gametic selection, on the one hand, and the case of zygotic selection combined with non-Mendelian segregation, on the other hand. Scudo (1967) has studied a deterministic model of zygotic selection combined with gametic selection; Hiraizumi, Sandler and Crow (1960) have analyzed a deterministic model of zygotic selection combined with non-Mendelian segregation. Although there is no experimental precedent for combining the two sorts of models, they are conveniently combined for theoretical purposes and the results of the analysis specialized to whichever particular model may be of interest. The whole question of the effects of stochastic variation in selection intensity has come under close scrutiny recently, especially in diploid, one-locus, diallelic models with variation in the intensity of zygotic selection (Hartl and Cook 1973, 1976; Gillespie 1973; Karlin and Lieberman 1974, 1975; Karlin and Levikson 1974; Levikson 1974; Levikson and Karlin 1975; Cook and Hartl 1975; Norman * Work supported by N. I. H. grants GM21732 and GM21623. The author is supported by Research Career Award GM2301.
264
D.L. Hartl:
1975). A principal motivation has been to determine whether the nuisance caused by variation in selection intensity renders the results of classical constant fitness models totally inappropriate in the sense that variation in selection intensity generates qualitatively new phenomena that are not encompassed by the deterministic theory9 This is an important question because selective values in nature are almost certainly not constants, at least not for any appreciable length of time. (How much selective values may vary is, of course, unknown; autocorrelations between the fitnesses in successive generations are also unknown.) The most optimistic hope for the classical deterministic models has been that their results, suitably re-interpreted, would be valid for the limiting behaviour of stochastic models9 So far, in infinite population models, this has been the case, although some care must be exercised in the re-interpretation. It seems likely, however, that the deterministic results will not carry over straightforwardly in the cases of finite population models9 The results of the present heuristic analysis have been briefly summarized elsewhere (Hartl 1975), and the details of the analysis will be presented here. Throughout the analysis I will use the methods employed by Levikson (1974) and Levikson and Karlin (1975), methods which are originally due largely to Feller (1954)9 The Deterministic Model
Consider a locus having two alleles, A and a, in a large, randomly-mating, diploid population. Let the frequency of A among newly-formed zygotes be denoted p, and let q = 1 - p. The frequencies of A A, A a and a a genotypes among newly-formed zygotes will therefore be p2, 2 p q and q2, respectively9 Suppose the egg-to-adult viabilities of AA, A a and a a zygotes are w11, w12 and wz2, where wis>O. Among adults, the frequencies of AA, Aa and aa will be p2w11/~,, 2 p q w 12/# and q2 w22/w, respectively, where # = p2 Wll + 2 p q w 12 + q2 w22. Suppose now that the fraction of functional A-bearing gametes produced by an A a heterozygote is k, where 0 < k < 1. (This is the main assumption in the model of Hiraizumi et al. 1960). Then, at the moment of gamete formation, the frequency of functional A-bearing gametes will be (p2 Wl1 "at-2 p q k Wl2)/W.Suppose that these gametes are subject to further selection, and that A- and a-bearing gametes survive with probabilities vl and v2, respectively9 (This is the main assumption in the model of Scudo 1967.) After the haplophase selection, let the frequency of functional A-bearing gametes be denoted by p'. This frequency will perforce be the frequency of A among newly-formed zygotes of the next generation and
\W12/')2,]
pr \
12 2 /
\ /)2 ,//
\Wl2J
It will be useful in the sequel to have an expression for A p = p' - p :
Stochastic Selection in Both Haplophase and Diplophase
w11/)1
2(1-k)
PLwI2/)2
Pq Ap=
+q v2
PZ(W w12 ~ l /v2 ) l /) + p q [
265
w12
.
(1)
/)2 + 2 ( l - k ) ] +qz(w22]kw12 J
This deterministic model has a fixed point,/3, at 2 k W12 /)1 --W22 /)2
/3= 2 k w12/)1 + 2 (1 -k)
WI.2 I)2 - - W l l /)1 --W22 /)2 '
2 (1 - k) v z , or if (b) w22 > 2 k/)1 and 0 W12
then p--*l, whereas if (d) w2z > W12
2kvl /)2
and
wll
3. Thus the process converges weakly to the diffusion process with drift and diffusion coefficients in (3) and (4).
Analysis of the Stochastic Model A biologically pathological case arises if there exists a value of p, call it ~, such that # (~) = o-z (~) = 0. Inspection of (3) and (4) shows that this can happen only if the random variables are linearly dependent in such a way that $ 4 - $ 2 = ~ ( $ 1 - $ 2 - $ 3 + $ 4 ) , where e is some constant, 0 < c ~ < l . In this case, /3=e. Assume (for the moment) that the random variables are not linearly dependent in this manner. Then the diffusion process has only two boundaries, p = 0 and p = 1, and the limiting behaviour of the process is largely determined by the stochastic stability of these boundaries. Local stochastic stability of a boundary can be investigated using the function S (p)= ~f,os (z) dz, where po=p (0) and s (p)=exp ( - 2 ~o # (x)/0-2 (x) dx). The integral in S (p) converges near a boundary if and only if the boundary is locally stochastically stable. Partial fraction expansion of #(p)/o-2 (p) and
M2-M4-R24+
V4
integration yields s(p)=p-2a(1-p)2Bg(p) where A V2_2R24+V4 , M 1 - M 3 +R13 - V1 B, and where 9 (p) is a function that is continuous in
V1-2R13+V3
[0, 1] and bounded away from zero. Therefore, around p = 0, S (p) converges iff -2A>-1, which is to say iff M 2 - V2/2 - 1, which is to say M1 - V1/2 > M3 - V3/2; this is a necessary and sufficient condition for p-- 1 to be locally stochastically stable. A key theorem in working out the global behavior of the diffusion process is the following: Given p (0) = Po and two constants, ~ and 6, with 0 < e < Po -< 1 - 6 < 1, then
s (po)- s (~)
Prob {p (t) reaches 1 - 6 before reachin9 e I Po, 0 < e_< Po < 1 - 6 < 1} - S (1 - 6 ) - S (~) = L (e, 6) where S (p) is as defined above. We are interested in the behavior of L (e, 6) as ~, 6 4 0 . Consider first the case M 2 - V z / 2 < M 4 - V 4 / 2 and M 1 - V~/2 M 4 - V4/2 and M ~ - V I / 2 > M 3 - 1/3/2, the case in which p = 0 is stochastically unstable and p = 1 is locally stochastically stable. Then S ( 1 - 6) remains bounded as 6--.0 whereas S (e)~oo as e~0. Therefore L (e, 6)--,1 as e, 6--*0, which implies p ( 0 - ' 1 almost surely. Thus, in this case, p = 1 is globally stochastically stable. The third case to consider is M 2 - V2/2 > M 4 - V4/2 and M 1 - V1/2 < M 3 - V3/2, in which case both p = 0 and p = 1 are stochastically unstable. In this case, a necessary and sufficient condition for the gene frequency to converge to a unique stationary distribution is that 1/0-2 (p) s (p) be integrable on [0, 1]. Now, 1/0- 2 (p) s ( p ) = p -2+2A ( 1 - - p ) - Z - 2 ~ g * (p), where A and B are as defined above and g* (P) is continuous on [0, 1] and bounded away from zero. Obviously 1/o-2 (p)S (/9) is integrable on the open interval (0, 1). The function is integrable near 0 iff - 2 + 2 A > - 1, which is to say M 2 - - V 2 / 2 > M 4 - V 4 / 2 . The function is integrable near 1 iff - 2 - 2 B > - 1 , which is to say M 1 - V 1 / 2 < M 3 - V3/2. These are the identical conditions that make p = 0 and p = 1 stochastically unstable. Therefore, in this case, the gene frequency always converges to a unique stationary distribution. The final case to consider is M 2 - 1/2/2 < M 4 - V4/2 and m l - V~/2 > M 3 - V3/2 , in which case both p = 0 and p = 1 are locally stochastically stable. In this case both S ( 1 - 6 ) and S(e) remain finite as 5, 6--.0. Thus L (e, 6) remains finite as e, 6"~0. In this case, as a consequence, some sample paths will converge to p = 0 and the rest will converge to p = 1, with the probability of convergence to a given boundary being determined by the initial gene frequency. In summary, presenting the cases in the same order as the analogous cases in the deterministic model, we have if (a) M 2 - V2/2> M 4 - 1/4/2 and M 1 - VI/2 < M 3 - V3/2, then p converges to a unique stationary distribution; if (b) M 2 - V2/2 < M , - V4/2 and M 1 - I/1/2 > M 3 - V3/2 , then, with some probability dependent on the initial gene frequency, p will converge to 1, and with the remaining probability it will converge to 0; if (c) M 2 - V2/2 > M 4 - V4/2 and MI - V~/2 > M 3 - 1/3/2, then p ~ 1 almost surely; and if (d) M 2 - V2/2 < M 4 - I/-J2 and m 1 - V1/2 < M 3 - V3/2 , then p-~0 almost surely. Returning now to the one case excluded at the beginning of this section, assume that $ 4 - Sz = c~( S I - $ 2 - $3 + $4) where e is a constant and 0 < c~< 1. Since the case is of little interest biologically, local stochastic stability analysis will suffice. Partial fraction expansion of #(p)/0-2(p) by tedious algebra reveals that s ( p ) = p -2A ( 1 - p ) 2 B ( p - ~ ) - 2 c where A and B are as previously defined and C-M1-Mz-M3+M4-R14+R23 R 1 4 - R 1 2 - - R 3 4 q-R23 Journ. Math. Biol. 3 / 3 ~
. Thus the conditions for p = 0 and p = 1 to 19
268
D.L. Hartl:
be locally stochastically stable are as stated above, and the point p--~ is locally stochastically stable if and only if C < 1/2.
Discussion
To the order of the diffusion approximation, E In(1 + St) = M i - V~/2.Thus the four cases of interest may be restated in termsofthe original random variables as follows: if (a) Eln (W2g/W12)Eln [2 (1 - k ) v2/vi], then p-o0 or p ~ l with probabilities depending on the initial gene frequency; if (c) Eln (w22/w12)< Eln (2 k vi/v2) and Eln (wii/Wig)>Eln [2 (1 -k) v2/vl], then p--,1 almost surely; and if (d) Eln (w22/w12)>Eln (2 k vl/v2) and Eln (wli/wl2)
