Bug 16890 - nlme::intervals.lmList with seems not to work with unpooled SE estimates
Summary: nlme::intervals.lmList with seems not to work with unpooled SE estimates
Status: ASSIGNED
Alias: None
Product: R
Classification: Unclassified
Component: Models (show other bugs)
Version: R-devel (trunk)
Hardware: All All
: P5 enhancement
Assignee: R-core
URL:
Depends on:
Blocks:
 
Reported: 2016-05-07 00:54 UTC by Ben Bolker
Modified: 2016-05-07 16:36 UTC (History)
1 user (show)

See Also:


Attachments

Note You need to log in before you can comment on or make changes to this bug.
Description Ben Bolker 2016-05-07 00:54:36 UTC
?intervals.lmList suggests that we can get un-pooled
confidence intervals on the fly:
>    pool: an optional logical value indicating whether a pooled
>          estimate of the residual standard error should be used.
>          Default is ‘attr(object, "pool")’.

However, intervals.lmList doesn't seem to work for unpooled SE estimates,
whether the unpooling is specified in the original model or as an override
argument.

library(nlme)
fm1 <- lmList(height~age|Subject, Oxboys)  ## pool=TRUE (default)
fm2 <- lmList(height~age|Subject, Oxboys, pool=FALSE)
i1 <- intervals(fm1)  ## works (uses pooled SE)
attr(fm1,"pool") ## TRUE
attr(fm2,"pool") ## FALSE

intervals(fm1,pool=FALSE)
## Error in qf(level, 1, smry$df.residual) :
##   Non-numeric argument to mathematical function
intervals(fm2)  ## same error

The proximal problem seems to be that summary.lme gives
a NULL df.residual when pool=FALSE ...

summary(fm1)$df.residual  ## TRUE
summary(fm1, pool=FALSE)$df.residual  ## NULL
Comment 1 Martin Maechler 2016-05-07 16:13:12 UTC
Thank you, Ben.

At least this is a *very* old bug. I've checked a 2005 version of R's nlme...

Further suggestions {what are the correct 'df' ???}  are welcome.

Or should we (R-core) just update the documentation with its "suggestion" ?
Comment 2 Ben Bolker 2016-05-07 16:36:02 UTC
(In reply to Martin Maechler from comment #1)
> Thank you, Ben.
> 
> At least this is a *very* old bug. I've checked a 2005 version of R's nlme...
> 
> Further suggestions {what are the correct 'df' ???}  are welcome.
> 
> Or should we (R-core) just update the documentation with its "suggestion" ?

I would think that the correct df.residuals would be a vector, equivalent to whatever would have been estimated for the *individual* group fits (i.e. nobs per group - npar per group).  If having length(df.residuals)>1 would mess things up elsewhere, then intervals.lmList should probably just compute these values itself (rather than expecting it to be in the summary object)