使用 matlab 实现 SEIRS 模型,模拟传染病流行过程
该模型基于疾病的自然史将人群分为四类:易感者、暴露者、传染者、恢复者。易感者与传染者发生有效接触之后会成为暴露者;暴露者经过潜伏期之后开始具有传染性,成为传染者;传染者在传染期内可以感染易感者,传染期在康复、住院、隔离等事件发生之后结束,此时传染者成为恢复者;恢复者具有对疾病的免疫力,不会被感染,但经过一段恢复期之后,免疫力会消失,重新成为易感者。
参数 | 含义 |
---|---|
S | 易感者仓室,表示易感者的数量 |
E | 暴露者仓室,表示暴露者的数量 |
I | 传染者仓室,表示传染者的数量 |
R | 恢复者仓室,表示恢复者的数量 |
$\beta$ | 每个个体在单位时间内发生的足以传染疾病的平均接触次数 |
$\frac{1}{\omega}$ | 平均潜伏期 |
$\frac{1}{\gamma}$ | 平均传染期 |
$\frac{1}{\alpha}$ | 平均恢复期 |
$$
\begin{aligned}
\frac{dS}{dt} &= \alpha R - \beta S \frac{I}{N}\\
\frac{dE}{dt} &= \beta S \frac{I}{N} - \omega E\\
\frac{dI}{dt} &= \omega E - \gamma I\\
\frac{dR}{dt} &= \gamma I - \alpha R
\end{aligned}
$$
该模型适用于存在潜伏期、感染之后会获得免疫力但免疫力在一段时间后会消失的短期暴发传染病。
适用的疾病类型有普通感冒、甲型流感、手足口病等。
适用的疾病暴发条件有人群在一个相对封闭的场所聚集,发生的接触相对比较密切。或者是在一个更广泛的空间内,有多个传染源相对均匀地分布在人群之中。
参数 | 设定值 | 单位 |
---|---|---|
S0 | 9990 | 人 |
E0 | 0 | 人 |
I0 | 10 | 人 |
R0 | 0 | 人 |
$\beta$ | $\frac{2}{7}$ | /天 |
$\omega$ | $\frac{1}{5}$ | /天 |
$\gamma$ | $\frac{1}{7}$ | /天 |
$\alpha$ | $\frac{1}{3}$ | /天 |
下面展示的代码分为两部分,一部分是函数代码,一部分是脚本代码。函数文件用于创建基于常微分方程的SEIRS模型,函数的输入是模型中所需要的一些参数,返回的是一阶常微分方程组等号右边的内容。脚本文件主要是用来测试模型结果,通过设置特定的参数值,以及初始值,利用求解器ode45求模型的数值解,并且将结果可视化。通过比较不同参数组合以及初始值下的模型结果,我们可以更好地理解疾病传播的动力学特性。
%函数
function dydt = SEIRSmodel(t, y, beta, omega, gamma, alpha)
%SEIRAMODEL 创建SEIRS模型
% input
%t 时间,可以处理常微分方程中的时间依赖项
%y 未知的一元函数组
%beta 传播速率系数
%omega 平均潜伏期的倒数
%gamma 平均传染期的倒数
%alpha 平均恢复期的倒数
% output
%dydt 数据类型为 single 或 double 的列向量 dydt,包含一阶方程组等号右边的内容
dydt = zeros(4,1);
dydt(1) = alpha.*y(4) - beta.*y(1).*y(3)./sum(y);
dydt(2) = beta.*y(1).*y(3)./sum(y) - omega.*y(2);
dydt(3) = omega.*y(2) - gamma.*y(3);
dydt(4) = gamma.*y(3) - alpha.*y(4);
end
%测试脚本
clc
clear
close all
%设置参数
beta = 2./7; %beta = R.*gamma,单位是$天_{-1}$
omega = 1./5; %单位是$天_{-1}$
gamma = 1./7; %单位是$天_{-1}$
alpha = 1./3; %单位是$天_{-1}$
tspan = [0, 200]; %积分区间,指定初始时间和最终时间,单位是天
y0 = [9990, 0, 10, 0]; %时间tspan(1)下的初始条件,1*4的行向量
%求解ODE
[t,y] = ode45(@(t,y) SEIRSmodel(t, y, beta, omega, gamma, alpha), tspan, y0);
y = interp1(t, y, 0:200); %使用插值来寻找指定时刻t下的解
%绘图
plot(0:200, y)
legend(["S", "E", "I", "R"])
从四类人群随时间 t 变化的趋势图我们可以看出在这个人口总数为10000的人群中,以及在4.1条件设定部分中所设置的条件下,该疾病不断蔓延,最终在 t = 120 以后趋向平稳,形成稳定的地方病平衡点,这一结果与模型计算出的基本再生数 $R_{0} = \frac{\beta}{\gamma} = 2$ 相一致。因为恢复者的免疫力会逐渐消失,导致人群中的易感者会得到补充,所以当基本再生数 $R_{0}$ 大于1时,在足够长的一段时间过后,人群中的患者不会消失,而是会维持在一定的水平。
尚无。