第 2 章：单物种离散种群模型（Discrete Population Models for a Single Species）

2.1 引言：简单模型（Introduction: Simple Models）

微分方程模型（无论常微分、延迟微分、偏微分或随机形式）默认种群世代之间存在连续重叠。但许多物种世代之间完全不重叠，种群增长以离散步骤发生。对原始生物而言步骤可以非常短，此时连续时间模型是合理近似；然而依物种不同，步骤长度差异显著——一年是常见周期，果蝇从蛹羽化为一天，细胞为数小时，细菌和病毒更短。本章及第 5 章的模型中时间步长统一缩放为 1，需要研究差分方程（即离散模型）$N_{t+1} = N_t F(N_t) = f(N_t)$（式 2.1），其中 $f(N_t)$ 一般为 $N_t$ 的非线性函数。第一种写法强调零稳态的存在。这样的方程通常无法解析求解，但无需解析解亦能提取大量种群动力学信息。差分方程的数学正在被深入研究并应用于众多领域：癌症生长（Cross 和 Cotton 1994）、老化（Lipsitz 和 Goldberger 1992）、细胞增殖（Hall 和 Levinson 1990）、遗传学（Hoppensteadt 和 Peskin 1992、Roughgarden 1996）；最近在婚姻互动和离婚预测的动力学建模中显示出惊人用途，将在第 5 章讨论。迄今最大量的应用在生态学领域，Hassell（1978）和 Kot（2001）列举大量实例。实际应用中只要知道 $f(N_t)$ 的形式，便可通过递归使用式 2.1 简单地计算 $N_{t+1}$ 及后续世代。

无论 $f(N_t)$ 取何形式，我们只对非负种群值感兴趣。建模某一特定种群增长动力学的技巧在于选取合适的 $f(N_t)$ 形式以反映已知观察事实。这要求我们必须理解 $f(N_t)$ 形式和参数改变对解的主要影响，也要了解若干有实际兴趣的典型例子下式 2.1 的解看起来是什么样。数学问题是一个映射问题——给定初值 $N_0 > 0$，求非线性映射的轨道或轨迹。值得指出的是：差分方程模型与看似连续微分方程类比物之间没有简单联系，即便有限差分近似会导出离散方程，下文即可看清这一点。

设 $F(N_t) = r > 0$（即种群下一步正比于当前种群），由式 2.1 得 $N_{t+1} = r N_t \Rightarrow N_t = r^t N_0$（式 2.2），种群按几何级数增长或衰减，分别取决于 $r>1$ 还是 $r<1$；这里 $r$ 为净生殖率。这一特别简单的模型对大多数种群不切实际，也不适用于长时间，但对某些细菌生长的早期阶段有过一定程度的正当使用——是第 1 章 Malthus 模型的离散版本。一个引入拥挤效应的简单修正可写作 $N_{t+1} = r N_t S$, $N_t^S = N_t^{1-b}$，$b$ 为常数，$N_t^S$ 是存活到繁殖期的种群，当然对 $b$ 有限制要求 $N_t^S \le N_t$，否则存活到繁殖期的数目将超过它们所属的种群。

Fibonacci 数列：Leonardo of Pisa 直到 18 世纪才被赋予 Fibonacci 的绰号，他在 1202 年的算术书中设置了一道关于虚构增长兔子种群的建模练习：从繁殖季之初一对未成熟雌雄兔开始，一对兔子在经过一个繁殖季后产下两对未成熟雌雄幼兔，停止繁殖；它们的后代继续完全相同的繁衍。问题是要确定每一繁殖期兔子的对数。以 $N_t$ 表示雄雌配对的兔子数，将繁殖期归一为 1，则第 $t$ 繁殖阶段满足 $N_{t+1} = N_t + N_{t-1}$, $t = 2, 3, \dots$（式 2.3）。取 $N_0 = 1$, $N_1 = 1$ 即得 Fibonacci 数列 $1, 1, 2, 3, 5, 8, 13, \dots$。每一项是前两项之和。式 2.3 是线性差分方程，可以寻找形如 $N_t \propto \lambda^t$ 的解，代入得特征方程 $\lambda^2 - \lambda - 1 = 0$，解为 $\lambda_{1,2} = \frac{1}{2}(1 \pm \sqrt{5})$。以 $N_0 = 1, N_1 = 1$ 的初始条件得式 2.4 的解；大 $t$ 时因 $\lambda_1 > \lambda_2$，$N_t \approx \frac{1}{2}(1 + 1/\sqrt{5})\lambda_1^t$。式 2.3 是一个更新方程，若将年龄与繁殖期及繁殖期后不再繁殖联系起来考虑，可从这个模型直观地看出年龄结构，这催生了更新矩阵与 Leslie 矩阵。取连续 Fibonacci 数之比 $N_t/N_{t+1} \approx (\sqrt{5} - 1)/2$，即所谓黄金比或黄金数。古典绘画中以之为天空与陆地之比的追求目标。这一数列及其极限值在令人惊奇的众多场合出现：松果、向日葵花盘、雏菊花瓣、许多植物连续分支间的角度等。在向日葵花盘上可看到从中心发出的交错螺旋（松果也可由底部观察），螺旋的条数变化但总是 Fibonacci 数。图 2.1 给出了两个自然出现的交错对数螺旋的例子，图 2.1(b) 中每个鳞片同时属于顺时针与逆时针螺旋：仔细计数得 8 条顺时针、13 条逆时针螺旋，是 Fibonacci 数列的连续两项；雏菊花盘上为 21 条顺时针、34 条逆时针，仍是连续两项。在分枝叶序情形，若将许多植物和树木的分支投影到平面上，相邻分支间的角度几乎恒定，约为 $137.5^\circ$。与 Fibonacci 数列的联系是：用 $360^\circ$ 乘以上述 Fibonacci 比的极限值 $(\sqrt{5} - 1)/2$ 得 $222.5^\circ$，大于 $180^\circ$，需从 $360^\circ$ 中减去得 $137.5^\circ$，即 Fibonacci 角。

已有若干尝试建模植物形态的图案生成过程以产生 Fibonacci 角与 Fibonacci 数列，但目前问题仍未解决。方法从反应扩散机制的操控（Thornley 1976），到考察前 $n$ 个自然数的置换代数关系（Kunz 和 Rothen 1992），再到磁场中磁性液滴的实验（Douady 和 Couder 1992）。本书后面将详细讨论生成空间图案的各种可能机制，包括反应扩散系统。作者坚信该过程是机制性的而非遗传性的，Douady 和 Couder（1992, 1993a, 1993b）的工作虽然是物理而非生物性质的，但支持了这一信念。他们的工作巧妙而有趣：考虑叶序分枝中原始体顺序出现在生长顶点上，并以等时间间隔向外移动到生长尖周围的一个圆上；这些原始体被假定为向外移动时相互排斥，从而最大化彼此距离，以高度规则的空间图案自组织。若如此，模拟此情景的实验应当给出元素分布，元素之间的角度应为 Fibonacci 角。他们用直径 8 cm 的圆形盘装满硅油，放入垂直磁场（场强向盘周缘递增），每隔相等时间向盘中心的一个小截锥（模拟植物顶点）滴下少量铁磁流体；液滴被磁场极化，形成相互排斥的小磁偶极子，并因磁场梯度向外移动至周缘，最终规则分布；与之前液滴的相互作用使新液滴在能量最小方向落下；为防止液滴在周缘堆积，它们最终落入周缘的沟槽。铁磁流体滴落的时间间隔影响螺旋的生成以及到达周缘时最终角度；大量实验中角度本质上就是 Fibonacci 角、螺旋条数也是 Fibonacci 数；并以计算机模拟得到确认。

通常由于拥挤和自调节，式 2.1 中的 $f(N_t)$ 在 $N_m$ 处有极大（对 $N_t > N_m$ 时递减），图 2.2 给出典型形式；实际中已使用过多种 $f(N_t)$ 形式（如 Kot 2001）。一种称为 Verhulst 过程的模型是 $N_{t+1} = r N_t (1 - N_t/K)$, $r>0$, $K>0$（式 2.5），看起来像是连续 logistic 增长模型的离散类比但其实不然：稳态不是 $N = K$。该模型的一个明显缺点是若 $N_t > K$ 则 $N_{t+1} < 0$。从连续 Verhulst 方程导出它的更恰当方式（图 2.11 图注中亦见）是将导数 $dN/dt$ 替换为步长 1 的差分形式得 $N(t+1) - N(t) = r N(t)(1 - N(t)/K) \Rightarrow N(t+1) = (1 + r - r N(t)/K) N(t)$（式 2.6）；以 $N(t) = ((1+r)/r) K x(t)$ 重缩放并设 $1+r = r'$，最后方程化为与式 2.2 同形的 $x(t+1) = r' x(t)(1 - x(t))$（式 2.7）。更现实的模型应满足大 $N_t$ 时增长率降低但 $N_{t+1}$ 仍非负；图 2.2 中 $f(N_t)$ 的定性形式即一例。一种常用模型是 Ricker 曲线 $N_{t+1} = N_t \exp[r(1 - N_t/K)]$, $r>0$, $K>0$（式 2.8），可看作式 2.2 的修正，加入一个对大 $N_t$ 更严厉的死亡因子 $\exp(-r N_t/K)$；只要 $N_0 > 0$ 即对所有 $t$ 有 $N_t > 0$。由于 $t$ 以离散步长增加，种群登记变化本身就有某种内在延迟，因此这些差分方程与第 1 章讨论的延迟微分方程之间存在某种启发式联系——视延迟长短可能产生振荡解。由于将通式 (2.1) 中时间步长缩放为 1，其他参数就是决定解是否周期的控制因素。对式 2.5 和式 2.8 而言，决定性参数是 $r$（$K$ 可通过 $N_t \to N_t/K$ 缩放出去）。

2.2 蛛网作图：图解求解法（Cobwebbing: A Graphical Procedure of Solution）

可以通过简单的图解法获取种群增长行为的丰富信息。考虑式 2.1 中 $f$ 如图 2.2 所示的情况。稳态 $N^*$ 是 $N^* = f(N^*) = N^* F(N^*)$ 的解，即 $N^* = 0$ 或 $F(N^*) = 1$（式 2.9）。一般只用式 2.9 中第一个解，第二个主要用以强调 $N^* = 0$ 总是稳态。图解中稳态是曲线 $N_{t+1} = f(N_t)$ 与直线 $N_{t+1} = N_t$ 的交点；图 2.3(a) 给出曲线最大值 $N_m > N^*$ 时的情况。解 $N_t$ 的动态演化可用下述图解方式得到：从图 2.3(a) 中 $N_0$ 出发，沿 $N_{t+1}$ 轴移动至与曲线 $N_{t+1} = f(N_t)$ 相交即得 $N_1 = f(N_0)$；再用直线 $N_{t+1} = N_t$ 把 $N_1$ 作为新的 $N_0$；依次得 $N_2, N_3, \dots$，箭头显示路径。路径只是直线 $N_{t+1} = N_t$ 上的一系列反射。可见 $t \to \infty$ 时 $N_t \to N^*$，且单调趋于 $N^*$，如图 2.3(b) 所示。若从 $N_0' > N^*$ 出发，再次 $N_t \to N^*$ 且第一步之后仍单调。若起始足够接近稳态 $N^*$，只要曲线以适当方式与 $N_{t+1} = N_t$ 相交，趋近就是单调的；这意味着 $0 < |df(N_t)/dN_t|_{N_t=N^*} = f'(N^*) < 1$（式 2.10）。$f'(N^*)$ 是重要的参数（即系统在稳态 $N^*$ 处的特征值）；任何关于 $N^*$ 的小扰动简单衰减到零，故 $N^*$ 是线性稳定平衡态。

设 $f(N_t)$ 使平衡态 $N^* > N_m$（如图 2.4）。种群动力学行为关键取决于 $N^*$ 处曲线相交的几何形状（图 2.4(a)(b)(c) 内嵌放大显示）：分别为 $-1 < f'(N^*) < 0$, $f'(N^*) = -1$, $f'(N^*) < -1$。解 $N_t$ 在 $N^*$ 附近振荡。若振荡幅度递减而 $N_t \to N^*$，则 $N^*$ 稳定（图 2.4(a)）；若振荡幅度增大则不稳定（图 2.4(c)）。图 2.4(b) 中振荡是周期的，提示 $N_{t+1} = f(N_t)$ 存在周期解。稳态在受到小扰动不趋于零时严格说是不稳定的。三种情形下的种群动力学行为如图 2.5 所示。

参数 $\lambda = f'(N^*)$（即平衡态 $N^*$ 在 $N_{t+1} = f(N_t)$ 下的特征值）对于确定稳态附近局部行为至关重要。行为清晰而决定性的情形是 $0 < \lambda < 1$（图 2.3(a)）以及 $-1 < \lambda < 0$ 和 $\lambda < -1$（分别为图 2.4(b)(c)）。平衡态在 $-1 < \lambda < 1$ 时稳定，称为吸引平衡态。临界分岔值 $\lambda = \pm 1$ 处解的行为性质发生变化。$\lambda = 1$ 时曲线 $N_{t+1} = f(N_t)$ 在稳态处与 $N_{t+1} = N_t$ 相切，显然称为切分岔。$\lambda = -1$ 时原先称为 pitchfork bifurcation，但现在称为 period-doubling bifurcation（周期倍化分岔）。图 2.3, 2.4, 2.6 中可清楚看出此图解法之所以命名为"蛛网作图"（cobwebbing）。对于类型 (2.1) 的单方程，这是用于提示 $N_t$ 动力学行为的极为有用的程序；虽主要集中在稳态附近的局部行为，它也给出定量全局行为。若稳态不稳定，可预兆这类方程解表现出的特有行为。例：$\lambda = f'(N^*) < -1$（即 $N^*$ 附近局部行为如图 2.4(c)）；对此情形作蛛网图（图 2.6），由于没有方式产生更大的 $N_t$（虽然可以由更大的初值开始），解不能趋于 $N^*$，种群必被 $N_{\max}$ 所界；故解全局有界但不趋于任何稳态，事实上它看上去以表面随机的方式游走（图 2.6(b)）。具有这种行为的解称为混沌解。配合蛛网作图所示的多种类型解以及特征值 $\lambda$ 在特殊临界值处暗示的敏感性，现在必须以解析方式考察这些方程。图解方法所提示的结果对分析非常有用。

2.3 离散 logistic 型模型：混沌（Discrete Logistic-Type Model: Chaos）

作为具体实例考虑非线性 logistic 型模型 $u_{t+1} = r u_t (1 - u_t)$, $r > 0$（式 2.11），设 $0 < u_0 < 1$，关注 $u_t \ge 0$ 的解。相对连续微分 logistic 模型，这里的 $r$ 严格说是 $1 + r$。稳态和对应特征值 $\lambda$ 为 $u^* = 0$, $\lambda = f'(0) = r$；以及 $u^* = (r-1)/r$, $\lambda = f'(u^*) = 2 - r$（式 2.12）。$r$ 从零开始增大但 $0 < r < 1$ 时，唯一的现实（即非负）平衡态是 $u^* = 0$，因 $0 < \lambda < 1$ 故稳定。这从式 2.11 在 $0 < r < 1$ 时的蛛网图或解析地由 $u_1 < u_0 < 1$ 及 $u_{t+1} < u_t$ 对所有 $t$ 成立可清楚看出，故 $t \to \infty$ 时 $u_t \to 0$。

第一次分岔出现在 $r = 1$：因 $r > 1$ 时 $u^* = 0$ 的特征值 $\lambda > 1$ 而变为不稳定；同时正稳态 $u^* = (r-1)/r > 0$ 在 $1 < r < 3$ 时 $-1 < \lambda < 1$ 故稳定。第二次分岔在 $r = 3$ 时，$\lambda = -1$。此时 $f'(u^*) = -1$，$u^*$ 附近局部处于图 2.4(b) 的情况，呈现周期解。

为看清 $r$ 通过分岔值 $r = 3$ 时发生什么，先引入迭代记号（式 2.13）：$u_1 = f(u_0)$, $u_2 = f(f(u_0)) = f^2(u_0)$, $\dots$, $u_t = f^t(u_0)$。对式 2.11 而言，第一次迭代即式 2.11 本身，第二次迭代为 $u_{t+2} = f^2(u_t) = r[r u_t(1 - u_t)][1 - r u_t(1 - u_t)]$（式 2.14）。图 2.7(a) 显示第一次迭代随 $r$ 的变化效果；特征值 $\lambda = f'(u^*)$ 随 $r$ 增大而减小，$r = 3$ 时 $\lambda = -1$。现在考察第二次迭代 (2.14)，问它是否有平衡态（即 $u_{t+2} = u_t = u^*_2$ 的解）。稍作代数得 $u^*_2 [r u^*_2 - (r-1)] [r^2 (u^*_2)^2 - r(r+1) u^*_2 + (r+1)] = 0$（式 2.15），其解为 $u^*_2 = 0$；或 $u^*_2 = (r-1)/r > 0$（若 $r > 1$）；或 $u^*_2 = (r+1 \pm [(r+1)(r-3)]^{1/2})/(2r) > 0$（若 $r > 3$）（式 2.16）。可见若 $r > 3$，式 2.11 的 $f(u_t)$ 给出的 $u_{t+2} = f^2(u_t)$ 还有两个新的实稳态。对应图 2.7(b) 中 A、B、C 是正平衡态 $u^*_2$，B 等于 $(r-1)/r$，位于式 2.16 中出现于 $r > 3$ 的两个新解之间。

可以将式 2.14 看作时间步长为 2 的模型中的第一次迭代。可计算 A、B、C 点处的特征值 $\lambda$。由图 2.7(b) 显然 $\lambda_B = f'(u^*_B) > 1$（$u^*_B$ 指 B 处的 $u^*_2$），A 和 C 类似。对刚大于 3 的 $r$，$-1 < \lambda_A < 1$ 且 $-1 < \lambda_C < 1$（由图可见或通过在 A 和 C 处对式 2.14 计算 $\partial f^2(u_t)/\partial u_t$ 得）。故第二次迭代 (2.14) 的稳态 $u^*_A$ 和 $u^*_C$ 稳定。这意味着第二次迭代 (2.14) 有稳定平衡态，等价于式 2.11 存在稳定的周期 2 解。换言之，若从 A 出发，2 次迭代后回到 A，即 $u^*_{A+2} = f^2(u^*_A)$ 但 $u^*_{A+1} = f(u^*_A) \neq u^*_A$；事实上 $u^*_{A+1} = u^*_C$ 而 $u^*_{C+1} = u^*_A$。

随着 $r$ 继续增大，A 和 C 处的特征值 $\lambda$ 通过 $\lambda = -1$，于是这两个 2 周期解变为不稳定。此时考察第 4 次迭代，正如可预期的那样，$u_{t+4}$ 作为 $u_t$ 的函数有 4 个峰（与图 2.7(b) 中的 2 个峰相比），出现 4 周期解。因此 $r$ 通过一系列分岔值时，解 $u_t$ 的特性经历一系列分岔——这里是周期解的倍周期。分岔情形如图 2.8(a) 所示。$\lambda = -1$ 处的这些分岔原本因图 2.8(a) 中产生的形状而称为 pitchfork bifurcation，但因为仅从 2 周期视角看才是叉形，故现在称为 period-doubling bifurcation。例如 $3 < r < r_4$（$r_4$ 是通向 4 周期解的分岔值）时，周期解位于图 2.8(a) 中两个 $u^*$ 之间（即竖直线通过 $r$ 值与平衡态曲线的两个交点）。图 2.8(b) 是 4 周期解的例子（$r_4 < r < r_8$），$u_t$ 值由竖直线与平衡态曲线的 4 个交点给出。

随着 $r$ 通过连续分岔值，每个偶数 $p$ 周期解分支进入 $2p$ 周期解，发生于 $p$ 周期解的特征值通过 $-1$ 时。$r$ 空间中分岔间距越来越小：启发地看这是合理的，因为更高阶迭代意味着更多峰（图 2.7(b)），所有这些峰都装入同一区间 $(0, 1)$。故对每个 $n$，存在周期 $2^n$ 解的层级，每个对应一个参数区间使该周期解稳定。在所有周期 $2^n$ 解都失稳的极限值 $r_c$ 处，对 $r > r_c$ 所有的原始 $2^n$ 周期都不稳定。行为相当复杂。对 $r > r_c$，奇数周期解开始出现，简单 3 周期解在 $r \approx 3.828$ 时最终出现，并出现局部吸引的 $k, 2k, 4k, \dots$ 周期循环（但 $k$ 现在为奇数）。例如另一个稳定 4 周期解在 $r \approx 3.96$ 出现。

$r_c$ 是模型 (2.11) 中当奇周期解恰可能时的临界参数值。当第三次迭代有 3 个与直线 $u_{t+3} = u_t$ 相切的稳态（即这些稳态处特征值 $\lambda = 1$）时，即出现 3 周期。该情形示意如图 2.9。对模型 (2.11)，$r \approx 3.828$。

Sarkovskii 定理（1964）：关于一维映射的一篇重要论文，有重大实际意义，与图 2.9 直接相关。他证明（其中包括）：若对某个 $r_3$ 值存在奇（$\ge 3$）周期解，则对 $r > r_3$ 存在非周期或混沌解。此处的分岔（$r_3$ 处）称为切分岔：名称提示了图 2.9 中所示情形。图 2.10 展示模型方程 (2.11) 对若干 $r$ 值的解，包括图 2.10(d)(f) 中的混沌例子。注意图 2.10(f) 中的行为：种群爆发、崩溃、缓慢恢复。

Sarkovskii 定理被 Stefan（1977）进一步推广。Li 和 Yorke（1975）的结果——若存在周期 3 解则对所有 $n \ge 1$ 都存在周期 $n$ 解——是 Sarkovskii 定理的特例。

虽然这里集中于 logistic 模型 (2.11)，但这类行为是形如 (2.1) 的差分方程模型的典型；它们都展示分岔到更高周期解最终导致混沌。

图 2.10(d)–(f) 展示通向混沌路径的一个有趣侧面：$r$ 从产生图 2.10(d) 非周期解的值进一步增大时，重新得到周期解（图 2.10(e)）；更大的 $r$ 时非周期解再次出现（图 2.10(f)）。故 $r$ 增大到首次出现混沌之后，参数值窗口处解行为是周期的。故存在周期性与非周期性交织的参数窗口。图 2.11 是迭代映射运行长时间后的典型图形：先运行几千次迭代稳定化，然后继续运行很多次迭代记录 $u_t$ 值。考察图 2.11 并考虑增大 $r$ 时对解的影响：$r_2 < r < r_4$ 时解 $u_t$ 仅在两点（如 A、B）之间振荡，对应竖直线通过 $r$ 值的交点；$r_4 < r < r_8$ 时 $u_t$ 展示 4 周期解，值仍为竖直线与曲线的交点；$r_c < r < r_p$ 区间解是混沌的；$r_p$ 之后小窗口再次出现规则周期解，然后又是非周期的——非周期–周期–非周期序列重复。考察插图（对 $r$ 方向放大倍数大于 $u_t$ 方向的小矩形放大）可见同样的分岔序列以分形方式重复。第 14 章给出了分形的简要介绍，第 3 章 3.9 节讨论了其在生物背景中的应用。Peitgen 和 Richter（1986）的著作以二维模型展示了令人眼花缭乱的分形序列和图形；第 5 章将讨论一个实际应用。

混沌相关行为研究兴趣日增、研究量大、应用广泛。流行媒体现在称之为混沌理论或新的非线性理论（这种容易识别的名称最能吸引公众关注，还有突变理论和分形理论）。兴趣不限于离散模型：它最初由常微分方程组——Lorenz 系统（Lorenz 1963；参见 Sparrow 1982, 1986）——所展示。此类混沌研究产生了许多有趣而意外的结果，例如 Feigenbaum（1978）证明：若 $r_2, r_4, \dots, r_{2^n}, \dots$ 是周期倍化分岔值序列，则 $\lim_{n \to \infty} |r_{2(n+1)} - r_{2^n}|/|r_{2(n+2)} - r_{2(n+1)}| = \delta = 4.66920\dots$。他证明 $\delta$ 是普适常数；对形式 $u_{t+1} = f(u_t)$（$f$ 有如图 2.2 所示极大）的通用迭代映射，只要呈现周期倍化，$\delta$ 取相同值。

一个实用而快捷的判别混沌存在性的方法由 Li 等（1982）给出：他们证明若对某 $u_t$ 及任意 $f(u_t)$，存在奇整数 $n$ 使 $f^n(u_t;r) < u_t < f(u_t;r)$，则存在奇周期解，从而混沌存在。例如对 $u_{t+1} = f(u_t;r) = u_t \exp[r(1 - u_t)]$，$r = 3.0$, $u_0 = 0.1$，前几项计算给出 $u_7 = f^5(u_2) < u_2 < f(u_2) = u_3$，即上述不等式要求中 $n = 5$；故此 $f(u_t;r)$ $r = 3$ 时是混沌的。

2.4 稳定性、周期解与分岔（Stability, Periodic Solutions and Bifurcations）

所有相关种群模型至少含一个参数 $r$。由以上讨论，随着该参数变化，通用模型方程 $u_{t+1} = f(u_t;r)$（式 2.17）的解通常在 $r$ 的特定值处分岔。这类分岔可以到达具逐次更高周期的周期解，最终对 $r$ 大于某有限临界值 $r_c$ 时产生混沌解。由图解分析可见，这类分岔出现于相应特征值 $\lambda$ 通过 $\lambda = 1$ 或 $\lambda = -1$ 时。这里讨论与这些分岔相关的一些解析结果。为代数简洁常略去 $f(u_t;r)$ 中的 $r$ 写作 $f(u_t)$，但对参数的依赖始终理解。心中设想的函数 $f$ 定性上类似于图 2.2。

方程 (2.17) 的平衡点或不动点是 $u^* = f(u^*;r)$ 的解，即 $u^*(r)$（式 2.18）。为研究 $u^*$ 的线性稳定性，通常写 $u_t = u^* + v_t$, $|v_t| \ll 1$（式 2.19）。代入式 2.17 并对小 $v_t$ 展开（Taylor 展开）得 $u^* + v_{t+1} = f(u^* + v_t) = f(u^*) + v_t f'(u^*) + O(v_t^2)$, $|v_t| \ll 1$。因 $u^* = f(u^*)$，决定 $u^*$ 线性稳定性的线性（$v_t$ 一次）方程为 $v_{t+1} = v_t f'(u^*) = \lambda v_t$, $\lambda = f'(u^*)$，其中 $\lambda$ 是首次迭代 (2.17) 在不动点 $u^*$ 处的特征值。解为 $v_t = \lambda^t v_0$，当 $|\lambda| < 1$ 时 $t \to \infty$ 时趋于 0，$|\lambda| > 1$ 时趋于 $\pm \infty$。故 $u^*$ 在 $-1 < f'(u^*) < 1$ 时稳定，$|f'(u^*)| > 1$ 时不稳定（式 2.20）。若 $u^*$ 稳定，从该平衡态的任何小扰动衰减到零——若 $0 < f'(u^*) < 1$ 单调衰减，若 $-1 < f'(u^*) < 0$ 振荡衰减。若 $u^*$ 不稳定，任何扰动若 $f'(u^*) > 1$ 单调增长，若 $f'(u^*) < -1$ 振荡增长。这与之前图解论证完全一致。

作为例子，重缩放模型 (2.8) 为 $u_{t+1} = u_t \exp[r(1 - u_t)]$, $r > 0$（式 2.21）。稳态为 $u^* = 0$ 或 $1 = \exp[r(1 - u^*)] \Rightarrow u^* = 1$（式 2.22）。对应特征值 $\lambda_{u^*=0} = f'(0) = e^r > 1$（$r > 0$），故 $u^* = 0$ 不稳定（单调地）；$\lambda_{u^*=1} = f'(1) = 1 - r$（式 2.23）。故 $u^* = 1$ 在 $0 < r < 2$ 时稳定（若 $1 < r < 2$ 则振荡回到平衡）；$r > 2$ 时因振荡增长而不稳定。故 $r = 2$ 是首次分岔值。基于上述可预期 $r$ 通过分岔值 $r = 2$ 时，从 $u^* = 1$ 分岔出周期解。对 $|1 - u_t|$ 小，式 2.21 近似为 $u_{t+1} \approx u_t[1 + r(1 - u_t)]$，正是图 2.11 模拟的形式。写作 $U_{t+1} = (1 + r) U_t[1 - U_t]$，其中 $U_t = r u_t/(1 + r)$，就得到与 logistic 模型 (2.11) 同样形式（$r+1$ 替换 $r$）。那里首次分岔处出现稳定周期 2 解。对式 2.21 而言，下一分岔（到 4 周期解）发生在 $r = r_4 \approx 2.45$，6 周期解发生在 $r = r_6 \approx 2.54$，$r > r_c \approx 2.57$ 时出现非周期或混沌行为。$r$ 的连续分岔值越来越接近。$r > 2$ 时解对小 $r$ 变化极敏感——多数此类模型均如此，至少对 $r$ 超过前几个分岔值后。

经过 $t$ 次 $u_0$ 迭代后，$u_t = f^t(u_0)$，使用式 2.13 中记号。由 $u_0$ 生成的轨道是点集 $\{u_0, u_1, u_2, \dots\}$，其中 $u_{i+1} = f(u_i) = f^{i+1}(u_0)$, $i = 0, 1, 2, \dots$。称一点周期 $m$ 或 $m$ 周期若 $f^m(u_0;r) = u_0$ 且 $f^i(u_0;r) \neq u_0$（$i = 1, 2, \dots, m-1$）（式 2.24），$u_0$ 是映射 $f^m$ (2.24) 的不动点，是映射 $f$ (2.17) 的周期 $m$ 不动点。$u_0, u_1, \dots, u_{m-1}$ 形成一个 $m$ 周期循环。

对不动点（解）的稳定性要求特征值；对平衡态 $u^*$ 即为 $f'(u^*)$。现将该定义推广到 $m$ 周期点 $u_0, u_1, \dots, u_{m-1}$。方便起见引入 $F(u;r) = f^m(u;r)$, $G(u;r) = f^{m-1}(u;r)$。则 $m$ 周期的特征值 $\lambda_m$ 定义为 $\lambda_m = \partial f^m(u;r)/\partial u|_{u=u_i}$, $i = 0$ 或 $1$ 或 $\dots$ 或 $m-1$（式 2.25），$= F'(u_i;r) = f'(G(u_i;r)) G'(u_i;r) = f'(u_{i-1};r) G'(u_i;r) = f'(u_{i-1};r) \cdot \partial f^{m-1}(u_i;r)/\partial u|_{u=u_i}$，故 $\lambda_m = \prod_{i=0}^{m-1} f'(u_i;r)$（式 2.26），与 $i$ 无关。

综上，分岔出现在某参数值 $r_0$ 处若 $r < r_0$ 与 $r > r_0$ 时解的动力学性质发生定性变化。由以上讨论现在可预期这是从一个周期解到另一个具不同周期的周期解的分岔。当偶数周期序列分岔到奇周期解时，Sarkovskii（1964）定理说所有整数周期循环都存在，意味着混沌。$\lambda = -1$ 处的分岔是周期倍化分岔，$\lambda = 1$ 处是切分岔。利用现成代数操作的计算机程序包，容易计算各次迭代的特征值并生成 $r$ 的分岔值序列 (2.25) 或 (2.26)。也有系统解析方法，如 Gumowski 和 Mira（1980）；还有若干近似方法，如 Hoppensteadt 和 Hyman（1977）。推荐入门书 Strogatz（1994）。Cvitanović（1984）编辑的早期论文集展示了早期对混沌研究的兴趣；Holden（1986）的综述文章书；化学方面 Scott（1991）是良好起点。混沌还可用于通过叠加混沌掩码加密秘密消息——掩码仅发收双方知道，接收方解掩恢复原消息，Strogatz（1994）详细讨论。这仅展示了混沌被发现和研究的少数领域。

2.5 离散延迟模型（Discrete Delay Models）

迄今讨论的所有离散模型基于的假设是：物种在时间 $t$ 的每个成员都对 $t+1$ 时刻的种群有贡献（通式 (2.1) 或缩放版 (2.17) 蕴含此意）。这适用于多数昆虫，但对许多其他动物（如性成熟前有相当成熟期）并不成立。因此这些情形下种群动力学模型必须包含延迟效应——某种意义上即包含年龄结构。若该延迟（到成熟期）为 $T$ 个时间步，则研究差分延迟模型 $u_{t+1} = f(u_t, u_{t-T})$（式 2.27）。下文讨论的须鲸模型中延迟 $T$ 是数年量级。

为说明这类模型线性稳定性分析的问题并获取对延迟方程的预期，考虑如下简单模型（即便如此亦有实际意义）：$u_{t+1} = u_t \exp[r(1 - u_{t-1})]$, $r > 0$（式 2.28）。这是式 2.21 的延迟版本。平衡态仍为 $u^* = 0$ 和 $u^* = 1$。$u^* = 0$ 几乎可看即不稳定；在 $u^* = 0$ 附近线性化即知。

在 $u^* = 1$ 附近线性化，按通常方式令 $u_t = 1 + v_t$, $|v_t| \ll 1$，式 2.28 给出 $1 + v_{t+1} = (1 + v_t) \exp[-r v_{t-1}] \approx (1 + v_t)(1 - r v_{t-1})$，故 $v_{t+1} - v_t + r v_{t-1} = 0$（式 2.29）。寻找此差分方程形如 $v_t = z^t$ 的解，得 $z^2 - z + r = 0$，对 $z$ 给出两个值 $z_1$ 和 $z_2$，其中若 $r < 1/4$, $z_{1,2} = \frac{1}{2}[1 \pm (1 - 4r)^{1/2}]$；若 $r > 1/4$, $z_{1,2} = \rho e^{\pm i\theta}$（式 2.30），其中 $\rho = r^{1/2}$, $\theta = \tan^{-1}((4r-1)^{1/2})$（$r > 1/4$）。式 2.29 的解（其特征方程是 $z$ 的二次方程）为 $v_t = A z_1^t + B z_2^t$（式 2.31），其中 $A, B$ 为任意常数。

若 $0 < r < 1/4$，$z_1, z_2$ 为实数，$0 < z_1 < 1$, $0 < z_2 < 1$，故由式 2.31，$t \to \infty$ 时 $v_t \to 0$，因此 $u^* = 1$ 是线性稳定平衡态；此外小扰动后回到平衡是单调的。若 $r > 1/4$，$z_1, z_2$ 为复数且 $z_2 = \bar{z}_1$；又 $z_1 z_2 = |z_1|^2 = \rho^2 = r$。故对 $1/4 < r < 1$，$|z_1||z_2| < 1$。此时解为 $v_t = A z_1^t + B \bar{z}_1^t$，因实数解必须 $B = \bar{A}$，结合式 2.30，实解为 $v_t = 2|A| \rho^t \cos(t \theta + \gamma)$, $\gamma = \arg A$, $\theta = \tan^{-1}((4r-1)^{1/2})$（式 2.32）。$r \to 1$ 时 $\theta \to \tan^{-1}\sqrt{3} = \pi/3$。当 $r$ 通过临界值 $r_c = 1$ 时，$|z_1| > 1$，故 $t \to \infty$ 时 $v_t$ 无界增长，$u^*$ 不稳定。因 $r \approx 1$ 时 $\theta \approx \pi/3$, $v_t \approx 2|A| \cos(t\pi/3 + \gamma)$（周期为 6），可预期式 2.28 的解在 $r$ 刚大于 $r_c (=1)$ 时呈现 6 周期。图 2.12 显示三个 $r > 1$ 值的计算解。图 2.12(b) 中仍存在 6 周期元素但已不规则；图 2.12(c) 中 6 周期元素消失，解变得更尖刺状（spikelike）——这是混沌的早期征兆。

上章中可见延迟有失稳效应且随延迟增大而加剧。离散模型中延迟也有类似失稳效应，比较图 2.10 和图 2.12 的 $r$ 值即可看出：前者临界 $r_c = 2$，解分岔到 2 周期解；后者（延迟情形）临界 $r_c = 1$，分岔到 6 周期解。同样，延迟越长失稳效应越大。这无疑是下例建模和分析引发关注的原因。更高周期解常以大幅度种群摆动为特征，若从前一次极高的种群崩溃回落到低水平足够剧烈，灭绝是完全可能的。2.7 节简要讨论了通向灭绝的可能路径。

作为本节结论，简短描述国际捕鲸委员会（IWC）用于须鲸的实际模型。IWC 的目标是管理鲸群实现持续产量、防止灭绝等；IWC 面临巨大的商业和文化压力。要现实地履行其章程要求，必须理解鲸群增长动力学和生态学。

IWC 使用的（现受保护的）须鲸模型基于成熟鲸群在时间 $t$ 的离散延迟模型 $N_{t+1} = (1 - \mu) N_t + R(N_{t-T})$（式 2.33）。这里 $(1 - \mu) N_t$（$0 < \mu < 1$）是存活到下一年对种群有贡献的鲸群份额；$R(N_{t-T})$ 是 $T$ 年前出生的、补充到成年种群的数目。延迟 $T$ 是达到性成熟的时间，量级 5–10 年。该模型假设性别比为 1 且各性别死亡率相同。模型关键在于补充项 $R(N_{t-T})$ 的形式，在 IWC 模型（IWC 1979）中为 $R(N) = \frac{1}{2}(1 - \mu)^T N [P + Q (1 - (N/K)^z)]$（式 2.34）。其中 $K$ 是未捕捞时的平衡密度，$P$ 是 $N = K$ 时雌性的每头生育力，$Q$ 是种群密度降到低水平时生育力的最大可能增量，$z$ 是该密度被感知的严重程度度量。$1 - \mu$ 是新生鲸每年存活概率，故 $(1 - \mu)^T$ 是经过所需 $T$ 年存活到成年的比例；$1/2$ 是因为一半鲸为雌性，故雌性生育力须乘以 $N/2$。该具体模型由 Clark（1976a）详细研究。Getz 和 Haight（1989）讨论了渔业管理的一般模型。

式 2.33 和 2.34 中的参数 $\mu, T, P$ 并非独立。平衡态为 $N^* = N_{t+1} = N_t = N_{t-T} = K$，由之得 $\mu = \frac{1}{2}(1 - \mu)^T P = h$（式 2.35），除定义 $h$ 之外，将生育力 $P$ 与死亡率 $\mu$ 和延迟 $T$ 联系起来。独立测量这些量给出大致一致性检验。现以 $u_t = N_t/K$ 重缩放模型，式 2.33 配合式 2.34 变为 $u_{t+1} = (1 - \mu) u_t + h u_{t-T} [1 + q(1 - u_{t-T}^z)]$（式 2.36），其中 $h$ 由式 2.35 定义，$q = Q/P$。在稳态 $u^* = 1$ 附近线性化，写 $u_t = 1 + v_t$ 得扰动方程 $v_{t+1} = (1 - \mu) v_t + h(1 - qz) v_{t-T}$（式 2.37）。设 $v_t \propto s^t$，得 $s^{T+1} - (1 - \mu) s^T + h(qz - 1) = 0$（式 2.38），这是特征方程。稳态在 $|s| > 1$ 时失稳。此处有 4 个参数 $\mu, T, h, qz$，分析围绕式 2.38 根的研究展开（参见 Clark 1976b）。虽然根的形式复杂，但可以通过 Jury 条件（见附录 B）确定参数使 $|s| < 1 的条件——Jury 条件是实多项式系数须满足的不等式，使根的模小于 1。对阶数大于约 4 的多项式，条件繁琐得难以处理。$|s| > 1$ 时（正如可预期的）式 2.33 的解展示到逐次更高周期解的分岔，最终导致混沌；响应参数 $z$ 是关键的。

混沌与数据：混沌并不是描述确定性方程（如本章所讨论）解的看似随机混沌行为的特别好的名称。当我们观察复杂实验数据并寻求用简单模型建模时，是在暗示其潜在机制实际上相当简单。故在面对真实数据时，重要的是知道该随机性质究竟是真正随机的还是本章所讨论的确定性意义上的混沌的。不出意外，这是一个困难且有争议的问题。虽然我们可能对生成数据和过程的机制有某种生物学洞见，但不太可能了解得足够确定以写下机制的精确模型。已经发展出若干方法试图判定数据究竟是否随机的还是确定性混沌的，但无一万无一失。

为体会这种困难，设有数据点 $N_t$（度量离散时间 $t$ 的某种种群）。若画 $N_t$ 对 $N_{t+1}$ 的图得到相对光滑的曲线（如定性如图 2.2），那么合理地可对生成机制提出确定性模型，即本章讨论的、能产生确定性混沌的模型。换言之我们在为式 2.1 中的 $f(N_t)$ 寻找定性形式。然而若不能给出任何合理的曲线，不能断定潜在机制不是确定性的。例如本节所见延迟可以自然地出现在更新过程中；那种情况下或许可以做 $N_{t-1}$ 与 $N_t$ 对 $N_{t+1}$ 的三维图。若得到相对光滑的曲面，可能是确定性机制。同样若在该空间仍给出随机分布的点，不必然指向非确定性模型，因为 $N_t$ 与 $N_{t+1}$（或 $N_{t-1}$ 或其他更早时间种群值）之间的关系可能只是更复杂的离散模型或涉及多个延迟。从数据中寻求确定关系时选择几乎无限。坚实的生物学知识当然能大幅减少可能性。故一种方法是先验地尝试确定一个合理模型，并且若似乎在任一时间步只涉及 $N_t$ 和 $N_{t+1}$，则有时可用数据确定 $N_t$ 与 $N_{t+1}$ 函数关系的定量细节。一个惊人成功的例子出现在婚姻互动和离婚预测这一不太可能的领域，将在第 5 章讨论；参见 Cook 等（1995）和 Gottman 等（2002）关于一般婚姻理论的著作。这里离散耦合方程构成预备模型。混沌解给出生物学过程洞见的另一个完全不同例子由 Cross 和 Cotton（1994）给出，将在 2.8 节讨论其问题与模型分析。

2.6 渔业管理模型（Fishery Management Model）

离散模型在渔业管理中已使用相当长时间。它们常被证明可用于评估各种捕捞策略以优化经济产量并维持产量。但第 1 章 1.6 节末的评论应当谨记。Clark（1976b, 1985, 1990）、Goh（1982）、Getz 和 Haight（1989）、Hilborn 和 Mangel（1997）、Cohen（1987）编辑的论文集、Levin（1994）的相关章节是若干相关著作。以下模型原则上适用于任何被收获的可再生资源；详细分析适用于任何动力学可由离散模型描述的种群。

设无捕捞时种群密度由 $N_{t+1} = f(N_t)$ 支配。令 $h_t$ 为时间 $t$ 从种群取出（产生 $t+1$ 时刻下一代）的收获，则种群动力学模型为 $N_{t+1} = f(N_t) - h_t$（式 2.39）。这里讨论的问题是：(i) 最大持续生物产量是多少？(ii) 最大经济产量是多少？平衡态下 $N_t = N^* = N_{t+1}$, $h_t = h^*$，由式 2.39 得 $h^* = f(N^*) - N^*$（式 2.40）。最大持续稳态产量 $Y_M$ 在 $N^* = N_M$ 处，由 $\partial h^*/\partial N^* = 0 \Rightarrow f'(N^*) = 1$ 且 $Y_M = f(N_M) - N_M$（式 2.41）。当然唯一值得关注的情形是 $Y_M \ge 0$。

一种管理策略可以简单地维持种群以获得最大产量 $Y_M$。由于实际鱼群难以掌握，做到这点很困难。已知的是实际产量及投入的努力量。故更好的做法是用产量和努力量来表述优化问题。假设单位努力捕鱼从种群 $N$ 中捕获 $cN$，常数 $c$ 为"可捕性"参数，与种群密度无关。则将 $N$ 减少 1 单位的努力量为 $1/(cN)$，将 $f(N)$ 减少 1 单位的努力量为 $1/(cf(N))$。故提供产量 $Y_M = f(N_M) - N_M$ 的努力量 $E_M$ 为 $E_M = \sum_{N_i = N_M}^{f(N_M)} (c N_i)^{-1}$。若 $cN$ 远大于 1 单位，可以将上式中的求和近似为积分得 $E_M \approx (1/c) \int_{N_M}^{f(N_M)} N^{-1} dN = (1/c) \ln[f(N_M)/N_M]$（式 2.42）。两式 (2.41) 和 (2.42) 给出 $E_M$ 与 $Y_M$ 通过 $N_M$ 的参数关系。

举例，设未捕捞动力学由 $N_{t+1} = f(N_t) = b N_t/(a + N_t)$（$0 < a < b$）支配，则 $N_M: 1 = f'(N_M) = ab/(a + N_M)^2 \Rightarrow N_M = a^{1/2}(b^{1/2} - a^{1/2})$。代入式 2.41 和 2.42 得 $Y_M = [b N_M/(a + N_M) - N_M]$, $E_M = (1/c) \ln[b/(a + N_M)]$（式 2.43）。本例中消去 $N_M$ 可得 $Y_M$ 与 $E_M$ 之间的显式关系 $Y_M = [b \exp(-c E_M) - a][\exp(c E_M) - 1]$（式 2.44）。图 2.13(a) 说明 $Y_M$-$E_M$ 关系。利用这一点，管理策略的一个关键之处是注意到：若增加努力量反而减少产量，则最大持续产量已被超越，必须减小努力量使种群恢复。努力量随后可重新调整以尝试获得 $Y_c$（对应 $E_c$），两者都可由式 2.44 计算。该分析是关于最大持续生物产量的。最大经济产量必须包括收获的价格和努力量的成本。作为初步模型可将二者纳入经济回报表达式 $R = p Y_M - k E_M$，其中 $p$ 是单位产量的价格，$k$ 是单位努力量的成本。利用式 2.43 中 $Y_M(N_M)$ 和 $E_M(N_M)$，得到 $R(N_M)$ 并对其求极大。故最大回报 $R$ 作为努力量 $E$ 的函数得曲线；如图 2.13(b) 所示。

这样的"模型"结果在未经实验观察佐证前不应被过于严肃地对待；但它们可以给出重要定性提示。本分析基于被收获的种群具有稳态的事实。鱼类尤其具高每头增长率，在详细模型中这与参数 $r$ 相关。故可预期鱼群展示周期波动，这已知是事实。可能增长率足够高使行为在某些情形下处于混沌区域。由于收获从某种意义上说是繁殖率的有效降低，它有可能起稳定作用——例如从混沌变为周期甚至稳态。

2.7 生态学含义与注意事项（Ecological Implications and Caveats）

建模种群动力学的主要原因是理解控制特征并能在环境参数变化时预测可能的发展模式。在建立模型时可以不同程度地拥有物种的生物学知识及观察数据，用以与分析结果比较。总结本章所考虑类型模型能告诉我们种群动力学的哪些方面、并指出若干困难与限制可能有益。

一旦得到种群增长动力学的合理模型，就能确定全局动力学。利用图解法可看出当主要环境参数变化时解的变化。例如从图 2.4 可见，若从低种群开始，会先增长一段时间，然后看似准规则振荡，然后稳定到常数状态，或展现周期行为，或看似随机地大幅振荡——某阶段种群很高、下一时间步崩溃到很低密度。无论何种模型，只要它有图 2.6 的一般形式，种群密度最终总有界。

这种看似随机的动力学从建模角度提出严重问题：表现出这种行为的数据究竟由确定性模型生成还是随机情形生成？这是决定哪个更恰当的问题，且可能不能在具体情形中真正解决。但建模能做的，是指出种群动力学对环境参数变化的敏感程度（环境参数估计常常困难且通常重要）。

以 $f(N_t)$ 如图 2.6 形式所展现的动力学类型可见：长时间后种群总被某最大值 $N_{\max}$ 和最小值 $N_{\min}$ 所界；若 $N_0$ 足够小，前几次迭代可低于 $N_{\min}$。以图 2.6 为例，$N_{\max}$ 由 $df/dN_t = 0 \Rightarrow N_m$, $N_{\max} = f(N_m)$ 给出；$N_{\min}$ 是 $N_{\max}$ 的下一迭代，即 $N_{\min} = f(N_{\max}) = f(f(N_m)) = f^2(N_m)$（式 2.45）。这些最终极限种群大小对给定模型易于算出。例如对 $N_{t+1} = f(N_t) = N_t \exp[r(1 - N_t/K)]$, $f'(N_t) = 0 \Rightarrow N_m = K/r$, $N_{\max} = f(N_m) = K/r \cdot e^{r-1}$, $N_{\min} = f(f(N_m)) = K/r \cdot \exp[2r - 1 - e^{r-1}]$（式 2.46）。

随着 $f(N_t)$ 动力学曲线在 $N_t > N_m$ 后的陡降行为，种群剧烈下降到接近 $N_{\min}$ 的低值提出物种灭绝的问题。若种群下降到 $N_t < 1$ 则显然灭绝。事实上若 $N_t$ 下降到低值灭绝几乎不可避免。此阶段需要随机模型。然而可以对何时种群降到 1 或以下（即灭绝）作出估计，方法是对给定模型求 $N_{\min}$。条件即式 2.45 的 $N_{\min} = f^2(N_m) \le 1$（式 2.47）。对式 2.46 的例子，条件为 $K/r \cdot \exp[2r - 1 - e^{r-1}] \le 1$。故若 $r = 3.5$ 且 $K < 1600$（约）时，种群最终将灭绝。该模型 (2.46) 的分析指出了一个重要现象：繁殖参数 $r$ 越大，$N_{\min}$ 越小，种群崩溃使物种灭绝的可能性越大。还需注意灾难性下降前的种群大小通常很大。对上述例子若 $r = 3.5$ 则几乎为 3500（由式 2.46）。灭绝概念的一个有趣且潜在实际的应用是引入害虫不育物种以控制数量；见下面习题 6。但该过程的高成本通常令人望而却步。

一类迄今未具体讨论但属于通类 (2.1) 的重要模型是展示 Allee 效应的模型。展示该效应的生物种群当种群低于某阈值水平 $N_c$ 时规模减小。一个典型密度依赖种群模型如图 2.14 所示。若初始种群 $N_0$ 使 $f^2(N_0) < N_c$ 则 $N_t \to 0$。这类模型常由捕食引起。第 1 章云杉蚜虫方程 (1.6) 的连续时间模型就具有这种行为。区域 $N_t < N_c$ 有时称为捕食坑。这里 $N_t = 0, N_c, N^*$ 都是稳态，$N_t = 0$ 稳定，$N_c$ 不稳定，$N^*$ 按 $f'(N^*)$ 符号可为稳定或不稳定。在此类动力学下，若 $N_t < N_c$ 则灭绝不可避免，无论 $N_c$ 多大。展示 Allee 效应的模型比上文考虑的模型表现出更丰富的行为谱——所有上述奇异振荡行为加上对任意 $m$ 若某迭代 $f^m(N_t) < N_c$ 时的灭绝可能性。

本章考虑的非线性离散模型的含义关键依赖于从观察数据分析得到的生物学参数。Southwood（1981）除其他外讨论了这些种群参数并给出若干物种的具体事实。Hassell 等（1976）分析了大量物种生活数据并拟合模型 $N_{t+1} = f(N_t) = r N_t/(1 + a N_t)^b$（$r, a, b$ 为正参数）；另见 Kot（2001）。$b > 1$ 时该 $f(N_t)$ 有图 2.2 中的单峰形式。例如 Colorado beetle 显然处于稳定周期区，而 Nicholson（1954）的 blowflies 可能处于混沌区。

最后应强调，这些解行为的丰富性是非线性模型的结果。同样有趣的是，许多定性特征可以由相当初等的方法找到，尽管它们呈现若干精致且具挑战性的数学问题。

2.8 肿瘤细胞生长（Tumour Cell Growth）

Cross 和 Cotton（1994）讨论了病理学中的一个数据问题：给定数据为一种种群（记为 $N_t$），由肿瘤细胞组成。分析中他们选取式 2.5 的简单 logistic 形式（$K = 1$），即 $N_{t+1} = r N_t (1 - N_t)$（式 2.48），其中 $r$ 反映肿瘤细胞的增长率。$N_t$ 归一化为 1 意味着 $N_t$ 是细胞培养容器中可维持的总细胞种群的比例。由 2.2–2.4 节分析知，$r < 3$ 时种群 $N_t$ 简单增长到稳态 $(r-1)/r$，若 $N_0$ 不太小则相对快速——例如 $N_0 = 0.001, r = 2$ 时种群每个时间步大致翻倍。$r > 3$ 时周期解出现，最终 $r > r_c$ 时产生混沌。处于混沌区时，任何时间 $t$ 的细胞种群本质上取决于初始条件。图 2.15 显示不同 $r$ 值的典型种群增长。图 2.15(b) 中 $N_t$ 趋近周期解但早期阶段也展示准 S 形增长曲线；图 2.15(c) 中解是混沌的。

假设模型中包含多克隆性，各细胞克隆有不同的初始种群。再设其增长率不同但都有 $r > r_c$，因此都展示混沌行为。病理学的主要兴趣在于肿瘤总大小，即细胞总数。Cross 和 Cotton（1994）首先考虑 5 个克隆并对种群求和得到总种群。一个典型结果如图 2.16(a)。开始看到混沌行为的某种平滑化以及图 2.15(a) 中种群 S 形状特征的初步出现。当包含 200 个克隆时，平滑效应显著得多（图 2.16(b)）。多克隆性在肿瘤生长中常见，数据展示图 2.16(b) 所示的增长模式。通过该简单例显然多克隆性能掩盖潜在的确定性混沌。该模型有诸多粗略假设，如假设每个克隆的增长率 $r$ 对所有时间不变。建模细胞分裂如何随时间变化本身就是一个有趣的问题，因为初始新生细胞群体从离散分裂到实质连续分裂的转变。它由 Murray 和 Frenzen（1986）讨论。多克隆情形下变化增长率提示年龄结构模型可能更恰当。研究多克隆系统（各克隆处于混沌增长区）的增长特征及可变增长率与年龄结构如何在实验观察中显现，是值得研究的。许多生物过程是混沌的，或即便不严格混沌、至少是随机的——但在神经学、病理学、生理学等领域观察到时（例如）仍表现出相当秩序。最近对癫痫的可能联系的综述由 Iasemidis 和 Sackellares（1996）给出。除其他外，他们认为其相关的理由是许多混沌系统在混沌区域——高度无序态——与高度有序区域之间有急剧间歇过渡；图 2.11 即典型例子。他们假设癫痫可能是混沌的一个例子，并以这种思路细致分析脑电图，对整个癫痫发生过程提供了新洞见。他们认为这可能在癫痫的诊断和治疗两方面都可用。

Larter 等（1999）关于混沌在脑活动中作用的有趣研究也提示混沌是常态，但在癫痫发作时活动变得异常规则。在部分癫痫发作中（这种患者最难响应药物），仅部分脑开始展示规则性，并相应地扩散和发作扩散。他们研究一千个互连神经元并对系统施加扰动以试图理解通信如何发生。除其他外，他们感兴趣的是影响规则行为区域向邻近混沌区域转移速率的因素。其目的也是将结果应用于治疗部分癫痫患者。关于心脏中的波活动（见第 1 卷第 1 章）则相反。若规则活动变得混沌（称为心颤的疾病），除非心脏能被电击回规则状态，否则是致命的——通常方法是使用大电流电击。

Cross 和 Cotton（1994）提到的另一个例子是人类毛发的正常生长是异步的，但某些情况下会被各种（通常是疾病）刺激同步。一个常见例子是妊娠和分娩期间所有毛发被同步到休止期（即细胞周期的休止阶段），后果是暂时性脱发——该症状称为休止期脱发（telogen effluvium，Benedict 等 1991）。从图 2.16 来看，随着克隆数目增加秩序增加，提示存在一种"反混沌"机制——该词在文献中越来越常出现。它与各种生物情形中发生的同步现象有相似之处。一个例子是培养中的某些细胞最初具不同细胞周期，但可用适当刺激诱导其同步。另一个是萤火虫的同步，第 9 章讨论生物振荡器时再论。一种完全不同的反混沌由 Benchetrit 等（1987）和 Demongeot 等（1987, 1996）在混沌呼吸模式分析中发现。后者使用 Aubin（1991）发展的可变性理论新领域的概念，展示了如何从潜在混沌（本章意义上的）中提取某种相干秩序。Demongeot 和 Jacob（1989）以及 Cosnard 和 Demongeot（1985）发展了关于吸引子和限制子的有趣方法。

本章个人批注

本章的核心主题是"差分方程如何从最简单的非线性映射开始产生极其丰富的动力学——稳态、周期解、周期倍化、混沌窗口与分形结构"。读完最深的印象是 logistic map $u_{t+1} = r u_t(1-u_t)$ 是少数几个以"显式初等方程"身份承载整套现代动力学语言（分岔、Lyapunov 指数、Feigenbaum 常数、Sarkovskii 排序）的模型。Murray 把它与生态学、渔业、肿瘤、癫痫、婚姻互动一连串真实或潜在的生物学应用串联起来，反映出 80–90 年代混沌研究跨学科的爆发式扩展。

我对几个具体技术点有几点疑问和评论： - 2.1 节末对 Fibonacci 数列与植物叶序的讨论虽然华丽，但严格来说它属于应用而非第 2 章核心离散模型。作者的论点（即"机制性而非遗传性"）是信念而非定理——Douady-Couder 的磁液滴实验是有力的支持但并不是证明。这一点放到批注是因为我不想把它误读为"数学证明"。 - 2.2 节 cobwebbing 的几何论证非常清晰，$|f'(N^*)| < 1$ 与 $-1 < f'(N^*) < 1$ 这两个等价条件在 2.4 节再次出现——值得把这两节当一组读。 - 2.3 节 Feigenbaum 常数 $\delta = 4.66920\dots$ 的"普适性"是一个惊人结论：所有具有单峰极大、经历周期倍化通向混沌的一维映射都具有相同的极限比例。这一点与统计物理中的重整化群不动点有同构关系，是少数几个从纯数学/计算中浮出来、后来在实验流体力学中也被验证的"普适类"之一。 - 2.4 节 m-cycle 特征值的链式法则 $\lambda_m = \prod_i f'(u_i)$ 推导漂亮，是后续讨论周期解稳定性的基础工具。 - 2.5 节离散延迟模型 $u_{t+1} = u_t \exp[r(1 - u_{t-1})]$ 在 $r = 1$ 分岔到 6 周期（而非 2 周期）——这是连续延迟模型在 2.4 章 1.3 节中没有的"新现象"：延迟不只延迟稳定化，更改变了分岔路径。 - IWC 须鲸模型的生物合理性被作者带过——参数 $P, Q, z$ 的实际测量难度和模型的拟合数据能力没有充分讨论，Clark（1976a）的更细致研究是必读项。 - 2.8 节肿瘤多克隆"混沌平滑"是一个反直觉的结果：总种群曲线变 S 形，看起来像确定性的 logistic 增长，但每个克隆单独看都是混沌。这意味着从实验数据判定"机制是确定性的"远比想象困难——这与本章末尾"混沌与数据"的讨论形成了完整的因果链。 - "反混沌"（antichaos）这个词在 1990 年代末曾经时髦，后来很快被吸收进"涌现同步"和"临界性"等更成熟的框架，今天的文献中已不常单列。 - 全章的反复出现的主题：简单非线性 + 离散化 = 复杂行为。它与第 1 章（连续 ODE）、第 6 章以后（化学反应动力学、振荡器）、第 8 章（Belousov-Zhabotinsky 反应）等形成完整的"模型离散性如何放大动力学复杂度"叙事线。

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

第 2 章在全书结构中位于"连续 → 离散"的关键转折点：第 1 章用连续 ODE 和延迟 ODE 处理单种群增长，并在 1.2 节（蚜虫模型）首次给出多稳态与滞后的几何图像、在 1.5 节（Cheyne-Stokes 与造血）展示延迟如何破坏稳定性产生极限环。本章则把时间步长设为 1，把"延迟"以"多步延迟"形式纳入 2.5 节，而把多稳态、滞后、极限环的故事推到一个新的高度——周期倍化、混沌、Sarkovskii 排序、Feigenbaum 普适性，并以离散 logistic map 为统一载体；这些主题将在后续章节中以更复杂的差分方程组（多物种相互作用第 3 章、捕食者–食饵第 4 章、婚姻互动第 5 章）和反应扩散系统（卷二第 2–6 章）继续展开。从方法论上看，本章是从"用图解法窥探一维动力学"过渡到"用解析方法（高次迭代、Jury 条件、特征值乘积）分析高周期和混沌"的桥梁；从生物学应用上看，本章末 2.6–2.8 节（渔业管理、生态学警示、肿瘤生长）建立了"模型行为—数据—管理决策"的完整反馈环，为后续多物种、空间模型的政策含义讨论奠定方法论基础。