Error in Model Fitting
48 views (last 30 days)
Show older comments
Abdullahi
on 7 Aug 2024 at 23:31
Hello! i have been to run this model fitting code but keeps on getting these error messages as:
- FUN must have two input arguments.
- userFcn_ME = initEvalErrorHandler(userFcn_ME,funfcn_x_xdata{3}, ...
- [params_fit, ~] = lsqcurvefit(obj_fun, params0, time, observed_I, [], [], options);
Find attached, the code,, Please help
function dydt = sir_model(t, y, beta, gamma)
S = y(1);
I = y(2);
R = y(3);
N = S + I + R; % Total population
dSdt = -beta * S * I / N;
dIdt = beta * S * I / N - gamma * I;
dRdt = gamma * I;
dydt = [dSdt; dIdt; dRdt];
end
function error = model_error(params, time, observed_I, y0)
beta = params(1);
gamma = params(2);
% Solve the differential equations
[~, y] = ode45(@(t, y) sir_model(t, y, beta, gamma), time, y0);
% Extract the infectious data from the solution
model_I = y(:, 2);
% Compute the error as the difference between model and observed data
error = model_I - observed_I';
end
% Define the time vector and observed data
time = [0, 1, 2, 3, 4, 5]; % Example time points
observed_I = [10, 15, 20, 25, 30, 35]; % Example observed infectious data
% Initial guess for parameters
beta_guess = 0.3;
gamma_guess = 0.1;
params0 = [beta_guess, gamma_guess];
% Initial conditions
S0 = 1000; % Initial number of susceptible individuals
I0 = observed_I(1); % Initial number of infectious individuals
R0 = 0; % Initial number of recovered individuals
y0 = [S0; I0; R0];
% Define the objective function for fitting
obj_fun = @(params) model_error(params, time, observed_I, y0);
% Perform the curve fitting
options = optimoptions('lsqcurvefit', 'Display', 'off');
[params_fit, ~] = lsqcurvefit(obj_fun, params0, time, observed_I, [], [], options);
% Display the fitted parameters
beta_fit = params_fit(1);
gamma_fit = params_fit(2);
disp(['Fitted beta: ', num2str(beta_fit)]);
disp(['Fitted gamma: ', num2str(gamma_fit)]);
0 Comments
Accepted Answer
Torsten
on 7 Aug 2024 at 23:46
Edited: Torsten
on 7 Aug 2024 at 23:55
% Define the time vector and observed data
time = [0, 1, 2, 3, 4, 5]; % Example time points
observed_I = [10, 15, 20, 25, 30, 35]; % Example observed infectious data
% Initial guess for parameters
beta_guess = 0.3;
gamma_guess = 0.1;
params0 = [beta_guess, gamma_guess];
% Initial conditions
S0 = 1000; % Initial number of susceptible individuals
I0 = observed_I(1); % Initial number of infectious individuals
R0 = 0; % Initial number of recovered individuals
y0 = [S0; I0; R0];
% Define the objective function for fitting
obj_fun = @(params,time) model_error(params, time, y0);
% Perform the curve fitting
options = optimoptions('lsqcurvefit', 'Display', 'off');
[params_fit, ~] = lsqcurvefit(obj_fun, params0, time, observed_I, [], [], options);
% Display the fitted parameters
beta_fit = params_fit(1);
gamma_fit = params_fit(2);
disp(['Fitted beta: ', num2str(beta_fit)]);
disp(['Fitted gamma: ', num2str(gamma_fit)]);
hold on
plot(time,observed_I,'o')
plot(time,model_error(params_fit,time,y0))
hold off
grid on
function error = model_error(params, time, y0)
beta = params(1);
gamma = params(2);
% Solve the differential equations
[~, y] = ode45(@(t, y) sir_model(t, y, beta, gamma), time, y0);
% Extract the infectious data from the solution
model_I = y(:, 2);
% Compute the error as the difference between model and observed data
error = model_I.';
end
function dydt = sir_model(t, y, beta, gamma)
S = y(1);
I = y(2);
R = y(3);
N = S + I + R; % Total population
dSdt = -beta * S * I / N;
dIdt = beta * S * I / N - gamma * I;
dRdt = gamma * I;
dydt = [dSdt; dIdt; dRdt];
end
6 Comments
More Answers (0)
See Also
Categories
Find more on Surface and Mesh Plots in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!