[R-sig-dyn-mod] Model giving negative values

[Thomas Petzoldt]

Dear Caroline,

first of all: the mailing list address is: 
r-sig-dynamic-models using r-project.org whereas the -owner using ... - address 
writes only to the moderators.

Thank you for using R for your problems and for writing to the forum. 
Which package have you used? In the meantime, several ODE solver 
packages are available in R. Non-negativity is indeed a common request 
and there are different ways how software can handle this. Some other 
systems have an option to force non-negative solutions, but we (i.e. the 
deSolve developers) decided not to offer such an option. It can be 
possible, that another solver package offers such an option.

The reason why deSolve has none, is that negative values can have 
different reasons. One reason is indeed numerical error of the solver 
that may be handled with a non-negativity constraint. However, the in 
our experience most frequent reason is an incomplete model 
specification. The typical case is an export from an empty pool.

Example: in the following equation, X can become negative by model design:

dX/dt = -k

whereas the following it cannot:

dX/dt = -k * X

If you need a linear decrease "most of the time" that reduces to zero 
shortly before X becomes zero, you can use a Monod-type safeguard (or 
something similar):

dX/dt = -k*X / (k2 + X)

The selection of k2 is relatively uncritical. It should be small enough 
not to influence the overall behavior and not too small, compared to the 
numerical accuracy of the solver.

Another method to avoid negative values is to work in log-transformed 
space.

Here are some related threads:

https://stat.ethz.ch/pipermail/r-sig-dynamic-models/2010q2/000028.html

https://stat.ethz.ch/pipermail/r-sig-dynamic-models/2013q3/000222.html

https://stat.ethz.ch/pipermail/r-sig-dynamic-models/2016/000437.html

https://stackoverflow.com/questions/56614347/non-negative-ode-solutions-with-functools-in-r/56692927#56692927

If your model is small enough or (even better and highly recommended) if 
you can reduce it to a "minimal reproducible example", you can include 
it in your mail. If it is a little bit bigger, you may consider a public 
web page or a link to your cloud store (e.g. dropbox). But, if it is a 
really complex application, you may not get helpful answers because only 
very few people will have enough time and motivation for this.

Hope it helps

Thomas

On 03.07.2019 at 20:34 Caroline Mburu wrote:
> Hello,
>
> I have a compartmental model with 9 compartments.I coded it in r and 
> used runsteady to solve the ODEs. The problem is I am getting very 
> small values and when I try estimating parameters the model is not 
> converging.
>
> Where should I share my code?
>
> I will be very grateful for your assistance.
>
> Regards,
>
> Caroline


-- 
Dr. Thomas Petzoldt
wiss. Mitarbeiter                     member of scientific staff

Technische Universität Dresden        Technische Universitaet Dresden
Fakultät Umweltwissenschaften         Faculty of Environmental Sciences
Institut für Hydrobiologie            Institute of Hydrobiology
Professur für Limnologie              Chair of Limnology
01062 Dresden                         01062 Dresden, Germany

Tel.: +49 351 463 34954
Fax:  +49 351 463 37108
E-Mail: thomas.petzoldt using tu-dresden.de
http://tu-dresden.de/Members/thomas.petzoldt


	[[alternative HTML version deleted]]