新冠病毒的SEIR模型仿真---matlab

    技术2023-08-30  129

    ## SEIR模型

    关于什么是SEIR模型我这里也不做过多解释了,我也看到过很多关于它的文章,都写的很好。我这里主要讲一下分阶段模型的建立,这个需要根据疫情防控的不同阶段来进行各个参数的设置,使模型的预测值更接近真实值,具体参考大家可以看这篇论文《多阶段动态时滞动力学模型的COVID⁃19 传播分析》。 下面上代码

    %SEIR模型 clear;clc; %参数设置 [data,~]=xlsread('data.xls');%读入实际疫情数据可以和预测的画图做比较 I=41;%114号初始传染者 R=7;%康复者 D=1;%死亡患者数量 E=0;%潜伏者 N=1400000000;%全国人口 N1=10890000;%武汉市人口 S=N1-I;%易感染者 r1=12;%开始发病初期,接触病患的人数 r2=8;%疫情发生后还没采取更严格的措施接触病患的人数 r3=1.5;%疫情强制管控 r4=1; a=0.125;%潜伏者患病概率 B=0.6;%感染概率 %B2=0.05;%感染概率 y1=0.23;%医疗物资缺乏时康复概率 y2=0.95;%医疗物资等充足时康复概率 k1=0.065373;%疫情前期病死率 k2=0.05373;%武汉基本病死率 k3=0.035373;%湖北省病死率 k4=0.0173;%全国病死率 t_1=2; T=1:75; %% 时间 t0=datetime(2020,01,14);%武汉 data1=[]; for i=1:length(T) data1=[data1,t0+caldays(i)]; % caldays自增,获取数组内的日期格式数据,便于绘图,下同 end %% for idx =1:length(T)-1 if idx<14%第一阶段 S(idx+1)=S(idx)-r1*B*I(idx)*S(idx)/N1;%易感人数迭代 E(idx+1)=E(idx)+r1*B*S(idx)*I(idx)/N1-a*E(idx)%潜伏者人数迭代 I(idx+1)=I(idx)+a*E(idx)-(y1+k1)*I(idx);%患病人数迭代 R(idx+1)=R(idx)+y1*I(idx);%康复人数迭代 D(idx+1)=D(idx)+k1*I(idx);%死亡患者人数迭代 elseif idx>=14&idx<28%2阶段 S(idx+1)=S(idx)-r1*B*I(idx)*S(idx)/N;%易感人数迭代 E(idx+1)=E(idx)+r1*B*S(idx)*I(idx)/N-a*E(idx-t_1)%潜伏者人数迭代 I(idx+1)=I(idx)+a*E(idx)-(y1+k2)*I(idx);%患病人数迭代 R(idx+1)=R(idx)+y1*I(idx);%康复人数迭代 D(idx+1)=D(idx)+k2*I(idx);%死亡患者人数迭代 elseif idx>=28&idx<42%3阶段 S(idx+1)=S(idx)-r2*B*I(idx)*S(idx)/N;%易感人数迭代 E(idx+1)=E(idx)+r2*B*S(idx)*I(idx)/N-a*E(idx-t_1)%潜伏者人数迭代 I(idx+1)=I(idx)+a*E(idx)-(y1+k3)*I(idx);%患病人数迭代 R(idx+1)=R(idx)+y1*I(idx);%康复人数迭代 D(idx+1)=D(idx)+k3*I(idx);%死亡患者人数迭代 elseif idx>=42&idx<56%4阶段 S(idx+1)=S(idx)-r3*B*I(idx)*S(idx)/N;%易感人数迭代 E(idx+1)=E(idx)+r3*B*S(idx)*I(idx)/N-a*E(idx-t_1)%潜伏者人数迭代 I(idx+1)=I(idx)+a*E(idx)-(y2+k3)*I(idx);%患病人数迭代 R(idx+1)=R(idx)+y2*I(idx);%康复人数迭代 D(idx+1)=D(idx)+k3*I(idx);%死亡患者人数迭代 elseif idx>=56&idx<70%5阶段 S(idx+1)=S(idx)-r4*B*I(idx)*S(idx)/N;%易感人数迭代 E(idx+1)=E(idx)+r4*B*S(idx)*I(idx)/N-a*E(idx-t_1)%潜伏者人数迭代 I(idx+1)=I(idx)+a*E(idx)-(y2+k4)*I(idx);%患病人数迭代 R(idx+1)=R(idx)+y2*I(idx);%康复人数迭代 D(idx+1)=D(idx)+k4*I(idx);%死亡患者人数迭代 elseif idx>=70%6阶段 S(idx+1)=S(idx)-r4*B*I(idx)*S(idx)/N;%易感人数迭代 E(idx+1)=E(idx)+r4*B*S(idx)*I(idx)/N-a*E(idx-t_1)%潜伏者人数迭代 I(idx+1)=I(idx)+a*E(idx)-(y2+k4)*I(idx);%患病人数迭代 R(idx+1)=R(idx)+y2*I(idx);%康复人数迭代 D(idx+1)=D(idx)+k4*I(idx);%死亡患者人数迭代 end end figure plot(data1,E,data1,I,data1,R,data1,D); grid on; xlabel('日期'); ylabel('人数'); legend('潜伏者','传染者','康复者','死亡者'); title('预测模型');

    可以看到模型预测出的患病人数在8万左右,和我们国内疫情患病人数是非常接近的。这里的峰值拐点比实际时间提前到,这和实际的疫情防控、检测能力等外在因素是密切相关的。 现在疫情在国外大多数国家的发展差不多只有前两个阶段,没用强有力的隔离防控措施的话是很难抑制它的传播的。

    Processed: 0.019, SQL: 9