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

**References**:**Re: Fitting parameters of nonlinear diff equation system***From:*Szabolcs Horvát <szhorvat@gmail.com>