How to apply Genetic Algorithm solver to do parameter estimation for ODE based model?

I have a ODE based model and a set of experimental data for my project. There are seven parameters (Constant) inside the model need to do parameter estimation through Genetic Algorithm solver. My fitness function/Objective function will be sum of quadratic difference between simulated and experimental result.

pm(..) is the parameter I need to do estimation to best fit my curve. The function cgtase_cont is the mathematical model used. I tried to use GA solver to run it but it always appeared "Optimization running. Error running optimization. Too many output arguments.". Function square_error_ga is my objective function.

Using GA to solve this problem is the requirement of the project. Thanks in advance!

----------------------------------------------------------------------

Function cgtase_cont code:

function dx = cgtase_cont(t,x,ux,Ks,K1,H,KI,Kox,Ag,Ea,R,Ad,Ed,dH,dS,dCp,Temp0,rol,N,Di,V,Vs,...
  K2,Kd,Ysx,mx,Clx,qO2,YH,Ah,Cp,U,uxp,I,Ki,K,D,Tempj)

dx = zeros(6,1);

 %state variables
u = (ux*x(2))/((Ks+x(2))*(1+K1/H+H/K2)*(x(2)+KI+(x(2)^2)/Ks)*(Kox+1.4e-3))*(Ag*exp(-Ea/(R*x(4)))-Ad*exp(-Ed/(R*x(4))));
dGTS = dH-x(4)*dS+dCp*(x(4)-Temp0-x(4)*log(x(4)/Temp0));
Pg = 0.4*6*rol*N^3*Di^5;
Kla = 0.0195*((Pg/V)^0.55)*((Vs)^0.64)*(1+2.12*x(1)+0.2*x(1)^2)^-0.25;
Kn = 111.3*exp(-((x(4)-313.3)/7.42)^2);

%models
dx(1) = (u-Kd)*x(1);
dx(2) = (u/Ysx+mx)*(-x(1));
dx(3) = (Kla*(Clx-x(3))-qO2*x(1));
dx(4) = ((u*x(1))/(YH*rol*Cp)-(U*Ah)/(V*rol*Cp)*(x(4)-Tempj));               
dx(5) = ((uxp*I)/(I+Ki)*u*x(1)-K*x(5));
dx(6) = (D*exp(-dGTS/(R*x(4))))*x(5)-Kn*x(6);

------------------------------------------------------------

Function square_error_ga code:

function error = squared_error_ga(data_time, data_points, fn_time, fn_points)

    % Intepolate the fit to match the time points of the data
    fn_values = interp1(fn_time, fn_points, data_time);

    % Square the error values and sum them up
    error = sum((fn_values - data_points).^2);

--------------------------------------------------------

Function ga_error code:

function error = ga_error(pm)

 %Estimated Parameter set
Ks = pm(1); %half saturation constant (g/l)  pm(1)
Ki = pm(2); %saturation constant for inducer (g/l) pm(2)
Ag=10^pm(3);    %Ag = 10^7.9; %Arrhenius constant for growth pm(3)
Ad=10^pm(4);    %Ad = 10^10; %Arrhenius constant for death pm(4)
KI= pm(5); %inhibition constant due to excessive subtrate (g/l)pm(5)
mx = pm(6); %cell maintenance coefficient on substrate pm(6)
Ysx = pm(7); %yield coefficient (g/g)pm(7)
D = pm(8); %preexponential factor (per h) pm(8)
Kox = pm(9); %oxygen limitation constant (g/l) pm(9)

%Actual Constant
%Ks = 7.5; %half saturation constant (g/l)
%Ki = 0.001; %saturation constant for inducer (g/l)
%Ag = 10^7.9; %Arrhenius constant for growth
%Ad = 10^10; %Arrhenius constant for death
%KI= 0.001; %inhibition constant due to excessive subtrate (g/l)
%mx = 0.003; %cell maintenance coefficient
%Ysx = 0.0105; %yield coefficient (g/g)
%D = 1.65*10^14; %preexponential factor (per h)
%Kox = 7.104*10^-7; %oxygen limitation constant (g/l

ux = 0.1027; %maximum specific growth rate (per h)
K1 = 7.05*10^-10; %constant (Molar)
K2 = 6.86*10^-8; %constant (Molar)
H = 10^(-7.1); %ion hydrogen concentration (M)
Ea = 10000; %activation energy for growth (cal/mol)
Ed = 80000; %activation energy for death (cal/mol)
R = 1.987; %gas constant (cal/mol K)
uxp = 0.5872; %maximum specific production rate (per h)
I = 0.162; %inducer concentration (g/l)
K = 0.0625; %protein degradation rate
dH = 9.4*10^3; %change in enthalpy (cal/mol)
dS = -0.01*10^3; %change in entropy (cal/mol)
dCp = -0.32*10^3; %change in the heat capacity of the protein between fold and unfold state
Temp0 = 295; %reference temperature (K)
rol = 998.23; %density medium (g/l)
Di = 46*10^-3; %impeller diameter (m), Rushton
V = 1.5; %volume culture (l)
Vs = 1.1372; %oxygen superficial velocity (m/s)
Clx = 7*10^-3; %saturated DO concentration (g/l)
qO2 = 384*10^-3; %specific uptake rate of oxygen (g/g h)
U = 3.6482*10^4; %overall heat transfer coefficient (cal/h m2 K)
Ah = 0.606; %contacted surface area (m2)
YH = 0.42*10^-3; %1/YH is metabolic heat evolved per g biomass (g/cal)
Cp = 1; %heat capacity of medium (cal/ g K)
N = 10/3; %agitation speed (rev/s)
Kd = 0.02; %specific death rate (per h)
Tempj=293;

%-------------------------------------------

ode45(@(t,x)cgtase_cont(t,x,ux,Ks,K1,H,KI,Kox,Ag,Ea,R,Ad,Ed,dH,dS,dCp,Temp0,rol,N,Di,V,Vs,...
  K2,Kd,Ysx,mx,Clx,qO2,YH,Ah,Cp,U,uxp,I,Ki,K,D,Tempj), ...
  [0 60], [0.12 50 1.4e-3 310 0 0]);

error = squared_error_ga([0 4 8 12 16 20 24 28 32 36 40 44 48] , ...
                     [0 0 0.7586 1.2944 4.9136 6.9158 10.3277...
                      10.0373 10.8607 12.9424 12.6819 11.7554 10.6795], ...
                     t, ...
                     y(:,6));

------------------------------------------------------

The graph of my project is something like:

%

[Joker]

Hello, Bendo! have you address this problem. i have the same question.

[Jan]

@Joker: Really exactly the same question? In the original question, the output of ode45 was not caught, so I cannot guess, why the code runs at all. As far as I can see, the GA-part has not been posted and the error message was shown partially only. So please be so kind and open a new thread providing the required details to help you.