modular | R Documentation |
Modular functions for mixed model fits
lFormula(formula, data = NULL, REML = TRUE,
subset, weights, na.action, offset, contrasts = NULL,
control = lmerControl(), ...)
mkLmerDevfun(fr, X, reTrms, REML = TRUE, start = NULL,
verbose = 0, control = lmerControl(), ...)
optimizeLmer(devfun,
optimizer = formals(lmerControl)$optimizer,
restart_edge = formals(lmerControl)$restart_edge,
boundary.tol = formals(lmerControl)$boundary.tol,
start = NULL, verbose = 0L,
control = list(), ...)
glFormula(formula, data = NULL, family = gaussian,
subset, weights, na.action, offset, contrasts = NULL,
start, mustart, etastart, control = glmerControl(), ...)
mkGlmerDevfun(fr, X, reTrms, family, nAGQ = 1L,
verbose = 0L, maxit = 100L, control = glmerControl(), ...)
optimizeGlmer(devfun,
optimizer = if(stage == 1) "bobyqa" else "Nelder_Mead",
restart_edge = FALSE,
boundary.tol = formals(glmerControl)$boundary.tol,
verbose = 0L, control = list(),
nAGQ = 1L, stage = 1, start = NULL, ...)
updateGlmerDevfun(devfun, reTrms, nAGQ = 1L)
formula |
a two-sided linear formula object
describing both the fixed-effects and random-effects parts
of the model, with the response on the left of a |
data |
an optional data frame containing the
variables named in |
REML |
(logical) indicating to fit restricted maximum likelihood model. |
subset |
an optional expression indicating the
subset of the rows of |
weights |
an optional vector of ‘prior weights’ to be used
in the fitting process. Should be |
na.action |
a function that indicates what should
happen when the data contain |
offset |
this can be used to specify an a priori known
component to be included in the linear predictor during
fitting. This should be |
contrasts |
an optional |
control |
a list giving
|
fr |
A model frame containing the variables needed to create an
|
X |
fixed-effects design matrix |
reTrms |
information on random effects structure (see
|
start |
starting values (see |
verbose |
print output? |
maxit |
maximal number of Pwrss update iterations. |
devfun |
a deviance function, as generated by |
nAGQ |
number of Gauss-Hermite quadrature points |
stage |
optimization stage (1: nAGQ=0, optimize over theta only; 2: nAGQ possibly >0, optimize over theta and beta) |
optimizer |
character - name of optimizing
function(s). A character vector or list of functions:
length 1 for
Special provisions are made for For |
restart_edge |
see |
boundary.tol |
see |
family |
a GLM family; see |
mustart |
optional starting values on the scale of
the conditional mean; see |
etastart |
optional starting values on the scale of
the unbounded predictor; see |
... |
other potential arguments; for |
These functions make up the internal components of an [gn]lmer fit.
[g]lFormula
takes the arguments that would normally be
passed to [g]lmer
, checking for errors and processing the
formula and data input to create a list of objects required to fit a
mixed model.
mk(Gl|L)merDevfun
takes the output of the previous
step (minus the formula
component) and creates a
deviance function
optimize(Gl|L)mer
takes a
deviance function and optimizes over theta
(or
over theta
and beta
, if stage
is set
to 2 for optimizeGlmer
updateGlmerDevfun
takes the first stage of a GLMM
optimization (with nAGQ=0
, optimizing over theta
only)
and produces a second-stage deviance function
mkMerMod
takes the environment of a
deviance function, the results of an optimization, a list of
random-effect terms, a model frame, and a model all and produces a
[g]lmerMod
object.
lFormula
and glFormula
return a list containing
components:
model frame
fixed-effect design matrix
list containing information on random effects structure:
result of mkReTrms
(lFormula only): logical indicating if restricted maximum likelihood was used (Copy of argument.)
mkLmerDevfun
and mkGlmerDevfun
return a function to
calculate deviance (or restricted deviance) as a function of the
theta (random-effect) parameters. updateGlmerDevfun
returns a function to calculate the deviance as a function of a
concatenation of theta and beta (fixed-effect) parameters. These
deviance functions have an environment containing objects required
for their evaluation. CAUTION: The environment
of
functions returned by mk(Gl|L)merDevfun
contains reference
class objects (see ReferenceClasses
,
merPredD-class
, lmResp-class
), which
behave in ways that may surprise many users. For example, if the
output of mk(Gl|L)merDevfun
is naively copied, then
modifications to the original will also appear in the copy (and
vice versa). To avoid this behavior one must make a deep copy (see
ReferenceClasses
for details).
optimizeLmer
and optimizeGlmer
return the results of an
optimization.
### Fitting a linear mixed model in 4 modularized steps
## 1. Parse the data and formula:
lmod <- lFormula(Reaction ~ Days + (Days|Subject), sleepstudy)
names(lmod)
## 2. Create the deviance function to be optimized:
(devfun <- do.call(mkLmerDevfun, lmod))
ls(environment(devfun)) # the environment of 'devfun' contains objects
# required for its evaluation
## 3. Optimize the deviance function:
opt <- optimizeLmer(devfun)
opt[1:3]
## 4. Package up the results:
mkMerMod(environment(devfun), opt, lmod$reTrms, fr = lmod$fr)
### Same model in one line
lmer(Reaction ~ Days + (Days|Subject), sleepstudy)
### Fitting a generalized linear mixed model in six modularized steps
## 1. Parse the data and formula:
glmod <- glFormula(cbind(incidence, size - incidence) ~ period + (1 | herd),
data = cbpp, family = binomial)
#.... see what've got :
str(glmod, max=1, give.attr=FALSE)
## 2. Create the deviance function for optimizing over theta:
(devfun <- do.call(mkGlmerDevfun, glmod))
ls(environment(devfun)) # the environment of devfun contains lots of info
## 3. Optimize over theta using a rough approximation (i.e. nAGQ = 0):
(opt <- optimizeGlmer(devfun))
## 4. Update the deviance function for optimizing over theta and beta:
(devfun <- updateGlmerDevfun(devfun, glmod$reTrms))
## 5. Optimize over theta and beta:
opt <- optimizeGlmer(devfun, stage=2)
str(opt, max=1) # seeing what we'got
## 6. Package up the results:
(fMod <- mkMerMod(environment(devfun), opt, glmod$reTrms, fr = glmod$fr))
### Same model in one line
fM <- glmer(cbind(incidence, size - incidence) ~ period + (1 | herd),
data = cbpp, family = binomial)
all.equal(fMod, fM, check.attributes=FALSE, tolerance = 1e-12)
# ---- -- even tolerance = 0 may work
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.