改进SEIR模型的matlab代码.zip

%SEIR模型修正 %假设1月15日开始出现第一例确诊，；1月23号武汉政府开始封城，此时其他省市也相应作出响应隔离措施，大约距离15号后的11天政府管控发挥明显作用 % N=13.95*100000000;%人口数 N=1395380000; load quanguo.mat %第一列为累计确诊人数，第二列为累计死亡人数，第三列为累计治愈人数 E=0;%潜伏者 D=0;%死亡患者人数 I=1;%感染人数 S=N-I;%易感人数 R=0;%康复者人数 r=17;%感染者接触数量 % r=19; B=0.602;%传染概率 % a=0.17;%潜伏者转化为感染者的概率 % a=0.175; a=0.198;%潜伏者转化为感染者的概率 % r2=8;%潜伏者接触人数 r2=15;%潜伏者接触人数 % B2=0.03;%潜伏者传染正常人的概率 B2=0.05; y=0.05;%康复概率 k=0.0001;%日致死率 B3=0.1;%转阴率 % T=1:200; T=1:180; for idx=1:length(T)-1 %若以1月18日为疫情起点，在政府发出管控号召时间以及各个地方响应延迟时间，此处采用11天后为临界点， %相当于11天后，感染者与潜伏者流动性和医疗配置发生明显变化，具体为接触人数 if idx>=11 r=0.20;%感染者接触数量 r2=1.8;%感染者接触人数量 y=0.15;%康复率上升为0.16 a=0.12;%潜伏者转化为感染者的概率 % a=0.10; k=0.0001;%日致死率暂无改变 end if idx<11 S(idx+1)=S(idx)-r*B*S(idx)*I(idx)/N-r2*B2*S(idx)*E(idx)/N;%易感人群迭代 E(idx+1)=E(idx)+r*B*S(idx)*I(idx)/N-a*E(idx)+r2*B2*S(idx)*E(idx)/N;%潜伏者迭代 I(idx+1)=I(idx)+a*E(idx)-(k+y)*I(idx);%感染人数迭代 R(idx+1)=R(idx)+0.05*I(idx);%康复人数迭代 D(idx+1)=R(idx)+k*I(idx);%死亡患者人数迭代 else S(idx+1)=S(idx)-r*B*S(idx)*I(idx)/N-r2*B2*S(idx)*E(idx)/N+B3*E(idx-10);%易感人群迭代 E(idx+1)=E(idx)+r*B*S(idx)*I(idx)/N-a*E(idx)+r2*B2*S(idx)*E(idx)/N-B3*E(idx-10);%潜伏者迭代 I(idx+1)=I(idx)+a*E(idx)-(k+y)*I(idx);%感染人数迭代 % Y参数有问题 R(idx+1)=R(idx)+0.045*I(idx-9);%康复人数迭代 D(idx+1)=R(idx)+k*I(idx);%死亡患者人数迭代 end end B={'01-19','02-08','02-28','03-19','04-08','04-28','05-18','06-07','06-27','07-17','08-06'}; % plot(1:1:102,huibei(:,1)-huibei(:,2)-huibei(:,3),'r*');hold on plot(1:1:102,quanguo(:,1)-quanguo(:,2)-quanguo(:,3),'g-');hold on plot(1:1:102,quanguo(:,3),'k-');hold on % legend('实际患病','实际康复') % xlabel('天数'); % ylabel('人数'); % legend('实际患病') plot(T,R,'b',T,I,'r'); grid on; hold on; plot([7 7],[0 1000]); set(gca,'XTickLabel',B) xlabel('日期'); ylabel('人数'); legend('实际患病','实际康复','预测康复者','预测患者'); title('采取隔离措施的SEIR模型');