🍠Python流感常微分方程房室数学模型

关键词

Python | 数学 | 流感 | 病情 | 感染 | 病毒 | 常微分方程 | 模型 | 求解器 | 推理 | 人口 | 分布 | 感染率 | 免疫力 | 表征 | 相图 | 流图 | 时间 | 衰减率 | 携带者 | 宿主 | 异质 | 矩阵 | 接触网络 | 年龄 | 种群 | 病原体 | 传播 | 动力学 | 交叉免疫 | 平衡态 | 稳定性 | 晶格 | 欧拉

🏈指点迷津 | Brief

🎯要点

🎯常微分方程-用例：Python机器人动力学和细胞酶常微分方程。

🍇Python欧拉法求解一般流感常微分方程模型

这些模型通常使用常微分方程（确定性）运行，但也可以与随机框架一起使用，这种框架更现实，但分析起来要复杂得多。

此类模型可以预测诸如疾病如何传播、感染总数或流行病持续时间之类的事情，并估计各种流行病学参数，如再生数。还可以显示不同的公共卫生干预措施如何影响流感的结果，例如，在特定人群中发放有限数量的疫苗的最有效技术是什么。

模型参数假设：

平均而言，群体中的 $S$ 个体每单位时间会遇到 $\beta$ 个体。
每单位时间离开房室 $I$ 的感染者的比率是 $\gamma$ $I$ （一旦一个人被感染，他就会对该疾病产生免疫力）。
人口规模 $N=S+I+R$ 是恒定的

这是模型的方程组：

\left\{\begin{array}{l} \frac{ d S}{ d t}=-\frac{\beta S}{N} I \\ \frac{ d I}{ d t}=\frac{\beta S}{N} I-\gamma I, \\ \frac{ d R}{ d t}=\gamma I . \end{array}\right.

模型图示：

问题公式化

X=\left(\begin{array}{l} S \\ I \\ R \end{array}\right)

微分系统写为：

\dot{X}=\left(\begin{array}{c} -\frac{\beta S}{N} I \\ \frac{\beta S}{N} I-\gamma I \\ \gamma I \end{array}\right)=f(X)

💦通用方法求解

 import numpy as np
 import matplotlib.pyplot as plt
 import matplotlib as mpl
 import pandas as pd
 import numba
 from scipy import integrate

💦给定相应数据：

 N = 350. 
 I0, R0 = 1., 0 
 S0 = N - I0 - R0 
 beta, gamma = 0.4, 0.1 
 tmax = 160 
 Nt = 160
 t = np.linspace(0, tmax, Nt+1)

💦定义微分函数：

 def derivative(X, t):
     S, I, R = X
     dotS = -beta * S * I / N
     dotI = beta * S * I / N - gamma * I
     dotR = gamma * I
     return np.array([dotS, dotI, dotR])

 X0 = S0, I0, R0 #Initial conditions vector
 res = integrate.odeint(derivative, X0, t)
 S, I, R = res.T
 Seuil = 1 - 1 / (beta/gamma)
 Seuil

 0.75

绘制个体疑似感染、感染和恢复趋势图：

 plt.figure()
 plt.grid()
 plt.title("odeint method")
 plt.plot(t, S, 'orange', label='Susceptible')
 plt.plot(t, I, 'r', label='Infected')
 plt.plot(t, R, 'g', label='Recovered with immunity')
 plt.xlabel('Time t, [days]')
 plt.ylabel('Numbers of individuals')
 plt.ylim([0,N])
 plt.legend()
 
 plt.show();

💦欧拉法求解：

 def Euler(func, X0, t):
     dt = t[1] - t[0]
     nt = len(t)
     X  = np.zeros([nt, len(X0)])
     X[0] = X0
     for i in range(nt-1):
         X[i+1] = X[i] + func(X[i], t[i]) * dt
     return X

 Nt = 100
 Xe = Euler(derivative, X0, t)

 plt.figure()
 
 plt.title("Euler method")
 plt.plot(t, Xe[:,0], color = 'orange', linestyle = '-.', label='Susceptible')
 plt.plot(t, Xe[:,1], 'r-.', label='Infected')
 plt.plot(t, Xe[:,2], 'g-.', label='Recovered with immunity')
 plt.grid()
 plt.xlabel("Time, $t$ [s]")
 plt.ylabel("Numbers of individuals")
 plt.legend(loc = "best")
 
 plt.show();

PreviousPython流感传播感染康复图模型计算和算法 NextPython脑溶质扩散生理几何模型计算

Last updated 6 months ago

plt.figure() plt.grid() plt.title("odeint method") plt.plot(t, S, 'orange', label='Susceptible') plt.plot(t, I, 'r', label='Infected') plt.plot(t, R, 'g', label='Recovered with immunity') plt.xlabel('Time t, [days]') plt.ylabel('Numbers of individuals') plt.ylim([0,N]) plt.legend() plt.show();

def Euler(func, X0, t): dt = t[1] - t[0] nt = len(t) X = np.zeros([nt, len(X0)]) X[0] = X0 for i in range(nt-1): X[i+1] = X[i] + func(X[i], t[i]) * dt return X

plt.figure() plt.title("Euler method") plt.plot(t, Xe[:,0], color = 'orange', linestyle = '-.', label='Susceptible') plt.plot(t, Xe[:,1], 'r-.', label='Infected') plt.plot(t, Xe[:,2], 'g-.', label='Recovered with immunity') plt.grid() plt.xlabel("Time, $t$ [s]") plt.ylabel("Numbers of individuals") plt.legend(loc = "best") plt.show();