求解微分方程参数拟合的问题并模拟混沌现象
新人金币终于够了，可以发帖了。
现在写论文遇到了一个障碍，是希望通过实际数据进行一个微分方程的参数辨识，然后利用拟合出来的参数进行混沌现象的模拟以及时间序列图的模拟等。
之前看过一篇参考文献对实际数据的处理如下图所示，参考这篇文献我也想对我的方程进行实际数据的参数辨识，我的方程如下所示：
用到的原始数据如下：
x y z
2008 3923 797.24 6589.35
2009 5843 481.14 7295.21
2010 8194 574.47 9387.21
2011 12781 737.65 13217.37
2012 15247 833.77 15486.97
2013 18052 923.19 16735.09
2014 22168 1022.37 19691.06
2015 20820 1076.14 20126.82
2016 25212 1133.57 22801.22
2017 24211 1246.67 22653.54
归一化插值法以后得到的数据如下：
X Y Z
0 0.412916541 0
0.045093711 0.206458271 0.021769852
0.090187421 0 0.043539703
0.14540373 0.060957768 0.108060329
0.200620039 0.121915536 0.172580955
0.308351731 0.228495291 0.290709215
0.416083423 0.335075046 0.408837475
0.474000658 0.397855081 0.478835569
0.531917892 0.460635116 0.548833663
0.597796985 0.519039098 0.58732768
0.663676077 0.577443079 0.625821697
0.760345719 0.642221729 0.71698854
0.85701536 0.707000379 0.808155382
0.825355818 0.742119839 0.821594918
0.793696275 0.777239298 0.835034453
0.896848138 0.814749259 0.917517227
1 0.852259219 1
0.976490206 0.92612961 0.995445313
0.952980412 1 0.990890625
求大神帮忙看看怎么用代码实现，可有偿！
微信图片_20200109201816.png
图片1.png@月只蓝
返回小木虫查看更多
京公网安备 11010802022153号
微分方程组参数拟合 ，参考
http://blog.sina.com.cn/s/blog_c0cb8ce60102ysqu.html
function ODEfunctiongroup
clear all;
clc;
k0 = [0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01]; % 参数初值
lb = [0 0 0 0 0 0 0 0 0 0 0 0 0 0]; % 参数下限
ub = [+inf +inf +inf +inf +inf +inf +inf +inf +inf +inf +inf +inf +inf +inf]; % 参数上限
x0 = [0.015 0.758 1.83]; %初始状态
% 动力学数据:
% t y1 y2 y3
ExpData = ...
[ 0 0 0.412916541 0
1 0.045093711 0.206458271 0.021769852
2 0.090187421 0 0.043539703
3 0.14540373 0.060957768 0.108060329
4 0.200620039 0.121915536 0.172580955
5 0.308351731 0.228495291 0.290709215
6 0.416083423 0.335075046 0.408837475
7 0.474000658 0.397855081 0.478835569
8 0.531917892 0.460635116 0.548833663
9 0.597796985 0.519039098 0.58732768
10 0.663676077 0.577443079 0.625821697
11 0.760345719 0.642221729 0.71698854
12 0.85701536 0.707000379 0.808155382
13 0.825355818 0.742119839 0.821594918
14 0.793696275 0.777239298 0.835034453
15 0.896848138 0.814749259 0.917517227
16 1 0.852259219 1
17 0.976490206 0.92612961 0.995445313
18 0.952980412 1 0.990890625
];
yexp = ExpData(:,2:4); % yexp: 对应实验数据[x1 x4 x5 x6]
% 使用函数fmincon()进行参数估计
[k,fval,flag] = fmincon(@ObjFunc4Fmincon,k0,[],[],[],[],[],[],[],[],[],[],[],[],[],lb,ub,[],[],x0,yexp);
fprintf('\n使用函数fmincon()估计得到的参数值为:\n')
fprintf('\tk1 = %.4f\n',k(1))
fprintf('\tk2 = %.4f\n',k(2))
fprintf('\tk3 = %.4f\n',k(3))
fprintf('\tk4 = %.4f\n',k(4))
fprintf('\tk5 = %.4f\n',k(5))
fprintf('\tk5 = %.4f\n',k(6))
fprintf('\tk5 = %.4f\n',k(7))
fprintf('\tk5 = %.4f\n',k(8))
fprintf('\tk5 = %.4f\n',k(9))
fprintf('\tk5 = %.4f\n',k(10))
fprintf('\tk5 = %.4f\n',k(11))
fprintf('\tk5 = %.4f\n',k(12))
fprintf('\tk5 = %.4f\n',k(13))
fprintf('\tk5 = %.4f\n',k(14))
fprintf('The sum of the squares is: %.1e\n\n',fval)
k_fmincon = k;
% 这一步通常被省略，通过反复迭代初始值得到最优解，加上后可以降低对初始值的依赖。
% 以函数fmincon()估计得到的结果为初值，使用函数lsqnonlin()进行参数估计
% 需要指出，这种方法并非在所有场合均有效，但有时确实可以改善求解效果。
% 使用函数lsqnonlin()进行参数估计
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
lsqnonlin(@ObjFunc4LNL,k0,lb,ub,[],x0,yexp);
ci = nlparci(k,residual,jacobian);
fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
Output(k,ci,resnorm)
% 以函数fmincon()估计得到的结果为初值，使用函数lsqnonlin()进行参数估计
k0 = k_fmincon;
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
lsqnonlin(@ObjFunc4LNL,k0,lb,ub,[],x0,yexp);
ci = nlparci(k,residual,jacobian);
fprintf('\n\n以fmincon()的结果为初值，使用函数lsqnonlin()估计得到的参数值为:\n')
Output(k,ci,resnorm)
function dxdt = KineticEqs(t,x,k)
dxdt =[k(1)*x(1)*(1-x(2)/k(10))+k(2)*x(2)+k(3)*x(3)
k(4)*x(1)*(1-x(1)/k(14))+k(5)*x(2)*(1-x(2)/k(12))+k(6)*x(3)*(1-x(3)/k(13))
k(8)*x(1)*(1-x(1)/k(14))+k(8)*x(2)*(x(2)/k(12)-1)+k(9)*x(3)*(x(2)/k(11)-1)];
这是我根据您的那篇文章修改的代码，但是跑不出来，希望您指导一下，MATLAB小白太痛苦了
，
数据归一化后，模型也需要归一化处理。
另外，你的模型只有9个参数。