Re: Data fitting from coupled differential equations
- To: mathgroup at smc.vnet.net
- Subject: [mg93402] Re: Data fitting from coupled differential equations
- From: margeorias at gmail.com
- Date: Thu, 6 Nov 2008 04:06:18 -0500 (EST)
- References: <gemjhd$4s9$1@smc.vnet.net> <gepb3g$si9$1@smc.vnet.net>
On 4 Nov, 06:19, dh <d... at metrohm.ch> wrote:
> Hi,
>
> draw your example solution. Is an exponential growning example really
>
> what you want? Your measured data will probably not be of the order of
>
> 10^50.
>
> hope this helps, Daniel
>
>
>
> margeor... at gmail.com wrote:
> > Hi,
>
> > I'm working on a predator-prey model for 2-species:
> > n1'[t] == a1 n1[t] - b1 n1[t]^2 + c1 n1[t] n2[t]
> > n2'[t] == a2 n2[t] - b2 n2[t]^2 + c2 n1[t] n2[t]
> > n1[0] == n10,
> > n2[0] == n20
>
> > I have real data for (t,n1[t],n2[t]), i.e. data={t,n1[t],n2[t]} and my
> > objective is to estimate the coefficients a1,a2,b1,b2,c1,c2.
>
> > To do this, I generate synthetic data by solving the above equations
> > using {n10,n20,a1,a2,b1,b2,c1,c2}={0.01,0.001,2,3,1,1,1,2} and sample
> > the solutions,i.e.
> > data = Table[{t, n1[t], n2[t]} /. sol, {t, 0, 20, 0.25}];
>
> > I then set-up the function below:
>
> > sse[a1_?NumericQ,a2_?NumericQ,b1_?NumericQ,b2_?NumericQ,c1_?
> > NumericQ,c2_?NumericQ]:=
> > Block[{sol,n1,n2},
> > sol=NDSolve[{
> > n1'[t]==a1 n1[t]-b1 n1[t]-c1 n1[t] n2[t],
> > n2'[t]==a2 n2[t]-b2 n2[t]-c2 n1[t] n2[t],
> > n1[0]==0.01,n2[0]==0.001},
> > {n1,n2},{t,0,20}][[1]];
> > Plus@@Apply[((n1[#1]-#2)^2+(n2[#1]-#3)^2)&,data,{1}]/.sol]
>
> > and use FindMinimum to retrieve the coefficients:
> > FindMinimum[sse[a1,a2,b1,b2,c1,c2],{{a1,1.},{a2,1.},{b1,1.},{b2,1.},
> > {c1,1.},{c2,1.}},AccuracyGoal->20,PrecisionGoal->18,WorkingPrecision-
> >> 40]
>
> > My problem however is that the resulting coefficients are not the same
> > as the original one, ie:
> > {a1->-9.49992942438124754467310,a2->-9.50000018840483650972573,b1-
> >> -9.50007117410767898402213,b2->-9.50000041008883222559689,c1-
> >> -9.50000036535088621114653,c2->-9.50000030036113407572174}
>
> > I can improve the results by making starting points very close to the
> > actual value of the coefficients i.e. {a1,2}. Well this somehow
> > defeats the purpose and in reality, I don't know whats the actual
> > value of the coefficients - they range from -Infinity to Infinity.How
> > do I go about obtaining the correct coefficients? Is there something
> > wrong with the function above?
>
> > Another thing, wouldn't NMinimize be appropriate for this as well
> > since I want to find the global minimum of sse? I have tried it but it
> > gave the wrong values for the coefficients as well :(
>
> > I'm trying to make FindFit work on the coupled differential equations
> > but no luck yet. Any guidance on this?
>
> > Appreciate any feedback you have.
>
> > Thanks again.
>
> --
>
> Daniel Huber
>
> Metrohm Ltd.
>
> Oberdorfstr. 68
>
> CH-9100 Herisau
>
> Tel. +41 71 353 8585, Fax +41 71 353 8907
>
> E-Mail:<mailto:d... at metrohm.com>
>
> Internet:<http://www.metrohm.com>
Thanks for the feedback Daniel,
Sorry I made a mistake in the writing the differential equations. (I
"cleaned up" the equations after copying from Mathematica). It should
be:
n1'[t] == a1 n1[t] - b1 n1[t]^2 - c1 n1[t] n2[t]
n2'[t] == a2 n2[t] - b2 n2[t]^2 - c2 n1[t] n2[t]
n1[0] == n10,
n2[0] == n20
What are the differences of using FindMinimum and NMinimize with
regards to extracting coefficients from differential equations model?
Wouldn't NMinimize be more appropriate since it is global in nature?
If I were to increase the number of species from the current 2 with
two differential equations to say 100 coupled differential equations,
what are the things I need to worry about? How would it affect the
computed coefficients extracted from data?
Thanks again for the insights.