Python实战：利用莱斯利模型预测种群动态变化

1. 莱斯利模型入门：从动物种群到Python代码

第一次听说莱斯利模型时，我正在研究一个有趣的生态问题：如何预测动物园里某种濒危动物的数量变化。传统的指数增长模型显然不够用，因为不同年龄段的动物繁殖能力和死亡率差异巨大。这时我发现了莱斯利模型这个神器，它完美解决了年龄结构对种群增长的影响问题。

莱斯利模型本质上是个离散时间的矩阵模型，由人口生态学家Patrick Leslie在1945年提出。它的精妙之处在于用矩阵运算来刻画种群动态：矩阵的每一行代表一个年龄组的生育率，每一列则对应存活率。举个例子，假设我们要研究某种昆虫，可以将其生命周期划分为幼虫、成虫和老虫三个阶段。通过这个模型，我们不仅能预测未来总数量，还能知道各年龄组的比例变化。

在实际应用中，莱斯利模型特别适合处理以下场景：

• 野生动物保护：预测濒危物种数量变化
• 农业养殖：优化鱼群或家畜的年龄结构
• 流行病研究：分析不同年龄段人群的疾病传播

# 一个简单的莱斯利矩阵示例 import numpy as np L = np.array([[0, 2, 1], # 生育率 [0.6, 0, 0], # 幼虫存活率 [0, 0.3, 0]]) # 成虫存活率

这个3×3矩阵就完整描述了三阶段种群的动态规律。第一行[0,2,1]表示成虫(第二列)每个体能产生2个后代，老虫(第三列)能产生1个。而次对角线上的0.6和0.3则代表年龄组间的存活转换率。

2. 模型搭建五步法：从数据到预测

2.1 数据收集与年龄分组

去年帮朋友分析池塘鱼群时，我深刻体会到数据收集的重要性。首先要确定最大寿命，比如某种鱼最长活5年，那么可以按年划分为5组。关键要获取两类数据：

1. 生育率(f)：每个年龄组平均产生的后代数
2. 存活率(s)：存活到下一个年龄组的概率

记得当时我们用了三个月跟踪记录，最终得到这样一组数据：

# 鱼类种群示例数据 age_groups = ["0-1年", "1-2年", "2-3年", "3-4年", "4-5年"] fertility = [0, 0.5, 1.2, 0.8, 0] # 生育率 survival = [0.2, 0.4, 0.6, 0.1] # 存活率(注意比年龄组少一个)

2.2 构建莱斯利矩阵

将上述数据转换为莱斯利矩阵时有几个技巧：

1. 第一行放生育率
2. 次对角线放存活率
3. 其余位置补零

def build_leslie_matrix(f, s): n = len(f) L = np.zeros((n, n)) L[0] = f # 设置生育率 for i in range(n-1): L[i+1,i] = s[i] # 设置存活率 return L fish_L = build_leslie_matrix(fertility, survival)

这个转换过程看似简单，但实际建模时我经常犯两个错误：一是搞混生育率和存活率的位置，二是忘记存活率数量比年龄组少一个。建议大家把这个函数保存为工具函数。

2.3 初始种群设置

初始种群向量的设置直接影响预测结果。在保护东北虎的项目中，我们就遇到初始数据不准确导致预测偏差的问题。好的做法是：

• 进行多次实地调查取平均值
• 对缺失数据使用相邻年龄组插值
• 记录数据采集时间和环境条件

# 初始种群向量的三种常见设置方式 X0_single = np.array([0,0,100,0,0]) # 只有特定年龄组有个体 X0_uniform = np.array([20]*5) # 均匀分布 X0_real = np.array([35,28,15,8,2]) # 实际调查数据

2.4 时间步长选择

时间步长的选择很有讲究。在研究昆虫种群时，开始我按天计算，结果矩阵太大计算困难。后来发现按发育阶段划分更合理：

• 寿命短的生物(如果蝇)：按天或周
• 中等寿命(如小鼠)：按月
• 长寿物种(如大象)：按年

# 不同时间步长的比较 daily_L = build_daily_matrix() # 每日变化 weekly_L = build_weekly_matrix() # 每周变化

2.5 模型验证技巧

模型建好后必须验证。我常用的三种方法：

1. 回溯验证：用历史数据检验预测准确性
2. 灵敏度分析：微调参数看结果稳定性
3. 对比现实：咨询领域专家评估合理性

def validate_model(L, X0, historical_data): predictions = [] current = X0 for _ in range(len(historical_data)): current = L @ current predictions.append(current) return np.mean(np.abs(predictions - historical_data))

3. Python实现详解：从基础到进阶

3.1 基础实现：NumPy版

使用NumPy实现莱斯利模型是最直接的方式。记得第一次实现时，我花了半天调试矩阵维度问题。下面是经过多次优化的版本：

import numpy as np class LeslieModel: def __init__(self, fertility, survival): self.L = build_leslie_matrix(fertility, survival) def predict(self, X0, steps): """预测未来多期的种群分布""" result = [X0] for _ in range(steps): result.append(self.L @ result[-1]) return np.array(result) def long_term_ratio(self): """计算长期年龄结构比例""" eigvals, eigvecs = np.linalg.eig(self.L) max_idx = np.argmax(eigvals) stable_dist = eigvecs[:,max_idx].real return stable_dist / stable_dist.sum()

使用时只需要：

fertility = [0, 0.8, 1.5] survival = [0.7, 0.2] model = LeslieModel(fertility, survival) predictions = model.predict([100,50,30], 10)

3.2 可视化分析

好的可视化能让结果一目了然。我常用的三种图形：

1. 种群趋势图：看总量变化
2. 年龄结构图：看比例变化
3. 热力图：看矩阵参数敏感性

import matplotlib.pyplot as plt def plot_population_trend(predictions): plt.figure(figsize=(10,6)) total = predictions.sum(axis=1) plt.plot(total, label='总数量') for i in range(predictions.shape[1]): plt.plot(predictions[:,i], label=f'年龄组{i+1}') plt.legend() plt.xlabel('时间周期') plt.ylabel('种群数量')

3.3 使用SymPy进行符号计算

当需要精确计算时，SymPy是更好的选择。在研究某个珍稀物种时，浮点误差导致结果异常，改用SymPy后问题解决：

import sympy as sp def symbolic_leslie(f, s): """创建符号计算的莱斯利矩阵""" n = len(f) L = sp.zeros(n) for j in range(n): L[0,j] = f[j] for i in range(n-1): L[i+1,i] = s[i] return L # 示例：精确计算特征值 L_sym = symbolic_leslie([0,sp.Rational(4),3], [sp.Rational('1/2'),sp.Rational('1/4')]) eigenvals = L_sym.eigenvals() # 精确特征值

3.4 性能优化技巧

处理大规模年龄组时，我总结了这些优化方法：

1. 使用稀疏矩阵(scipy.sparse)
2. 并行计算多场景
3. 预计算特征分解

from scipy.sparse import dia_matrix def build_sparse_leslie(f, s): """创建稀疏莱斯利矩阵""" n = len(f) diagonals = [f, [0]*(n-1) + [0], s + [0]] return dia_matrix((diagonals, [0,-1,-2]), shape=(n,n))

4. 实战案例：从建模到决策

4.1 案例背景：草原鼠害防治

去年参与的一个实际项目：某草原鼠害严重，需要预测未来5年的种群变化以制定防治计划。我们收集到以下数据：

• 最大寿命：3年
• 年龄组划分：幼鼠(0-1年)、成年鼠(1-2年)、老年鼠(2-3年)
• 生育率：[0, 5, 3]
• 存活率：[0.2, 0.1]

# 建立鼠害模型 rat_L = np.array([[0,5,3], [0.2,0,0], [0,0.1,0]]) rat_X0 = np.array([1000,500,200])

4.2 模型预测与分析

运行10个周期的预测后，发现了几个关键现象：

1. 前两年种群数量会下降(由于初始老年鼠较多)
2. 第三年开始指数增长
3. 长期年龄结构稳定在82:15:3

model = LeslieModel([0,5,3], [0.2,0.1]) predictions = model.predict(rat_X0, 10) # 计算关键指标 growth_rate = predictions[1:].sum(axis=1) / predictions[:-1].sum(axis=1) stable_dist = model.long_term_ratio()

4.3 防治方案模拟

我们测试了三种防治方案：

1. 灭幼方案：降低幼鼠存活率50%
2. 节育方案：降低成年鼠生育率50%
3. 综合方案：同时实施上述两种

# 测试不同防治方案 scenarios = { '基线': ([0,5,3], [0.2,0.1]), '灭幼': ([0,5,3], [0.1,0.1]), '节育': ([0,2.5,3], [0.2,0.1]), '综合': ([0,2.5,3], [0.1,0.1]) } results = {} for name, (f,s) in scenarios.items(): m = LeslieModel(f, s) results[name] = m.predict(rat_X0, 5).sum(axis=1)

4.4 结果可视化与决策

通过对比发现，虽然综合方案效果最好，但节育方案性价比更高。最终建议采用节育方案结合局部灭幼：

# 绘制不同场景对比图 plt.figure(figsize=(10,6)) for name, data in results.items(): plt.plot(data, label=name) plt.legend() plt.title('不同防治方案效果对比') plt.xlabel('年份') plt.ylabel('鼠群总量')

这个案例让我深刻体会到，好的数学模型不仅要准确预测，更要为决策提供清晰依据。莱斯利模型的价值不仅在于计算那几个数字，更在于帮助我们理解种群变化的驱动因素。