Yingfu Xie
2004-Aug-11 14:38 UTC
[R] always NaN after some running in R, but all fine in S-plus
Hello, S-plus and R helpers,(sorry for cross-post) I wrote some simple C code for one likelihood to be optimized (using optim(MASS)). I use same function, same data, same starting points and same DLL in R and S-plus for comparison. (I compiled it with 'Rcmd SHLIB likelihood.c' and the header files of it include only R.h and math.h). While it works quite fine in S-plus, it forever returns NaN for the value of log-likelihood after several normal running in R. What is worse, it can never return a normal value once NaN appears. I just cannot understand why. Has anyone had similar experience? I searched R and S-plus archive for nothing. I attach more details below and really appreciate any help from you! R:1.9.0 S-plus:6.2 OS: Windows 2000 Best Regards, Yingfu The output in R are:> symbol.C("like")[1] "like"> likelihood(data=data03) # The likelihood function calling .C("like",...)[1] 5850.12> myoptim(data=data03) # Optimization using optim with fn=likelihoodInitial parameters are 1 10 0.2 0.2 0.2 0.1 0.1 0.1 0.1 0.2 0.3 0.5 1 10 0.2 0.2 0.2[1] 5850.154 #I print first 5 parameters and the value 1.001 10 0.2 0.2 0.2[1] 5850.308 #of the minus log-likelihood every run 0.999 10 0.2 0.2 0.2[1] 5849.997 ........... #similar outputs are omitted 1 10 0.2 0.2 0.2[1] 5848.578 0 0 0 0 0[1] NaN #It seems that here it jumps to the edge of 0.001 0 0 0 0[1] NaN #box constrain, all zeros in my case. ............. #And it returns NaN. From now on, it is 0 0 0 0 0[1] NaN # forever NaN, 0 0 0 0 0[1] NaN # even if it jumps back near to the initial 0.995734 9.95734 0.1991468 0.1991468 0.1991468[1] NaN #parameters. 0.996734 9.95734 0.1991468 0.1991468 0.1991468[1] NaN ............. 1 10 0.2 0.2 0.201[1] NaN 1 10 0.2 0.2 0.199[1] NaN $par [1] 1.0 10.0 0.2 0.2 0.2 0.1 0.1 0.1 0.1 0.2 0.3 0.5 $value # Since NaN is not allowed with method L-BFGS-B, I set it [1] 1e+05 # returns a bigger 10000 if the value is NaN or NA. $counts function gradient 7 7 $convergence [1] 0 What is more, it became worse now:> likelihood(data=data03) #It returns a NaN also!![1] NaN> myoptim(data=data03)Initial parameters are 1 10 0.2 0.2 0.2 0.1 0.1 0.1 0.1 0.2 0.3 0.5 1 10 0.2 0.2 0.2[1] NaN 1.001 10 0.2 0.2 0.2[1] NaN ................ 1 10 0.2 0.2 0.2[1] NaN $par [1] 1.0 10.0 0.2 0.2 0.2 0.1 0.1 0.1 0.1 0.2 0.3 0.5 $value [1] 1e+05 $counts function gradient 1 1 $convergence [1] 0 $message [1] "CONVERGENCE: NORM OF PROJECTED GRADIENT <= PGTOL" Now, the problem seems lie in my C code and algorithm. However, in S-plus, it works fine:> is.loaded("like")[1] T> myoptim(data=data03)Initial parameters are 1 10 0.2 0.2 0.2 0.2 0.1 0.1 0.1 0.1 0.3 0.5 [1] 1.0 10.0 0.2 0.2 0.2 [1] 5834.419 [1] 1.001 10.000 0.200 0.200 [1] 5834.579 .............. #it never jumps to zeros [1] 0.7701213 4.3442599 0.1424875 0.3072238 0.1274840 [1] 5106.236 $par: [1] 0.77012134 4.34425988 0.14248754 0.30722383 0.12748400 [6] 0.48420116 0.00000000 0.02095689 0.00000000 0.61156935 [11] 0.77179635 0.00000000 $value: [1] 5106.049 $counts: function gradient 21 21 $convergence: [1] 0 $message: [1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH" I tried the 'bad' parameters in S-plus:> likelihood(par=rep(0,12),data=data03) #It is NA.[1] NA> likelihood(data=data03) #It is normal again.[1] 5834.421 Thank you very much for your time! ########################################### This message has been scanned by F-Secure\ Anti-Virus for Mi...{{dropped}}