Evolutionary branching in a finite population: deterministic branching vs. stochastic branching

Abstract

Adaptive dynamics formalism demonstrates that, in a constant environment, a continuous trait may first converge to a singular point followed by spontaneous transition from a unimodal trait distribution into a bimodal one, which is called "evolutionary branching." Most previous analyses of evolutionary branching have been conducted in an infinitely large population. Here, we study the effect of stochasticity caused by the finiteness of the population size on evolutionary branching. By analyzing the dynamics of trait variance, we obtain the condition for evolutionary branching as the one under which trait variance explodes. Genetic drift reduces the trait variance and causes stochastic fluctuation. In a very small population, evolutionary branching does not occur. In larger populations, evolutionary branching may occur, but it occurs in two different manners: in deterministic branching, branching occurs quickly when the population reaches the singular point, while in stochastic branching, the population stays at singularity for a period before branching out. The conditions for these cases and the mean branching-out times are calculated in terms of population size, mutational effects, and selection intensity and are confirmed by direct computer simulations of the individual-based model.

Figure 1

Evolutionary dynamics. (A) A simulation run with N = 8000. (B and C) Two different simulation runs with N = 1000. (D) A simulation run with N = 200. Parameters: $\mu =0.01,\hspace{0.17em}\sigma =0.02,\hspace{0.17em}(b_1,b_2,c_1,c_2)=(6.0,-1.4,\hspace{0.17em}4.56,-1.6)$.

Figure 2

A schematic illustration of phenotype distribution $\varphi (z)$ and fitness function $w(z)$ when the mean trait value $\stackrel{-}{z}$ has reached the singular strategy z*. The trait variance is $\text{Var}[\varphi (z)]=m_2$ so the standard deviation is $\sqrt{m_2}$.

Figure 3

“Force” plotted against the trait variance based on Equation 10. Dotted, solid, and thick solid curves represent stabilizing selection ($w_2$ = −0.2), weak disruptive selection ($w_2$= 0.2), and strong disruptive selection ($w_2$ = 0.6), respectively. They also correspond to cases i, ii, and iii in the main text. Parameters: N = 400, $\mu =0.01,\hspace{0.17em}\sigma =0.02$.

Figure 4

Stable (solid) and unstable (dashed) branches of the trait variance $m_2$ with disruptive selection intensity $w_2$ as a bifurcation parameter. Arrows indicate the direction of systematic change in $m_2$. Parameters: N = 1000, $\mu =0.01,\hspace{0.17em}\sigma =0.02$.

Figure 5

Stable (solid) and unstable (dashed) branches of the trait variance $m_2$ with population size N as a bifurcation parameter. A Log plot is shown. Arrows indicate the direction of systematic change in $m_2$. Parameters: $\mu =0.01,\hspace{0.17em}\sigma =0.02,\hspace{0.17em}(b_1,b_2,c_1,c_2)=(6.0,-1.4,\hspace{0.17em}4.56,-1.6)$.

Figure 6

The stationary distribution $\tilde{p}(m_2)$ of the trait variance $m_2$ (Equation 14). A Log-Log plot is shown, and thus the probability density is much higher at the peak than it appears. The locally stable value of the variance according to Equation 11 is $m_2*\approx$ 0.00081, which is close to $N\mu {\sigma }^2$ = 0.0008 (dashed line). Note that a local peak lies at a value <0.0008 because of geometric Brownian motion. Parameters: N = 200, $\mu =0.01,\hspace{0.17em}\sigma =0.02,\hspace{0.17em}(b_1,b_2,c_1,c_2)=(6.0,-1.4,\hspace{0.17em}4.56,-1.6)$.

Figure 7

The stationary distribution of the trait variance $m_2$. Curves represent an analytic result (Equation 14), while bars are a histogram of realized trait variances in a simulation run for ${10}^6$ time steps. The dashed line represents $N\mu {\sigma }^2$ = 0.0008. In the simulation, the initial trait distribution was set as monomorphic at z = 0.2. The transition from this initial state to the stationary state was negligible as it took only ∼5000 time steps. Parameters: N = 200, $\mu =0.01,\hspace{0.17em}\sigma =0.02,\hspace{0.17em}(b_1,b_2,c_1,c_2)=(6.0,-1.4,\hspace{0.17em}4.56,-1.6)$.

Figure 8

A Log plot of waiting times for evolutionary branching for different population sizes (N). The thick and thin curves represent analytic results based on a SDE model (Equation 16) and those based on an ODE model (Equation 11), respectively. Each open circle represents the result of a single run in individual-based simulation. The first time the trait variance exceeded 0.1 is shown. In the simulation, the initial trait distribution was set as monomorphic at z = 0.2 and each run consisted of ${10}^6$ time steps. Simulation runs in which the variance did not explode are plotted at ${10}^6$. Time required to reach a singular point in the simulation was ∼5000 when N = 200 and ∼2000 when N = 1800, which is included to show open circles. Thus, open circles are expected to lie slightly above the curve. Parameters: $\mu =0.01,\hspace{0.17em}\sigma =0.02,\hspace{0.17em}(b_1,b_2,c_1,c_2)=(6.0,\hspace{0.17em}-1.4,\hspace{0.17em}4.56,\hspace{0.17em}-1.6)$.

Figure 9

The distribution of waiting times. Curves represent an analytic result (see Appendix E), while bars represent a histogram of the waiting time obtained in 1000 different runs of individual-based simulation. Parameters: N = 1000, $\mu =0.01,\hspace{0.17em}\sigma =0.02,\hspace{0.17em}(b_1,b_2,c_1,c_2)=(6.0,\hspace{0.17em}-1.4,\hspace{0.17em}4.56,\hspace{0.17em}-1.6)$.

Figure 10

A Log plot of waiting times for evolutionary branching for different population sizes and mutation parameters. The top panel ($\mu =0.0025,\sigma =0.02$) depicts a result for a smaller mutation rate than in Figure 8, while the bottom panel ($\mu =0.01,\sigma =0.01$) depicts a result for a smaller mutation step size. As $\mu {\sigma }^2={10}^{-6}$ is the same, the analytic prediction curves are the same. The two simulation results are also undistinguishable. Parameters: $(b_1,b_2,c_1,c_2)=(6.0,-1.4,\hspace{0.17em}4.56,-1.6)$.