I am currently in the final phase of developing my first R package, which should conform to the multidimensional processing (MPT) models (see its homepage for the current version). Model fitting is achieved by the R optim function.
Today I first spoke with him on a Windows 7 machine and noticed something really strange: optim does not succeed when using the 64-bit version of R. This looks like an error for me (especially how nlminb converges for both versions of R). Since optim is at the core of my package, any help on this is greatly appreciated.
Here is a minimally reproducible example (usually the model is specified using expressions and not specified in the objective function, but for simplicity I put everything in the objective function here):
# The objective function: llk.tree <- function (Q, data) { p <- Q[1] q <- Q[2] r <- Q[3] e <- rep(NA,4) e[1] <- p * q * r e[2] <- p * q * (1-r) e[3] <- p * (1-q) * r e[4] <- p * (1-q) * (1-r) + (1-p) llk <- sum(data * log(e)) if (is.na(llk)) llk <- -1e+19 if (llk == -Inf) llk <- -1e+19 return(-llk) }
This example reproduces the MPT seed example from Riefer and Batchelder , namely row 1 in table 1 s. 327 (the expected parameter values ββwill be p = 1, q = .49 and r = .30).
Running it on a 32-bit R always gives the correct result (verified with versions 2.12.2 and 2.13.0):
$par [1] 1.0000000 0.4944449 0.3000001 $value [1] 234.7110 $counts function gradient 11 11 $convergence [1] 0 $message [1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"
(Please note that the quantity may vary due to random starting values.)
Running it on 64 bit R, on the other hand, can lead to this (wrong) result:
$par [1] 0.8668081 0.6326655 0.1433857 $value [1] 257.7328 $counts function gradient 3 3 $convergence [1] 0 $message [1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"
The return value of the objective function and the return values ββof the parameters are different for each run, but the score is always 3!
Note that running nlminb gives the correct results for 32-bit and 64-bit R:
nlminb(runif(3, 0.1, 0.9), llk.tree, data = c(24, 65, 30, 61), lower = 0, upper = 1) $par [1] 1.0000000 0.4944445 0.3000000 $objective [1] 234.711 $convergence [1] 0 $iterations [1] 14 $evaluations function gradient 19 55 $message [1] "relative convergence (4)"
Finally, we have examples (this is our simple sample model) that worked on 64-bit R and optim , but more examples (for example, shown here) did not work.
And the score is always 3 ...
EDIT:
When fixing the initial values ββ(thanks to Joshua Ulrich) optim does not optim from these fixed values ββto 64 bits of R:
optim(c(0.5, 0.5, 0.5), llk.tree, data = c(24, 65, 30, 61), method = "L-BFGS-B", lower = rep(0, 3), upper = rep(1, 3)) $par [1] 0.5 0.5 0.5 $value [1] 276.1238 $counts function gradient 3 3 $convergence [1] 0 $message [1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"