problems with using genetic algorithm to fit kinetic parameters

2 次查看（过去 30 天）
Wendy
Wendy 2014-6-19
评论： Star Strider 2014-6-19
Hi all, I am learning how to fit kinetic parameters of a system of ODE using the genetic algorithm in the optimization toolbox. I set up a testing case but got the following error. Any advise on how to solve the problem is appreciated. I have also attached by testing code here. Thank you in advance.
Best, Wendy
----
Error using sir
Too many input arguments.
Error in odearguments (line 88)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode15s (line 149)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
Error in fitness (line 5)
[t, y] = ode15s('sir', [0 14], [50000, 1, 0], [], param);
Error in createAnonymousFcn>@(x)fcn(x,FcnArgs{:}) (line 11)
fcn_handle = @(x) fcn(x,FcnArgs{:});
Error in fcnvectorizer (line 14)
y(i,:) = feval(fun,(pop(i,:)));
Error in makeState (line 47)
Score =
fcnvectorizer(state.Population(initScoreProvided+1:end,:),FitnessFcn,1,options.SerialUserFcn);
Error in gaunc (line 41)
state = makeState(GenomeLength,FitnessFcn,Iterate,output.problemtype,options);
Error in ga (line 351)
[x,fval,exitFlag,output,population,scores] = gaunc(FitnessFcn,nvars,
...
Error in driver_fitness (line 5)
[val, fval] = ga(@fitness, 2);
Caused by:
Failure in user-supplied fitness function evaluation. GA cannot continue.
28 rethrow(userFcn_ME)
%%%%main file
clear all; close all; clc
global param
[val, fval] = ga(@fitness, 2);
% compare to the input data
param = val;
time_obs = [2 4 5 6 7 8 10 12 14]; %observation time
y_obs = [1 5 35 150 250 260 160 90 45];
[tx, yx] = ode15s('sir', [0 14], [50000, 1, 0],[],param);
yxx = interp1(tx, yx(:,2), time_obs);
figure
plot(tx, yx(:,2), '-b', time_obs, y_obs, 'rs')
%%%fitness.m - objective function
function fit = fitness(param)
% Get the model output for the given values of beta and delta
[t, y] = ode15s('sir', [0 14], [50000, 1, 0], [], param);
tint = [2 4 5 6 7 8 10 12 14]; % time point
yint = interp1(t, y(:,2), tint); % simulated data
yobs = [1 5 35 150 250 260 160 90 45]; % observed data
error = sum((yint - yobs).^2);
fit = error^-1;
%%%sir.m - ODE system
function dx = sir(t, x, param)
beta = param(1);
delta = param(2);
dx = [0; 0; 0];
dx(1) = -beta * x(1) * x(2);
dx(2) = beta * x(1) * x(2) - delta * x(2);
dx(3) = delta * x(2);

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

回答（1 个）

Star Strider
Star Strider 2014-6-19
I didn’t run your code, but I believe I see the problem.
Try changing your ode15s call to:
[t, y] = ode15s(@(t,x) sir(t,x,param), [0 14], [50000, 1, 0], [], param);
Your ‘param’ variable gets passed to ‘sir’ but ode15s doesn’t see the extra argument.
You can also specify the time vector argument to the ODE functions as your ‘tint’ vector. The ODE functions will evaluate your ODE at those points (or very close to them), so you can eliminate the interpolation step.
  3 个评论
显示 1更早的评论
Star Strider
Star Strider 2014-6-19
When in doubt, I start from scratch. I have nothing against Genetic Algorithms, but I have more experience fitting ODEs with curve-fitting functions. That said, I didn’t get a good fit with this:
tint = [2 4 5 6 7 8 10 12 14]; % time point
yobs = [1 5 35 150 250 260 160 90 45]; % observed data
B = lsqcurvefit(@sirobj, rand(2,1)*0.0001, tint', yobs')
[t, y] = ode15s(@(t,x) sir(t,x,B), [tint], [50000, 1, 0]);
figure(1)
plot(tint, yobs, 'pb')
hold on
plot(tint, y(:,2), '-r')
hold off
grid
% legend('Data', 'y_1', 'y_2', 'y_3', 'Location', 'NE')
legend('Data', 'y_2', 'Location', 'NE')
axis([0 15 0 300])
% --------------- SEPARATE FUNCTION FILE ---------------
function [ y2 ] = sirobj(param, t)
%UNTITLED2 Summary of this function goes here
% Detailed explanation goes here
beta = param(1);
delta = param(2);
sir = @(t,x,param) [-beta .* x(1) .* x(2); beta .* x(1) .* x(2) - delta .* x(2); delta .* x(2)];
[t, y] = ode15s(@(t,x) sir(t,x,param), t, [50000, 1, 0]);
y2 = y(:,2);
end
The values I got for the parameters are:
beta = 211.7341e-006
delta = 9.4255e+000
I don’t know what you’re doing, so I don’t know whether these are appropriate. The only thing I can say for it is that it does run.
Star Strider
Star Strider 2014-6-19
One problem I didn’t see last night is this line:
fit = error^-1;
The GA is going to search for the minimum of the fitness function, so you’re having it search for the maximum (in this instance, the maximum error) with this definition of your ‘fit’ variable.
You’ll likely have better luck defining ‘fit’ as ‘error’ instead.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by