文章索引
显示
SEIR模型
\( \frac{\mathrm{d}}{\mathrm{d}t}E=r\beta I\frac{S}{N}-\alpha E+r_2 \beta_2 E\frac{S}{N} \\ \frac{\mathrm{d}}{\mathrm{d}t}S=-r\beta I\frac{S}{N}-r_2 \beta_2 E\frac{S}{N} \\ \frac{\mathrm{d}}{\mathrm{d}t}R=\gamma I \\ \frac{\mathrm{d}}{\mathrm{d}t}I=\alpha E-\gamma I\)
查看所有模型及说明:点击病毒传染的动态扩散模拟及数学模型
clear
figure
N = 100; %人口总数
E = 0; %潜伏者
I = 1; %传染者
S = N - I; %易感者
R = 0; %康复者
r = 2; %感染者接触易感者的人数
B = 0.005; %传染概率
a = 0.1; %潜伏者转化为感染者概率
r2 = 15; %潜伏者接触易感者的人数
B2 = 0.015; %潜伏者传染正常人的概率
y = 0.025; %康复概率
D = 20;
days = 400;
Dynamic = false; %设置为true显示动态模拟过程
state=zeros(1,N);
index=randperm(N,I);
state(index)=1;
if Dynamic
axis([0 N 0 N]);
hold on
end
for i=1:days
xpos=randperm(N);
ypos=randperm(N);
if i>=D
r=5;
r2=5;
end
for j=1:N
if state(j)==1
if rand()<=y
state(j)=2;
continue
end
dis=sqrt((xpos-xpos(j)).^2+(ypos-ypos(j)).^2);
l=0;
for k=1:N
if state(k)==1 || state(k)==2 || state(k)==3
continue
end
l=l+1;
peo(l).dis=dis(k);
peo(l).num=k;
end
T = struct2table(peo);
sortedT = sortrows(T,'dis');
sortedS = table2struct(sortedT);
tp=sortedT{1:r,2};
for k=1:length(tp)
if rand()<=B
state(tp(k))=3;
end
end
elseif state(j)==3
if rand()<=a
state(j)=1;
end
dis=sqrt((xpos-xpos(j)).^2+(ypos-ypos(j)).^2);
l=0;
for k=1:N
if state(k)==1 || state(k)==2 || state(k)==3
continue
end
l=l+1;
peo(l).dis=dis(k);
peo(l).num=k;
end
T = struct2table(peo);
sortedT = sortrows(T,'dis');
sortedS = table2struct(sortedT);
tp=sortedT{1:r2,2};
for k=1:length(tp)
if rand()<=B2
state(tp(k))=3;
end
end
end
end
Idata(i)=length(find(state==1));
Sdata(i)=length(find(state==0));
Rdata(i)=length(find(state==2));
Edata(i)=length(find(state==3));
%a=-8;
%b=8;
%xpos=xpos+ a + (b-a).*rand(1,N);
%ypos=ypos+ a + (b-a).*rand(1,N);
if Dynamic && rem(i,10)==0
phd=scatter(xpos,ypos,[],state,"filled");
title(['第' num2str(i) '天'])
drawnow
delete(phd)
end
end
figure
set(gcf,'visible',true)
days=length(Sdata);
plot(1:days,Sdata,1:days,Idata,1:days,Rdata,1:days,Edata)
hold on
T = 1:days+100;
for idx = 1:length(T)-1
if S(idx)<0
S(idx)=0;
elseif S(idx)>N
S(idx)=N;
end
if E(idx)<0
E(idx)=0;
elseif E(idx)>N
E(idx)=N;
end
if I(idx)<0
I(idx)=0;
elseif I(idx)>N
I(idx)=N;
end
if R(idx)<0
R(idx)=0;
elseif R(idx)>N
R(idx)=N;
end
if idx>=D
r=5;
r2=5;
end
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);
end
plot(T,S,T,E,T,I,T,R);
plot([D D],[0 N])
legend('模拟易感者','模拟传染者','模拟康复者','模拟潜伏者','易感者','潜伏者','传染者','康复者','执行戒严措施')
xlabel('天');ylabel('人数')
title('SEIR模型2-戒严措施')
Code language: Matlab (matlab)
模拟最终截图如下(动态过程未给出):
未实施戒严措施:
从20天开始实施戒严措施: