Derek, The problem is that your function is poorly scaled. You can see that the parameters vary over 10 orders of magnitude (from 1e-04 to 1e06). You can get good convergence once you properly scale your function. Here is how you do it:
par.scale <- c(1.e06, 1.e06, 1.e-06, 1.0) SPoptim <- optim(pars, SPsse, B=d$catch, CPE=d$cpe, control=list(maxit=1500, parscale=par.scale)) > SPoptim $par B0 K q r 7.329553e+05 1.160097e+06 1.484375e-04 4.050476e-01 $value [1] 1619.487 $counts function gradient 1401 NA $convergence [1] 0 $message NULL Hope this helps, Ravi. ____________________________________________________________________ Ravi Varadhan, Ph.D. Assistant Professor, Division of Geriatric Medicine and Gerontology School of Medicine Johns Hopkins University Ph. (410) 502-2619 email: rvarad...@jhmi.edu ----- Original Message ----- From: Derek Ogle <do...@northland.edu> Date: Saturday, June 26, 2010 4:28 pm Subject: [R] optim() not finding optimal values To: "R (r-help@R-project.org)" <r-help@r-project.org> > I am trying to use optim() to minimize a sum-of-squared deviations > function based upon four parameters. The basic function is defined as > ... > > SPsse <- function(par,B,CPE,SSE.only=TRUE) { > n <- length(B) # get number of years of > data > B0 <- par["B0"] # isolate B0 parameter > K <- par["K"] # isolate K parameter > q <- par["q"] # isolate q parameter > r <- par["r"] # isolate r parameter > predB <- numeric(n) > predB[1] <- B0 > for (i in 2:n) predB[i] <- predB[i-1]+r*predB[i-1]*(1-predB[i-1]/K)-B[i-1] > predCPE <- q*predB > sse <- sum((CPE-predCPE)^2) > if (SSE.only) sse > else list(sse=sse,predB=predB,predCPE=predCPE) > } > > My call to optim() looks like this > > # the data > d <- data.frame(catch= > c(90000,113300,155860,181128,198584,198395,139040,109969,71896,59314,62300,65343,76990,88606,118016,108250,108674), > > cpe=c(109.1,112.4,110.5,99.1,84.5,95.7,74.1,70.2,63.1,66.4,60.5,89.9,117.0,93.0,116.6,90.0,105.1)) > > pars <- c(800000,1000000,0.0001,0.17) # put all > parameters into one vector > names(pars) <- c("B0","K","q","r") # name the parameters > ( SPoptim <- optim(pars,SPsse,B=d$catch,CPE=d$cpe) ) # run optim() > > > This produces parameter estimates, however, that are not at the > minimum value of the SPsse function. For example, these parameter > estimates produce a smaller SPsse, > > parsbox <- c(732506,1160771,0.0001484,0.4049) > names(parsbox) <- c("B0","K","q","r") > ( res2 <- SPsse(parsbox,d$catch,d$cpe,SSE.only=FALSE) ) > > Setting the starting values near the parameters shown in parsbox even > resulted in a movement away from (to a larger SSE) those parameter values. > > ( SPoptim2 <- optim(parsbox,SPsse,B=d$catch,CPE=d$cpe) ) # run optim() > > > This "issue" most likely has to do with my lack of understanding of > optimization routines but I'm thinking that it may have to do with the > optimization method used, tolerance levels in the optim algorithm, or > the shape of the surface being minimized. > > Ultimately I was hoping to provide an alternative method to fisheries > biologists who use Excel's solver routine. > > If anyone can offer any help or insight into my problem here I would > be greatly appreciative. Thank you in advance. > > ______________________________________________ > R-help@r-project.org mailing list > > PLEASE do read the posting guide > and provide commented, minimal, self-contained, reproducible code. ______________________________________________ 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.