头歌实践教学平台 Python答案 Python 计算思维训练——常微分新冠预测

任务描述

本关任务：使用常微分方程来预测新冠疫情。

相关知识

为了完成本关任务，你需要掌握： 1.常微分方程的基础知识； 2.前向欧拉法，后向欧拉法与梯形法； 3.SIR 模型。

其中 1，2 知识点见 Python 计算思维训练——常微分方程关卡。

SIR模型

SIR 模型是一种传播模型，会对信息传播过程进行抽象描述，是传染病模型中最经典的模型。模型中把传染病流行范围内的人群分成三类：

• S 类，易感者（Susceptible），指未得病者，但缺乏免疫能力，与感病者接触后容易受到感染；
• I 类，感病者（Infective），指染上传染病的人，它可以传播给 S 类成员；
• R 类，移出者（Removal），指被隔离，或因病愈而具有免疫力的人。

感染机制和 SIR 模型如下所示：

任务描述

本关任务：使用常微分方程来预测新冠疫情。

相关知识

为了完成本关任务，你需要掌握： 1.常微分方程的基础知识； 2.前向欧拉法，后向欧拉法与梯形法； 3.SIR 模型。

其中 1，2 知识点见 Python 计算思维训练——常微分方程关卡。

SIR模型

SIR 模型是一种传播模型，会对信息传播过程进行抽象描述，是传染病模型中最经典的模型。模型中把传染病流行范围内的人群分成三类：

• S 类，易感者（Susceptible），指未得病者，但缺乏免疫能力，与感病者接触后容易受到感染；
• I 类，感病者（Infective），指染上传染病的人，它可以传播给 S 类成员；
• R 类，移出者（Removal），指被隔离，或因病愈而具有免疫力的人。

感染机制和 SIR 模型如下所示：

任务描述

本关任务：使用常微分方程来预测新冠疫情。

相关知识

为了完成本关任务，你需要掌握： 1.常微分方程的基础知识； 2.前向欧拉法，后向欧拉法与梯形法； 3.SIR 模型。

其中 1，2 知识点见 Python 计算思维训练——常微分方程关卡。

SIR模型

SIR 模型是一种传播模型，会对信息传播过程进行抽象描述，是传染病模型中最经典的模型。模型中把传染病流行范围内的人群分成三类：

• S 类，易感者（Susceptible），指未得病者，但缺乏免疫能力，与感病者接触后容易受到感染；
• I 类，感病者（Infective），指染上传染病的人，它可以传播给 S 类成员；
• R 类，移出者（Removal），指被隔离，或因病愈而具有免疫力的人。

感染机制和 SIR 模型如下所示：

基于微分方程组对 SIR 模型进行求解，即可有效地拟合曲线。

编程要求

根据提示，在右侧编辑器 Begin-End 区间补充代码，实现新冠疫情的预测。

测试说明

平台会对你编写的代码进行测试。

---

开始你的任务吧，祝你成功！

import ode_7 as ode
import numpy as np
import matplotlib.pyplot as plt
# 采用matplotlib作图时默认设置下是无法显示中文的，凡是汉字都会显示成方块。
# 实际上，matplotlib是支持unicode编码的，不能正常显示汉字主要是没有找到合适的中文字体。
from pylab import mpl
mpl.rcParams['font.sans-serif']= ['SimHei']
# 解决负号显示问题
import matplotlib
matplotlib.rcParams['axes.unicode_minus']=False

'''
根据传染病传播的SIR模型，微分方程组可以转换为向量运算：
y = [S, I, R]
f(t, y) = [-beta*y[0]*y[1], beta*y[0]*y[1]-gamma*y[1], gamma*y[1]]
前向欧拉公式：
    y[k+1] = y[k] + dt * f(t[k], y[k])
后向欧拉公式（两步迭代）：
    y_p = y[k] + dt * f(t[k], y[k])
    y[k+1] = y[k] + dt * f(x[k+1], y_p)
改进的欧拉法：
    y_p = y[k] + dt * f(t[k], y[k])
    y_c = y[k] + dt * f(t[k+1], y_p)
    y[k+1] = 1/2 * (y_p + y_c)
'''
class SIR(ode.ODE):
    def f(self, t, y): # Y为向量[S,I,R]，t为标量
        #请在此处编程，实现f函数
        ########### Begin ###########
        return np.array([-self.beta*y[0]*y[1], self.beta*y[0]*y[1]-self.gamma*y[1], self.gamma*y[1]])

		########### End ###########

    def __init__(self, beta, gamma, y0):
        #请在此处编程，调用ODE类对SIR进行初始化
        ########### Begin ###########
        self.beta = beta; self.gamma = gamma; self.y0 = y0  # 参数属性
        ode.ODE.__init__(self, self.f, self.y0)

		########### End ###########

N     = 1e8             # 武汉总人数：1000万人
beta  = 1.0/N           # 1名病毒携带者平均每天感染1.0人
gamma = 1/25            # 假设肺炎患者平均25天治愈（15天潜伏 + 10天治疗）
y0    = [N-1, 1, 0]     # 初始发病1人，其他人员正常 [S0, I0, R0]

t = np.arange(0, 61, 1) # 模拟60天的发展情况，单位时间为1天
simulation = SIR(beta=beta, gamma=gamma, y0=y0) 
y = simulation.solve(t, 'T')  # 'T' 代表改进的欧拉法
print(y)
plt.plot(t, y[:, 0], '-g', label='易感者')
plt.plot(t, y[:, 1], '-r', label='感染者')
plt.plot(t, y[:, 2], '-k', label='恢复（隔离）者')

plt.legend(loc=0)
plt.title('继承ODE类的新冠病毒传播分析（SIR模型）')
plt.xlabel('时间（天）')
plt.ylabel('人数')
Vx = [1.0e6, 5.0e6] + [i * 1.e7 for i in range(11)]
plt.yticks(Vx, ['%d'%e for e in Vx])

plt.show()
plt.savefig('test1/test1.png')

关注更新后续答案