传染病模型SIR及其变体(python版本)

已于 2024-05-14 08:29:01 修改
阅读量1.2k 收藏 35
点赞数 10
分类专栏： python 知识分享 文章标签： python 开发语言
于 2024-04-24 21:03:56 首次发布
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/qq_36714950/article/details/138058867
版权

文章目录

传染病模型及其变体

1. SI模型

在该模型里面，群体中只有两种人：易感者和感染者。
感染者每天会感染一定的数量的易感者。
S表示易感者（尚未感染但是容易被感染的人） 的数量
I表示感染者（已经感染的人）的数量
N表示人口总数
r表示一个感染者平均每天感染易感者的人数
那么每天易感者和感染者的数量变化为
d S = − r S I N d I = r S I N \begin{align} dS &= \frac{-rSI}{N} \\ dI &= \frac{rSI}{N} \end{align} dSdI​=N−rSI​=NrSI​​​

1.1代码

import numpy as np   # 科学计算工具包
import matplotlib.pyplot as plt  # 画图工具包

plt.rcParams["font.family"] = 'Arial Unicode MS'  # 用来正常显示中文
# 定义SI函数——根据SI模型计算每天新增的易感人数和感染人数，返回累计人数
def SI(si,dt):
    S,I = si          # 每天初始SI的数值
    dS = -(r*I*S)/N   # 易感者微分方程
    dI = r*I*S/N      # 感染者微分方程
    S = 0 if S+dS*dt<=0 else S+dS*dt      # 当天易感者人数    
    I = N if I+dI*dt>=N else I+dI*dt      # 当天感染者人数
    return [S, I]

def calculate(func,si,days):  
    dt = 1
    t = np.arange(0,days,dt) # 设置时间步
    res = []
    for itm in t:
        si=func(si,dt)      # 运行SI模型函数
        res.append(si)      # 存储每天人数结果
    return np.array(res)

# 画图函数
def plot_graph(np_res): 
    plt.figure(figsize=(10,6),dpi=300)
    plt.plot(np_res[:,0])
    plt.plot(np_res[:,1])
    plt.title("SI模型")
    plt.xlabel("天数")
    plt.ylabel("人数")
    plt.legend(['易感者','感染者'])
    plt.show()
    

N = 10000
I = 1
r = 1
days = 50
si= [N-I,    # 易感人数
       I]    # 感染人数

result = calculate(SI,si,days)
plot_graph(result)

在这里插入图片描述

2. SIS模型

在该模型里面，群体中依然只有两种人：易感者和感染者。

感染者每天会感染一定的数量的易感者，同时每天会有一定数量的感染者康复，但是他们康复之后依然有可能被感染。

S S S表示易感者（尚未感染但是容易被感染的人） 的数量

I I I表示感染者（已经感染的人）的数量

N N N表示人口总数

r r r表示一个感染者平均每天感染易感者的人数

μ \mu μ 表示感染者每天康复的比例

那么每天易感者和感染者的数量变化为

d S = − r S I N + μ I d I = r S I N − μ I \begin{align} dS &= \frac{-rSI}{N} + \mu I\\ dI &= \frac{rSI}{N} - \mu I \\ \end{align} dSdI​=N−rSI​+μI=NrSI​−μI​​

2.1 代码

# 定义SIS函数——根据SIS模型计算每天新增的易感人数和感染人数，返回累计人数
def SIS(sis,dt):
    S,I = sis          # 每天初始SI的数值
    dS = -(r*I*S)/N + u*I   # 易感者微分方程
    dI = r*I*S/N - u*I     # 感染者微分方程
    S = 0 if S+dS*dt<=0 else S+dS*dt    # 当天易感者人数
    I = N if I+dI*dt>=N else I+dI*dt    # 当天感染者人数
    return [S, I]

def calculate(func,sis,days):  
    dt = 1
    t = np.arange(0,days,dt) # 设置时间步
    res=[]
    for itm in t:
        sis=func(sis,dt)      # 运行SI模型函数
        res.append(sis)      # 存储每天人数结果
    return np.array(res)

# 画图函数
def plot_graph(np_res): 
    plt.figure(figsize=(10,6),dpi=300)
    plt.plot(np_res[:,0])
    plt.plot(np_res[:,1])
    plt.title("SIS模型")
    plt.xlabel("天数")
    plt.ylabel("人数")
    plt.legend(['易感者','感染者'])
    plt.text(10,4500,'N=%d \nI=%d \nr=%2.1f \nu=%2.1f' % (N,I,r,u), bbox={'facecolor': 'red', 'alpha': 0.4, 'pad': 8})
    final = [round(x) for x in result[-1]]  # 取整表示
    plt.text(90,final[0],final[0])
    plt.text(93,final[1],final[1])
    plt.show()

# 赋值绘图
N = 10000
I = 1
r = 0.6
u = 0.2
days = 100

sis= [N-I,    # 易感人数
       I]    # 感染人数

result = calculate(SIS,sis,days)
plot_graph(result)
print(result[-1])

在这里插入图片描述
可以看出，在SIS模型中，并非所有的人都会被感染，最终感染人数为：

I = N ( 1 − u r ) I = N(1-\frac{u}{r}) I=N(1−ru​)

3. 基本再生数 basic reproductive number

这里，从SIS模型中引出了一个概念，就是基本再生数，其定义为：

R 0 = r u R_0 = \frac{r}{u} R0​=ur​

r 表示感染的速率，u 表示治愈的速率

如果 r>u 也就是 R 0 > 1 R_0>1 R0​>1 时，感染的速度快于治愈的速度，那么疾病就会传播开来，

而且 R 0 R_0 R0​越大，传播速度越快，感染人数越多；

如果 r<u 也就是 R 0 < 1 R_0<1 R0​<1 时，感染的速度慢于治愈的速度，那么疾病就不会传播开来。

这里我们可以看看其他传染病的基本再生数
%%html

4. SIR模型

在该模型里面，群体中有三种人：易感者，感染者，康复者（包括死亡者）。

感染者每天会感染一定的数量的易感者，同时每天会有一定数量的感染者康复（或者死亡），而且他们康复之后不可能再次被感染（拥有了抗体）。

S表示易感者（尚未感染但是容易被感染的人） 的数量

I表示感染者（已经感染的人）的数量

R表示康复者（包括治愈和死亡的人）的数量

N表示人口总数

r表示一个感染者平均每天感染易感者的人数

μ \mu μ 表示感染者每天康复（含死亡）的比例

那么每天易感者/感染者/康复者的数量变化为

d S = − r S I N d I = r S I N − μ I d R = μ I \begin{align} dS &= \frac{-rSI}{N} \\ dI &= \frac{rSI}{N} - \mu I \\ dR &= \mu I \\ \end{align} dSdIdR​=N−rSI​=NrSI​−μI=μI​​

4.1 代码

# 定义SIR函数——根据SIR模型计算每天新增的易感人数，感染人数，康复人数，返回累计人数
def SIR(sir,dt):
    S,I,R = sir             # 每天初始SIR的数值
    dS = -(r*I*S)/N         # 易感者微分方程
    dI = r*I*S/N - u*I      # 感染者微分方程
    dR = u*I                # 康复者微分方程
    S = 0 if S+dS*dt<=0 else S+dS*dt                        # 当天易感者人数
    I = N if I+dI*dt>=N else (0 if I+dI*dt<=0 else I+dI*dt) # 当天感染者人数
    R = N if R+dR*dt>=N else (0 if R+dR*dt<=0 else R+dR*dt) # 当天康复者人数
    return [S, I, R]

def calculate(func,sir,days):  
    dt = 1
    t = np.arange(0,days,dt) # 设置时间步
    res=[]
    for itm in t:
        sir=func(sir,dt)      # 运行SI模型函数
        res.append(sir)      # 存储每天人数结果
    return np.array(res)

# 画图函数
def plot_graph(np_res): 
    plt.figure(figsize=(10,6),dpi=300)
    plt.plot(np_res[:,0])
    plt.plot(np_res[:,1])
    plt.plot(np_res[:,2])
    plt.title("SIR模型")
    plt.xlabel("天数")
    plt.ylabel("人数")
    plt.legend(['易感者','感染者','康复者'])
    plt.text(1,4500,'N=%d \nI=%d \nr=%2.1f \nu=%2.1f' % (N,I,r,u), bbox={'facecolor': 'red', 'alpha': 0.4, 'pad': 8})
    final = [round(x) for x in result[-1]]  # 取整表示
    plt.text(45,final[0],final[0])
    plt.text(46,final[1],final[1])
    plt.text(45,final[2],final[2])
    plt.show()

# 赋初值、绘图
N = 10000
I = 1
r = 1
u = 0.8
days = 100

sir= [N-I,    # 易感人数
       I,     # 感染人数
       0]     # 康复人数

result = calculate(SIR,sir,days)
plot_graph(result)
final = [round(x) for x in result[-1]]  # 取整表示
print(final)

在这里插入图片描述

5. SEIR模型

考虑到感染病毒之后存在一定的潜伏期，而且病毒携带者未必一定会出现明显的症状（即患病），SEIR模型就被提出来了。

在SEIR模型中存在四种人：易感者，潜伏者，感染者（特指患病的人），康复者（含死亡者）。

感染者每天会感染一定的数量的易感者，同时每天会有一定数量的感染者康复（或者死亡），而且他们康复之后不可能再次被感染（拥有了抗体）。

S表示易感者（尚未感染但是容易被感染的人） 的数量

E表示潜伏者（已经感染尚未发病的人） 的数量

I表示感染者（已经感染的人）的数量

R表示康复者（包括治愈和死亡的人）的数量

N表示人口总数

r表示一个感染者平均每天感染易感者的人数

α \alpha α 表示潜伏者每天发展感染者的比例

μ \mu μ 表示感染者每天康复（含死亡）的比例

那么每天易感者和感染者的数量变化为

d S = − r S I N d E = r S I N − α E d I = α E − μ I d R = μ I \begin{align} dS &= \frac{-rSI}{N} \\ dE &= \frac{rSI}{N} - \alpha E \\ dI &= \alpha E - \mu I \\ dR &= \mu I \\ \end{align} dSdEdIdR​=N−rSI​=NrSI​−αE=αE−μI=μI​​

5.1 代码

import numpy as np
import matplotlib.pyplot as plt

plt.rcParams["font.family"] = 'Arial Unicode MS'  # 用来正常显示中文


# 定义SEIR函数——根据SEIR模型计算每天新增的易感人数，潜伏人数，感染人数，康复人数，返回累计人数
def SEIR(seir,dt):
    S,E,I,R = seir             # 每天初始SIR的数值
    dS = -(r*I*S)/N            # 易感者微分方程
    dE = (r*I*S)/N - a*E       # 潜伏者微分方程
    dI = a*E - u*I             # 感染者微分方程
    dR = u*I                   # 康复者微分方程
    S = 0 if S+dS*dt<=0 else S+dS*dt                        # 当天易感者人数
    E = 0 if E+dE*dt<=0 else (N if E+dE*dt>=N else E+dE*dt) # 当天潜伏者人数
    I = N if I+dI*dt>=N else (0 if I+dI*dt<=0 else I+dI*dt) # 当天感染者人数
    R = N if R+dR*dt>=N else (0 if R+dR*dt<=0 else R+dR*dt) # 当天康复者人数
    return [S, E, I, R]

def calculate(func,seir,days):  
    dt = 1
    t = np.arange(0,days,dt) # 设置时间步
    res=[]
    for itm in t:
        seir=func(seir,dt)      # 运行SEIR模型函数
        res.append(seir)       # 存储每天人数结果
    return np.array(res)

# 画图函数
def plot_graph(np_res): 
    plt.figure(figsize=(10,6),dpi=300)
    plt.plot(np_res[:,0])
    plt.plot(np_res[:,1])
    plt.plot(np_res[:,2])
    plt.plot(np_res[:,3])
    plt.title("SEIR模型")
    plt.xlabel("天数")
    plt.ylabel("人数")
    plt.legend(['易感者','潜伏者','感染者','康复者'])
    plt.text(1,4500,'N=%d \nE=%d  \nI=%d \nr=%2.1f \na=%2.1f \nu=%2.1f' % (N,E,I,r,a,u), bbox={'facecolor': 'red', 'alpha': 0.4, 'pad': 8})
    final = [round(x) for x in result[-1]]  # 取整表示
    plt.text(90,650,'%d' % final[0])
    plt.text(92,100,'%d' % final[2])
    plt.text(90,final[3],'%d' % final[3])
    plt.show()

# 赋初值、绘图
N = 10000
E = 0 
I = 10
r = 0.6
a = 0.3
u = 0.3
days = 100

seir= [N-I-E,    # 易感人数
       E,      # 潜伏人数
       I,      # 感染人数
       0]      # 康复人数

result = calculate(SEIR,seir,days)
plot_graph(result)
final = [round(x) for x in result[-1]]  # 取整表示
print(final)

在这里插入图片描述

6. SEIJR模型

在SEIR模型的基础上考虑确诊情况，提出SEIJR模型。确诊病例 confirmed case.

在SEIJR模型中存在五种人：易感者，潜伏者，感染者（指患病的但未被确诊的人），确诊者，康复者（含死亡者）。

潜伏者有可能继续感染易感者，未确认感染者自我防护意识强不会感染其他人，确诊患者由于被强制隔离页不会感染其他人。

S表示易感者（尚未感染但是容易被感染的人） 的数量

E表示潜伏者（已经感染尚未发病的人）的数量

I表示感染者（已经感染的人）的数量

J表示确诊者（感染疾病并确诊的人）的数量

R表示康复者（包括治愈和死亡的人）的数量

N表示人口总数

r表示一个潜伏者平均每天感染易感者的人数

α \alpha α 表示潜伏者每天发展为感染者的比例

π \pi π 表示潜伏者每天自愈的比例

μ \mu μ 表示感染者每天康复（含死亡）的比例

τ \tau τ 表示感染者每天确诊的比例

那么每天易感者和感染者的数量变化为

d S = − r S E N d E = r S E N − α E − π E d I = α E − μ I d J = τ I d R = π E + μ I \begin{align} dS &= \frac{-rSE}{N} \\ dE &= \frac{rSE}{N} - \alpha E - \pi E \\ dI &= \alpha E - \mu I \\ dJ &= \tau I \\ dR &= \pi E + \mu I \\ \end{align} dSdEdIdJdR​=N−rSE​=NrSE​−αE−πE=αE−μI=τI=πE+μI​​

6.1 代码

def SEIJR(seijr,dt):
    S,E,I,J,R = seijr
    dS = -(r*E*S)/N                  # 易感者微分方程
    dE = (r*E*S)/N - a*E -p*E        # 潜伏者微分方程
    dI = a*E - u*I                   # 感染者微分方程
    dJ = t*I                         # 确诊者微分方程
    dR = p*E + u*I                   # 康复者微分方程
    S = 0 if S+dS*dt<=0 else S+dS*dt                        # 当天易感者人数
    E = 0 if E+dE*dt<=0 else (N if E+dE*dt>=N else E+dE*dt) # 当天潜伏者人数
    I = N if I+dI*dt>=N else (0 if I+dI*dt<=0 else I+dI*dt) # 当天感染者人数
    J = N if J+dJ*dt>=N else (0 if J+dJ*dt<=0 else J+dJ*dt) # 当天感染者人数
    R = N if R+dR*dt>=N else (0 if R+dR*dt<=0 else R+dR*dt) # 当天康复者人数
    return [S, E, I, J, R]

def calculate(func,seijr,days):  
    dt = 1
    t = np.arange(0,days,dt) # 设置时间步
    res=[]
    for itm in t:
        seijr=func(seijr,dt)      # 运行SEIR模型函数
        res.append(seijr)       # 存储每天人数结果
    return np.array(res)

# 画图函数
def plot_graph(np_res): 
    plt.figure(figsize=(10,6),dpi=300)
#     plt.plot(np_res[:,0])
    plt.plot(np_res[:,1])
    plt.plot(np_res[:,2])
    plt.plot(np_res[:,3])
    plt.plot(np_res[:,4])
    plt.title("SEIJR模型")
    plt.xlabel("天数")
    plt.ylabel("人数")
    plt.legend(['潜伏者','感染者','确诊者','康复者']) # '易感者',
    plt.text(1,500,'N=%d \nE=%d \nr=%3.2f \na=%3.2f \nu=%3.2f' % (N,E,r,a,u), bbox={'facecolor': 'red', 'alpha': 0.4, 'pad': 8})
    final = [round(x) for x in result[-1]]  # 取整表示
#     plt.text(90,final[0]+100,'%d' % final[0])  # 易感人数
#     plt.text(90,final[1]+100,'%d' % final[1])  # 潜伏人数
    plt.text(days-10,final[2]-50,'%d' % final[2])  # 感染人数
    plt.text(days-10,final[3]-50,'%d' % final[3])  # 确诊人数
#     plt.text(90,final[4]+100,'%d' % final[4])  # 康复人数
    plt.show()

# 赋初值、绘图
N = 14000000
E = 1000
I = 100
J = 41
r = 0.25
a = 0.2
u = 0.1
p = 0.08
t = 0.2
days = 200

seijr = [N-E,    # 易感人数
         E,      # 潜伏人数
         I,      # 感染人数
         J,      # 确诊人数
         0]      # 康复人数

result = calculate(SEIJR,seijr,days)
plot_graph(result)
final = [round(x) for x in result[-1]]  # 取整表示
print(final)

在这里插入图片描述

7. SEIJRD模型

考虑到全面防控等严格隔离措施，将感染率设为参数可调，提出SEIJRD模型。

该模型分为两个阶段，采取隔离措施前r较大；采取隔离措施后r较大；增加死亡者微分方程。

S表示易感者（尚未感染但是容易被感染的人） 的数量

E表示潜伏者（已经感染尚未发病的人）的数量

I表示感染者（已经感染的人）的数量

J表示确诊者（感染疾病并确诊的人）的数量

R表示康复者（治愈患者）的数量

D表示康复者（死亡患者）的数量

N表示人口总数

r表示一个潜伏者平均每天感染易感者的人数

α \alpha α 表示潜伏者每天发展为感染者的比例

π \pi π 表示潜伏者每天自愈的比例

μ \mu μ 表示感染者每天康复的比例

τ \tau τ 表示感染者每天确诊的比例

ω \omega ω 表示感染者每天死亡的比例

那么每天易感者和感染者的数量变化为

d S = − r S E N d E = r S E N − α E − π E d I = α E − μ I − ω I d J = τ I d R = π E + μ I d D = ω I \begin{align} dS &= \frac{-rSE}{N} \\ dE &= \frac{rSE}{N} - \alpha E - \pi E \\ dI &= \alpha E - \mu I - \omega I\\ dJ &= \tau I \\ dR &= \pi E + \mu I \\ dD &= \omega I \\ \end{align} dSdEdIdJdRdD​=N−rSE​=NrSE​−αE−πE=αE−μI−ωI=τI=πE+μI=ωI​​

7.1 代码

def SEIJRD(seijrd,para,dt):
    S,E,I,J,R,D = seijrd
    r,a,u,p,t,w = para   # 赋值
    dS = -(r*E*S)/N                  # 易感者微分方程
    dE = (r*E*S)/N - a*E -p*E        # 潜伏者微分方程
    dI = a*E - u*I - w*I             # 感染者微分方程
    dJ = t*I                         # 确诊者微分方程
    dR = p*E + u*I                   # 康复者微分方程
    dD = w*I                         # 死亡者微分方程
    
    S = 0 if S+dS*dt<=0 else S+dS*dt                        # 当天易感者人数
    E = 0 if E+dE*dt<=0 else (N if E+dE*dt>=N else E+dE*dt) # 当天潜伏者人数
    I = N if I+dI*dt>=N else (0 if I+dI*dt<=0 else I+dI*dt) # 当天感染者人数
    J = N if J+dJ*dt>=N else (0 if J+dJ*dt<=0 else J+dJ*dt) # 当天感染者人数
    R = N if R+dR*dt>=N else (0 if R+dR*dt<=0 else R+dR*dt) # 当天康复者人数
    D = N if D+dD*dt>=N else (0 if D+dD*dt<=0 else D+dD*dt) # 当天死亡者人数
    return [S, E, I, J, R, D]

def calculate(func,seijrd,days):  
    dt = 1
    t = np.arange(0,days,dt) # 设置时间步
    res=[]
    for itm in t:
        if itm < 13:
            seijrd=func(seijrd,para1,dt)      # 运行SEIJR模型函数
        else:
            seijrd=func(seijrd,para2,dt)      # 运行SEIJR模型函数
        res.append(seijrd)                   # 存储每天人数结果
    return np.array(res)

# 画图函数
def plot_graph(np_res): 
    plt.figure(figsize=(10,6),dpi=300)
    ax = plt.subplot(111)
#     plt.plot(np_res[:,0])
    plt.plot(np_res[:,1])
    plt.plot(np_res[:,2])
    plt.plot(np_res[:,3])
    plt.plot(np_res[:,4])
    plt.plot(np_res[:,5])
    ax.set_xticks([0,20,40,60,80,100])
    ax.set_xticklabels(['0110','0130','0209','0229','0320','0409'])
    plt.title("武汉市COVID-19疫情预测")
    plt.xlabel('日期')
    plt.ylabel("人数")
    plt.legend(['潜伏者','感染者','确诊者','康复者','死亡者']) # '易感者',
    plt.text(60,20000,'N=%d, r1=%3.2f, r2=%3.2f' \
             % (14000000,0.52,0.2), \
             bbox={'facecolor': 'red', 'alpha': 0.4, 'pad': 8})
    plt.text(60,15000,r'$\alpha$=%3.2f, $\mu$=%3.2f, $\pi$=%3.2f, $\tau$=%3.2f, $\omega$=%3.2f' \
             % (0.2,0.1,0.1,0.18,0.009), \
             bbox={'facecolor': 'red', 'alpha': 0.4, 'pad': 8})
    final = [round(x) for x in result[-1]]  # 取整表示
#     plt.text(90,final[0]+100,'%d' % final[0])  # 易感人数
#     plt.text(90,final[1]+100,'%d' % final[1])  # 潜伏人数
    plt.text(days-20,final[2]+100,'感染人数为 %d' % final[2])   # 感染人数
    plt.text(days-20,final[3]+100,'确诊人数 %d' % final[3])     # 确诊人数
    plt.text(days-20,final[4]+100,'康复人数为 %d' % final[4])  # 康复人数
    plt.text(days-20,final[5]+100,'死亡人数为 %d' % final[5])  # 死亡人数
    
    plt.show()

# 武汉
N = 14000000
E = 800
I = 100
J = 41
days = 100
para1 = [0.52,0.2,0.1,0.1,0.18,0.009]
para2 = [0.2,0.2,0.1,0.1,0.18,0.009]

seijrd = [N-E,   # 易感人数
         E,      # 潜伏人数
         I,      # 感染人数
         J,      # 确诊人数
         0,      # 死亡人数
         0]      # 康复人数


result = calculate(SEIJRD,seijrd,days)
plot_graph(result)
final = [round(x) for x in result[-1]]  # 取整表示
print(final)

在这里插入图片描述

# 美国 3.3亿人口

def SEIJRD(seijrd,para,dt):
    S,E,I,J,R,D = seijrd
    r,a,u,p,t,w = para   # 赋值
    dS = -(r*E*S)/N                  # 易感者微分方程
    dE = (r*E*S)/N - a*E -p*E        # 潜伏者微分方程
    dI = a*E - u*I - w*I             # 感染者微分方程
    dJ = t*I                         # 确诊者微分方程
    dR = p*E + u*I                   # 康复者微分方程
    dD = w*I                         # 死亡者微分方程
    
    S = 0 if S+dS*dt<=0 else S+dS*dt                        # 当天易感者人数
    E = 0 if E+dE*dt<=0 else (N if E+dE*dt>=N else E+dE*dt) # 当天潜伏者人数
    I = N if I+dI*dt>=N else (0 if I+dI*dt<=0 else I+dI*dt) # 当天感染者人数
    J = N if J+dJ*dt>=N else (0 if J+dJ*dt<=0 else J+dJ*dt) # 当天感染者人数
    R = N if R+dR*dt>=N else (0 if R+dR*dt<=0 else R+dR*dt) # 当天康复者人数
    D = N if D+dD*dt>=N else (0 if D+dD*dt<=0 else D+dD*dt) # 当天死亡者人数
    return [S, E, I, J, R, D]

def calculate(func,seijrd,days):  
    dt = 1
    t = np.arange(0,days,dt) # 设置时间步
    res=[]
    for itm in t:
        if itm < 22:
            seijrd=func(seijrd,para1,dt)      # 运行SEIJR模型函数
        else:
            seijrd=func(seijrd,para2,dt)      # 运行SEIJR模型函数
        res.append(seijrd)                   # 存储每天人数结果
    return np.array(res)

# 画图函数
def plot_graph(np_res): 
    plt.figure(figsize=(10,6),dpi=300)
    ax = plt.subplot(111)
#     plt.plot(np_res[:,0])
    plt.plot(np_res[:,1])
    plt.plot(np_res[:,2])
    plt.plot(np_res[:,3])
    plt.plot(np_res[:,4])
    plt.plot(np_res[:,5])
    ax.set_xticks([0,20,40,60,80,100])
    ax.set_xticklabels(['0301','0321','0410','0430','0520','0610'])
    plt.title("美国COVID-19疫情预测")
    plt.xlabel('日期')
    plt.ylabel("人数")
    plt.legend(['潜伏者','感染者','确诊者','康复者','死亡者']) # '易感者',
    plt.text(60,150000,'N=%d, r1=%3.2f, r2=%3.2f' \
             % (330000000,0.52,0.2), \
             bbox={'facecolor': 'red', 'alpha': 0.4, 'pad': 8})
    plt.text(60,100000,r'$\alpha$=%3.2f, $\mu$=%3.2f, $\pi$=%3.2f, $\tau$=%3.2f, $\omega$=%3.2f' \
             % (0.2,0.1,0.1,0.18,0.009), \
             bbox={'facecolor': 'red', 'alpha': 0.4, 'pad': 8})
    final = [round(x) for x in result[-1]]  # 取整表示
#     plt.text(90,final[0]+100,'%d' % final[0])  # 易感人数
#     plt.text(90,final[1]+100,'%d' % final[1])  # 潜伏人数
    plt.text(days-20,final[2]+100,'感染人数为 %d' % final[2])   # 感染人数
    plt.text(days-20,final[3]+100,'确诊人数 %d' % final[3])     # 确诊人数
    plt.text(days-20,final[4]+100,'康复人数为 %d' % final[4])  # 康复人数
    plt.text(days-20,final[5]+100,'死亡人数为 %d' % final[5])  # 死亡人数
    
    plt.show()

# 赋初值、绘图
N = 330000000
E = 1000
I = 200
J = 69
days = 100
para1 = [0.52,0.2,0.1,0.1,0.18,0.009]
para2 = [0.2,0.2,0.1,0.1,0.18,0.009]

seijrd = [N-E,   # 易感人数
         E,      # 潜伏人数
         I,      # 感染人数
         J,      # 确诊人数
         0,      # 死亡人数
         0]      # 康复人数


result = calculate(SEIJRD,seijrd,days)
plot_graph(result)
final = [round(x) for x in result[-1]]  # 取整表示
print(final)

在这里插入图片描述

yanxiaoyu110
关注
专栏目录
SIR模型python实现
Yygod_p的博客
12-08 1347
python实现SIR传播模型
BasicSIRModel:Python中的简单SIR模型
05-02
基本SIR模型 该项目包含一个基本的SIR模型，并带有绘图功能。 可以在找到该项目的配套博客文章。 运行模型 >> > import model >> > m = model . SIR () >> > m . run () >> > m . plot () 更改模型参数 创建model.SIR类的新实例时，必须指出特定的参数。 例如，可以如下更改beta（S到I）率： >> > import model >> > m = model . SIR ( rateSI = 0.05 ) 可变参数包括： 'eons'（模型的时间点数，默认为1000） “易感”（时间0时易感个体的数量，默认为950） “已感染”（时间0处的已感染人数，默认值为50） “抗性”（时间0的抗性个体数，默认为0） 'rateSI'（从S到I的基本速率'beta'，默认值为0.05） 'rateIR'（从
SIR模型
02-08
SIR模型.CDC constructed the EbolaResponse modeling tool in spreadsheet (available at http://dx.doi.org/10.15620/ dc.24900) using Microsoft Excel 2010 and used the model o estimate the increase in Ebola cases in Liberia and Sierra eone (see Appendix for additional results and technical notes).
python传染病模型 供参考
02-11
python传染病模型 供参考
【保姆级教程】使用python实现SIR模型（包含数据集的制作与导入及最终结果的可视化）
m0_65238732的博客
08-20 6354
使用python实现SIR传染病模型，详细介绍了数据集的制作与导入，整体模型框架，和最后以图片形式对结果进行呈现，并额外介绍了如何利用python制作散点图。
用python实现传染病模型传染病模型
qq_42839893的博客
01-31 8317
用python实现传染病模型1.SI模型1.1 代码实现1.2 模型的结果2.SIS （治愈后仍然还是易感者）2.1 代码实现2.2模型的结果3 SIR模型（治愈后直接移除）3.2代码实现3.2绘制图像：4.SEIR 模型（新增一个人群，叫潜伏者E）4.1代码实现4.2模型的结果 参考文章为：https://www.kesci.com/mw/project/60161bf6ac79f40016b7d7d9 1.SI模型 示意图： 我们假设城市有一千万（N=10的7次方）人，每个患者每天接触感染每天0.
传染病SEIR传播动力模型python代码
04-17
代码建立了传染病SEIR传播动力模型，可以通过修改其中参数从而能够有效模拟传染病传播情况，并做了可视化，能够直观观察。
SIR模型和Python实现
weixin_41168304的博客
01-10 1万+
一、SIR模型介绍 SIR模型时传染病中最基础最核心的模型，研究的是某个封闭地区的疫情传播规律。 SIR模型的动力学关系如下图： 健康人数S的变化与 健康人数S和正感人数I的乘积（代表健康人数和正感人数的接触）成正比，其中α代表交叉感染率 移出人数的变化与正感人数的数量成正比，其中β代表回复率。 基于上面的是自，SIR模型可以表示成一个常微分方程组如下图： 当s(t)=β/α时就是病毒最严重的时候； 表示S的反函数，代表S越大，t越小，越早结束疫情 封城就是为了较低交叉感染率α，并.
SIR传染病模型（微分方程系列1）
qq_46211790的博客
07-25 2万+
一：基本参数 SIR模型是常见的一种描述传染病传播的数学模型，其基本假设是将人群分为以下三类： S:(Susceptible):易感人群，指未得病者，但缺乏免疫能力，与感病者接触后容易受到感染。 I:(Infective):患病人群，指染上传染病的人，他可以传播给易感人群。 R:(Removed）：移除人群，被移出系统的人，因病愈（具有免疫力）或死亡的人。这部分人不再参与感染和被感染过程。 二：问题假设 1：没有天生携带抗体的人 2：总人数N不变且仅包含三种健康状态，健康人、病人和移出者的比例
Python学习：SIR传染病模型
weixin_45237889的博客
02-07 3550
SIR模型是传染病模型中最经典的一个，类似的还有SI和SIS两种。SIR是三个单词首字母的缩写，其中S是Susceptible的缩写，表示易感者；I是Infective的缩写，表示感染者，R是Removal的缩写，表示移除者。 这个模型本身是在研究这三者的关系。在病毒最开始的时候，所有人都是易感者，也就是所有人都有可能中病毒；当一部分人在接触到病毒以后中病毒了，变成了感染者；感染者会接受各种治疗，...
改进的SIR模型评估k个重要点算法(Python)
06-27
本资源针对一种边权重存在重尾分布复杂网络，改进原本的SIR模型对TopK重要节点进行性能评估。并将传播过程绘制成可视化图。本资源使用networkx工具包。
模型python执行代码
05-31
模型python执行代码模型python执行代码模型python执行代码模型python执行代码
传染病模型
07-19
传染病模型的算法支持，侧重于系统动力学的学习和运转。
【数学建模】传染病SIR模型
热门推荐
SuperSources的博客
08-20 3万+
MATLAB实现 转自：SIR模型实现（matlab） function dydt = odefun(t,y,A,B) dydt = zeros(2,1); dydt(1) = A*y(1)*y(2) - B*y(1); dydt(2) = -A*y(1)*y(2); end clc clear close all; A = 0.4; B = 0.1; I = 0.4; S = 0.5; %ode tspan = [0 50]; y0 = [I S]; [t, y] = ode45(@(t,y)ode
传染病模型python_用matplotlib绘制python传染病SIR模型
weixin_34125336的博客
02-05 547
我试图用Python中的matplotlib绘制一个SIR模型，它具有：a)人口2200b)并显示了30天内受影响的人口过程c)δt=0.01-Delta是变化的，t是时间这是我的代码：import matplotlib.pyplot as pyplotdef SIR(population, days, dt):susceptible = population - 1infected = 1.0r...
Python：SEIR传染病模型
FightingBob的博客
04-12 1万+
Hello，大家好！好久没有更新了，今天给大家介绍一个非常经典的模型——SEIR传染病模型！ SEIR传染病模型模型介绍S、E、I、R代表的含义模型建立变量解释建立微分方程转化迭代形式Python实现总结获得源代码 模型介绍 传染病的基本数学模型，研究传染病的传播速度、空间范围、传播途径、动力学机理等问题，以指导对传染病的有效地预防和控制。常见的传染病模型按照传染病类型分为 SI、SIR、SIRS、SEIR 模型等，按照传播机理又分为基于常微分方程、偏微分方程、网络动力学的不同类型。 S、E、I、R代表的含
传染病SIR模型
千里的博客
02-06 1万+
SIR模型简介 在典型的传染病模型中，种群（Population）内的N个个体的状态可分为如下几类： 易染状态S(Susceptible)：即健康状态，可被感染的个体。 感染状态I(Infected)：处于感染状态的个体还能够感染将康状态的个体。 移除状态R(Removed,Refractory or Recovered)：也称为免疫状态和恢复状态。一个个体经历过一个完整的感染周期后，...
python代码实现， SIR 模型进行信息传播模拟
最新发布
weixin_42602368的博客
12-20 487
SIR模型是一种常用的信息传播模型，它用来模拟传染病(如冠状病毒)的流行情况。SIR模型包括三种状态： S(Susceptible)：易感者，还没有患上传染病的人。 I(Infected)：感染者，已经患上传染病的人。 R(Recovered)：恢复者，已经患过传染病并康复的人。 SIR模型假设传染病是通过直接接触传播的，并且每个人只有在感染后才能被传染。因此，SIR模型可以用来模拟传染病的流...
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

yanxiaoyu110

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付：¥1
获取中
扫码支付

您的余额不足，请更换扫码支付或充值

打赏作者

实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值