SEIR 传染病模型 Python 实战：基于 2020 年新冠数据拟合与预测

传染病动力学模型是理解疫情传播规律的重要工具。在众多模型中，SEIR（Susceptible-Exposed-Infectious-Recovered）模型因其考虑了潜伏期人群而特别适合新冠疫情的建模分析。本文将带您用 Python 实现完整的 SEIR 模型构建、参数拟合和预测流程，并提供可直接运行的 Jupyter Notebook 代码。

1. SEIR 模型基础与数学原理

SEIR 模型将人群划分为四个互斥类别：

• 易感者 (S)：可能被感染但尚未感染的人群
• 暴露者 (E)：已被感染但尚未发病且无传染性的人群
• 感染者 (I)：具有传染性的患病个体
• 康复者 (R)：痊愈并获得免疫力的个体

模型的核心是一组非线性微分方程：

def seir_model(y, t, beta, sigma, gamma): S, E, I, R = y dSdt = -beta * S * I / N dEdt = beta * S * I / N - sigma * E dIdt = sigma * E - gamma * I dRdt = gamma * I return [dSdt, dEdt, dIdt, dRdt]

其中关键参数包括：

• β：有效接触率（感染者每天传染易感者的概率）
• σ：潜伏期倒数（1/σ 为平均潜伏期）
• γ：恢复率（1/γ 为平均感染期）

基本再生数 R₀是衡量传播能力的关键指标：

R₀ = β / γ

当 R₀ > 1 时疫情会扩散，R₀ < 1 时疫情将逐渐消退。

2. 数据准备与预处理

我们使用约翰霍普金斯大学系统科学与工程中心（JHU CSSE）提供的全球新冠疫情数据：

import pandas as pd # 下载原始数据 url = "https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_covid19_confirmed_global.csv" df = pd.read_csv(url) # 处理中国地区数据 china_data = df[df['Country/Region'] == 'China'].groupby('Country/Region').sum().iloc[:, 4:] daily_cases = china_data.diff(axis=1).fillna(0).iloc[0]

典型的数据清洗步骤包括：

1. 处理缺失值和异常值
2. 计算每日新增病例
3. 数据平滑处理（如7日移动平均）

日期 确诊病例 2020-01-22 548 2020-01-23 643 2020-01-24 920 ... ...

3. 模型参数估计与拟合

参数估计是建模中最具挑战性的环节。我们使用最小二乘法拟合模型参数：

from scipy.optimize import minimize from scipy.integrate import odeint def fit_seir(params, t, data): beta, sigma, gamma = params solution = odeint(seir_model, [N-1,1,0,0], t, args=(beta, sigma, gamma)) return np.sum((solution[:,2] - data)**2) initial_guess = [0.4, 0.2, 0.1] # β, σ, γ 的初始猜测值 bounds = ((0,1), (0,1), (0,1)) # 参数范围约束 result = minimize(fit_seir, initial_guess, args=(t, observed_cases), bounds=bounds)

参数估计的注意事项：

1. 使用滑动窗口方法处理参数时变特性
2. 考虑不同防控阶段参数的变化
3. 验证参数的生物合理性（如潜伏期应在2-14天范围内）

4. 完整建模流程实现

以下是完整的 Python 实现代码：

import numpy as np import matplotlib.pyplot as plt from scipy.integrate import odeint from scipy.optimize import minimize # 1. 定义SEIR模型 def seir_model(y, t, beta, sigma, gamma): S, E, I, R = y N = S + E + I + R dSdt = -beta * S * I / N dEdt = beta * S * I / N - sigma * E dIdt = sigma * E - gamma * I dRdt = gamma * I return [dSdt, dEdt, dIdt, dRdt] # 2. 加载和处理数据 # [此处添加数据加载和预处理代码] # 3. 定义损失函数 def loss_function(params, t, data): beta, sigma, gamma = params solution = odeint(seir_model, [N-E0-I0, E0, I0, 0], t, args=(beta, sigma, gamma)) predicted = solution[:,2] # 预测的感染人数 return np.sum((predicted - data)**2) # 4. 参数估计 initial_guess = [0.4, 1/5.2, 1/14] # R0≈2.8 bounds = ((0,1), (1/14,1/2), (1/21,1/7)) # 合理范围约束 result = minimize(loss_function, initial_guess, args=(time_points, observed_cases), bounds=bounds, method='L-BFGS-B') # 5. 模型验证与预测 best_params = result.x t_predict = np.arange(0, 180, 1) solution = odeint(seir_model, [N-E0-I0, E0, I0, 0], t_predict, args=tuple(best_params)) # 6. 可视化结果 plt.figure(figsize=(12,6)) plt.plot(t_predict, solution[:,0], label='Susceptible') plt.plot(t_predict, solution[:,1], label='Exposed') plt.plot(t_predict, solution[:,2], label='Infectious') plt.plot(t_predict, solution[:,3], label='Recovered') plt.scatter(observed_days, observed_cases, color='red', label='Observed Cases') plt.legend() plt.xlabel('Days') plt.ylabel('Population') plt.title('SEIR Model Fit to COVID-19 Data') plt.show()

5. 模型改进与扩展

基础 SEIR 模型可通过以下方式增强：

1. 考虑潜伏期传染性：

dSdt = -beta*S*I/N - beta2*S*E/N dEdt = beta*S*I/N + beta2*S*E/N - sigma*E

2. 加入时变参数：

def time_varying_beta(t): return beta0 * np.exp(-alpha*t) # 反映防控措施加强

3. 空间异质性扩展：

# 多地区耦合模型 def multi_regional_seir(y, t, params): # 包含地区间迁移项 dS1dt = -beta1*S1*I1/N1 + m*(S2-S1) # ...其他方程类似

6. 结果分析与应用

模型输出可支持多种决策场景：

1. 疫情趋势预测：

   • 峰值时间估计
   • 最终感染规模预测
   • 不同防控场景模拟

2. 防控措施评估：

   • 社交距离效果量化
   • 隔离政策优化
   • 医疗资源需求预测

3. 参数敏感性分析：

7. 局限性与挑战

1. 数据质量依赖：

   • 检测能力变化影响病例报告数
   • 无症状感染者难以统计

2. 模型简化假设：

   • 均质混合假设不现实
   • 忽略年龄结构和空间异质性

3. 长期预测不确定性：

   • 病毒变异影响参数
   • 人群行为动态变化

# 不确定性分析示例 num_simulations = 100 params_distribution = { 'beta': np.random.normal(beta, 0.1*beta, num_simulations), 'sigma': np.random.uniform(1/7, 1/3, num_simulations), 'gamma': np.random.uniform(1/21, 1/10, num_simulations) } plt.figure(figsize=(12,6)) for i in range(num_simulations): sol = odeint(seir_model, [N-E0-I0,E0,I0,0], t_predict, args=(params_distribution['beta'][i], params_distribution['sigma'][i], params_distribution['gamma'][i])) plt.plot(t_predict, sol[:,2], color='gray', alpha=0.1) plt.plot(observed_days, observed_cases, 'ro', label='Observed') plt.title('Uncertainty Analysis of SEIR Model') plt.xlabel('Days') plt.ylabel('Infectious Cases') plt.legend() plt.show()

通过完整的代码实现和系统分析，我们建立了一个可解释、可验证的疫情分析工具。这种数据驱动的方法不仅适用于历史数据分析，也能为未来公共卫生决策提供科学参考。