第 3 章：相互作用种群的模型（Models for Interacting Populations）

物种相互作用时，每一种群的动态都会受到其他种群的影响。自然界中通常存在一整套相互作用的物种网络（即食物网 trophic web），结构非常复杂。本章考虑涉及两个或多个物种的系统，重点放在两物种情形。三类主要相互作用为：(i) 一个种群增长率下降、另一个上升——捕食者–猎物关系；(ii) 两个种群增长率都下降——竞争；(iii) 两个种群增长率都上升——互惠共生或共生。本章的数学方法与分析技术可以直接应用于第 6 章反应动力学，形式相同的方程描述的是化学浓度而非生物种群。

3.1 捕食者–猎物模型：Lotka–Volterra 系统（Predator–Prey Models: Lotka–Volterra Systems）

Volterra（1926）最早提出一个简单的捕食–被捕食模型，用以解释亚得里亚海某些鱼类捕获量的周期性振荡。设 $N(t)$ 为猎物、$P(t)$ 为捕食者在 $t$ 时刻的数量，模型为 $\frac{dN}{dt} = N(a - bP)$ $\frac{dP}{dt} = P(cN - d)$（式 3.1、3.2），其中 $a,b,c,d>0$ 为常数。模型假设：(i) 猎物在无捕食时按 Malthus 方式无限增长——即式 3.1 的 $aN$ 项；(ii) 捕食作用使猎物的每员增长率按与两群数量成正比的量减少，即 $-bNP$；(iii) 缺乏猎物时捕食者指数衰减，即式 3.2 的 $-dP$ 项；(iv) 猎物对捕食者增长率的贡献为 $cNP$，既与可用猎物量成正比，也与捕食者种群大小成正比。$NP$ 项可视为能量的种间转化——$bNP$ 取自猎物，$cNP$ 流向捕食者。这一模型有严重缺陷，但作为"提出正确问题"的起点以及更现实模型的出发点，仍具有不可替代的价值。

该模型又称 Lotka–Volterra 模型：Lotka（1920；1925）从一假想化学反应独立推导出同一方程，并指出该反应可在化学浓度上展现周期行为。这条动机在第 6 章还会再次出现。

分析的第一步是将系统无量纲化。令 $u(\tau) = \frac{cN(t)}{d},\quad v(\tau) = \frac{bP(t)}{a},\quad \tau = at,\quad \alpha = d/a$（式 3.3），则模型化为 $\frac{du}{d\tau} = u(1 - v),\qquad \frac{dv}{d\tau} = \alpha v(u - 1)$（式 3.4）。在 $(u,v)$ 相平面上给出 $\frac{dv}{du} = \frac{\alpha v(u-1)}{u(1-v)}$（式 3.5），其奇点在 $u=v=0$ 与 $u=v=1$。积分（3.5）可得相轨迹 $\alpha u + v - \ln u^{\alpha v} = H$（式 3.6），其中 $H > H_{\min} = 1 + \alpha$ 为常数，$H_{\min}$ 在 $u=v=1$ 处取到。对任意 $H > 1 + \alpha$，相轨迹是如图 3.1 所示的闭合曲线。闭合轨迹意味着式 3.4 关于 $\tau$ 有周期解，初始条件 $u(0),v(0)$ 决定常数 $H$ 继而确定具体轨迹；典型周期解 $u(\tau), v(\tau)$ 见图 3.2。由式 3.4 立即可见 $u$ 在 $v=1$ 处取极值，$v$ 在 $u=1$ 处取极值。

Lotka–Volterra 模型的一个主要缺陷是解不具备结构稳定性。若初始条件落在靠近 $u$ 或 $v$ 轴的 $H_4$ 轨迹上，任何微小扰动都会把解推到另一条并不处处接近 $H_4$ 的轨迹上，从而对小扰动产生极大影响（至少对振幅如此）。这一现象是任何具有如式 3.6 这种"闭合相轨迹的首积分"的系统所共有的，这类系统称为保守系统，其"守恒律"使得它们通常不适合作为真实种群相互作用的模型。

回到式 3.4，对奇点进行线性化可判断奇点类型与稳态稳定性，这套做法可推广到任意方程个数的等价系统。先考虑 $(u,v)=(0,0)$：令 $x,y$ 为 $(0,0)$ 的小扰动，仅保留线性项得式 3.7 的矩阵形式，其中 $\mathbf{A} = \begin{pmatrix} 1 & 0 \\ 0 & -\alpha \end{pmatrix}$，特征方程给出 $\lambda_1 = 1, \lambda_2 = -\alpha$。因为 $\lambda_1>0$，$u=v=0$ 线性不稳定；又因为 $\lambda_1>0$、$\lambda_2<0$，该奇点为鞍点。再考虑 $u=v=1$：令 $u=1+x$, $v=1+y$，$|x|,|y|$ 很小，式 3.4 化为式 3.8，矩阵 $\mathbf{A} = \begin{pmatrix} 0 & -1 \\ \alpha & 0 \end{pmatrix}$，特征值为 $\lambda_{1,2} = \pm i\sqrt{\alpha}$（式 3.9）。特征值为纯虚数，故 $u=v=1$ 是中心奇点，$\mathrm{Re}\,\lambda=0$，稳态为中性稳定。式 3.8 的解形如 $\mathbf{x}(\tau) = l e^{i\sqrt{\alpha}\tau} + m e^{-i\sqrt{\alpha}\tau}$，故在 $u=v=1$ 附近解关于 $\tau$ 周期为 $2\pi/\sqrt{\alpha}$。回到原始量纲由式 3.3 知周期为 $T = 2\pi\sqrt{a/d}$，即与猎物线性增长率 $a$ 与捕食者死亡率 $d$ 之比的平方根成正比。即便我们只在小扰动下观察这一关系，周期如何依赖内禀增长与死亡速率也已显现：例如，猎物增长率升高会使周期变长，捕食者死亡率下降亦然——这是否符合直觉？

在生态学语境下，线性化方程（式 3.7、3.8）的矩阵 $\mathbf{A}$ 称为群落矩阵（community matrix），其特征值 $\lambda$ 决定稳态稳定性：若 $\mathrm{Re}\,\lambda > 0$ 则不稳定；若所有 $\mathrm{Re}\,\lambda < 0$ 则稳定；临界情形 $\mathrm{Re}\,\lambda=0$ 称为中性稳定。

Lotka–Volterra 模型被多次尝试用于解释真实世界的振荡现象，但由于其结构不稳定性，这些尝试几乎都无法在定量上实用。它的价值在于引发"应问哪些正确问题"的讨论。一个特别有意思的例子是将模型用于解释哈德逊湾公司 1845 至 1930 年代山猫–雪兔毛皮捕获记录中的长期周期性数据，假设捕获数与种群总数成固定比例（这一假设未必准确）。Williamson（1996）的书是含周期或准周期行为的种群数据的好来源。图 3.3 显示了较为规律的周期振荡，相平面图近似闭合——正如时间周期行为所预期。Leigh（1968）用标准 Lotka–Volterra 模型尝试解释数据；Gilpin（1973）用修正版做同样尝试。仔细看图 3.3(c) 中相轨迹方向为顺时针，与图 3.1 中的相反；这反映了图 3.3(a)(b) 中山猫（捕食者）振荡先于雪兔发生——与图 3.2 所示捕食者–猎物情形相反。这意味着按数据表现雪兔在"吃"山猫，造成严重的解释困难。Gilpin（1973）猜想雪兔或许会将所携疾病传给山猫，于是把流行病效应纳入模型，数值结果确实重现了图 3.3(c) 的轨迹方向；可惜目前并无此类已知疾病。Gilpin 同时提出或许更对的解释——毛皮捕猎者本身就是"病原"：种群密度低的年份猎人会另谋生计，等到雪兔再次繁衍才回到猎捕线；山猫利润更高，猎人花更多时间捕山猫，这正好产生图 3.3(b)(c) 所示的现象。Schaffer（1984）进一步猜想山猫–雪兔数据可能是自然界中奇怪吸引子（即混沌行为）的证据。这一案例的教训是：仅给出能复现振荡的模型并不够，必须提供在生态与生物学层面经得起检验的解释。

3.2 复杂度与稳定性（Complexity and Stability）

为说明复杂度对稳定性的影响，简要考察一般化的 Lotka–Volterra 捕食者–猎物系统，其中存在 $k$ 个猎物物种与 $k$ 个捕食者物种，每个捕食者以不同强度捕食所有猎物。替代式 3.1、3.2 的方程组为 $\frac{dN_i}{dt} = N_i\left(a_i - \sum_{j=1}^k b_{ij}P_j\right),\quad i=1,\dots,k$

[ \frac{dP_i}{dt} = P_i\left(\sum_{j=1}^k c_{ij}N_j - d_i\right)$（式 3.10），所有 $a_i, b_{ij}, c_{ij}, d_i > 0$。平凡稳态为 $N_i=P_i=0$，对应的群落矩阵是对角矩阵 $\mathbf{A} = \mathrm{diag}(a_1,\dots,a_k,-d_1,\dots,-d_k) \] ，$2k$ 个特征值为 $\lambda_i = a_i > 0$（$i=1,\dots,k$）与 $\lambda_{k+i} = -d_i < 0$，故该稳态不稳定。

非平凡稳态为满足 $\sum_{j=1}^k b_{ij}P_j^* = a_i,\qquad \sum_{j=1}^k c_{ij}N_j^* = d_i,\qquad i=1,\dots,k$ 的列向量 $(N^*,P^*)$，或用矩阵记号写为 $\mathbf{B}P^* = \mathbf{a}$，$\mathbf{C}N^* = \mathbf{d}$（式 3.11），其中 $\mathbf{B}, \mathbf{C}$ 为 $k\times k$ 矩阵 $[b_{ij}]$、$[c_{ij}]$。

式 3.10 可写为 $\frac{dN}{dt} = N^T \cdot [\mathbf{a} - \mathbf{B}P],\qquad \frac{dP}{dt} = P^T \cdot [\mathbf{C}N - \mathbf{d}]$。在 $(N^*,P^*)$ 附近线性化：令 $N = N^*+u$, $P = P^*+v$，其中 $|u|,|v|$ 远小于 $|N^*|,|P^*|$，得 $\frac{du}{dt} \approx -N^{*T}\cdot \mathbf{B}v,\qquad \frac{dv}{dt} \approx P^{*T}\cdot \mathbf{C}u$ 由此得到式 3.12 的群落矩阵 $\mathbf{A} = \begin{pmatrix} \mathbf{0} & -N^{*T}\mathbf{B} \\ P^{*T}\mathbf{C} & \mathbf{0} \end{pmatrix}$，是 $2k\times 2k$ 的分块矩阵，对角块为零。$\mathbf{A}$ 的 $2k$ 个特征值 $\lambda_i$ 满足 $|\mathbf{A}-\lambda \mathbf{I}|=0$，由迹的线性性得 $\sum_{i=1}^{2k} \lambda_i = \mathrm{tr}\,\mathbf{A} = 0$（式 3.13）。由于 $\mathbf{A}$ 的元素为实数，复特征值成共轭对出现。结合式 3.13 有两种情况：所有特征值纯虚，或存在非零实部。若全部 $\mathrm{Re}\,\lambda_i = 0$，则 $(N^*,P^*)$ 中性稳定，与两物种情形一致；若存在 $\lambda_i$ 满足 $\mathrm{Re}\,\lambda_i \neq 0$，则因共轭对称必同时存在 $\mathrm{Re}\,\lambda > 0$ 与 $\mathrm{Re}\,\lambda < 0$，$(N^*,P^*)$ 不稳定。

分析表明，种群相互作用网络的复杂度引入了不稳定的可能性：若模型偶然给出仅含纯虚特征值（扰动关于时间为周期）的情况，那么群落矩阵任一参数稍有变化就会产生 $\mathrm{Re}\,\lambda \neq 0$ 的特征值从而引发不稳定。这一结论虽只针对式 3.12 形式的群落矩阵，但它指向了一个相当普遍且重要的结论：复杂度通常带来不稳定而非稳定。

3.3 更现实的捕食者–猎物模型（Realistic Predator–Prey Models）

尽管 Lotka–Volterra 模型并不现实，但它提示了简单捕食–被捕食相互作用能够产生种群周期行为。直观上这并不意外：若猎物增多，会刺激其捕食者增长；更多捕食者消耗更多猎物，使猎物下降；猎物减少后捕食者因食物短缺而下降；当捕食者低到一定程度，猎物得以再次增长，循环重启。视具体系统而定，这类振荡可放大或衰减，可形成稳定的极限环振荡，甚至出现混沌行为——但要出现混沌，至少需要 3 个相互作用物种，或模型含有延迟项。

极限环解是相平面上的闭合轨迹，但不属于像 Lotka–Volterra 那样一族连续闭合轨迹。稳定极限环是指对该轨迹的任何小扰动都会随时间衰减至零。两物种捕食 $(P)$–被捕食 $(N)$ 相互作用的典型极限环示意图见图 3.4。其存在条件在附录 A 给出。

式 3.1、3.2（以及一般化的式 3.10）一个不现实之处是：猎物在无捕食时无限增长。在式 3.1、3.2 当前的形式下，括号内的项为密度依赖的每员增长率。更现实地，这些增长率应同时依赖于猎物与捕食者的密度，如 $\frac{dN}{dt} = N F(N,P),\qquad \frac{dP}{dt} = P G(N,P)$（式 3.14），其中 $F,G$ 的具体形式取决于相互作用、物种等。

合理的第一步是让猎物种群在无捕食者时满足 logistic 增长（类似第 1 章式 1.2），即有最大承载量 $K$。一个更现实的猎物方程形如 $\frac{dN}{dt} = N F(N,P),\quad F(N,P) = r\left(1 - \frac{N}{K}\right) - P R(N)$（式 3.15），其中 $R(N)$ 为下文讨论的捕食项之一（图 3.5），$K$ 为 $P\equiv 0$ 时猎物的承载量。

捕食项即捕食者对猎物种群密度变化的功能性反应，通常表现出饱和效应。替代 Lotka–Volterra 模型中的 $bNP$，考虑 $PN R(N)$ 形式，其中 $N R(N)$ 在 $N$ 较大时饱和。一些例子（图 3.5(b)–(d)）为 $R(N) = \frac{A}{N+B},\qquad R(N) = \frac{AN}{N^2+B^2},\qquad R(N) = \frac{A(1-e^{-aN})}{N}$（式 3.16），其中 $A,B,a>0$。第二个（图 3.5(c)）类似于第 1 章式 1.6 的芽虫模型，也是蚜虫（Aphidicus zbeckistamicus）捕食的典型形式。式 3.16(b)(c) 在低密度时近似线性，$N$ 大时的饱和反映出捕食者能力的上限。

式 3.14 的第二个方程（捕食者方程）也应比 $G = -d + cN$（式 3.2）更现实。可能的形式为 $G(N,P) = k\left(1 - \frac{hP}{N}\right),\qquad G(N,P) = -d + e R(N)$（式 3.17），其中 $k,h,d,e>0$，$R(N)$ 如式 3.16。式 3.17 的第一个表明捕食者的承载量与猎物种群密度成正比。

式 3.14–3.17 仅是众多已提出并研究过的模型中的一些例子，都比经典 Lotka–Volterra 模型更现实。其他例子可参见 Nisbet 与 Gurney（1982）的书以及 Levin（1994）编辑的文集。

3.4 带极限环周期行为的捕食者–猎物模型分析：参数稳定域（Analysis of a Predator–Prey Model with Limit Cycle Periodic Behaviour: Parameter Domains of Stability）

作为分析此类现实两物种模型的范例，考察以下系统 $\frac{dN}{dt} = N\left[r\left(1 - \frac{N}{K}\right) - \frac{kP}{N+D}\right]$

[ \frac{dP}{dt} = P\left[s\left(1 - \frac{hP}{N}\right)\right] ] （式 3.18），其中 $r,K,k,D,s,h>0$ 共 6 个参数。

一如既往，无量纲化非常有用。无量纲化的方式不唯一，但通常与某个关键参数对齐。这里把 $N,P$ 表达为无捕食者承载量 $K$ 的比例。令 $u(\tau) = \frac{N(t)}{K},\quad v(\tau) = \frac{hP(t)}{K},\quad \tau = rt,\quad a = \frac{k}{hr},\quad b = \frac{s}{r},\quad d = \frac{D}{K}$（式 3.19），式 3.18 化为 $\frac{du}{d\tau} = u(1-u) - \frac{a u v}{u+d} = f(u,v)$

[ \frac{dv}{d\tau} = b v\left(1 - \frac{v}{u}\right) = g(u,v) ] （式 3.20），仅含 3 个无量纲参数 $a,b,d$。无量纲化以有意义的方式把参数分组，使参数个数减少；无量纲分组一般给出量纲参数的相对效应。例如 $b$ 是捕食者与猎物线性增长率之比，$b>1$ 与 $b<1$ 具有明确的生态意义：后者表示猎物比捕食者繁殖更快。

平衡态 $(u^*, v^*)$ 是 $du/d\tau = 0$、$dv/d\tau = 0$ 的解，即 $f(u^*,v^*)=0, g(u^*,v^*)=0$（式 3.21）。只看正解，可推出 $v^* = u^*,\qquad u^{*2} + (a+d-1)u^* - d = 0$ 唯一正解为 $u^* = \frac{(1-a-d) + \{(1-a-d)^2 + 4d\}^{1/2}}{2},\qquad v^* = u^*$（式 3.22）。

对稳态（式 3.20 相平面上的奇点）做稳定性分析需做线性化。令 $x(\tau) = u(\tau) - u^*,\quad y(\tau) = v(\tau) - v^*$（式 3.23）代入式 3.20 并保留线性项（利用式 3.21）得式 3.24 的群落矩阵 $\mathbf{A} = \begin{pmatrix} \frac{\partial f}{\partial u} & \frac{\partial f}{\partial v} \\ \frac{\partial g}{\partial u} & \frac{\partial g}{\partial v} \end{pmatrix}_{(u^*,v^*)} = \begin{pmatrix} u^*\left[\frac{a u^*}{(u^*+d)^2} - 1\right] & -\frac{a u^*}{u^*+d} \\ b & -b \end{pmatrix}$ 特征值 $\lambda$ 满足 $|\mathbf{A}-\lambda\mathbf{I}|=0$，即 $\lambda^2 - (\mathrm{tr}\,\mathbf{A})\lambda + \det \mathbf{A} = 0$（式 3.25）。

线性稳定要求 $\mathrm{Re}\,\lambda < 0$，充要条件（由式 3.25）为 $\mathrm{tr}\,\mathbf{A} < 0 \Rightarrow u^*\left[\frac{a u^*}{(u^*+d)^2} - 1\right] < b$

[ \det \mathbf{A} > 0 \Rightarrow 1 + \frac{a}{u^+d} - \frac{a u^}{(u^*+d)^2} > 0 ] （式 3.26）。

把 $u^*$ 用式 3.22 代入，稳定性条件即用 $a,b,d$ 表达，进而也可用原始参数 $r,K,k,D,s,h$ 表达。一般而言 $(a,b,d)$ 空间中有一域：当参数位于域内时 $(u^*,v^*)$ 稳定（$\mathrm{Re}\,\lambda < 0$）；位于域外则不稳定，至少违反式 3.26 之一。

利用 $u^*$ 的表达式（3.22）以及 $f(u^*,v^*)=0$、$v^*=u^*$，可证 $\det \mathbf{A} = \left[1 + \frac{a d}{(u^*+d)^2}\right] b u^* > 0$（式 3.27），对所有 $a,b,d>0$ 成立——故式 3.26 的第二式总成立。不稳定域仅由第一式 $\mathrm{tr}\,\mathbf{A} < 0$ 决定。代入 $u^*$（式 3.22）以及 $f=0$、$v^*=u^*$ 后，第一式变为 $b > \frac{[a - \{(1-a-d)^2 + 4d\}^{1/2}]\,[1 + a + d - \{(1-a-d)^2 + 4d\}^{1/2}]}{2a}$（式 3.28），定义 $(a,b,d)$ 参数空间中的一个三维曲面。

仅考虑 $a,b,d>0$。式 3.28 的第二个方括号是 $d$ 的单调递减函数且总为正。第一个方括号也是 $d$ 的单调递减函数，在 $d=0$ 处取极大值。因此式 3.28 在 $d=0$ 给出 $b > (2a-1)/a$（若 $0<a\leq 1$）或 $b > 1/a$（若 $a\geq 1$）。对 $0 < a < 1/2$ 及所有 $d>0$，稳定性条件（3.28）对任意 $b>0$ 都满足——即稳态 $(u^*,v^*)$ 对 $0<a<1/2$、所有 $b>0$、所有 $d>0$ 线性稳定。另一方面若 $a > 1/2$，$(a,b,d)$ 空间中 $(b>0, d>0)$ 存在使（3.28）不满足的域，从而第一式被违反，式 3.25 中至少一个特征值 $\mathrm{Re}\,\lambda > 0$，稳态 $(u^*,v^*)$ 对小扰动不稳定。分界曲面由式 3.28 给出，它与 $b=0$ 平面的交点 $d = d_m(a)$ 满足 $a = \{(1-a-d_m)^2 + 4d_m\}^{1/2}$，解为 $d_m(a) = d_{b=0} = (a^2 + 4a)^{1/2} - (1+a)$ $d_m(a)$ 是 $a$ 的单调递增函数，以 $d=1$ 为上界。注意 $a>1/2$ 时总有 $d < a$。图 3.6 示意了 $(a,b,d)$ 空间中的稳定性/不稳定性域。

$\mathrm{Re}\,\lambda < 0$ 时稳态稳定：若式 3.25 两特征值为实，则式 3.21 的奇点 $(u^*,v^*)$ 是式 3.20 相平面上的稳定结点；若为复，则为稳定螺线。$\mathrm{Re}\,\lambda > 0$ 时奇点为不稳定结点或螺线，此时为应用 Poincaré–Bendixson 定理以判定极限环振荡存在（见附录 A），需要在 $(u,v)$ 相平面正象限内找到一个简单闭合边界曲线，使该曲线上的相轨迹都指向域内。即对曲线上每一点，若 $\mathbf{n}$ 是该边界的外向单位法向，则需 $\mathbf{n}\cdot \left(\frac{du}{d\tau}, \frac{dv}{d\tau}\right) < 0$（参见式 3.20 后的讨论）。该不等式在边界上某点成立意味着速度向量 $(du/d\tau, dv/d\tau)$ 指向域内；直观上若一条解轨迹在边界上的速度都指向内，则轨迹不可能离开该域。

寻找有界集必要且总有信息量的做法是画出零倾线，即相平面上 $du/d\tau=0$ 与 $dv/d\tau=0$ 的曲线。式 3.20 的零倾线即 $f(u,v)=0$ 与 $g(u,v)=0$，示于图 3.7。向量分量 $(f,g)$（即 $(du/d\tau, dv/d\tau)$）的符号给出轨迹方向：若在某一域内 $f>0$，则 $u$ 在该处递增。在 $DE, EA, AB, BC$ 各段上，轨迹明显由 $f$ 与 $g$ 的符号决定而指向内。可以（虽然繁琐地）证明存在一条直线 $DC$ 使得 $\mathbf{n}\cdot(du/d\tau, dv/d\tau)<0$，其中 $\mathbf{n}$ 是 $DC$ 的单位法向。

如此当 $(u^*,v^*)$ 不稳定时即得到 Poincaré–Bendixson 定理所需的"有界集"。该定理说明：从不稳定稳态 $(u^*,v^*)$ 出来的任何轨迹都不能穿越约束边界 $ABCDEA$，因此必将演化为一条闭合极限环，定性上与图 3.4 类似。对模型（3.20）而言，图 3.8(a) 展示了这样的闭合轨迹，图 3.8(b) 给出对应的种群时间变化。该图所用参数下稳态在相平面上是不稳定结点（两个特征值均为正实数）；任何对极限环的扰动都迅速衰减。

该模型与大多数容许极限环行为的模型一样，随参数变化展示分岔性质，不过复杂度不及第 2、5 章的离散模型与第 1 章的延迟模型。立即可从图 3.6 看出：取定 $a>1/2$，存在有限的不稳定域（图 3.9）；固定 $0<d<d_m$，对应线 $DEF$。让 $b$ 从 $D$ 处的值连续减小，越过分岔线于 $E$ 时，稳态变为不稳定并出现周期极限环解——均匀稳态分岔为振荡解。任何从稳定到不稳定域的参数变化都发生类似情况。无量纲变量经过分岔值的事实为等价的量纲参数效应提供了有用信息。例如 $b = s/r$ 为捕食者与猎物线性增长率之比：稳态稳定时，捕食者增长率 $s$ 减小使 $b$ 减小，从而更可能进入不稳定域引发周期行为；反之 $r$ 减小使 $b$ 增大，可能减弱振荡倾向。但后者并不干脆——由式 3.19 知 $r$ 减小同时使 $a$ 增大，由图 3.6(c) 又趋于增加振荡可能。量纲分岔空间是 6 维，难以图形表达；无量纲化将其降为简洁的 3 维空间，由式 3.19 给出不同量纲参数变化的等价效应。例如承载量 $K$ 加倍恰等价于捕食者响应参数 $D$ 减半。无量纲参数是决定分岔的关键。

3.5 竞争模型：竞争排斥原理（Competition Models: Principle of Competitive Exclusion）

两个或更多物种竞争同一有限食物源，或在某些方面相互抑制增长。例如为领地竞争，领地直接关系到食物资源。实际竞争模型的研究已经发现若干有趣现象（见 Hsu 等 1979）。本节讨论一个非常简单的竞争模型，展示一个相当一般的、且在自然界中确实观察到的原理：两个物种竞争同一有限资源时，其中一个通常会灭绝。

考察两物种 Lotka–Volterra 竞争模型，每个物种 $N_1,N_2$ 在缺少另一物种时按 logistic 增长。引入 logistic 项使 Lotka–Volterra 系统更现实，但为凸显原理，考察较简单的模型 $\frac{dN_1}{dt} = r_1 N_1\left[1 - \frac{N_1}{K_1} - b_{12}\frac{N_2}{K_1}\right]$

[ \frac{dN_2}{dt} = r_2 N_2\left[1 - \frac{N_2}{K_2} - b_{21}\frac{N_1}{K_2}\right] ] （式 3.29、3.30），其中 $r_1,K_1,r_2,K_2,b_{12},b_{21}>0$，$r$ 为线性出生率，$K$ 为承载量。$b_{12}$ 与 $b_{21}$ 分别度量 $N_2$ 对 $N_1$、$N_1$ 对 $N_2$ 的竞争效应，一般不相等。注意竞争模型（3.29）（3.30）不是 Lotka–Volterra 捕食模型那样的保守系统。

按 $u_1 = \frac{N_1}{K_1},\quad u_2 = \frac{N_2}{K_2},\quad \tau = r_1 t,\quad \rho = \frac{r_2}{r_1},\quad a_{12} = b_{12}\frac{K_2}{K_1},\quad a_{21} = b_{21}\frac{K_1}{K_2}$（式 3.31）无量纲化，式 3.29、3.30 化为 $\frac{du_1}{d\tau} = u_1(1 - u_1 - a_{12}u_2) = f_1(u_1, u_2)$

[ \frac{du_2}{d\tau} = \rho u_2(1 - u_2 - a_{21}u_1) = f_2(u_1, u_2) ] （式 3.32）。

稳态（相平面奇点）$(u_1^*, u_2^*)$ 为 $f_1(u_1,u_2)=0$, $f_2(u_1,u_2)=0$ 的解，即 $u_1^* = 0, u_2^* = 0;\quad u_1^* = 1, u_2^* = 0;\quad u_1^* = 0, u_2^* = 1;\quad u_1^* = \frac{1-a_{12}}{1-a_{12}a_{21}},\ u_2^* = \frac{1-a_{21}}{1-a_{12}a_{21}}$（式 3.33），最后一个仅在 $u_1^*\geq 0, u_2^*\geq 0$ 且有限时相关，需 $a_{12}a_{21}\neq 1$。

在 $(u_1,u_2)$ 相平面上画出零倾线 $f_1=0$、$f_2=0$（图 3.10），四种情形立即可见。零倾线的关键部分（由式 3.32）为直线 $1 - u_1 - a_{12}u_2 = 0,\qquad 1 - u_2 - a_{21}u_1 = 0$ 第一条连同 $u_2$ 轴构成 $f_1=0$，第二条连同 $u_1$ 轴构成 $f_2=0$。

稳态稳定性仍由群落矩阵（式 3.32 下）决定 $\mathbf{A} = \begin{pmatrix} 1 - 2u_1 - a_{12}u_2 & -a_{12}u_1 \\ -\rho a_{21}u_2 & \rho(1 - 2u_2 - a_{21}u_1) \end{pmatrix}_{(u_1^*, u_2^*)}$（式 3.34）。

对第一稳态 $(0,0)$，特征方程 $|\mathbf{A}-\lambda\mathbf{I}|=0$ 给出 $\lambda_1=1, \lambda_2=\rho$（皆正），故 $(0,0)$ 不稳定。对第二稳态 $(1,0)$，$\lambda_1 = -1, \lambda_2 = \rho(1-a_{21})$，于是 $u_1^* = 1, u_2^* = 0 \text{ 是 } \begin{cases} \text{稳定，若 } a_{21}>1 \\ \text{不稳定，若 } a_{21}<1 \end{cases}$（式 3.35）类似地第三稳态 $(0,1)$ 的特征值为 $\lambda_1 = -\rho, \lambda_2 = (1-a_{12})$，于是 $u_1^* = 0, u_2^* = 1 \text{ 是 } \begin{cases} \text{稳定，若 } a_{12}>1 \\ \text{不稳定，若 } a_{12}<1 \end{cases}$（式 3.36）

对第四稳态（当其存在于正象限时），由式 3.34 得 $\mathbf{A} = (1-a_{12}a_{21})^{-1}\begin{pmatrix} a_{12}-1 & a_{12}(a_{12}-1) \\ \rho a_{21}(a_{21}-1) & \rho(a_{21}-1) \end{pmatrix}$ 特征值为 $\lambda_{1,2} = \frac{(a_{12}-1) + \rho(a_{21}-1) \pm \{[(a_{12}-1)+\rho(a_{21}-1)]^2 - 4\rho(1-a_{12}a_{21})(a_{12}-1)(a_{21}-1)\}^{1/2}}{2(1-a_{12}a_{21})}$（式 3.37）$\lambda$ 的符号（若复数则为 $\mathrm{Re}\,\lambda$）从而稳态稳定性取决于 $\rho,a_{12},a_{21}$ 的大小，存在多个需分别考察的情形，都有生态学含义。

在讨论各情形前注意：存在一个有界集，其边界上的导数向量 $(du_1/d\tau, du_2/d\tau)$ 沿边界或指向内——这里是 $(u_1,u_2)$ 平面中的一个矩形。由式 3.32，该条件在 $u_1$ 与 $u_2$ 轴上成立。矩形外边例如 $u_1=U_1$ 满足 $1-U_1-a_{12}u_2<0$、$u_2=U_2$ 满足 $1-U_2-a_{21}u_1<0$，任取 $U_1>1, U_2>1$ 即可。故系统全局稳定。

各种情形为：(i) $a_{12}<1, a_{21}<1$；(ii) $a_{12}>1, a_{21}>1$；(iii) $a_{12}<1, a_{21}>1$；(iv) $a_{12}>1, a_{21}<1$。分析方式类似，图 3.10(a)–(d) 与图 3.11(a)–(d) 对应 (i)–(iv)。作为示例，仅考察情形 (ii)：由式 3.35、3.36 知 $(1,0)$ 与 $(0,1)$ 稳定；因 $1-a_{12}a_{21}<0$，式 3.33 的第四稳态 $(u_1^*, u_2^*)$ 落在正象限；由式 3.37 知其特征值满足 $\lambda_2<0<\lambda_1$，故对小幅扰动不稳定——是鞍点。此时相轨迹可趋向两个稳态中的任一个（如图 3.11(b)），每个稳态有吸引域。一条分隔线（separatrix）将正象限分为不重叠的 I、II 两域，分隔线穿过 $(u_1^*,u_2^*)$，是鞍点轨迹之一。

以下看这些结果的生态学含义。情形 (i)：$a_{12}<1$ 且 $a_{21}<1$，两物种可共存的稳定稳态存在（图 3.10(a)）。用原始参数（3.31）表达即 $b_{12}K_2/K_1<1$ 与 $b_{21}K_1/K_2<1$。例如 $K_1,K_2$ 大致相等且种间竞争 $b_{12}, b_{21}$ 不太强时，这些条件意味着两个物种共存于比无竞争时更低但仍可维持的种群水平——竞争不激烈。若 $b_{12},b_{21}$ 大致相等而 $K_1,K_2$ 不同，则不易直接判断，必须借助无量纲分组 $a_{12}, a_{21}$ 进行比较。

情形 (ii)：$a_{12}>1$ 且 $a_{21}>1$。$K$ 大致相等时 $b_{12}, b_{21}$ 不小。分析表明三个非平凡稳态可同时存在，但只有 $(1,0)$ 与 $(0,1)$ 稳定（图 3.11(b)）。哪个物种最终存活相当微妙，取决于各物种的起始优势：若初始条件在 I 区，则 $u_2\to 0, u_1\to 1$，$N_2$ 灭绝（$N_1\to K_1$）；若 $N_2$ 有起始优势使初始点在 II 区，则 $u_1\to 0, u_2\to 1$，$N_1$ 灭绝（$N_2\to K_2$）。即便初始种群非常接近分隔线，自然涨落也会不可避免地把某 $u_i$ 推向零，因此仍可预期一个物种灭绝。

情形 (iii)、(iv)：某物种的种间竞争远强于另一者，或承载量差异大到使 $a_{12}=b_{12}K_2/K_1<1$ 与 $a_{21}=b_{21}K_1/K_2>1$（或反之），结果较为确定：(iii) 中（图 3.11(c)）$u_1$ 物种无量纲意义下的强竞争压制另一者，$u_2$ 灭绝；(iv) 中情形相反，$u_1$ 灭绝。

虽然不是所有情形都导致物种灭绝，但 (iii)、(iv) 必然如此，(ii) 也因自然涨落而不可避免。这项工作奠定了上文提到的竞争排斥原理。注意此原理成立的条件取决于无量纲参数分组 $a_{12}, a_{21}$：增长率比 $\rho$ 不影响大尺度稳定性结果，只影响动力学细节。由 $a_{12}=b_{12}K_2/K_1, a_{21}=b_{21}K_1/K_2$，竞争排斥的条件关键取决于竞争与承载量的相互作用，以及情形 (ii) 中的初始条件。

举例：两物种由大型动物与小型动物组成，竞争同一片固定面积草地。设二者竞争能力相等 $b_{12}=b_{21}$。以 $N_1$ 为大型动物、$N_2$ 为小型，则 $K_1<K_2$，从而 $a_{12}=b_{12}K_2/K_1<b_{21}K_2/K_1=a_{21}$。例如 $b_{12}=1=b_{21}$ 时 $a_{12}<1, a_{21}>1$，于是 $N_1\to 0, N_2\to K_2$——大型动物灭绝。

$a_{12}=1=a_{21}$ 的情形特殊，自然界中由于随机变异几乎不会严格成立；该情形下某一物种也仍会被竞争排斥。

物种竞争在自然界中的重要性不言而喻。本节仅讨论了一个特别简单的模型，但分析方法具有相当的一般性。Pianka（1981）的综述介绍了竞争的实际方面，Waltman（1984）的讲义集也是。Flores（1998）将略简化的竞争模型（见练习 2）应用于尼安德特人因早期现代人到来而灭绝的问题，其模型基于两物种略不同的死亡率，证明共存不可能。他从独立来源估计相关参数，所得到的灭绝时段与古生物学公认的 5000–10000 年吻合。第 II 卷第 1 章与第 14 章将讨论一些空间竞争的实际案例——松鼠、狼–鹿生存、转基因生物释放等。

3.6 互惠共生（Mutualism or Symbiosis）

许多情形中两个或多个物种的相互作用对各方都有利。互惠共生或共生常在促进乃至维持物种生存中起关键作用——植物与种子传播即是一例。即便生存不成问题，互惠或共生的互利也可能极其重要。作为理论生态学的主题，即使仅两个物种的情形，这一方向的研究远不及捕食–被捕食与竞争那样广泛，虽然其重要性与后两者相当。部分原因是 Lotka–Volterra 式的简单模型给出荒谬结果。最简单的互惠模型对应经典 Lotka–Volterra 捕食模型：$\frac{dN_1}{dt} = r_1 N_1 + a_1 N_1 N_2,\qquad \frac{dN_2}{dt} = r_2 N_2 + a_2 N_2 N_1$ 其中 $r_1,r_2,a_1,a_2>0$。因 $dN_1/dt>0, dN_2/dt>0$，$N_1,N_2$ 无限增长——正如 May（1981）所贴切形容的"一场互惠盛宴"。

现实模型至少须对各物种提供互惠，并具有正稳态或极限环型振荡。Whittaker（1975）描述了一些这样的模型，May（1975）讨论了一个实际例子。

构造合理两物种模型的第一步是为两物种引入有限的承载量，考察 $\frac{dN_1}{dt} = r_1 N_1\left[1 - \frac{N_1}{K_1} + b_{12}\frac{N_2}{K_1}\right]$

[ \frac{dN_2}{dt} = r_2 N_2\left[1 - \frac{N_2}{K_2} + b_{21}\frac{N_1}{K_2}\right]$（式 3.38），其中 $r_1,r_2,K_1,K_2,b_{12},b_{21}>0$。采用与竞争模型相同的无量纲化（注意此处 $b$ 前的符号为正，竞争模型为负），即式 3.31、3.40：$u_1 = \frac{N_1}{K_1},\quad u_2 = \frac{N_2}{K_2},\quad \tau = r_1 t,\quad \rho = \frac{r_2}{r_1},\quad a_{12} = b_{12}\frac{K_2}{K_1},\quad a_{21} = b_{21}\frac{K_1}{K_2}$ 式 3.38 化为 $\frac{du_1}{d\tau} = u_1(1 - u_1 + a_{12}u_2) = f_1(u_1, u_2)$ $\frac{du_2}{d\tau} = \rho u_2(1 - u_2 + a_{21}u_1) = f_2(u_1, u_2) ] （式 3.39）

按常规分析稳态 $(u_1^*,u_2^*)$，由式 3.39 得 $(0,0),\ (1,0),\ (0,1),\ \left(\frac{1+a_{12}}{\delta}, \frac{1+a_{21}}{\delta}\right),\ \text{其中 } \delta = 1 - a_{12}a_{21} > 0 \text{ 时正}$（式 3.41）

计算式 3.39 的群落矩阵并求各态特征值后可直接证 $(0,0)$、$(1,0)$、$(0,1)$ 均不稳定：$(0,0)$ 是不稳定结点，$(1,0)$ 与 $(0,1)$ 是鞍点。若 $1-a_{12}a_{21}<0$ 则只有前三个稳态，种群无界增长。在相平面上画零倾线 $f_1=0$, $f_2=0$，可见相轨迹在某些域内向无穷远移，$u_1\to\infty$, $u_2\to\infty$（图 3.12(a)）。

当 $1-a_{12}a_{21}>0$，式 3.41 的第四稳态存在于正象限。群落矩阵特征值显示它是稳定平衡点，是相平面上的结点。该情形如图 3.12(b)，正象限所有轨迹趋于 $u_1^*>1$ 与 $u_2^*>1$——即 $N_1>K_1$, $N_2>K_2$，每个物种的稳态数量都超过了孤立时的最大值。

此模型有一定缺陷。一是"无界增长"与"有限正稳态"之间的敏感性完全取决于不等式 $a_{12}a_{21}<1$；按式 3.40 用原始参数表达即 $b_{12}b_{21}<1$（$b$ 已无量纲）。即若任一物种的共生效应过大，条件被违反，两种群都无界增长。

3.7 一般模型与若干一般性及警示性评注（General Models and Some General and Cautionary Remarks）

本章讨论的所有模型都化为如下非线性微分方程组 $\frac{dN_i}{dt} = N_i F_i(N_1, N_2, \dots, N_n),\quad i = 1, 2, \dots$（式 3.42）这凸显了种群向量 $\mathbf{N}$ 必有 $\mathbf{N}=0$ 作为稳态。两物种版本有时称为 Kolmogorov 模型或 Kolmogorov 方程。

虽然本章主要讨论两物种相互作用，自然界中（尤其海洋中）存在许多物种或营养级，能量（以食物形式）从一种群流向另一种群。某物种的总质量常称为其生物量（biomass），即种群数量乘以单位质量。最终能量来源是太阳；例如在海洋中食物网经过浮游生物、鱼、鲨鱼直至鲸鱼再至人类，中间物种繁多。同一营养级的物种可捕食其下层若干物种。一般而言模型涉及多种群相互作用。

多种群模型为 $\frac{d\mathbf{u}}{dt} = \mathbf{f}(\mathbf{u}),\ \text{即}\ \frac{du_i}{dt} = f_i(u_1, \dots, u_n),\ i=1,\dots,n$（式 3.43）其中 $\mathbf{u}(t)$ 是种群密度的 $n$ 维向量，$\mathbf{f}(\mathbf{u})$ 描述种间的非线性相互作用。$\mathbf{f}(\mathbf{u})$ 含表征所研究系统各种生长与相互作用特征的参数，$f_i$ 给出第 $i$ 物种的总增长率。稳态稳定性分析方法与前相同：在稳态 $\mathbf{u}^*$（$\mathbf{f}(\mathbf{u}^*)=0$）附近线性化，检查群落矩阵或稳定矩阵 $\mathbf{A} = (a_{ij}) = \left.\left(\frac{\partial f_i}{\partial u_j}\right)\right|_{\mathbf{u}=\mathbf{u}^*}$（式 3.44）的特征值 $\lambda$。特征值满足 $|\mathbf{A}-\lambda\mathbf{I}|=0$，要求所有 $\mathrm{Re}\,\lambda < 0$ 的充要条件由 Routh–Hurwitz 条件给出（见附录 B）。

若稳态不稳定，解 $\mathbf{u}$ 可能无界增长，可能演化到另一稳态，或演化为如极限环的稳定振荡形态。对两物种模型此类方程的理论基本完整：它们是相平面系统，分析方法见附录 A 简要回顾。对三物种或更多物种相互作用系统，一般理论少得多。在稳态变为不稳定时，至少对稳态附近的解，常可用 Hopf 分岔理论找到部分结果（Strogatz 1994 一书对此有很好的教学讨论）。该理论的最简形式说：若系统某参数 $p$ 有临界值 $p_c$，使 $p<p_c$ 时最大实部特征值 $\mathrm{Re}\,\lambda < 0$；$p=p_c$ 时 $\mathrm{Re}\,\lambda = 0$、$\mathrm{Im}\,\lambda \neq 0$；$p>p_c$ 时 $\mathrm{Re}\,\lambda > 0$、$\mathrm{Im}\,\lambda \neq 0$；则当 $p-p_c > 0$ 且很小时，解 $\mathbf{u}$ 在 $\mathbf{u}^*$ 附近表现小振幅极限环行为。Smith（1993）发展了新方法来研究三（或更多）竞争或合作物种，其方法通过把流映射为拓扑等价的二维流，使 Poincaré–Bendixson 定理可应用于三物种系统。

群落矩阵 $\mathbf{A}$（式 3.44）在决定稳态线性稳定性中至关重要，并具有直接的生物学意义：元素 $a_{ij}$ 衡量第 $j$ 物种在稳态附近对第 $i$ 物种的效应。设 $U_i$ 为偏离稳态 $u_i^*$ 的扰动，则 $U_i$ 的方程为 $\frac{dU_i}{dt} = \sum_{j=1}^n a_{ij} U_j$（式 3.45），$a_{ij}U_j$ 是物种 $U_j$ 对 $U_i$ 生长的效应：若 $a_{ij}>0$，$U_j$ 直接促进 $U_i$ 增长；若 $a_{ij}<0$，则抑制。若 $a_{ij}>0$ 且 $a_{ji}>0$，则 $U_i$ 与 $U_j$ 互促增长，即共生；若 $a_{ij}<0$ 且 $a_{ji}<0$，则竞争。May（1975）综述了一些一般化模型，并在稳定性与复杂性的讨论中给出了基于群落矩阵性质的稳定性结论。

已有相当多的研究关注群落矩阵具对角对称、反对称或其他特殊性质的系统，对这类系统可给出特征值与稳态稳定性的一般结论。但实际情形下的模型很少具有如此简单的性质，且参数评估的随机性使此类模型的近似意义有限。然而正如经典 Lotka–Volterra 系统与现实世界不相关一样，这些特殊模型也常常促使人们提出正确问题。即便如此，若基本目标是理解现实世界，则须避免过分执着于这些模型或其推广。

本章未讨论的另一重要模型类是含延迟的相互作用模型。若物种存在不同或分布式的延迟，这类模型开启了真正的"潘多拉盒子"——其解行为在很大程度上仍是相对未探索的领域。

考察三物种及以上时，可出现非周期行为。Lorenz（1963）首次以下列模型系统展示此现象 $\frac{du}{dt} = a(v - u),\quad \frac{dv}{dt} = -uw + bu - v,\quad \frac{dw}{dt} = uv - cw$ 其中 $a,b,c>0$（该方程组实际来源于流体模型）。随参数变化解展现倍周期并最终出现混沌或非周期行为。Rössler（1976a,b; 1979; 1983）、Sparrow（1982; 1986）、Strogatz（1994）特别研究了此类系统并发现了若干其他具类似性质的基本例子；亦见 Holden（1986）编辑的书。若某些三物种及以上的种群相互作用模型不能展现类似性质反倒令人意外。三物种及以上的竞争模型会产生若干意想不到的结果。

自然界种群中混沌（甚至复杂振荡）的证据难以发现，也难以确证。因此有人提出进化倾向于使种群免于此类混沌行为。Ferrière 与 Gatto（1993）的有趣文章研究了由自然选择导致的可能混沌种群动力学，结果假设进化可能支持比此前所认为的更多混沌种群动力学，并颇具争议地提出混沌或许是种群动力学中的最优行为——这也与 Schaffer 与 Kot（1986）关于流行病的观点一致。

尽管如此，复杂种群相互作用的进化发育似乎总体上产生了相当稳定的系统。从本章对相互作用模型的研究可知，当某些参数（越过相应的分岔值）作适当变化时，系统可被驱动到不稳定。任何系统在被外部操纵改变之前，都应进行相当多的科学研究。利用模型研究人为干预此类食物网的效应是必要的且极具启发性。若此前做过此类模型评估，下述灾难或可避免。虽然使用现实的动力学模型不能给出完整答案（以引入新物种或消除一物种的预测形式），但它们能指出必须认真考虑的各种危险征兆。同理，不能过分信赖模型——相互作用常常极其复杂，建模者未必能构建足够好的模型。为结束本节，描述一项由试图操纵复杂食物网导致的主要生态灾难——发生在东非。

维多利亚湖与尼罗河鲈鱼灾难（1960）：1960 年尼罗河鲈（Lates niloticus）被引入维多利亚湖（东非最大湖），该湖被肯尼亚、坦桑尼亚、乌干达环抱，长期支撑着沿岸数百个小型渔业社区。引入这种大型肉食鱼类（可达 100 千克以上）被认为是高产且有价值的蛋白源，得到联合国粮农组织支持，但部分科学家的反对意见被忽视。35 年间，这种大型肉食鲈几乎消灭了湖中数百种较小的丽鱼；这些丽鱼曾是沿岸渔村经济的主要支柱。市场上充斥鲈鱼。据估计 1984 年湖泊总生产力比 1960 年前下降约 80%。湖内这一计划外的、引入主要且不合适物种的失误造成了生态灾难。除渔村的经济灾难外，还存在湖外的连锁效应，本应预料到：例如大型鲈鱼肉多油脂，无法日晒干燥，必须熏制保存，导致大量珍贵树木被砍伐用作燃料。更严重的是，许多几乎消失的丽鱼物种曾帮助控制某种特定蜗牛的数量（这种蜗牛生活在湖中及周边）。

这些淡水螺生活在许多大坝的大型水库中，是寄生虫病血吸虫病（schistosomiasis，又称 bilharzia）传播链中的关键一环。该病重要性仅次于疟疾，不受人体免疫系统攻击，若不治疗对人类致命。1990 年代另一灾难袭击维多利亚湖：水葫芦（Eichhornia crassipes）在湖周形成宽厚毡层，破坏鱼类繁殖场，堵塞水电站等。乌干达约 80% 海岸线已被感染。该植物原产巴西，可部分通过锈菌真菌天敌控制，但学界对外来天敌的使用仍持谨慎态度。然而生物防治与化学除草剂都在认真考虑之中。

3.8 阈值现象（Threshold Phenomena）

除 Lotka–Volterra 捕食模型外，本章已考察或提及的两物种模型要么具有稳定稳态（小扰动衰减），要么稳态不稳定（扰动无界增长或产生极限环周期解）。存在另一类有趣的模型具有非零稳定稳态，但若扰动足够大或为某种特定形式，种群密度会发生大幅变化后回到稳态。这类模型被称为展现阈值效应。本节研究其中一组。

考察如下捕食者–猎物模型 $\frac{dN}{dt} = N[F(N) - P] = f(N,P)$（式 3.46） $\frac{dP}{dt} = P[N - G(P)] = g(N,P)$（式 3.47）其中为方便起见所有参数已通过合适重标并入 $F,G$ 中；$F(N)$, $G(P)$ 的定性形式如图 3.13 所示。$F(N)$ 的特定形式体现了 Allee 效应，即猎物的每员增长率随猎物密度先增、在某 $N_m$ 处达到最大、随后随密度增大而下降。

式 3.46、3.47 的稳态 $(N^*, P^*)$ 满足 $N^*=0=P^*$ 以及非负解 $P^* = F(N^*),\quad N^* = G(P^*)$（式 3.48）。和往常一样，画零倾线 $f=0$, $g=0$（图 3.14）很有帮助。视 $F(N),G(P)$ 的具体参数，稳态可典型地处于 $S$ 或 $S'$ 处。这里考察 $N^*>N_m$ 的情形，即稳态处于图 3.14 的 $S$。

由式 3.46、3.47，零稳态 $N^*=0=P^*$ 的群落矩阵为 $\mathbf{A} = \begin{pmatrix} F(0) & 0 \\ 0 & -G(0) \end{pmatrix}$ 特征值 $\lambda = F(0) > 0$ 与 $\lambda = -G(0) < 0$，故 $(0,0)$ 不稳定，是 $(N,P)$ 相平面上的鞍点奇点。

对正稳态 $(N^*,P^*)$，群落矩阵为 $\mathbf{A} = \begin{pmatrix} N^* F'(N^*) & -N^* \\ P^* & -P^* G'(P^*) \end{pmatrix}$（' 表示求导）其中由图 3.14 知 $(N^*,P^*)$ 处于 $S$ 时 $G'(P^*)>0$ 且 $F'(N^*)<0$；处于 $S'$ 时 $G'(P^*)>0$ 且 $F'(N^*)>0$。特征值 $\lambda$ 满足式 3.49：$\lambda^2 - (\mathrm{tr}\,\mathbf{A})\lambda + \det\mathbf{A} = 0$，其中 $\mathrm{tr}\,\mathbf{A} = N^* F'(N^*) - P^* G'(P^*)$ $\det\mathbf{A} = N^* P^* [1 - F'(N^*) G'(P^*)]$（式 3.50）

当稳态处于图 3.14 的 $S$ 时 $\mathrm{tr}\,\mathbf{A}<0$、$\det\mathbf{A}>0$，由式 3.49 知 $\mathrm{Re}\,\lambda < 0$，对任意 $F(N), G(P)$ 都稳定。若稳态处于 $S'$，$\mathrm{tr}\,\mathbf{A}$ 与 $\det\mathbf{A}$ 可正可负（因 $F'(N^*)>0$），从而 $S'$ 可能稳定或不稳定，视具体 $F(N), G(P)$ 而定；若 $S'$ 不稳定，由于系统存在有界集，会产生极限环解（参见 3.4 节类似问题与图 3.8）。

本节关注的重点是稳态处于 $S$（永远稳定）的情形。假设把系统扰动到图 3.15(a) 中的 $X$ 点。该处 $f<0$、$g<0$，由式 3.46、3.47 知 $dN/dt<0$、$dP/dt<0$，轨迹定性按图 3.15(a) 所示移动，最终回到 $S$，但在相平面中要经历大幅偏离。把 $(N,P)$ 扰动到 $Y$ 发生类似行为；但若扰动到 $Z$，扰动保持在 $S$ 附近。图 3.15(b)(c) 示意了 $N$ 与 $P$ 的典型时间演化。

显然存在一个大致阈值扰动——低于它时扰动总保持在稳态附近，高于它时则不然，即便最终仍回到稳态。阈值扰动更像一条阈值曲线或阈值域：若扰动使轨迹越过图 3.15(a) 的最大 $N_m$，轨迹典型地如从 $X$、$Y$ 出发的情形；若轨迹在 $N>N_m$ 处穿越 $f=0$，则不发生大幅扰动。这种阈值性质的存在源于 $f=0$ 零倾线具最大值（如图所示）的形式——这是模型（3.46）动力学中 Allee 效应的体现。结合本章前述问题，仅从种群的时序行为上看，可能以为面对的是不稳定情形；仔细绘制零倾线显然是必要的。此时阈值的定义尚不精确。第 II 卷第 1 章将证明，若允许某物种扩散（例如通过扩散项），阈值行波解是可能的，它们具有重要的生物学意义。在该语境下阈值概念可以精确化。本书后续在生化背景下还会以重要方式出现此阈值行为——这些反应动力学方程与本章讨论的种群相互作用方程形式相同。

最后关于建模相互种群的评注：对给定情境不存在"正确"模型——许多模型都能给出定性相似的行为。得到正确的定性特征只是第一步，绝不应视为模型的正当性。一个好模型应当具备的条件包括：基于观察与真实事实的生长动力学合理性、各种参数的合理可估性，以及基于模型的预测是否能被后续实验与观察所验证。

3.9 相互作用种群的离散增长模型（Discrete Growth Models for Interacting Populations）

考察两个世代不重叠、互相影响种群动力学的物种。与连续增长模型相同，主要相互作用类型为捕食–被捕食、竞争与互惠。捕食–被捕食中一个种群增长率因另一者提高而下降；竞争中两者增长率均下降；互惠中两者增长率均提高。这些主题已被广泛研究，但远不及连续模型——后者在两物种情形已有完整的数学处理。Hassell（1978）的书讨论捕食–被捕食模型；Beddington 等（1975）给出耦合捕食系统动力学复杂性的若干结果；Gumowski 与 Mira（1980）的书更具数学味，讨论耦合系统数学并包含一些有趣的数值结果；亦见 Lauwerier（1986）的综述。May（1986）的综述与本章及前几章内容相关，核心议题是种群如何自我调节。他还讨论了叠加于确定性模型之上的不可预见环境因子以及资源管理的实际方面。鉴于单物种离散模型解行为的复杂性，耦合离散系统出现更复杂行为不足为奇。Peitgen 与 Richter（1986）以绚丽色彩展示了仅两个耦合方程系统就能生成的丰富模式——包括分形、Mandelbrot（1982）集、Julia 集、Hubbard 树等。第 14 章将给出分形的简要介绍。

本节关注捕食–被捕食模型。进化通过自然选择偏好高效的捕食者与机敏的猎物。在普遍情形中，我们主要考虑昆虫捕食–被捕食系统——既有相当规模的实验数据，昆虫的生活史也常适合用两物种离散模型。

设猎物 $(N)$ 与捕食者 $(P)$ 的相互作用由如下离散时间 $(t)$ 耦合方程组描述 $N_{t+1} = r N_t f(N_t, P_t)$（式 3.51） $P_{t+1} = N_t g(N_t, P_t)$（式 3.52）其中 $r>0$ 为猎物的净线性增长率，$f$ 与 $g$ 分别反映捕食压力下猎物的繁殖效率与捕食者的搜索效率。本节讨论的技术自然适用于其他种群相互作用。下一章讨论的"物种"是婚姻状态，理论与本章形式相似但有本质差别，关键差别在于下一章个体间相互作用是重叠的。

3.10 捕食者–猎物模型：详细分析（Predator–Prey Models: Detailed Analysis）

先考察一个简单模型：捕食者在固定面积上搜索且对猎物有无限消费能力。模型为 $N_{t+1} = r N_t \exp[-a P_t]$ $P_{t+1} = N_t \{1 - \exp[-a P_t]\}$ $a > 0$（式 3.53）

值得一提的是，看方程给出的定性行为对每物种意味着什么总是很有帮助——据此尝试预判稳定性分析的结果。一般若结果与直觉不符，这种事先的定性印象常有助于修正模型使之更现实。

式 3.53 的平衡值 $(N^*, P^*)$ 为 $N^* = 0,\ P^* = 0\ \text{或}\ 1 = r\exp[-aP^*],\ P^* = N^*(1 - \exp[-aP^*])$ 故正稳态种群为 $P^* = \frac{1}{a \ln r},\quad N^* = \frac{r}{a(r-1)\ln r},\qquad r > 1$（式 3.54）

用通常方法做线性稳定性分析：令 $N_t = N^* + n_t$, $P_t = P^* + p_t$，$|n_t/N^*| \ll 1$, $|p_t/P^*| \ll 1$（式 3.55），代入式 3.53 仅保留线性项。对零稳态 $(0,0)$ 特别简单：$n_{t+1} = r n_t$, $p_{t+1} = 0$。故 $r<1$ 时稳定（$N_t\to 0$），$r>1$ 时不稳定——这正是正稳态（3.54）存在的 $r$ 范围。

对正稳态，线性化系统为 $n_{t+1} = n_t - N^* a p_t$

[ p_{t+1} = n_t\left(1 - \frac{1}{r}\right) + \frac{N^* a}{r} p_t ] （式 3.56）其中用到 $1 = r\exp[-aP^*]$。

求解式 3.56 的一个直接方法是迭代第一式并利用第二式把 $n_t$ 化为单一方程：$n_{t+2} = n_{t+1} - N^* a p_{t+1} = \cdots = \left(1 + \frac{N^* a}{r}\right) n_{t+1} - N^* a n_t$ 即 $n_{t+2} - \left(1 + \frac{N^* a}{r}\right) n_{t+1} + N^* a n_t = 0$（式 3.57）设解形如 $n_t = A x^t$，代入得 $x^2 - \left(1 + \frac{N^* a}{r}\right) x + N^* a = 0$ 以式 3.54 的 $N^*$ 代入，特征多项式为 $x^2 - \left(1 + \frac{1}{r-1}\ln r\right) x + \frac{r}{r-1}\ln r = 0,\quad r > 1$（式 3.58）两根为 $x_{1,2} = \frac{1}{2}\left\{\left(1 + \frac{\ln r}{r-1}\right) \pm \left[\left(1 + \frac{\ln r}{r-1}\right)^2 - \frac{4r\ln r}{r-1}\right]^{1/2}\right\}$（式 3.59）于是 $n_t = A_1 x_1^t + A_2 x_2^t$（式 3.60）其中 $A_1,A_2$ 为任意常数；类似可得 $p_t = B_1 x_1^t + B_2 x_2^t$（式 3.61）

更优雅且易推广的方法是把线性化系统（3.56）写为矩阵形式 $\begin{pmatrix} n_{t+1} \\ p_{t+1} \end{pmatrix} = \mathbf{A}\begin{pmatrix} n_t \\ p_t \end{pmatrix},\quad \mathbf{A} = \begin{pmatrix} 1 & -N^* a \\ 1 - \frac{1}{r} & \frac{N^* a}{r} \end{pmatrix}$（式 3.62）并设解为 $\begin{pmatrix} n_t \\ p_t \end{pmatrix} = \mathbf{B}\begin{pmatrix} 1 \\ 1 \end{pmatrix} x^t$，代入得 $x\mathbf{B}\begin{pmatrix} 1 \\ 1 \end{pmatrix} = \mathbf{A}\mathbf{B}\begin{pmatrix} 1 \\ 1 \end{pmatrix}$。该方程有非平凡解当且仅当 $|\mathbf{A} - x\mathbf{I}| = 0$，即 $\begin{vmatrix} 1 - x & -N^* a \\ 1 - \frac{1}{r} & \frac{N^* a}{r} - x \end{vmatrix} = 0$ 重新给出特征方程（3.58）。$x_1,x_2$ 即矩阵 $\mathbf{A}$（式 3.62）的特征值。这一矩阵方法就是连续相互作用种群模型中所用方法的离散版本。推广到更高阶离散模型系统是显然的。

正稳态 $(N^*,P^*)$ 的稳定性由 $|x_1|, |x_2|$ 的大小决定：若任一 $|x_i| > 1$，则 $n_t, p_t$ 随 $t\to\infty$ 无界增长，从而 $(N^*,P^*)$ 不稳定（扰动随时间放大）。简单代数表明，式 3.59 中 $\left(1 + \frac{\ln r}{r-1}\right)^2 - \frac{4r\ln r}{r-1} < 0,\quad r>1$ 故 $x_1,x_2$ 为复共轭对。式 3.58（或 3.59）中两根之积为 $x_1 x_2 = |x_1|^2 = \frac{r\ln r}{r-1} > 1,\quad \forall r>1 \Rightarrow |x_1| > 1$（验证 $(r\ln r)/(r-1)>1$ 的简便方法：比较 $r>1$ 时 $\ln r$ 与 $(r-1)/r$ 的图像，可见 $d(\ln r)/dr > d[(r-1)/r]/dr$ 对所有 $r>1$ 成立。）因此式 3.60、3.61 的解 $(n_t,p_t)$ 随 $t\to\infty$ 无界增长，正平衡 $(N^*,P^*)$（式 3.54）不稳定，且因 $x_1,x_2$ 为复数，是增长型振荡不稳定。系统（3.53）的数值解表明对有限扰动也不稳定，解无界增长。因此该简单模型除在人造实验条件下且仅限短时间外，毫无实际应用价值。

密度依赖的捕食者–猎物模型：重新审视简单初始模型（3.53）背后的假设。方程形式意味着捕食者与猎物的相遇数随猎物种群密度无界增长，这相当不现实。更合理的设定是捕食者食欲存在上限。从另一等价角度看，若 $P_t=0$ 且 $r>1$，则 $N_t$ 无界增长；$0<r<1$ 时灭绝——即简单 Malthus 模型（2.2）。合理地修正 $N_t$ 方程（3.53）以纳入猎物种群的某种饱和或"猎物限制模型"——取更现实的模型为 $N_{t+1} = N_t \exp\left[r\left(1 - \frac{N_t}{K}\right) - a P_t\right]$ $P_{t+1} = N_t\{1 - \exp[-a P_t]\}$（式 3.63） $P_t=0$ 时该模型退化为单物种模型（2.8）（2.1 节），$0<r<2$ 时有稳定正平衡 $N^*=K$，$r>2$ 时出现振荡与周期解。可以合理预期该模型有类似分岔行为，但首次分岔点不必为 $r=2$，后续分岔的 $r$ 值也不同。Beddington 等（1975）对该模型作了详细研究。

式 3.63 的非平凡稳态满足 $1 = \exp\left[r\left(1 - \frac{N^*}{K}\right) - a P^*\right]$

[ P^ = N^(1 - \exp[-a P^])$（式 3.64）第一式给出 $P^* = \frac{r}{a}\left(1 - \frac{N^*}{K}\right)$（式 3.65）代入第二式得 $N^*$ 的超越方程 $\frac{r}{a}\left(1 - \frac{N^}{K}\right) \cdot \frac{1}{N^} = 1 - \exp\left[-r\left(1 - \frac{N^}{K}\right)\right] ] （式 3.66）显然 $N^*=K, P^*=0$ 是一个解。在 $N^*$ 轴上画出式 3.66 两边可见另一平衡 $0 < N_E^* < K$，是曲线的另一交点，依赖于 $r,a,K$。$N_E^*$ 确定后代入式 3.65 得 $P_E^*$。

此平衡的线性稳定性仍按前述方法处理，特征值 $x$ 仍由线性化系统的矩阵特征值给出，需数值求解。可证明对某些 $r>0$ 平衡稳定，更大 $r$ 时发生分岔。Beddington 等（1975）在 $r, N_E^*/K$ 参数空间中确定了从稳定到不稳定的分界——解展现周期行为并最终混沌。现实两物种模型的稳定性分析常需数值完成。对三物种及更高情形，Jury 条件（附录 B）可用来确定系数须满足的条件使线性解 $x$ 满足 $|x|<1$。但对更高阶系统，这些条件仅在数值方案内有用。

生物防治：一般性评注：利用天敌进行害虫防治的目的是以捕食者种群的相应增长抑制任何大幅的害虫增长，目标是将两种群都维持在可接受的低水平，目标是控制而非根除害虫。虽然许多真实捕食者–害虫的模型系统从稳定性角度看相当稳健，但有些可极为敏感——这正是现实模型分析如此重要的原因。当模型参数取自观察时，幸运的是许多参数给出稳态平衡或简单周期行为，混沌行为较少出现。因此在相当多实际情形下有效的参数操控更具可预测性。生物防治有许多显著成功，特别是对果树、森林作物等长期作物，由于持续的捕食者–猎物相互作用而成功。在多年生作物收获造成的重大生态变化情形下则较少成功。成功案例主要为捕食者作为寄生者的捕食–被捕食型，这在许多人类疾病中可能极其重要。Kot（2001）详细讨论了收获模型动力学，包括最优控制的重要方面。

在前述模型分析中，重点放在模型建立、参数范围内的稳定性研究以及稳态或周期行为的存在性。未讨论的是初始条件的影响——虽然通常不重要，但有时很重要。一个例子是控制红蜘蛛螨——温室番茄植物害虫，其中初始捕食者–猎物比至关重要。在那些振荡表现为爆发、崩溃与缓慢恢复的情形中，初始数据预计尤为重要。崩溃到低水平时物种可能接近灭绝临界从而真正灭绝。关于生物防治的书很多，例如 DeBach（1974）与 Huffaker（1971）。

一个相对较新且近乎未开垦的领域是捕食者与猎物时间步长不等的耦合系统。这在真实世界中显然存在。鉴于本章与第 2 章模型展现的有趣而意外的行为，若时间步长不等的模型不展现同样意外的行为反倒令人意外。

本章个人批注

本章给我印象最深的两点是：其一是 3.2 节那条"复杂度通常带来不稳定而非稳定"的论证——它从迹为零的群落矩阵出发，依赖复特征值成共轭对这一代数事实，立即得出"只要有一个 $\mathrm{Re}\,\lambda \neq 0$，必同时存在 $\mathrm{Re}\,\lambda > 0$"的结论，进而说明稍微变动参数就能把中性稳定推到不稳定。这种论证在生物学中极其少见——通常我们只能给出"经验上复杂系统更脆弱"的描述；这里能给出精确的代数机制，很值得记。其二是 3.5 节关于 $a_{12}a_{21}$ 与 1 的大小关系把竞争排斥原理条件清晰地写出来：竞争系数与承载量的乘积决定了哪个物种占优，而增长率比 $\rho$ 只影响动力学细节而不影响大尺度稳定性结论。这种"参数几何"的处理方法——把生态学问题化为几个无量纲分组之间的不等式——是 Murray 处理多物种问题的标志风格。

需要标出几点怀疑与待验证的地方：(1) 3.10 节式 3.63 的分岔分析被作者推到 Beddington 等（1975）的论文，具体首次分岔点未必等于 $r=2$ 但未给具体数值；以及"现实两物种模型的稳定性分析常需数值完成"——这条结论是限制性的，意味着我们面对真实生态场景时必须回到数值。(2) 3.8 节关于阈值扰动的论述，作者承认"此时阈值的定义尚不精确"；这一节读起来更像现象描述而非机制分析，后续第 II 卷第 1 章据说会给出严格化（threshold travelling waves）——这是后续章节需要回头核对的承诺。(3) 维多利亚湖与尼罗河鲈的灾难段属于生态学叙事而非建模分析，其作用是警示而非论证；学界对 Nile perch 引入的反思很多，作者立场鲜明（认为引入是错误的）但并未引用反对意见的具体文献，建议保留怀疑态度。(4) 关于 Flores（1998）尼安德特人灭绝模型——他说 s=0.995 是独立估计，但具体来源未给出；若读到该论文应核对原始数据来源。

与上下章的衔接（一段话）

本章在全书结构上承接第 1、2 章（单物种连续与离散种群模型）并直接为后续章节铺垫。从方法学看，第 1、2 章所用的无量纲化、相平面分析、群落矩阵与线性稳定性分析在第 3 章全部复用并扩展——单物种的 logistic 自限制与离散 logistic 的分岔结构在此变为两物种相互作用的稳定性与极限环分析；从主题看，第 3 章建立的三种典型相互作用（捕食–被捕食、竞争、互惠）以同样方程形式贯穿全书后续：作者明确指出第 6 章反应动力学的方程"数学上与本章相互作用的种群动力学方程同型"，那里的"物种"是化学浓度；第 4 章则专门讨论温度依赖性别决定（一种内源性单物种结构变化）；第 5 章虽然用了离散时间模型但"物种"代表婚姻状态，分析技术类似但概念核心完全不同；第 7、8 章的生物振荡与 BZ 反应则把第 3 章引入的极限环概念推广到更复杂的振荡系统；第 9 章讨论耦合振荡器的扰动，正是 3.10 节密度依赖捕食–猎物模型分析技术的多物种推广。第 3 章因此是"种群动力学"主线从单物种进入多物种的枢纽章节——它不仅是新内容，更是后续几乎所有章节共享的方法论基线。