带有捕食作用的随机CKTW模型维持物种多样性的研究

Maintenance of Species Diversity by Stochastic CKTW Models with Predation

摘要: 在海洋生态系统中，消耗细菌种群数量的因素包括原生动物的非选择性捕食与噬菌体的裂解。本文在具有共同进化机制的随机CKTW模型基础上，引入原生动物非选择性捕食作用，建立基于个体水平的带有非选择性捕食作用的随机CKTW模型。利用随机模拟算法对种群数量有限的随机模型进行数值模拟，捕捉到丧失物种多样性的三种灭绝路径与维持物种多样性的现象，并通过模拟实验计算出这些现象发生的概率。这里使用的随机模拟算法为Gillespie算法。最后以香农熵作为物种多样性的评价指标，对两类随机模型的物种多样性进行对比分析，说明共同进化可以在一定程度上维持物种多样性。

Abstract: In Marine ecosystems, factors that deplete bacterial populations include non-selective predation by protozoa and phage lysis. In this paper, based on the stochastic CKTW model with coevolutionary mechanism, we introduce the non-selective predation of protozoa and establish the stochastic CKTW model with non-selective predation based on individual level. The stochastic simulation algorithm was used to simulate the stochastic model with limited population to capture the three extinction paths that lost species diversity and the phenomena that maintained species diversity, and the probability of occurrence is calculated by a large number of simulation experiments. The stochastic simulation algorithm used here is Gillespie algorithm. Finally, Shannon entropy was used as the evaluation index of species diversity, and the species diversity of the two stochastic models was compared and analyzed, indicating that coevolution could maintain species diversity to a certain extent.

1. 引言

在生态系统中，物种多样性是一个重要的研究课题。Hutchinson等以海洋生态系统为例提出了浮游生物的悖论 [1] 。“Kill the Winner”(KTW)动力学可作为悖论的一种解释 [2] [3] [4] [5] 。KTW模型被Thingstad等最初表示为空间分布均匀、种群数量庞大的确定性Lotka-Volterra型方程 [6] [7] 。但是当种群数量较少时，确定性方程无法捕捉到因内部随机性导致物种灭绝的结果 [8] [9] [10] 。为找到维持随机KTW模型物种多样性的机制，Goldenfeld等根据细菌种群与噬菌体种群出现的变异现象，提出具有共同进化机制的随机CKTW模型，并说明共同进化可以在一定程度上维持系统的物种多样性 [10] - [15] 。

在具有共同进化机制的随机CKTW模型中，Goldenfeld等考虑了噬菌体种群对细菌种群特异性侵染的作用。在海洋生态系统中，控制细菌种群数量的因素还包括原生动物非选择性捕食作用 [16] [17] 。我们在Goldenfeld等工作的启发下，引入原生动物非选择性捕食作用，建立基于个体水平的带有非选择性捕食作用的随机CKTW模型。利用随机模拟算法对该模型进行多次数值模拟，捕捉到丧失物种多样性的三种灭绝路径与维持物种多样性的现象，并通过给定不同变异率进行大量数值模拟实验，计算出这些现象发生的概率 [18] [19] 。在数值模拟的过程中，观察到变异率值越高，物种多样性的现象越明显。最后，为了描述共同进化对模型物种多样性的维持程度，我们引入香农熵作为评价指标，利用随机模拟算法模拟出两类随机模型的变异率与香农熵的关系并进行对比分析，说明共同进化可以在一定程度上维持带有非选择性捕食作用的随机CKTW模型的物种多样性。

2. 模型建立

本文中我们考虑基于个体水平的带有非选择性捕食作用的随机CKTW模型：

$\begin{array}{l}{X}_i\stackrel{b_i}{\to }2X_i{X}_i+X_j\stackrel{e_{ij}}{\to }{X}_j\\ X_i+Y_i\stackrel{p_i}{\to }\left({\beta }_i+1\right)Y_i{Y}_i\stackrel{d_i}{\to }\varphi \\ X_i+Z\stackrel{g}{\to }2ZZ\stackrel{f}{\to }\varphi \\ X_i\stackrel{{\mu }_1/2}{\to }{X}_{i\pm 1}{Y}_i\stackrel{{\mu }_2/2}{\to }{Y}_{i\pm 1}\end{array}$ (1)

其中，所有的速率都为正数， $X_i$ 代表第i类细菌个体， $Y_i$ 代表第i类噬菌体个体，Z代表原生动物个体， $i,j=1,2,\cdots ,m$ 为角标。 $b_i$ 为第i类细菌种群的特异性生长速率，细菌种群之间以 $e_{ij}$ 的强度竞争隐式资源。第i类噬菌体以 $p_i$ 的速率感染特异性宿主细菌，裂解细菌后的病毒释放量为 ${\beta }_i$ ，并且第i类噬菌体种群以 $d_i$ 的速率衰减。g代表原生动物对细菌种群的非选择性捕食速率，f代表原生动物的衰减速率， ${\mu }_1/2$ 与 ${\mu }_2/2$ 分别代表细菌种群与噬菌体种群的变异速率。

接下来，考虑不带有共同进化机制的平均场速率方程：

$\begin{array}{l}{\stackrel{\dot{}}{B}}_i=b_i{B}_i-{\displaystyle \underset{j=1}{\overset{m}{\sum }}{e}_{ij}{B}_i{B}_j}-p_i{V}_i{B}_i-gP{B}_i\\ {\stackrel{\dot{}}{V}}_i={\beta }_i{p}_i{V}_i{B}_i-d_i{V}_i\\ \stackrel{\dot{}}{P}=gP{\displaystyle \underset{i=1}{\overset{m}{\sum }}{B}_i}-fP\end{array}$ (2)

其中， $B_i$ 与 $V_i$ 分别代表第i类细菌的种群密度和第i类噬菌体的种群密度，P代表原生动物的种群密度，参数与模型(1)相同。令细菌种群的特异性生长速率满足 [7] [17] ：

$b_1>b_2>b_3>\cdots >b_m$ (3)

经计算(取 $e_{ij}\equiv e$ )，模型(2)的平衡点为：

$\begin{array}{l}{B}^{*}=\left(\frac{d_1}{{\beta }_1{p}_1},\frac{d_2}{{\beta }_2{p}_2},\cdots ,\frac{d_{m-1}}{{\beta }_{m-1}{p}_{m-1}},\frac{f}{g}-{\displaystyle \underset{i=1}{\overset{m-1}{\sum }}\frac{d_i}{{\beta }_i{p}_i}}\right)\\ V^{*}=\left(\frac{1}{p_1}\left(b_1-b_m\right),\frac{1}{p_2}\left(b_2-b_m\right),\cdots ,\frac{1}{p_{m-1}}\left(b_{m-1}-b_m\right),0\right)\\ P^{*}=\frac{1}{g}\left(b_m-e\frac{f}{g}\right)\end{array}$ (4)

我们利用随机模拟算法对模型(1)进行数值模拟，将(4)中的平衡点作为初值点，研究当种群数量有限时，共同进化对模型(1)的影响。

3. 数值模拟

模型假设：

假设模型(1)有十类细菌种群、十类噬菌体种群与原生动物种群。十类细菌种群的特异性生长速率满足(3)式，取 $b=\left(1.2,1.15,1.1,1.05,1,0.95,0.9,0.85,0.8,0.75\right)$ ，细菌种群之间以 $e_{ij}=0.1$ 的强度竞争隐式资源。第i类噬菌体以 $p_i=0.1$ 的速率侵染第i类细菌，裂解细菌后的病毒释放量为 ${\beta }_i=10$ ，并且第i类噬菌体种群以 $d_i=0.5$ 的速率衰减。原生动物种群以 $g=0.23$ 的速率非选择性捕食细菌种群，并且原生动物种群以 $f=1$ 的速率衰减， $i,j=1,2,\cdots ,10$ 。

由上述参数可得平衡点：

$\begin{array}{l}{B}^{*}=\left(0.025,0.025,0.025,0.025,0.025,0.025,0.025,0.025,0.025,0.00\text{5}\right)\\ V^{*}=\left(0.225,0.2,0.175,0.15,0.125,0.1,0.075,0.05,0.025,0\right)\\ P^{*}=0.727\end{array}$ (5)

我们以(5)中的平衡点作为各类种群数值模拟的初值点，系统规模选为 $C=1000$ ，利用随机模拟算法对模型(1)进行数值模拟。由于算法的随机性，每次数值模拟得到的结果都不相同，在数值模拟一至数值模拟三中，将程序的终止时间设定为 $T_{\max }=200$ (若各类种群提前灭绝，直接结束程序)。

数值模拟一：

当变异率 ${\mu }_1={\mu }_2=0$ 时，共同进化作用无法发挥，种群数量有限的系统必然丧失物种多样性。利用随机模拟算法对模型(1)进行10⁵次数值模拟，得到丧失物种多样性的两种路径：

(a) 十类细菌种群与原生动物种群 (b) 十类噬菌体种群 (c) 十类细菌种群与原生动物种群 (d) 十类噬菌体种群

Figure 1. Two paths to lose species diversity

图1. 丧失物种多样性的两种路径

图1中图1(a)与图1(b)为丧失物种多样性的第一种路径，前期因噬菌体的特异性侵染作用，不同细菌种群的“赢家”交替出现，后期因大量细菌种群提前灭绝，原生动物种群与噬菌体种群因无法满足生存条件而全部灭绝。图1(c)与图1(d)为丧失物种多样性的第二种路径，细菌种群因捕食和竞争作用迅速衰减，使噬菌体种群无法满足生存条件而提前灭绝，最终只剩下原生动物种群对第十类细菌种群(不唯一)的非选择性捕食。

当物种的变异率 ${\mu }_1={\mu }_2=0$ 时，模型(1)丧失物种多样性是必然结果，我们利用随机模拟算法对模型(1)进行10⁵次数值模拟，计算两种路径发生的概率：

Table 1. The probability of two paths of lose species diversity

表1. 丧失物种多样性的两种路径概率

通过对模型(1)进行数值模拟，我们可以得到当物种的变异率 ${\mu }_1={\mu }_2=0$ 时，模型(1)必然丧失物种多样性，并且丧失物种多样性的第二种路径出现的概率较大。

数值模拟二：

当变异率 ${\mu }_1={\mu }_2=0.015$ 时，由于此时的变异率较低，利用随机模拟算法对模型(1)进行10⁵次数值模拟，可以得到丧失物种多样性的三种路径：

(a) 第一种路径的细菌与原生动物种群 (b) 第二种路径的细菌与原生动物种群 (c) 第三种路径的细菌与原生动物种群 (d) 第一种路径的噬菌体种群 (e) 第二种路径的噬菌体种群 (f) 第三种路径的噬菌体种群

Figure 2. Three paths to lose species diversity

图2. 丧失物种多样性的三种路径

图2中图2(a)与图2(d)为丧失物种多样性的第一种路径，十类细菌种群因随机捕食作用接连发生灭绝，噬菌体与原生动物种群因缺少猎物，继而发生灭绝的现象。图2(b)与图2(e)为丧失物种多样性的第二种路径，噬菌体种群因无法达到生存条件而提前灭绝。由于存在变异现象，细菌种群可以变异为相邻类别的细菌种群，使原生动物种群与几类细菌种群局部共存。图2(c)与图2(f)为丧失物种多样性的第三种路径，原生动物种群由于内部随机性迅速走向灭绝，细菌种群与噬菌体种群局部共存，并且细菌种群会受到噬菌体种群特异性捕食的影响，导致细菌的“赢家”交替产生。

同时，利用随机模拟算法对模型(1)进行数值模拟，也可以得到维持物种多样性的现象：

(a) 十类细菌种群与原生动物种群 ${\mu }_1={\mu }_2=0.015$ (b) 十类噬菌体种群 ${\mu }_1={\mu }_2=0.015$

Figure 3. Maintain species diversity

图3. 维持物种多样性

图3展示了维持物种多样性的现象。受到噬菌体种群特异性侵染作用的影响，细菌种群的“赢家”交替出现，没有单一的细菌种群占据主导地位，并且原生动物种群与中的噬菌体种群数量也受到细菌种群数量的影响，但没有出现丧失物种多样性的现象。

当物种的变异率 ${\mu }_1={\mu }_2=0.015$ 时，利用随机模拟算法对模型(1)进行10⁵次数值模拟，计算三种路径发生的概率与维持物种多样性的概率：

Table 2. The probability of the three paths to lose species diversity and the probability of maintaining species diversity

表2. 丧失物种多样性的三种路径与维持物种多样性的概率

通过利用随机模拟算法对模型(1)进行数值模拟，当物种的变异率 ${\mu }_1={\mu }_2=0.015$ 时，模型(1)可以出现维持物种多样性的现象。但由于此时变异率较低，共同进化作用发挥的不明显，丧失物种多样性现象出现的概率较大。

数值模拟三：

随着变异率的增加，模型(1)维持物种多样性的现象会逐渐明显。我们以变异率 ${\mu }_1={\mu }_2=0.1$ 与 ${\mu }_1={\mu }_2=1$ 为例，利用随机模拟算法对模型(1)分别进行10⁵次数值模拟可以得到维持物种多样性的现象：

(a) 十类细菌种群与原生动物种群 ${\mu }_1={\mu }_2=0.1$ (b) 十类噬菌体种群 ${\mu }_1={\mu }_2=0.1$ (c) 十类细菌种群与原生动物种群 ${\mu }_1={\mu }_2=1$ (d) 十类噬菌体种群 ${\mu }_1={\mu }_2=1$

Figure 4. Maintain species diversity with different variation rate

图4. 不同变异率下维持物种多样性的现象

图4中的模拟结果展示了在不同变异率下，维持物种多样性的现象。图4(a)与图4(b)为第一组模拟结果，此时的变异率相对较低。图4(a)中十类细菌种群在十类噬菌体种群的特异性侵染下，呈现出“赢家”交替产生的现象，没有单一细菌种群占据主导地位，所以物种都共存。并且不同细菌种群之间与不同噬菌体种群之间通过共同进化实现变异流通，变异的速度可以将濒临灭绝的种群恢复。图4(c)与图4(d)为第二组模拟结果，此时的变异率较高，不同细菌种群之间与不同噬菌体种群之间快速流通，变异的速度足以补偿较大的波动，时间序列走势趋于高度同步。

当变异率为 ${\mu }_1={\mu }_2=0.1$ 与 ${\mu }_1={\mu }_2=1$ 时，利用随机模拟算法对模型(1)在不同变异率下进行10⁵次数值模拟，计算三种路径发生的概率与维持物种多样性的概率：

Table 3. The probability of the three paths to lose species diversity and the probability of maintaining species diversity

表3. 丧失物种多样性的三种路径与维持物种多样性的概率

通过利用随机模拟算法对模型(1)进行数值模拟，随着变异率的增加，物种多样性的现象逐渐明显。

结合表1~3可以得出：随着变异率的增加，模型(1)维持物种多样性的现象逐渐明显。虽然存在丧失物种多样性的现象，但是发生这些现象的概率明显降低。说明共同进化在一定程度上可以维持模型(1)的物种多样性。

4. 物种多样性的评价指标

为了更直观地观察共同进化在模型(1)中维持物种多样性的程度，我们以香农熵作为物种多样性的评价指标。香农熵的公式为：

$S=-{\displaystyle \underset{i=1}{\overset{m}{\sum }}{f}_i\ln f_i}$ (6)

其中， $f_i$ 是第i类细菌菌株在全部细菌菌株中的占比。当所有细菌菌株都共存，并且处于相对稳定的状态时，香农熵为当前计算值；当只有一类细菌菌株存在时，香农熵 $S=0$ ；当细菌菌株全部灭绝时，香农熵记为 $S=-1$ 。在模型(1)中，虽然物种灭绝是必然结果，但我们只关注有限时间内的种群行为 [10] 。

数值模拟四：

将程序的终止时间设定为 $T_{\max }=400$ ，对每个变异率值循环1000次实验，求得香农熵的平均值。我们会得到模型(1)的变异率( ${\mu }_1={\mu }_2=\mu$ )与香农熵的关系如图5(a)所示。这里，我们与不带有非选择性捕食作用的随机CKTW模型的香农熵结果图5(b)进行对比 [10] 。

(a) 带有捕食作用的香农熵结果图 (b) 不带有捕食作用的香农熵结果图

Figure 5. The relationship between variation rate and Shannon entropy

图5. 变异率与香农熵的关系图

图5中图5(a)与图5(b)左上角的小图为标准坐标系下变异率与香农熵的关系，都可以看出共同进化在变异率较低时就开始发挥作用，但细节变化不是特别明显。于是我们转为主图所示的指数坐标系下变异率与香农熵的关系。从两幅图的主图中都可以看出，变异率将香农熵分为了三个不同区域，并且随着变异率的增加，香农熵值明显增大。

在图5(a)中，通过利用随机模拟算法对模型(1)进行数值模拟，可以看出变异率将香农熵分为了三个区域。在区域I中，由于变异率较低，细菌种群几乎是完全灭绝的，并且共同进化还没有发挥作用，物种多样性就已经丧失；在区域II中，所有物种对变异率都比较敏感，随着变异率的增加，共同进化机制驱动种群之间的流通，足以补偿物种的波动，使得物种多样性的现象逐渐明显。在此阶段主要发生的现象为原生动物种群与噬菌体种群对细菌种群的捕食与侵染作用，使细菌种群的“赢家”不断交替变化，没有任何一类细菌种群占据主导地位，物种多样性的现象逐渐明显；在区域III中，变异率相对较高，香农熵值都接近最大值，细菌种群、噬菌体种群、原生动物种群大量共存，物种多样性的现象十分明显。在图5(b)中，由于只考虑细菌种群与噬菌体种群之间的特异性侵染作用，因此，不存在丧失物种多样性的第三种路径。在区域I中，由于变异率较低，细菌种群几乎是完全灭绝的，并且共同进化还没有发挥作用，物种多样性就已经丧失；在区域II中，所有物种对变异率都比较敏感，随着变异率的增加，共同进化机制驱动种群之间的流通，足以补偿物种的波动，使得物种多样性的现象逐渐明显；在区域III中，变异率相对较高，香农熵值都接近最大值，细菌种群与噬菌体种群大量共存，物种多样性的现象十分明显。

通过图5(a)与图5(b)的对比，不同之处主要在区域II中，图5(a)中的香农熵结果明显低于图5(b)中的香农熵结果。这是由于模型(1)中引入原生动物非特异性捕食作用，当变异率较低时，丧失物种多样性的三种灭绝路径会占据一定比例，使得模型(1)丧失物种多样性的现象较为明显。但整体来说，在模型(1)中共同进化可以在一定程度上维持物种多样性。

5. 结论

本文在具有共同进化机制的随机CKTW模型基础上，引入原生动物非选择性捕食作用，建立基于个体水平的带有非选择性捕食作用的随机CKTW模型(1)。利用随机模拟算法对模型(1)进行数值模拟，捕捉到丧失物种多样性的三种灭绝路径与维持物种多样性的现象，并通过给定不同变异率进行十万次数值模拟实验，计算出这些现象发生的概率。根据数值模拟一可以得出：当变异率为0时，模型(1)必然丧失物种多样性。根据数值模拟二可以得出：当变异率较低时，模型(1)可以出现维持物种多样性的现象，但此时共同进化的作用发挥不明显，丧失物种多样性的现象占主导。根据数值模拟三可以得出：随着变异率的增加，细菌种群之间与噬菌体种群之间的时间系列走势趋于同步，物种共存的现象越明显。最后以香农熵作为物种多样性的评价指标，通过数值模拟四可以看出，引入原生动物非选择性捕食作用后，香农熵的计算结果会随着变异率的增加而增大，说明共同进化可以在一定程度上维持物种多样性。

基金项目

国家自然科学基金项目(批准号11901052)。