optim() and resultant hessian
On Fri, 16 Feb 2007, Brian Healey wrote:
R users; A question about optimization within R. I've been using both optim() and nlminb() to estimate parameters and all seems to be working fine. For context (but without getting into specifics - sorry), I'm working with a problem that is known to have correlated parameters, and parameter estimation can be difficult. I have a question on optim() - I'm using method="L-BFGS-B" to accommodate box constraints. For my dataset, I obtain parameter estimates using a few iterations of optim() - iterations in that I'm simply taking the results from a previous optim() call and using these as starting values in the next function call. The final call to optim() returns the following: $par [1] 0.2272361 0.8037642 26.8591998 3.0631280 0.2224566 $value [1] -46.13906 $counts function gradient 4 4 $convergence [1] 0 $message [1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH" $hessian [,1] [,2] [,3] [,4] [,5] [1,] 1.267070e+17 1.012691e+17 1.348054e+15 625551.58724 9.359559e+07 [2,] 1.012691e+17 8.189877e+16 1.144248e+15 569562.44945 8.699072e+07 [3,] 1.348054e+15 1.144248e+15 2.457323e+05 3426.60293 -2.297009e+03 [4,] 6.255516e+05 5.695624e+05 3.426603e+03 99.06880 -6.750806e+01 [5,] 9.359559e+07 8.699072e+07 -2.297009e+03 -67.50806 1.905247e+03 i.e. convergence and message report that things look "ok". However; if I take the hessian and compute eigenvalues and eigenvectors, the result is: $values [1] 2.080357e+17 5.889416e+14 1.907746e+03 9.648828e+01 -1.886641e+13 $vectors [,1] [,2] [,3] [,4] [,5] [1,] 7.797175e-01 6.246383e-01 -5.024932e-09 -2.214316e-10 -4.321721e-02 [2,] 6.260739e-01 -7.768497e-01 5.917821e-09 2.637892e-10 6.734993e-02 [3,] 8.496067e-03 -7.957109e-02 9.710924e-08 4.049691e-09 -9.967930e-01 [4,] 4.058779e-12 -8.828316e-11 3.729724e-02 -9.993042e-01 -4.192523e-10 [5,] 6.125907e-10 -1.547717e-08 -9.993042e-01 -3.729724e-02 -9.626469e-08 Optim() indicates convergence when one of the eigenvalues is negative !?!? Any reason to not be concerned?
I don't think you have scaled your problem as recommended on the help page, or the Hessian would not have elements differing by 15 orders of magnitude. So, yes, there is reason to be worried about the finite-difference approximations used, as well as the numerical accuracy of such extreme differences. (If you did scale, you don't want to do this on the rescaled Hessian.)
Possibly the Hessian is simply computed at the solution and not used to arrive at the estimate?)
Please do read what it says about the method on the help page and the references.
Brian D. Ripley, ripley at stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UK Fax: +44 1865 272595