Error using nlinfit (Check FunVals)

3 次查看（过去 30 天）
Abdullahi
Abdullahi 2024-8-27
编辑： Matt J 2024-8-27
Ran in:
Hi there.. i am trying to fit my odel equations to some data but i keep on getiing error.. please help
below is my code
% Data fitting for the covid-19 outbreak in Lagos
covid19_nlinfit
Norm of Norm of Iteration SSE Gradient Step ----------------------------------------------------------- 0 8.39013e+06
Warning: Rank deficient, rank = 2, tol = 1.269216e-11.
1 7.37178e+06 5.39372e+06 2.01582
Warning: Failure at t=1.591672e+01. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (2.842171e-14) at time t.
Error using nlinfit>checkFunVals (line 649)
The function you provided as the MODELFUN input has returned Inf or NaN values.

Error in nlinfit>LMfit (line 596)
if funValCheck && ~isfinite(sse), checkFunVals(r); end

Error in nlinfit (line 284)
[beta,J,~,cause,fullr] = LMfit(X,yw, modelw,beta,options,verbose,maxiter);

Error in solution>covid19_nlinfit (line 36)
[values] = nlinfit(day, Totalcases, @covidinfection, [a(1) a(2) a(3)], options);
function covid19_nlinfit
clear all
clc
global a x1 beta_c ep_g c_g sigma nu theta psi d_o d_d gamma_i alpha ...
delta ep_s ep_i ep_a E0 A0 I0
ep_g = 0.9; c_g = 0.9; delta = 0; ep = 0; ep_s = ep; ep_a=ep; ep_i=ep;
sigma = 1/5.2; nu = 0.5;
gamma_i = 1/15;
gamma_a = 0.13978; gamma_o = 0.13978;
d_o = 0.015; d_d = 0.015;
alpha = 0.5;
%==========================================================================
% Lagos data, starting from March 16, 2020 till May 3, 2020.
day = [1 3 4 6 7 8 9 10 11 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 ...
43 44 45 46 47 48];
Totalcases = [1 5 9 19 22 25 29 32 44 59 68 81 82 91 98 109 109 120 120 130 145 158 166 174 176 189 214 232 251 283 ...
306 376 376 430 504 582 657 689 731 764 844 931 976 1006 1068];
Activecases = [1 5 8 19 22 25 29 32 43 58 67 75 76 85 81 86 86 91 87 96 111 117 118 119 116 123 139 140 154 182 200 ...
266 266 309 383 460 525 548 584 602 682 718 756 753 791];
%==========================================================================
%Initial guesses for parameters to be estimated
beta_c=0; theta=0;psi=0;
% E0 = 16; A0 = 20; I0 = 10;
%======================================================================
% Using NLINFIT function
a = [beta_c theta psi]; % E0 A0 I0;
x1 = day;
options = statset('Display', 'iter');
[values] = nlinfit(day, Totalcases, @covidinfection, [a(1) a(2) a(3)], options);
%[values] = nlinfit(day, Activecases, @covidinfection, [a(1) a(2) a(3) a(4) a(5) a(6)], options);
beta_c = values(1)
theta = values(2)
psi = values(3)
% E0 = values(4)
% A0 = values(5)
% I0 = values(6)
%=====================================================================
%Control Reproduction Number
R_c = beta_c*(1-ep_g*c_g)*(1-delta)*(1-ep_s)*((alpha*nu*(1-ep_a)/...
(theta+gamma_a)) + ((1-nu)*(1-ep_i)/(psi+d_o+gamma_o)));
[t,y] = ode15s(@covidmodel1, [0,50],[14368332,16,20,10,1,0,0]); %E0 A0 I0;
plot(t,y(:,7),'-b',day,Totalcases,'rs')
xlabel('Time (days)'),ylabel('Cumulative number of reported cases')
legend('Model','Lagos COVID-19 Data')
%plot(t,y(:,5),'-b',day,Activecases,'rs')
%========================================================================
% Using NLINFIT function for data fitting
function y = covidinfection(a,x1)
beta_c = a(1);
theta = a(2);
psi = a(3);
% E0 = a(4);
% A0 = a(5);
% I0 = a(6);
[t,y] = ode15s(@covidmodel1, [0,x1(end)],[14368332,16,20,10,1,0,0]); %E0 A0 I0;
y = interp1(t,y(:,7),x1); % cumulative reported cases
%y = interp1(t,y(:,5),x1); % active cases
end
%=========================================================================
function dx = covidmodel1(t,x)
dx = [0;0;0;0;0;0;0];
N = x(1)+x(2)+x(3)+x(4)+x(5)+x(6);
beta = beta_c*(1-ep_g*c_g)*(1-delta)*(1-ep_s);
dx(1) = -beta*x(1)*(alpha*(1-ep_a)*x(3)+(1-ep_i)*x(4))/N;
dx(2) = beta*x(1)*(alpha*(1-ep_a)*x(3)+(1-ep_i)*x(4))/N - sigma*x(2);
dx(3) = nu*sigma*x(2) - theta*x(3) - gamma_a*x(3);
dx(4) = (1-nu)*sigma*x(2) - psi*x(4) - d_o*x(4) - gamma_o*x(4);
dx(5) = theta*x(3) + psi*x(4) - gamma_i*x(5) -d_d*x(5); % Active cases, detected (asymptomatic and symptomatic)
dx(6) = gamma_i*x(5) + gamma_a*x(3) + gamma_o*x(4);
dx(7) = theta*x(3) + psi*x(4); % Cumulative number of new cases reported (asymptomatic and symptomatic)
end
end

请先登录，再回答此问题。

回答（1 个）

Matt J
Matt J 2024-8-27
编辑：Matt J 2024-8-27
For one thing,
y = interp1(t,y(:,7),x1,'spline','extrap'); % cumulative reported cases
For another, some sort of bounds or constraints on the variables are needed to prevent the ODE solver from considering ill-behaved systmes of equations. nlinfit does not let you impose such constraints, so you may have to look to another solver, e.g. lsqcurvefit.

类别

在 Help Center 和 File Exchange 中查找有关 Biological and Health Sciences 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by