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

🎯要点

🎯房室数学模型：🖊常微分方程推理模型 | 🖊模型求解条件：静态人口、混合均匀分布、感染率相同、健康免疫力 | 🖊绘制模型表征接近平衡相图 | 🖊计算模型随时间变化有效传染数 | 🖊计算免疫力衰减率​ | 🖊计算反复感染和传染性关系模型 | 🖊计算病毒携带者模型​ | 🖊计算低传染性​。🎯宿主模型：🖊计算风险异质人群 | 🖊计算接触网络混合矩阵 | 🖊年龄异质性耦种群变化 | 🖊感染矩阵模型推理​。🎯病原体模型：🖊疟疾传播常微分方程模型 | 🖊绘制模型单相剖面流图 | 🖊媒介传播疾病推理 | 🖊人畜共患疾病模型。🎯病原体动力学：🖊无共感染无交叉免疫微分方程模型​。🎯疫苗建模 | 🎯流感时间动态平衡态和稳定性分析 | 🎯流感空间动态晶格模型。

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

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

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

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

模型参数假设：

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

这是模型的方程组：
$$\{\begin{array}{l}\frac{dS}{dt}=-\frac{\beta S}{N}I\\ \frac{dI}{dt}=\frac{\beta S}{N}I-\gamma I,\\ \frac{dR}{dt}=\gamma I.\end{array}$$
模型图示：

问题公式化
$$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();

🔗参阅一：计算思维

🔗参阅二：亚图跨际