estimation of parameters using lsqnonlin for system of differential equations

2 次查看（过去 30 天）

mallela ankamma rao 2022-3-29

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1683279-estimation-of-parameters-using-lsqnonlin-for-system-of-differential-equations

评论： mallela ankamma rao 2022-4-1

data 
1	36456
2	33761
3	38309
4	36638
5	36111
6	32272
7	27107
8	32067
9	26357
10	34666
11	30034
12	30354
13	27336
14	21941
15	26401
16	18164
17	26762
18	26991
19	26834
20	24589
21	19174
22	23886
23	24236
24	23924
25	22350
26	18574
27	20333
28	16072
29	20529
30	21957
31	19046
32	19522
33	19460
34	18938
35	18967
36	18591
37	18381
38	18252
39	18008
40	18105
41	18020
42	17478
43	17191
44	17297
45	16348
46	15877
47	15338
48	15034
49	14678
50	14376
51	14128
52	13959
53	13840
54	13818
55	13714
56	13588
57	13437
58	12905
59	13535
60	13354
61	13094
62	12850
63	12776
64	12529
65	12724
66	11794
67	11602
68	11459
69	11496
70	11515
71	11477
72	11422
73	10986
74	11047
75	11066
76	11054
77	11107
78	11230
79	11276
80	11831
81	12085
82	12331
83	12699
84	12900
85	13198
86	13780
87	14260
88	14638
89	15050
y = data(:,2);
t = data(:,1);
N = 1390000000;
e0 = 1217378;
a0 = 150;
i0 = 100;
q0 = 30;
h0 = 10;
r0 = 10;
d0 = 0;
s0 = N-e0-a0-i0-q0-h0-r0-d0;
y0 = [s0; e0; a0; i0; q0; h0; r0; d0];
beta0 = 0.5;
alpha0 = 0.004;
gamma0 = 0.1;
upsilon0 = 0.13;
epsilon0= 0.07;
lamda0= 0.1;
sigma0= 0.07;
kappa0= 0.03;
nu0 = 0.0001;	
xi0 = 0.0002;
p0 = [beta0; alpha0; gamma0;  upsilon0; epsilon0; lamda0; sigma0; kappa0; nu0; xi0 ];
options = optimoptions(@lsqnonlin,'Algorithm','trust-region-reflective');
y = lsqnonlin(Immanuel(t, y, p),y0,[],[],options);
function Immanuel = Immanuel(t, y, p)
% RHS of the ODE
N = 1390000000;
pi = 70000;
zeta = 0.1;
eta = 0.2;
theta = 0.3;
iota = 0.3;
delta = 0.1;
rho = 0.5;
mu = 0.0000425;
beta = p(1);
alpha =  p(2);
gamma =  p(3);
upsilon = p(4);
epsilon = p(5);
lamda = p(6);
sigma = p(7);
kappa = p(8);
nu = p(9);	
xi = p(10);
ydot = zeros(10,1);
S = y(1);
E = y(2);
A = y(3);
I = y(4);
Q = y(5);
H = y(6);
R = y(7);
D = y(8);
ydot(1) = pi -beta*(zeta*A+eta*I+theta*Q+iota*H)*(S/N) -mu*S;
ydot(2) = beta*(zeta*A+eta*I+theta*Q+iota*H)*(S/N) -(delta+mu)*E; 
ydot(3) = rho*delta*E-(lamda+gamma+nu+mu)*A;
ydot(4) = (1-rho)*delta*E-(sigma+kappa+mu)*I;
ydot(5) = lamda*A + sigma*I-(alpha+upsilon+mu)*Q;
ydot(6) = alpha*H  + kappa*I- (epsilon+xi+mu)*H;
ydot(7) =  gamma*A + upsilon*Q + epsilon*H;
ydot(8) = nu*A  + xi*H;
end