Dear R users,
I want to modify how the degrees of freedom are calculated from the
summary function of the package minpack.lm
My first thought was to replicate how minpack.lm:::summary.nls.lm
works, and modify the line of the code that is relevant for my task.
However, for some reason I am not able to use the code contained in
this function to perform this operation.
Here is a reproducible example.
First, I run a NLLS regression using minpack.lm::nls.lm
# From example 1 of the help page of ?minpack.lm::nls.lm
library(minpack.lm)
x <- seq(0,5,length=100)
getPred <- function(parS, xx) parS$a * exp(xx * parS$b) + parS$c
pp <- list(a=9,b=-1, c=6)
simDNoisy <- getPred(pp,x) + rnorm(length(x),sd=.1)
residFun <- function(p, observed, xx) observed - getPred(p,xx)
parStart <- list(a=3,b=-.001, c=1)
nls.out <- nls.lm(par=parStart, fn = residFun, observed = simDNoisy,
xx = x, control = nls.lm.control(nprint=1))
summary(nls.out)
Now, by running minpack.lm:::summary.nls.lm in the console, I get the following
> minpack.lm:::summary.nls.lm
function (object, ...)
{
param <- coef(object)
pnames <- names(param)
ibb <- chol(object$hessian)
ih <- chol2inv(ibb)
p <- length(param)
rdf <- length(object$fvec) - p
resvar <- deviance(object)/rdf
se <- sqrt(diag(ih) * resvar)
names(se) <- pnames
tval <- param/se
param <- cbind(param, se, tval, 2 * pt(abs(tval), rdf, lower.tail = FALSE))
dimnames(param) <- list(pnames, c("Estimate", "Std. Error",
"t value", "Pr(>|t|)"))
ans <- list(residuals = object$fvec, sigma = sqrt(object$deviance/rdf),
df = c(p, rdf), cov.unscaled = ih, info = object$info,
niter = object$niter, stopmess = object$message, coefficients = param)
class(ans) <- "summary.nls.lm"
ans
}
Specifically, my task requires to modify the object p of this
function, say I want to set p <- 2.
To this purpose, I replicate what the summary function does using my
object (nls.out), however an error is returned at line three:
> param <- coef(nls.out)
> pnames <- names(param)
> ibb <- chol(nls.out$hessian)
Error in array(x, c(length(x), 1L), if (!is.null(names(x))) list(names(x), :
'data' must be of a vector type, was 'NULL'
The reason of the error is that there is no nls.out$hessian in nls.out.
I guess that I am missing something about how summary functions work,
and this is the reason why I cannot use the code from
minpack.lm:::summary.nls.lm outside the minpack.lm environment.
Does anyone have any hints on how to proceed? Any other way of dealing
with this issue will be equally appreciated.
Thank you,
Valerio
______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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.
