[R] Obtaining SE from the hessian matrix
Timur Elzhov
Timur.Elzhov at jinr.ru
Thu Feb 19 14:38:03 CET 2004
Dear R experts,
In R-intro, under the 'Nonlinear least squares and maximum likelihood
models' there are ttwo examples considered how to use 'nlm' function.
In 'Least squares' the Standard Errors obtained as follows:
After the fitting, out$minimum is the SSE, and out$estimates are the
least squares estimates of the parameters. To obtain the approximate
standard errors (SE) of the estimates we do:
> sqrt(diag(2*out$minimum/(length(y) - 2) * solve(out$hessian)))
But under 'Maximum likelihood' section I've read:
After the fitting, out$minimum is the negative log-likelihood, and
out$estimates are the maximum likelihood estimates of the parameters.
To obtain the approximate SEs of the estimates we do:
> sqrt(diag(solve(out$hessian)))
As for me, I use MINPACK fortran library for NLS fitting in R, and there
I also get the hessian matrix. What formula should I use in _this_ case?
Well, some times ago I had a glance at gsl, GNU Scientific Library. It
use converted-to-C MINPACK for NLS fit too. And, in the GSL ref. manual
example
http://www.gnu.org/software/gsl/manual/html_node/gsl-ref_36.html#SEC475
SE calculated as
#define ERR(i) sqrt(gsl_matrix_get(covar,i,i))
where covar = (J^T * J)^-1 (i.e. how in the second formulae above).
So, what is the _right_ way for obtatining SE? Why two those formulas above
differ?
Thank you!
--
WBR,
Timur.
More information about the R-help
mailing list