Last chance! 50% off unlimited learning
Sale ends in
profile
() of [ng]lmer
fitted
models. The log()
method transforms a lmer profile
to the scale of the logarithm of the standard deviation of the random
effects, i.e. it returns a profile with all the .sigNN parameters
replaced by .lsigNN. The forward and backward splines for these
parameters are recalculated.
## S3 method for class 'merMod':
profile(fitted, which = 1:nptot, alphamax = 0.01,
maxpts = 100, delta = cutoff/8, verbose = 0, devtol = 1e-09,
maxmult = 10, startmethod = "prev", optimizer = "bobyqa",
signames = TRUE, ...)
## S3 method for class 'thpr':
as.data.frame(x, ...)
## S3 method for class 'thpr':
log(x, base = exp(1))
lmer(..)
.which
1 - alphamax
is the maximum alpha value for likelihood ratio confidence
regions; used to establish the range of values to be profiled.lmer
for details)..sigNN
should be used; otherwise, names are more meaningful
(but longer) of the form (sd|cor)_(effects)|(group)
. Note
that some code for profile transformations (e.thpr
(i.e., output of
profile
)profile()
returns
an object of S3 class "thpr"
, data.frame
-like. Methods for such a profile object
are notably confint()
, the three plots methods
xyplot
, densityplot
, and
splom
. Further, log()
(see above) and as.data.frame()
.
profile
(package
bootMer
.fm01ML <- lmer(Yield ~ 1|Batch, Dyestuff, REML = FALSE)
system.time(
tpr <- profile(fm01ML, optimizer="Nelder_Mead", which="beta_")
)## fast; as only *one* beta parameter is profiled over
## full profiling (default which means 'all) needs
## ~2.6s (on a 2010 Macbook Pro)
system.time( tpr <- profile(fm01ML))
## ~1s, + possible warning about bobyqa convergence
(confint(tpr) -> CIpr)
stopifnot(all.equal(CIpr,
array(c(12.1985292, 38.2299848, 1486.4515,
84.0630513, 67.6576964, 1568.54849), dim = 3:2,
dimnames = list(c(".sig01", ".sigma", "(Intercept)"),
c("2.5 %", "97.5 %"))),
tol= 1e-07))# 1.37e-9 {64b}
require(lattice)
xyplot(tpr)
xyplot(tpr, absVal=TRUE) # easier to see conf.int.s (and check symmetry)
densityplot(tpr, main="densityplot( profile(lmer(..)) )")
splom(tpr)
doMore <- lme4:::testLevel() > 1if(doMore) { ## not typically, for time constraint reasons
## Batch and residual variance only
system.time(tpr2 <- profile(fm01ML, which=1:2, optimizer="Nelder_Mead"))
print( xyplot(tpr2) )
print( xyplot(log(tpr2)) )# log(sigma) is better
## GLMM example
gm1 <- glmer(cbind(incidence, size - incidence) ~ period + (1 | herd),
data = cbpp, family = binomial)
## running time ~9 seconds on a modern machine:
print( system.time(pr4 <- profile(gm1)) )
print( xyplot(pr4,layout=c(5,1),as.table=TRUE) )
print( xyplot(log(pr4), absVal=TRUE) ) # log(sigma_1)
print( splom(pr4) )
print( system.time( # quicker: only sig01 and one fixed effect
pr2 <- profile(gm1, which=c("theta_", "period2"))))
print( confint(pr2) )
}
Run the code above in your browser using DataLab