Hi,
I try to fit parameters based on the comparison between a stochastic
model that uses these parameters and observed data.
The problem is that two runs with the same set of parameters will
produce a fitting criteria slightly different.
I found that optim() is highly sensitive to such a situation and
produces completly wrong result.
I reproduce this situation below.
v is a set of 50 "observed" values
par is the set of parameters to be fitted
fitNorm returns the -Ln L of the observations based on a Gaussian
distribution
fitNorm_error returns the -Ln L of the observations based on a Gaussian
distribution + a small error
Example of optim() result:
with fitNorm
final value 215.918650
converged
$par
[1] 156.27980 18.16392 (compared to 157.5, 18.5)
with fitNorm_error
final value 271.661187
converged
$par
[1] 169.98441 88.81387 (compared to 157.5, 18.5)
So the question: is there another solution for a fitting procedure that
is more resistant to change in fit criteria ?
Thanks a lot,
Marc
#The likelihood function that returns -Ln L of the data
fitNorm<-function(para) {
return(-sum(dnorm(v, para[1], para[2], log=TRUE)))
}
#The likelihood function that returns -Ln L of the data with a slight error
fitNorm_error<-function(para) {
return(-sum(dnorm(v, para[1], para[2], log=TRUE))+rnorm(1, 0, 0.1))
}
#Generate 50 random number from a Gaussian distribution
v<-rnorm(50, 157.5, 18.5)
par<-c(100, 10)
optim(par, fitNorm, method="BFGS", hessian=FALSE,
control=list(maxit=500, trace=1))
optim(par, fitNorm_error, method="BFGS", hessian=FALSE,
control=list(maxit=500, trace=1))
--
__________________________________________________________
Marc Girondot, Pr
Laboratoire Ecologie, Systématique et Evolution
Equipe de Conservation des Populations et des Communautés
CNRS, AgroParisTech et Université Paris-Sud 11 , UMR 8079
Bâtiment 362
91405 Orsay Cedex, France
Tel: 33 1 (0)1.69.15.72.30 Fax: 33 1 (0)1.69.15.73.53
e-mail: marc.giron...@u-psud.fr
Web: http://www.ese.u-psud.fr/epc/conservation/Marc.html
Skype: girondot
______________________________________________
R-help@r-project.org 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.
