用Python模拟兔子和羊的“地盘争夺战”：手把手教你实现Lotka-Volterra竞争模型
用Python模拟兔子和羊的“地盘争夺战”：手把手教你实现Lotka-Volterra竞争模型
生态学中的物种竞争关系一直是研究者关注的焦点。想象一片广袤的草原，兔子和羊作为主要的食草动物，它们之间存在着微妙的竞争关系——争夺有限的草资源。这种竞争关系可以通过Lotka-Volterra竞争模型进行数学描述和预测。本文将带你用Python一步步实现这个经典模型，并通过可视化直观展示不同参数下种群动态的变化。
1. 环境准备与模型基础
在开始编码之前，我们需要先理解Lotka-Volterra竞争模型的基本原理。这个模型由Alfred J. Lotka和Vito Volterra在20世纪早期独立提出，用于描述两个物种在共享资源环境中的竞争关系。
模型的核心是两个微分方程：
def lotka_volterra(t, y, r1, r2, K1, K2, alpha12, alpha21): N1, N2 = y dN1_dt = r1 * N1 * (1 - (N1 + alpha12 * N2) / K1) dN2_dt = r2 * N2 * (1 - (N2 + alpha21 * N1) / K2) return [dN1_dt, dN2_dt]其中：
N1和N2分别表示兔子和羊的种群数量r1和r2是各自的内禀增长率K1和K2是环境承载力alpha12表示羊对兔子的竞争系数alpha21表示兔子对羊的竞争系数
要运行这个模型，我们需要安装几个Python库：
pip install numpy matplotlib scipy2. 模型实现与参数设置
现在让我们用Python实现这个模型。我们将使用SciPy的solve_ivp函数来求解微分方程。
首先定义模型参数：
import numpy as np from scipy.integrate import solve_ivp import matplotlib.pyplot as plt # 基础参数设置 params = { 'r1': 0.1, # 兔子的增长率 'r2': 0.08, # 羊的增长率 'K1': 500, # 兔子的环境承载力 'K2': 600, # 羊的环境承载力 'alpha12': 0.5, # 羊对兔子的竞争影响 'alpha21': 0.6 # 兔子对羊的竞争影响 } # 初始种群数量 initial_populations = [100, 150] # [兔子, 羊] # 时间范围 (0到200天) t_span = (0, 200) t_eval = np.linspace(0, 200, 1000) # 评估点接下来，我们定义模型函数并求解：
def competition_model(t, y, r1, r2, K1, K2, alpha12, alpha21): N1, N2 = y dN1_dt = r1 * N1 * (1 - (N1 + alpha12 * N2)/K1) dN2_dt = r2 * N2 * (1 - (N2 + alpha21 * N1)/K2) return [dN1_dt, dN2_dt] # 求解微分方程 solution = solve_ivp( competition_model, t_span, initial_populations, args=tuple(params.values()), t_eval=t_eval, method='RK45' )3. 结果可视化与分析
有了计算结果后，我们可以绘制种群数量随时间变化的曲线：
plt.figure(figsize=(12, 6)) plt.plot(solution.t, solution.y[0], label='兔子', color='brown') plt.plot(solution.t, solution.y[1], label='羊', color='gray') plt.xlabel('时间 (天)') plt.ylabel('种群数量') plt.title('兔子和羊的种群动态变化') plt.legend() plt.grid(True) plt.show()为了更直观地理解两个物种的竞争关系，我们可以绘制相图：
plt.figure(figsize=(8, 8)) plt.plot(solution.y[0], solution.y[1], 'b-') plt.xlabel('兔子数量') plt.ylabel('羊数量') plt.title('种群竞争相图') plt.grid(True) # 添加零增长等斜线 N1_range = np.linspace(0, params['K1']*1.2, 100) N2_range = np.linspace(0, params['K2']*1.2, 100) # 兔子的零增长线 plt.plot(N1_range, (params['K1'] - N1_range)/params['alpha12'], 'r--', label='兔子零增长线') # 羊的零增长线 plt.plot((params['K2'] - N2_range)/params['alpha21'], N2_range, 'g--', label='羊零增长线') plt.legend() plt.xlim(0, params['K1']*1.2) plt.ylim(0, params['K2']*1.2) plt.show()4. 探索不同竞争情景
Lotka-Volterra模型可以预测四种不同的竞争结果，取决于参数设置。让我们创建一个函数来探索这些情景：
def explore_scenarios(scenario): scenarios = { 'rabbit_wins': {'alpha12': 0.5, 'alpha21': 1.2}, 'sheep_wins': {'alpha12': 1.2, 'alpha21': 0.5}, 'coexistence': {'alpha12': 0.8, 'alpha21': 0.7}, 'unpredictable': {'alpha12': 1.2, 'alpha21': 1.2} } params.update(scenarios[scenario]) solution = solve_ivp( competition_model, t_span, initial_populations, args=tuple(params.values()), t_eval=t_eval ) plt.figure(figsize=(12, 5)) plt.subplot(1, 2, 1) plt.plot(solution.t, solution.y[0], label='兔子') plt.plot(solution.t, solution.y[1], label='羊') plt.title(f'{scenario} - 种群动态') plt.legend() plt.subplot(1, 2, 2) plt.plot(solution.y[0], solution.y[1]) plt.title(f'{scenario} - 相图') plt.xlabel('兔子数量') plt.ylabel('羊数量') plt.tight_layout() plt.show() # 示例：探索兔子获胜的情景 explore_scenarios('rabbit_wins')四种情景的关键参数设置如下表所示：
| 情景类型 | α₁₂ (羊对兔) | α₂₁ (兔对羊) | 结果描述 |
|---|---|---|---|
| 兔子获胜 | <1 | >1 | 兔子最终占据主导地位 |
| 羊获胜 | >1 | <1 | 羊最终占据主导地位 |
| 稳定共存 | <1 | <1 | 两个物种达到平衡共存 |
| 结果不可预测 | >1 | >1 | 初始条件决定最终优势物种 |
5. 模型扩展与交互式探索
为了让学习体验更加丰富，我们可以创建一个交互式工具来实时调整参数并观察结果：
from ipywidgets import interact, FloatSlider def interactive_simulation(r1, r2, K1, K2, alpha12, alpha21): params = locals() solution = solve_ivp( competition_model, t_span, initial_populations, args=(r1, r2, K1, K2, alpha12, alpha21), t_eval=t_eval ) plt.figure(figsize=(12, 5)) plt.subplot(1, 2, 1) plt.plot(solution.t, solution.y[0], label='兔子') plt.plot(solution.t, solution.y[1], label='羊') plt.legend() plt.title('种群动态') plt.subplot(1, 2, 2) plt.plot(solution.y[0], solution.y[1]) plt.title('相图') plt.xlabel('兔子数量') plt.ylabel('羊数量') plt.tight_layout() plt.show() interact( interactive_simulation, r1=FloatSlider(min=0.01, max=0.2, step=0.01, value=0.1), r2=FloatSlider(min=0.01, max=0.2, step=0.01, value=0.08), K1=FloatSlider(min=100, max=1000, step=50, value=500), K2=FloatSlider(min=100, max=1000, step=50, value=600), alpha12=FloatSlider(min=0.1, max=2.0, step=0.1, value=0.5), alpha21=FloatSlider(min=0.1, max=2.0, step=0.1, value=0.6) )6. 实际应用与注意事项
在实际生态研究中应用Lotka-Volterra模型时，有几个关键点需要注意：
参数估计：
- 增长率(r)可以通过野外观察或实验室实验估计
- 环境承载力(K)通常取决于资源可用性
- 竞争系数(α)是最难估计的参数，需要精心设计的实验
模型局限性：
- 假设环境是均匀且恒定的
- 忽略空间异质性和个体差异
- 不考虑进化适应和长期变化
改进方向：
- 添加随机因素模拟环境波动
- 引入空间维度考虑栖息地分布
- 扩展为多物种交互模型
以下是一个考虑环境随机性的改进模型示例：
def stochastic_competition(t, y, r1, r2, K1, K2, alpha12, alpha21, noise_scale=0.05): N1, N2 = y noise1 = np.random.normal(0, noise_scale) noise2 = np.random.normal(0, noise_scale) dN1_dt = r1 * N1 * (1 - (N1 + alpha12 * N2)/K1) + noise1 dN2_dt = r2 * N2 * (1 - (N2 + alpha21 * N1)/K2) + noise2 return [dN1_dt, dN2_dt]在生态保护实践中，理解物种竞争关系对于制定合理的保护策略至关重要。例如，在引入新物种时，可以通过类似模型预测其对本地物种的潜在影响。