用Python模拟SIS模型从公式推导到可视化传播过程附完整代码在流行病学和信息传播研究中SISSusceptible-Infected-Susceptible模型是描述反复感染现象的经典工具。与一次性感染的SIR模型不同SIS模型更适用于模拟流感类疾病或社交媒体信息的传播模式——个体康复后仍可能再次被感染。本文将带您从零开始构建完整的Python模拟系统涵盖微分方程求解、参数调优和动态可视化三大核心环节。1. 环境准备与模型基础实现SIS模型模拟需要以下工具链科学计算栈NumPy用于数值运算SciPy处理微分方程可视化库Matplotlib基础绘图Plotly实现交互式图表辅助工具Jupyter Notebook作为实验环境可选安装依赖包pip install numpy scipy matplotlib plotlySIS模型的核心微分方程di/dt λ*i*(1-i) - μ*i其中λ感染率每个感染者每天传染他人的概率μ治愈率每天康复的概率i(t)t时刻感染人群比例关键参数对传播的影响参数组合传播特征稳态解λ μ持续传播i∞ 1 - μ/λλ ≤ μ自然消亡i∞ 02. 微分方程数值求解使用SciPy的odeint进行数值求解比解析解更适用于实际应用场景。以下是完整的求解实现import numpy as np from scipy.integrate import odeint def sis_model(i, t, λ, μ): 定义SIS微分方程 didt λ * i * (1 - i) - μ * i return didt # 参数设置 λ 0.3 # 感染率 μ 0.1 # 治愈率 i0 0.01 # 初始感染比例 t np.linspace(0, 100, 1000) # 100天内的1000个时间点 # 求解方程 solution odeint(sis_model, i0, t, args(λ, μ)) i solution.flatten()提示当λ≈μ时建议减小求解步长增加t的点数以获得更精确结果参数敏感性分析示例params [ (0.4, 0.1), # 高传播 (0.2, 0.15), # 低传播 (0.3, 0.3) # 临界状态 ] solutions [odeint(sis_model, i0, t, argsp) for p in params]3. 动态可视化实现3.1 Matplotlib静态图表基础感染曲线绘制import matplotlib.pyplot as plt plt.figure(figsize(10,6)) plt.plot(t, i, r-, labelInfected) plt.title(SIS Model Dynamics (λ0.3, μ0.1)) plt.xlabel(Days) plt.ylabel(Population Fraction) plt.legend() plt.grid(True) plt.show()多参数对比可视化技巧fig, ax plt.subplots(figsize(12,7)) for idx, (λ, μ) in enumerate(params): ax.plot(t, solutions[idx], labelfλ{λ}, μ{μ}, R0{λ/μ:.2f}) ax.set(xlabelTime (days), ylabelInfected Proportion, titleSIS Model Parameter Comparison) ax.legend() ax.grid(True)3.2 Plotly交互式可视化创建可缩放动态图表import plotly.graph_objects as go fig go.Figure() fig.add_trace(go.Scatter(xt, yi, modelines, nameInfected)) fig.update_layout( titleInteractive SIS Model, xaxis_titleDays, yaxis_titleProportion Infected, hovermodex unified ) fig.show()添加参数滑块控制from plotly.subplots import make_subplots fig make_subplots() buttons [] for λ_val in [0.2, 0.3, 0.4]: sol odeint(sis_model, i0, t, args(λ_val, μ)) fig.add_trace(go.Scatter( xt, ysol.flatten(), visibleFalse, namefλ{λ_val} )) # 设置按钮交互 for i, λ_val in enumerate([0.2, 0.3, 0.4]): buttons.append( dict(methodupdate, args[{visible: [False]*3}, {title: fSIS Model (λ{λ_val})}], labelfλ{λ_val}) ) buttons[i][args][0][visible][i] True fig.update_layout( updatemenus[dict(typedropdown, buttonsbuttons, x1.1, y1)], titleSIS Model with Parameter Selector )4. 构建完整模拟系统4.1 参数优化模块使用网格搜索寻找关键参数from itertools import product λ_range np.linspace(0.1, 0.5, 10) μ_range np.linspace(0.05, 0.3, 10) peak_infections [] for λ, μ in product(λ_range, μ_range): sol odeint(sis_model, i0, t, args(λ, μ)) peak_infections.append((λ, μ, max(sol)))4.2 实时模拟器基于IPython的交互式模拟from IPython.display import display import ipywidgets as widgets widgets.interact def sis_simulator(λ(0.1,0.5,0.01), μ(0.01,0.3,0.01)): sol odeint(sis_model, i0, t, args(λ, μ)) plt.figure(figsize(10,5)) plt.plot(t, sol) plt.title(fSIS Model (λ{λ}, μ{μ}, R0{λ/μ:.2f})) plt.ylim(0, 1) plt.show()4.3 网络传播模拟扩展对于复杂网络中的传播模拟需安装networkximport networkx as nx def network_sis_simulation(G, λ, μ, steps100): infected set(np.random.choice(G.nodes(), sizeint(len(G)*i0), replaceFalse)) history [] for _ in range(steps): new_infected set() # 感染传播 for node in infected: for neighbor in G.neighbors(node): if neighbor not in infected and np.random.rand() λ: new_infected.add(neighbor) # 康复过程 still_infected set() for node in infected: if np.random.rand() μ: still_infected.add(node) infected still_infected.union(new_infected) history.append(len(infected)/len(G)) return history5. 工程实践建议性能优化对于长期模拟建议使用Numba加速from numba import jit jit(nopythonTrue) def sis_model_numba(i, t, λ, μ): return λ * i * (1 - i) - μ * i结果验证检查稳态解是否满足理论值def verify_steady_state(λ, μ, tolerance1e-3): sol odeint(sis_model, i0, t, args(λ, μ)) theoretical max(0, 1 - μ/λ) if λ μ else 0 return abs(sol[-1] - theoretical) tolerance异常处理添加参数合法性检查def safe_sis_simulation(λ, μ, i0): assert 0 i0 1, Initial infection must be in [0,1] assert λ 0 and μ 0, Rates must be non-negative return odeint(sis_model, i0, t, args(λ, μ))实际项目中遇到的典型问题包括当μ0时模型退化为SI模型需特殊处理极小的初始感染比例(i00.001)可能导致数值不稳定多波次感染需要扩展为SIRS模型
