Skip to content
Snippets Groups Projects
Commit 93f2c210 authored by Rivnyák Tímea's avatar Rivnyák Tímea
Browse files

SIRC_RLSQ modificaticon

parent ab355113
Branches
No related tags found
No related merge requests found
...@@ -12,7 +12,7 @@ C0 = 0; ...@@ -12,7 +12,7 @@ C0 = 0;
S0 = N - I0 - R0 - C0; S0 = N - I0 - R0 - C0;
% Initial parameteres % Initial parameteres
nu = 14; nu = 14;
mu = 0.015; mu = 0.013;
epsilon = 1.2; epsilon = 1.2;
beta = 0.0002; beta = 0.0002;
gamma = 0.05; gamma = 0.05;
...@@ -31,31 +31,11 @@ xlabel('Time(days)'); ...@@ -31,31 +31,11 @@ xlabel('Time(days)');
ylabel('Number of individuals'); ylabel('Number of individuals');
legend('S', 'I', 'R', 'C'); legend('S', 'I', 'R', 'C');
%% Recursive Least Square without noise
%% Add noise
noise = wgn(size(t,1), 4, 0, 42);
y_noise = y + noise;
figure(2)
plot(t, y_noise);
xlabel('Time(days)');
ylabel('Number of individuals');
legend('S', 'I', 'R', 'C');
[t_noisy, y_noisy] = ode45(@(t, y) noisy_deriv(y, nu, mu, epsilon, beta, gamma, Gamma, q), [0 time], y0);
figure(3)
plot(t_noisy, y_noisy);
xlabel('Time(days)');
ylabel('Number of individuals');
legend('S', 'I', 'R', 'C');
%% Recursive Least Square
h = 1.0; h = 1.0;
lambda0 = 0.95; % 0.95-0.99 lambda0 = 0.95; % 0.95-0.99
alpha = 0.7; alpha = 0.5;
N = length(y); N = length(y);
% dSdt = nu - beta*I*S - epsilon*beta*C*S - mu*S; % dSdt = nu - beta*I*S - epsilon*beta*C*S - mu*S;
...@@ -67,22 +47,21 @@ phi_S(:,3) = y(1:end-1, 4) .* y(1:end-1, 1); ...@@ -67,22 +47,21 @@ phi_S(:,3) = y(1:end-1, 4) .* y(1:end-1, 1);
phi_S(:,4) = y(1:end-1, 1); phi_S(:,4) = y(1:end-1, 1);
[theta_S, p_S] = recursiveLSQ(phi_S, y_S, lambda0, alpha) [theta_S, p_S] = recursiveLSQ(phi_S, y_S, lambda0, alpha)
figure(4) figure(2)
plot(1:size(theta_S, 1), theta_S); plot(theta_S);
legend('nu', '-beta', '-epsilon*beta', '-mu') legend('nu', '-beta', '-epsilon*beta', '-mu')
% dIdt = beta*I*S + epsilon*beta*C*S - gamma*I - mu*I; % dIdt = beta*I*S + epsilon*beta*C*S - (gamma + mu)*I;
y_I = (y(2:end, 2) - y(1:end-1, 2)) / h; y_I = (y(2:end, 2) - y(1:end-1, 2)) / h;
phi_I = zeros(N-1, 4); phi_I = zeros(N-1, 3);
phi_I(:,1) = y(1:end-1, 2) .* y(1:end-1, 1); phi_I(:,1) = y(1:end-1, 2) .* y(1:end-1, 1);
phi_I(:,2) = y(1:end-1, 4) .* y(1:end-1, 1); phi_I(:,2) = y(1:end-1, 4) .* y(1:end-1, 1);
phi_I(:,3) = y(1:end-1, 2); phi_I(:,3) = y(1:end-1, 2);
phi_I(:,4) = y(1:end-1, 2);
[theta_I, p_I] = recursiveLSQ(phi_I, y_I, lambda0, alpha) [theta_I, p_I] = recursiveLSQ(phi_I, y_I, lambda0, alpha)
figure(5) figure(3)
plot(1:size(theta_I, 1), theta_I); plot(theta_I);
legend('beta', 'epsilon*beta', '-gamma', '-mu') legend('beta', 'epsilon*beta', '-(gamma+mu)')
% dRdt = gamma*(1-q)*I + Gamma*C - mu*R; % dRdt = gamma*(1-q)*I + Gamma*C - mu*R;
y_R = (y(2:end, 3) - y(1:end-1, 3)) / h; y_R = (y(2:end, 3) - y(1:end-1, 3)) / h;
...@@ -92,19 +71,18 @@ phi_R(:,2) = y(1:end-1, 4); ...@@ -92,19 +71,18 @@ phi_R(:,2) = y(1:end-1, 4);
phi_R(:,3) = y(1:end-1, 3); phi_R(:,3) = y(1:end-1, 3);
[theta_R, p_R] = recursiveLSQ(phi_R, y_R, lambda0, alpha) [theta_R, p_R] = recursiveLSQ(phi_R, y_R, lambda0, alpha)
figure(6) figure(4)
plot(1:size(theta_R, 1), theta_R); plot(theta_R);
legend('gamma*(1-q)', 'Gamma', '-mu') legend('gamma*(1-q)', 'Gamma', '-mu')
% dCdt = gamma*q*I - Gamma*C - mu*C; % dCdt = gamma*q*I - (Gamma + mu)*C;
y_C = (y(2:end, 4) - y(1:end-1, 4)) / h; y_C = (y(2:end, 4) - y(1:end-1, 4)) / h;
phi_C = zeros(N-1, 3); phi_C = zeros(N-1, 2);
phi_C(:,1) = y(1:end-1, 2); phi_C(:,1) = y(1:end-1, 2);
phi_C(:,2) = y(1:end-1, 4); phi_C(:,2) = y(1:end-1, 4);
phi_C(:,3) = y(1:end-1, 4);
[theta_C, p_C] = recursiveLSQ(phi_C, y_C, lambda0, alpha) [theta_C, p_C] = recursiveLSQ(phi_C, y_C, lambda0, alpha)
figure(7) figure(5)
plot(1:size(theta_C, 1), theta_C); plot(theta_C);
legend('gamma*q', '-Gamma', '-mu') legend('gamma*q', '-(Gamma+mu)')
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment