[R] Using odesolve to produce non-negative solutions
dave fournier
otter at otter-rsch.com
Tue Jul 3 19:34:10 CEST 2007
If you didn't get this solved.
I have done parameter estimation with models
defined by ODE's where negative solutions are a problem
and one can only avoid them with great difficulty if the
standard explicit methods for solving the ODE are used. I found that
using implicit methods could be a great help.
For example in the exponential case
dN/dt = -k*N
the simple finite difference approximation is
N_{t+1}-N+t
-------------- = -k*N_t , k>=0
h
or
N_{t+1} = N_t -k*h*N_t
and if k*h gets too large N_{t+1} goes negative and you are in trouble.
Consider instead the implicit formulation where the second
N_t on the RHS is replaced by N_{t+1} and one gets
N_{t+1} = N_t/(1+k*h)
which is correct for k*h=0 and as k*h--> infinity
For a more complicated example see
http://otter-rsch.com/admodel/cc4.html
for something I called "semi-implicit".
I hope these ideas will be useful for your problem.
Cheers,
Dave
--
David A. Fournier
P.O. Box 2040,
Sidney, B.C. V8l 3S3
Canada
Phone/FAX 250-655-3364
http://otter-rsch.com
More information about the R-help
mailing list