From: Spencer Graves <spencer.graves_at_pdf.com>

Date: Thu 10 Aug 2006 - 18:40:59 EST

}

R-help@stat.math.ethz.ch mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Thu Aug 10 19:55:28 2006

Date: Thu 10 Aug 2006 - 18:40:59 EST

<see in line>

Leslie Chavez wrote:

*> Hi,
**>
*

> I have a system of ODE's I can solve with lsoda.

*>
**> Model=function(t,x,parms)
**> {
**> #parameter definitions
**> lambda=parms[1]; beta=parms[2];
**> d = parms[3]; delta = parms[4];
**> p=parms[5]; c=parms[6]
**>
**> xdot[1] = lambda - (d*x[1])- (beta*x[3]*x[1])
**> xdot[2] = (beta*x[3]*x[1]) - (delta*x[2])
**> xdot[3] = (p*x[2]) - (c*x[3])
**>
**> return(list(xdot))
**> }
**>
**> I want to fit the output out[,4] to experimental data that is only
**> available on days 0, 7, 12, 14, 17, and 20. I don't know how to set up
**> optim or nls so that it takes out[,4] on the appropriate day, but still
**> runs lsoda on a time scale of 0.01 day.
**>
**> Below is the function I've been using to run 'optim', at the
**> course-grained time scale:
**>
*

SG: What about the following:

Modelfit=function(s) {

parms[1:4]=s[1:4]; times=c(0,7,12,14,17,20,25) lsodaTimes <- seq(min(times),max(times), by=0.01) out=lsoda(x0,lsodaTimes,Model,parms) obsTimes <- (100*times-1) mse=mean((log10(out[obsTimes,4])-log10(i(times)))^2) # cat(times) return(mse)

}

Your example is not self contained, so obviously I haven't tried this with it. However, something of this nature should work fine, I believe. Something similar but different should also work, I believe, with 'nls'; this would give you access to many helper functions (see "methods(class='nls')"). If 'nls' bombed on me, I'd then try 'optim' as it is less brittle. Then I might use the output of 'optim' as initial values for 'nls' to get confidence intervals etc.

hope this helps. Spencer Graves

> #x0=c(T0,I0,V0)

*> x0=c(2249,0,1)
**> #parms(lambda, beta, d, delta, p, c)
**> parms[5:6]=c(1.0,23)
**>
**> s0=c(49994,8456,6.16E-8,0.012) #initial values
**>
**> fit=optim(s0,Modelfit)
**>
**> Right now, lsoda is being run on too course-grained a time scale in the
**> function Modelfit. Most examples of optim and nls I have found compare
**> two data sets at the same times, and run lsoda on the time scale the
**> data is available at, but I would like to run lsoda at a finer scale, and
**> only compare the appropriate time points with the experiment. I have also
**> tried using nls, but I have the same problem. Does anyone have
**> suggestions?
**>
**> Thank you very much,
**>
**> Leslie
**>
**> ______________________________________________
**> R-help@stat.math.ethz.ch mailing list
**> https://stat.ethz.ch/mailman/listinfo/r-help
**> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
**> and provide commented, minimal, self-contained, reproducible code.
*

R-help@stat.math.ethz.ch mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Thu Aug 10 19:55:28 2006

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.1.8, at Thu 10 Aug 2006 - 22:22:31 EST.

*
Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help.
Please read the posting
guide before posting to the list.
*