MathGroup Archive 2007

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Re: Fitting parameters of nonlinear diff equation system

  • To: mathgroup at smc.vnet.net
  • Subject: [mg81567] Re: [mg81303] Re: Fitting parameters of nonlinear diff equation system
  • From: "Szabolcs Horvát" <szhorvat at gmail.com>
  • Date: Fri, 28 Sep 2007 02:05:42 -0400 (EDT)
  • References: <fcnl8h$s85$1@smc.vnet.net> <200709181002.GAA04068@smc.vnet.net>

On 26/09/2007, Paul Thomas <paul.pgt at gmail.com> wrote:
> Hi Szabolcs,
>
> A couple things, thanks for your  reply!
>
> Replacing f with f[t] in fun[] does seem to return a function. I tested this
> by just running fun with a given a,b,c and with the f[t] syntax it produces
> an interpolating function (as NDSolve should)--it doesn't with just f as I
> think it doesn't know what "f" is in that case--the output of just using "f"
> is literally just "f" with no function associated.

A few examples that demonstrate how NDSolve returns the result:

In[1]:= NDSolve[{f'[t]==f[t],f[0]==1},f,{t,0,1}]
Out[1]= {{f->InterpolatingFunction[{{0.,1.}},<>]}}

In[2]:= NDSolve[{f'[t]==f[t],f[0]==1},f[t],{t,0,1}]
Out[2]= {{f[t]->InterpolatingFunction[{{0.,1.}},<>][t]}}

In[3]:= f/.NDSolve[{f'[t]==f[t],f[0]==1},f,{t,0,1}]
Out[3]= {InterpolatingFunction[{{0.,1.}},<>]}

In[4]:= f[t]/.NDSolve[{f'[t]==f[t],f[0]==1},f[t],{t,0,1}]
Out[4]= {InterpolatingFunction[{{0.,1.}},<>][t]}

Using f[t] inside NDSolve gives the value of the
InpterpolatingFunction at the point 't'.

Using f (without the argument t) returns the function itself, not the
value of the function at a particular point.

> On the /@ priority issue, I actually changed one other thing, which maybe is
> a problem too.
> Instead of using "/@" there I used "/."--that I think does require the
> parentheses. I kept getting an error using "/@" and "/." seemed to fix it.

Here are two examples demonstrating /@ (Map) and /. (ReplaceAll).

Map applies f to each element in the list.

In[1]:= f/@{1,2,3}
Out[1]= {f[1],f[2],f[3]}

ReplaceAll can be used to replace parts of an expression with some
other expression. In this case we change 't' in f[t] to {1,2,3}.

In[2]:= f[t]/.{t->{1,2,3}}
Out[2]= f[{1,2,3}]

Both approaches may be used in the application we are discussing,
because InterpolatingFunctions behave as if they had the attribute
Listable, i.e. they accept a List as their argument, and "apply
themselves" to each element of the list.  (Actually I did not know
about this before, so thanks for the hint!)

In[3]:= SetAttributes[f,Listable]

In[4]:= f[{1,2,3}]
Out[4]= {f[1],f[2],f[3]}


> The error would say that essentially a non-numerical equation was resulting
> from the /@ and it would write it out--it did seem to be doing the wrong
> equation, but it looks like the right one with /.
> Now it seems to be doing exactly what you designed it to do.

Now I learned that it is a really bad idea to provide a not fully
tested, partial example.  It just leads to misunderstandings.  I am
sure that your approach is correct, just slightly different than mine
(e.g. you might have used variables such as tvals and fvals slightly
differently).

Here is a full example with the very simple differential equation
f'(t) = a*f(t):

We use the values of an exponential function for fitting:

In[1]:=
tvals = Range[0, 1, .1];
fvals = 2 Exp[3 tvals];

In[3]:= fun[a_,b_]:=First[f/.NDSolve[{f'[t]==a*f[t],f[0]==b},f,{t,0,1}]]

In[4]:= opt[a_?NumericQ,b_?NumericQ]:=Norm[fun[a,b]/@tvals-fvals]
(* The more compact fun[a,b][tvals] - fvals could be used instead! *)

In[5]:= FindMinimum[opt[a,b],{{a,1},{b,1}}]

During evaluation of In[5]:= FindMinimum::lstol: The line search
decreased the step size to within tolerance specified by AccuracyGoal
and PrecisionGoal but was unable to find a sufficient decrease in the
function.  You may need more than MachinePrecision digits of working
precision to meet these tolerances. >>

Out[5]= {0.0000218821,{a->3.,b->2.}}

> Once the parameters are fit I feed them back into "fun" and check the work
> essentially of the process and plot the optimization over time. For example,
> once I've done a particular fit:
>
> g[t_]:=fun[1,2,3] {pretending 1,2,3 were our best a,b,c fits}

This is where the misunderstanding came from:  In my example
fun[1,2,3] returned a function (e.g. f), in your example it returns
the value of the function at point 't' (e.g. f[t]).

> g[t] is now my best fit function and matches the findminimum minimum at the
> particular t points I have. The curves also look generally as they should.
>
> Do you think I have an error in this?

No, I think that everything is OK.  Here is the previous simplified
example converted to what I think you did:

In[1]:=
tvals=Range[0,1,.1];
fvals=2Exp[3tvals];

In[3]:= fun[a_,b_]:=First[f[t]/.NDSolve[{f'[t]==a*f[t],f[0]==b},f[t],{t,0,1}]]

In[4]:= opt[a_?NumericQ,b_?NumericQ]:=Norm[(fun[a,b]/.{t->tvals})-fvals]

In[5]:= FindMinimum[opt[a,b],{{a,1},{b,1}}]

During evaluation of In[5]:= FindMinimum::lstol: The line search
decreased the step size to within tolerance specified by AccuracyGoal
and PrecisionGoal but was unable to find a sufficient decrease in the
function.  You may need more than MachinePrecision digits of working
precision to meet these tolerances. >>

Out[5]= {0.0000218821,{a->3.,b->2.}}

> Thanks for all your help with it, like
> I said, it seems like a very powerful approach at least for the types of
> things that I've been trying to do--very flexible!

I'm glad you found it useful!

> Paul

Szabolcs


  • Prev by Date: Re: create a list with x,y,z coordinates
  • Next by Date: Re: CompleteCharacters Palette in v6?
  • Previous by thread: Re: Re: Fitting parameters of nonlinear diff equation system
  • Next by thread: Re: Fitting parameters of nonlinear diff equation system