数学建模学习笔记（十七）传染病模型（SIER)

传染病模型讲解比较清楚的是知乎这位博主,文章链接戳这在家宅着也能抵抗肺炎！玩一玩SEIR传染病模型本文基于这篇文章进行记录和整理

对于一般传染病来说，都具备潜伏者(E)，因此直接记录传统的SIER模型：

模型公式：

left{ begin{array}{l} frac{{dS}}{{dt}} = - frac{{rbeta IS}}{N}\ \ frac{{dE}}{{dt}} = frac{{rbeta IS}}{N} - sigma E\ \ frac{{dI}}{{dt}} = sigma E - gamma I\ \ frac{{dR}}{{dt}} = gamma I end{array} right.

迭代公式：

left{ begin{array}{l} {S_n} = {S_{n - 1}} - frac{{rbeta {I_{n - 1}}{S_{n - 1}}}}{N}\ \ {E_n} = {E_{n - 1}} frac{{rbeta {I_{n - 1}}{S_{n - 1}}}}{N} - sigma {E_{n - 1}}\ \ {I_n} = {I_{n - 1}} sigma {E_{n - 1}} - gamma {I_{n - 1}}\ \ {R_n} = {R_{n - 1}} gamma {I_{n - 1}} end{array} right.

引入潜伏者传染概率，改进SEIR模型，公式为

left{ begin{array}{l} {frac{{dS}}{{dt}} = - frac{{rbeta IS}}{N} - frac{{{r_2}{beta _2}ES}}{N}}\ {}\ {frac{{dE}}{{dt}} = frac{{rbeta IS}}{N} - sigma E frac{{{r_2}{beta _2}ES}}{N}}\ {}\ {frac{{dI}}{{dt}} = sigma E - gamma I}\ {}\ {frac{{dR}}{{dt}} = gamma I} end{array} right.

迭代公式为：

left{ begin{array}{l} {S_n} = {S_{n - 1}} - frac{{rbeta {I_{n - 1}}{S_{n - 1}}}}{N} - frac{{{r_2}{beta _2}{E_{n - 1}}{S_{n - 1}}}}{N}\ \ {E_n} = {E_{n - 1}} frac{{rbeta {I_{n - 1}}{S_{n - 1}}}}{N} - sigma {E_{n - 1}} frac{{{r_2}{beta _2}{E_{n - 1}}{S_{n - 1}}}}{N}\ \ {I_n} = {I_{n - 1}} sigma {E_{n - 1}} - gamma {I_{n - 1}}\ \ {R_n} = {R_{n - 1}} gamma {I_{n - 1}} end{array} right.

matlab代码：源代码：

代码语言：javascript复制

clear;clc;

%--------------------------------------------------------------------------
%   参数设置
%--------------------------------------------------------------------------
N = 12700000;                                                                  %人口总数
E = 0;                                                                      %潜伏者
I = 1;                                                                      %传染者
S = N - I;                                                                  %易感者
R = 0;                                                                      %康复者

r = 20;                                                                     %感染者接触易感者的人数
B = 0.03;                                                                   %传染概率
a = 0.1;                                                                    %潜伏者转化为感染者概率
y = 0.1;                                                                    %康复概率

T = 1:140;
for idx = 1:length(T)-1
    S(idx 1) = S(idx) - r*B*S(idx)*I(idx)/N;
    E(idx 1) = E(idx)   r*B*S(idx)*I(idx)/N-a*E(idx);
    I(idx 1) = I(idx)   a*E(idx) - y*I(idx);
    R(idx 1) = R(idx)   y*I(idx);
end

plot(T,S,T,E,T,I,T,R);grid on;
xlabel('天');ylabel('人数')
legend('易感者','潜伏者','传染者','康复者')

稍作改进，反应每日新增病例情况：

代码语言：javascript复制

%--------------------------------------------------------------------------
%   初始化
%--------------------------------------------------------------------------
clear;clc;

%--------------------------------------------------------------------------
%   参数设置
%--------------------------------------------------------------------------
N = 29000;                                                                  %人口总数
E = 0;                                                                      %潜伏者
I = 1;                                                                      %传染者
S = N - I;                                                                  %易感者
R = 0;                                                                      %康复者
m=1;

r = 25;                                                                     %感染者接触易感者的人数
B = 0.03;                                                                   %传染概率
a = 0.1;                                                                    %潜伏者转化为感染者概率
r2 = 3;                                                                     %潜伏者接触易感者的人数
B2 = 0.03;                                                                  %潜伏者传染正常人的概率
y = 0.1;                                                                    %康复概率

T = 1:182;
for idx = 1:length(T)-1
    S(idx 1) = S(idx) - r*B*S(idx)*I(idx)/N(1) - r2*B2*S(idx)*E(idx)/N;
    E(idx 1) = E(idx)   r*B*S(idx)*I(idx)/N(1)-a*E(idx)   r2*B2*S(idx)*E(idx)/N;
    I(idx 1) = I(idx)   a*E(idx) - y*I(idx);
    R(idx 1) = R(idx)   y*I(idx);
    m(idx 1) = E(idx 1)   I(idx 1);
end


x=1:182;
plot(x,m);grid on;
xlabel('day');ylabel('Demand for drugs')

0 人点赞