I am trying to simulate a model using ssa for a selection of initial values and
ID: 2902865 • Letter: I
Question
I am trying to simulate a model using ssa for a selection of initial values and be able to discuss the results. can anyone let me know if my code make sense?
Basically, it's a simple SIR model over 2 patches. Individual are able to move from patch A to patch B.the rate of being infected in patch A is different from that of patch B. similarly, the rate of recovering in patch A is different from that of patch B.
clear all
close all
clc
% part b
%a1 = (beta1/n)*S1*I1
%a2 = (beta2/n)*S2*I2
%a3 = (beta1/n)*I1
%a4 = (beta/n)*I2
%a5 = (omega/n)*S1
%a6 = (omega/n)*S2
%a7 = (omega/n)*I1
%a8 = (omega/n)*I2
% first set up the initial conditions
% the number of people that are there at the start of the simulation(both
% infected and susceptible)
%then when the simulation would start and where it will end
%initial conditions
S1 = 100; %Number of susceptible in patch 1
S2 = 100; %Number of susceptible in patch 2
I1 = 10; %Number of infected in patch 1
I2 = 10; %Number of infected in patch 2
b1 = 0.8; %Contact rate in Patch 1
b2 = 0.8; %Contact rate in patch 2
g1 = 0.4; %recovery rate in patch 1
g2 =0.4; % recovery rate in patch 2
w1 =0.08; % diffusion rate
t = 0;
t_end = 50;
k = [b1 b2 g1 g2 w1];
X = [S1;S2;I1;I2];
N = S1+S2+I1+I2;
%state change vector
v_1= [-1;1;0;0];
v_2= [1;-1;0;0];
v_3= [0;-1;0;0];
v_4= [0 ;0;0;-1];
v_5= [-1;0;1;0];
v_6= [1;0;-1;0];
v_7= [0;-1;0;1];
v_8 =[0;1;0;-1];
% Taking a time step
while t(end) a(1) =(k(1)/N)*X(1,end)*X(3,end);
a(2) =(k(2)/N)*X(2,end)*X(4,end);
a(3) = k(3)*X(3,end);
a(4) = k(4)*X(4,end);
a(5) = k(5)*X(1,end);
a(6) = k(5)*X(2,end);
a(7) = k(5)*X(3,end);
a(8) = k(5)*X(4,end);
a_0 = sum(a);
u = rand(2,1);
dt = -1/a_0* log(u(1));
t = [t, t(end) + dt];
if u(2) < a(1)/a_0;
v = v_1;
elseif u(2)< (a(1)+a(2))/a_0;
v = v_2;
elseif u(2)<(a(1)+a(2)+a(3))/a_0;
v = v_3;
elseif u(2)<(a(1)+a(2)+a(3)+ a(4))/a_0;
v = v_4;
elseif u(2)<(a(1)+a(2)+a(3)+a(4)+a(5))/a_0;
v = v_5;
elseif u(2)<(a(1)+a(2)+a(3)+a(4)+a(5)+a(6))/a_0;
v = v_6;
elseif u(2)<(a(1)+a(2)+a(3)+a(4)+a(5)+a(6)+a(7))/a_0;
v = v_7;
elseif u(2)<(a(1)+a(2)+a(3)+a(4)+a(5)+a(6)+a(7)+a(8))/a_0;
v = v_8;
end
X = [X, X(:,end) + v];
end
figure
plot(t,X(1,:)) % this plot the number of susceptible in patch 1
hold on
plot(t,X(2,:),'k') % this plot the number of susceptible in patch 2
hold on
plot(t,X(3,:),'c')% this plot the number of infected in patch 1
hold on
plot(t,X(4,:),'r'),% this plot the number of infected in patch 2
hold off
Explanation / Answer
clear all close all clc % part b %a1 = (beta1/n)*S1*I1 %a2 = (beta2/n)*S2*I2 %a3 = (beta1/n)*I1 %a4 = (beta/n)*I2 %a5 = (omega/n)*S1 %a6 = (omega/n)*S2 %a7 = (omega/n)*I1 %a8 = (omega/n)*I2 % first set up the initial conditions % the number of people that are there at the start of the simulation(both % infected and susceptible) %then when the simulation would start and where it will end %initial conditions S1 = 100; %Number of susceptible in patch 1 S2 = 100; %Number of susceptible in patch 2 I1 = 10; %Number of infected in patch 1 I2 = 10; %Number of infected in patch 2 b1 = 0.8; %Contact rate in Patch 1 b2 = 0.8; %Contact rate in patch 2 g1 = 0.4; %recovery rate in patch 1 g2 =0.4; % recovery rate in patch 2 w1 =0.08; % diffusion rate t = 0; t_end = 50; k = [b1 b2 g1 g2 w1]; X = [S1;S2;I1;I2]; N = S1+S2+I1+I2; %state change vector v_1= [-1;1;0;0]; v_2= [1;-1;0;0]; v_3= [0;-1;0;0]; v_4= [0 ;0;0;-1]; v_5= [-1;0;1;0]; v_6= [1;0;-1;0]; v_7= [0;-1;0;1]; v_8 =[0;1;0;-1]; % Taking a time step while t(end)Related Questions
drjack9650@gmail.com
Navigate
Integrity-first tutoring: explanations and feedback only — we do not complete graded work. Learn more.