我自己生了一些模擬的資料 2x2 ODE, 1x15 double
然後想用fminsearch找出最佳的參數 p(1) p(2) p(3)
但收到Error如下：
Index exceeds array bounds.
Error in TEST (line 5)
dydt(1) = -(p(1)+p(3))*y(1) + p(2)*y(2);
Error in odefit>@(t,y)TEST(t,y,p) (line 2)
[t,y] = ode45(@(t,y)TEST(t,y,p),exp_t,0) % i am using y0=0 you can choose
whatever.
Error in odearguments (line 90)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 115)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options,
varargin);
Error in odefit (line 2)
[t,y] = ode45(@(t,y)TEST(t,y,p),exp_t,0) % i am using y0=0 you can choose
whatever.
Error in @(p)odefit(exp_t,exp_y,p)
Error in fminsearch (line 200)
fv(:,1) = funfcn(x,varargin{:});
我的程式是由五個小.m組成的
其中的Main:
exp_y = transpose(call_ex);
exp_t = linspace(0.5,7.5,15);
p0 = [0.6, 0.1, 0.3];
p_estimate = fminsearch(@(p)odefit(exp_t,exp_y,p),p0);
求Residual的odefit:
function err = odefit(exp_t,exp_y,p)
[t,y] = ode45(@(t,y)TEST(t,y,p),exp_t,0)
err = sum((y-exp_y).^2);
end
但有三個未知參數的假設2x2 ODE:
function [dydt] = TEST(t,y,p)
dydt = zeros(2,1); % this creates an empty column
%vector that you can fill with your two derivatives:
dydt(1) = -(p(1)+p(3))*y(1) + p(2)*y(2);
dydt(2) = p(1)*y(1) - p(2)*y(2);
end
生產資料用的ex_ode和call_ex
function dydt = ex_ode(t,y)
dydt(1) = -0.95*y(1) + 0.15*y(2);
dydt(2) = 0.6*y(1) - 0.15*y(2);
end
function [exp_y] = call_ex()
tspan = linspace(0.5,7.5,15);
y1_0 = 1;
y2_0 = 0;
[T,Y] = ode15s(@osc,tspan,[y1_0 y2_0]);
plot(T,Y(:,2),'o')
exp_y = Y(:,2);
end