Matlab Code Assistance for Multistart or GA

I want to apply MultiStart global optimization technique for curve fitting and parameter estimation from data set and system of ordinary differential equation. I am new for that concept. I have Number of infected cases for time serious data and want to fit with my model. I am wondering if you could assist me how to estimate using that technique? I will be very happy if you send me any prototype matlabcode!

Any help is appreciated!

Answer 1

[Nihal Reddy]

I understand you require resources and reference examples for implementing MultiStart optimization technique in MATLAB.

Please refer to the following MATLAB documentation for more information and examples-

• How GlobalSearch and MultiStart Work - MATLAB & Simulink (mathworks.com)
• Find multiple local minima - MATLAB (mathworks.com)
• Using MultiStart for Optimization Problems - Video - MATLAB (mathworks.com)

Hope this helps!

[Temesgen]

Dear Nihal Readdy, I will refer these marials. Thank you very much!!

Answer 2

[Alan Weiss]

Perhaps more relevant documentation is here:

MultiStart with lsqnonlin, Problem-Based

MultiStart Using lsqcurvefit or lsqnonlin

Fit ODE Parameters Using Optimization Variables

Fit an Ordinary Differential Equation (ODE)

Good luck,

Alan Weiss

MATLAB mathematical toolbox documentation

[Temesgen]

Dear Alan Weiss, thank you for sharing this important information!!

[Temesgen]

Dear Alan, I have disease dynamics model with 10 compartments and 21 parameters. From these the inital values of the state variables(10) and 8 of model parameters have known values and want to estimate the remaining 13 parameters using MultiStart technique. Using the materials you sent me about MultiStart, I wrote Matlab code. But unable to run. could you please help me? the code is attached

Thank you

Temesgen

[Alan Weiss]

For more help, you need to give more information. I do not have your data so cannot try running your code. And I do not know what you see as output when you run your code, so have no idea what you see as a problem.

But I am suspicious about several lines:

err = norm(qdata - Ihh(tmeasure))
options = optimset('Display','iter','PlotFcns',@optimplotfval);
[x,fval,exitflag,output] = fmincon(@(par)ExpData(abs(par)),p0,[],[],[],[],[0 0],[],[],options) %Initiation of fminsearch

1. You use fmincon as the minimizer, and the norm of the difference as the objective function. I think that this is a mistake. You should use lsqnonlin as the minimizer and the vector of differences as the objective function.
2. Your lower bound is [0 0]. That implies that you have just two variables. Or maybe it means that only your first two variables have lower bounds. But it looks wrong to me.

Good luck,

Alan Weiss

MATLAB mathematical toolbox documentation

[Temesgen]

Sorry for not attaching my data, but here is. Sir the lower boud is expected to be 13 not 2 you are right. Please take a look with my data sir? I used also GA but it took a long time to finish more than 6 hours and I terminate the program.Here is also the GA code dear.

Thank you for helping the community!

[Alan Weiss]

I tried looking at your code but I got confused at what you are trying to do.

The idea should be as follows.

1. Read the data and get a bunch of known, unchanging parameters.
2. Put all the parameters that you want to fit into a vector x. How many components are in this vector? That should be exactly the same as the number of components in the vector x0. Also in the vector lb.
3. Write an ODE function. This needs to have the syntax dydt = odefun(t,y). Include the data in this function. See Passing Extra Parameters.
4. Write an ODE solver function. This is something like function fout = myfun(@odefun,tspan,y0), where myfun calls ode45 and where tspan is the set of times for your ODE, and y0 has the same number of elements as x0. The output fout is a matrix with one row per time step and one column for each component of the ODE. You might need to pass data into this function as well.
5. After you have debugged the ODE, so that it gives you the expected output for given inputs, you can try to optimize the ODE. Really, you need to take things one step at a time, and make sure that you know how to pass the data around. Get the ODE solver working before you try to optimize anything.

Good luck,

Alan Weiss

MATLAB mathematical toolbox documentation

[Temesgen]

Dear Alan, sorry for sending such a wrong code and confusing you, it is due to luck of exposure.

To make it more clear, I have a systems of differential equation with 10 equations (10 state variables with known initial values). The model have 21 parameters and from these 8 of them have fixed values. I want to estimate the remaining 13 parameters. These 13 parameters have inital gues. I have data in excel with time in the first column and infected cases in the second column.

My objective is to minimize the square of the total error which is the difference between the data in the second column with the solution of infected cases at fourth system and estimate model parameters by finding the global solution. I choose use MultiStart method.

Hence, I tryied to write the code as attached herewith. I try to put titles for parts of the code to be more clear. Could you please try to correct and send me the correct code my dear?

[Alan Weiss]

Here is what I would like for you to do.

Create an ODE solver that takes the parameter x in your problem, meaning a vector of 13 components that you are trying to fit, and gives a solution to the ODE. In other words, I am asking you to create a function that takes x and returns an ODE solution.

In order to do that you need to understand what an ODE solution is for your problem. How many components are in the solution? How many time steps does it represent? What is the size of the solution, 12 by 10 or what?

I am NOT asking you to optimize anything. I can show you how to do that next. What I want is a clear, unambiguous, and working piece of code that takes a set of parameters x and gives a matrix representing the ODE solution at a set of times.

If you need help for this step, that is OK, we can discuss the details as you work on it. But I want to see you do that first before we talk about optimization.

Alan Weiss

MATLAB mathematical toolbox documentation

[Temesgen]

Thank you for your imidiate reply sir. I have attached the code!

[Alan Weiss]

See if you like how this looks. I did not use MultiStart; before we do that you have to see how you like the plain lsqnonlin solution.

Alan Weiss

MATLAB mathematical toolbox documentation

[Temesgen]

Dear Alan, really I appreciate your effort, thank you so much. here are some thing I want to say.

1. I want to estimate the following (13) model parameters. Keep in mind, there are 8 fixed model paraeters in addition to another 10 state variables. I gave these initial values for simulation purpose.

x(1) = 0.000548 ; % x(1)=Lambdaa

x(2) = 0.004626; % x(2)= phia

x(3) = 0.00227; % x(3)=gammah

x(4) = 0.295; % x(4)=gammaa

x(5) = 0.0002857; % x(5)=mua

x(6) = 0.0014; % x(6)= deltaa

x(7) = 0.001; % x(7)=q

x(8) =0.003189 ; % x(8)=betah

x(9) = 0.000178; % x(9)=betaaa

x(10) = 0.09; % x(10)=betahe

x(11) = 0.0009; % x(11)=betaah

x(12) = 0.175;% x(12)=alphah

x(13) = 0.285; % x(13)=alphaa

2. the estimated values as indicated in ''yout'' shows for the ten state variables (y(1), y(2), ......y(10)). Am I right? But, the program should use those initial values and estimate those 13 parameters and since most of the values will be in between 0 and 1 since ther are rates except some parametes like x(1).

3. The lower and upper bounds should be for these 13 parameters not for the 10 state variables.

you can use these bounds for those parameters

%[Lama betah betaah betahe betaaa gammah gammaa alphah alphaa mua deltaa q phia];

lb=[1 0.14 0.01 0.01 0.01 0 0.01 0.01 0 0.5 0.5 0.01 0.1];

ub=[1000 0.428 0.95 0.95 0.95 0.95 0.75 0.95 0.5 0.85 0.85 0.75 0.9];

Dear, I am very greatful and I can't wait to see the result after you consider and incorporate these ideas into the program.

[Alan Weiss]

OK, I clearly misunderstood what are the parameters to vary in your problem. And I am not sure that I understood correctly what you are trying to match, meaning which data are you trying to fit and which output of the ODE you want to match.

Once more I ask you to perform the following.

1. Give a function F(x) that is the ODE solver, taking your parameters x to a solution y(t). Obviously, I am not sure what that function is in the code that you gave me.
2. After that we can try to identify the data that you want F(x) to fit. We will call an optimizer and minimize the sum of squares of some parts of y(t) minus some data d(t).

So I am asking you to give F(x). Not to optimize it, just give clear, unambiguous, working code that takes an input set of parameters x and returns an ODE solution y(t).

Alan Weiss

MATLAB mathematical toolbox documentation