[R] problems with optim, "for"-loops and machine precision
Setzer.Woodrow at epamail.epa.gov
Setzer.Woodrow at epamail.epa.gov
Wed Jan 10 17:57:25 CET 2007
Without more detail - a reproducible example - it is hard to give you
I wonder if the functions NLL23 and NLL21 depend on numerical solutions
of a system of ODEs, since you invoke the odesolve package? If so, try
switching to the Nelder-Mead optimizer, enforcing the parameter
constraints using transformation. Probably you are using the finite
difference derivatives calculated internally to optim for the gradient
information used in the L-BFGS-B optimizer. These can be unstable when
based on numerical solutions of odes, causing the optimizer to fail, or
sometimes to converge to a non-optimal point.
Some other points:
- you cannot change machine precision by changing values in .Machine.
To change the number of digits printed, use options(digits=8).
- use 'library()' instead of 'require()', unless you need to use the
return value from 'require()'
R. Woodrow Setzer, Ph. D.
National Center for Computational Toxicology
US Environmental Protection Agency
Mail Drop B205-01/US EPA/RTP, NC 27711
Ph: (919) 541-0128 Fax: (919) 541-1194
<s.ruegg at access.
Sent by: r-help at stat.math.ethz.ch
r-help-bounces at s cc
[R] problems with optim,
01/10/2007 07:18 "for"-loops and machine precision
Dear R experts,
I have been encountering problems with the "optim" routine using "for"
loops. I am determining the optimal parameters of several nested models
minimizing the negative Log-Likelihood (NLL) of a dataset.
The aim is to find the model which best describes the data. To this end,
am simulating artificial data sets based on the model with the least
of parameters (6) and the parameters determined with the field data. For
each artificial set I estimate the parameters of the model with 6
and the next more complex model with 7 parameters (two of these
are equal in the 6-parameter model) by minimizing the corresponding NLL
"optim". In theory the 7-parameter model should fit the data either
or better than the 6-parameter model. Therefore the difference of the
minimal NLLs should be 0 or larger.
For 500 data sets I use the following code:
for (s in 1:500)
nv=MyEnv) #reading a data set
For some strange reason the results do not correspond to what I expect:
about 10% of the solutions have a difference of NLL smaller than 0. I
verified the optimisation of these results manually and found that a
NLL was ignored and a higher NLL was returned at "convergence". To check
what was happening I inserted a printing line in the NLL function to
all parameters and the NLL as the procedure goes on. To my surprise
then stopped at the minimal NLL which had been ignored before. I have
reduced the machine precision to .Machine$double.digits=8 thinking, that
printing was slowing down the procedure and by reducing the machine
precision to speed up the calculations. For an individual calculation
solved my problem. However if I implemented the same procedure in the
above, the same impossible results occurred again.
Can anyone tell me where I should be looking for the problem, or what it
and how I could solve it?
Thanks a lot for your help
Dr.med.vet., PhD candidate
Institute for Parasitology
phone: +41 44 635 85 93
fax: +41 44 635 89 07
e-mail: s.ruegg at access.unizh.ch
[[alternative HTML version deleted]]
R-help at stat.math.ethz.ch mailing list
PLEASE do read the posting guide
and provide commented, minimal, self-contained, reproducible code.
More information about the R-help