Nothing
# Part of the lme4 package for estimating generalized linear
# mixed effect models. Copyright (C) 2003-2021 Douglas Bates,
# Martin Maechler, Ben Bolker, and Steven Walker.
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 2 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License along
# with this program; if not, write to the Free Software Foundation, Inc.,
# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
glFormula <- function(formula, data = NULL, subset, weights,
na.action, offset, contrasts = NULL,
start, mustart, etastart,
control = glmerControl(), verbose = FALSE, ...)
{
control <- control$checkControl
mf <- mc <- match.call()
ignoreArgs <- c("start", "verbose", "devFunOnly", "optimizer",
"control", "nAGQ")
l... <- list(...)
l... <- l...[names(l...) %not_in% ignoreArgs]
do.call(checkArgs, c(list("glmer"), l...))
cstr <- "check.formula.LHS"
checkCtrlLevels(cstr, control[[cstr]])
denv <- checkFormulaData(formula, data,
checkLHS = control$check.formula.LHS == "stop")
mc$formula <- formula <- as.formula(formula, env = denv)
m <- match(c("data", "subset", "weights", "na.action", "offset",
"mustart", "etastart"),
names(mf), 0L)
mf <- mf[c(1L, m)]
mf$drop.unused.levels <- TRUE
mf[[1L]] <- quote(stats::model.frame)
mf$na.action <- quote(stats::na.pass)
# fr and fr.form contains everything, build specialized
# frames below
fr.form <- subbart(subbars(formula))
environment(fr.form) <- environment(formula)
for (i in c("weights", "offset")) {
if (!eval(bquote(missing(x = .(i)))))
assign(i, get(i, parent.frame()), environment(fr.form))
}
mf$formula <- fr.form
fr <- eval(mf, parent.frame())
fr <- factorize(fr.form, fr, char.only = TRUE)
attr(fr, "formula") <- formula
attr(fr, "offset") <- mf$offset
if (!missing(start) && is.list(start)) {
attr(fr, "start") <- start$fixef
}
mf$na.action <- na.action
getVariableNames <- function(terms) {
# variables contains everything, even if there is a . - XX construction
# Subtracted out terms don't appear in the term labels.
factors <- attr(terms, "factors")
if (NROW(factors) > 0L) result <- rownames(factors)[rowSums(factors) > 0L]
else result <- attr(terms, "term.labels")
variables <- as.character(attr(terms, "variables"))[-1L]
response <- attr(terms, "response")
if (length(response) > 0L) {
response <- variables[response]
setdiff(result, response)
} else {
result
}
}
fixedform <- formula
RHSForm(fixedform) <- nobart(nobars(RHSForm(fixedform)))
RHSForm(mf$formula) <- RHSForm(fixedform)
fixedfr <- eval(mf, parent.frame())
fixedterms <- attr(fixedfr, "terms")
attr(attr(fr, "terms"), "predvars.fixed") <- attr(fixedterms, "predvars")
attr(attr(fr, "terms"), "varnames.fixed") <- getVariableNames(fixedterms)
attr(fr, "na.action.fixed") <- attr(fixedfr, "na.action")
bartform <- formula
RHSForm(bartform) <- allbart(nobars(RHSForm(bartform)))
RHSForm(mf$formula) <- RHSForm(bartform)
mf$drop.unused.levels <- FALSE
bartfr <- eval(mf, parent.frame())
bartterms <- attr(bartfr, "terms")
attr(attr(fr, "terms"), "predvars.bart") <- attr(bartterms, "predvars")
attr(attr(fr, "terms"), "varnames.bart") <- getVariableNames(bartterms)
bartlevels <-
lapply(colnames(attr(bartterms, "factors")), function(n) levels(bartfr[[n]]))
names(bartlevels) <- colnames(attr(bartterms, "factors"))
attr(attr(fr, "terms"), "levels.bart") <- bartlevels
attr(fr, "na.action.bart") <- attr(bartfr, "na.action")
ranform <- formula
RHSForm(ranform) <- subbars(RHSForm(reOnly(formula)))
mf$formula <- ranform
mf$drop.unused.levels <- TRUE
ranfr <- eval(mf, parent.frame())
ranterms <- attr(ranfr, "terms")
attr(attr(fr, "terms"), "predvars.random") <- attr(ranterms, "predvars")
attr(attr(fr, "terms"), "varnames.random") <- getVariableNames(ranterms)
attr(fr, "na.action.random") <- attr(ranfr, "na.action")
# We create a total frame with all of the data, and then subset the individual frames
# based on the total na situation
na.action.all <- c(attr(fixedfr, "na.action"), attr(bartfr, "na.action"), attr(ranfr, "na.action"))
if (length(na.action.all) > 0L) {
na.action.all <- sort(na.action.all[!duplicated(na.action.all)])
if (!is.null(class(attr(fr, "na.action.fixed"))))
class(na.action.all) <- class(attr(fr, "na.action.fixed"))
else if (!is.null(class(attr(fr, "na.action.bart"))))
class(na.action.all) <- class(attr(fr, "na.action.bart"))
else
class(na.action.all) <- class(attr(fr, "na.action.random"))
attr(fr, "na.action.all") <- na.action.all
all_rows <- seq_len(nrow(fr)) %not_in% na.action.all
if (!is.null(attr(fixedfr, "na.action"))) {
fixed_rows <- seq_len(nrow(fr)) %not_in% attr(fixedfr, "na.action")
fixedfr <- fixedfr[all_rows[fixed_rows],,drop = FALSE]
}
if (!is.null(attr(bartfr, "na.action"))) {
bart_rows <- seq_len(nrow(fr)) %not_in% attr(bartfr, "na.action")
bartfr <- bartfr[all_rows[bart_rows],,drop = FALSE]
}
if (!is.null(attr(ranfr, "na.action"))) {
random_rows <- seq_len(nrow(fr)) %not_in% attr(ranfr, "na.action")
ranfr <- ranfr[all_rows[random_rows],,drop = FALSE]
}
}
if (length(attr(ranterms, "term.labels")) > 0L) {
reTrms <- mkReTrms(findbars(RHSForm(formula)), ranfr)
wmsgNlev <- checkNlevels(reTrms$flist, n = nrow(ranfr), control, allow.n = TRUE)
wmsgZdims <- checkZdims(reTrms$Ztlist, n = nrow(ranfr), control, allow.n = TRUE)
wmsgZrank <- checkZrank(reTrms$Zt, n = nrow(ranfr), control, nonSmall = 1e+06,
allow.n = TRUE)
} else {
reTrms = list()
}
# dbartsData can't handle columns in the model frame not in the formula terms, so we
# pull out "(offset)" if it exists
bartprednames <- as.character(attr(bartterms, "predvars"))[-1L]
if (any(names(bartfr) %not_in% bartprednames)) {
keepcols <- names(bartfr) %in% bartprednames
bartfr <- bartfr[,keepcols, drop = FALSE]
attr(bartfr, "terms") <- bartterms
attr(attr(bartfr, "terms"), "dataClasses") <- attr(bartterms, "dataClasses")[keepcols]
bartterms <- attr(bartfr, "terms")
attr(attr(fr, "terms"), "varnames.bart") <- getVariableNames(bartterms)
}
bartData <- dbarts::dbartsData(bartform, bartfr)
if (ncol(bartData@x) == 0L)
stop("no bart component detected in formula; consider using rstanarm package instead")
X <- model.matrix(fixedform, fixedfr, contrasts)
# drop intercept, offset
keepcols <- colnames(X) %not_in% c("(Intercept)", "`(offset)`")
X <- X[,keepcols,drop=FALSE]
if (is.null(rankX.chk <- control[["check.rankX"]]))
rankX.chk <- eval(formals(lmerControl)[["check.rankX"]])[[1]]
X <- chkRank.drop.cols(X, kind = rankX.chk, tol = 1e-07)
if (is.null(scaleX.chk <- control[["check.scaleX"]]))
scaleX.chk <- eval(formals(lmerControl)[["check.scaleX"]])[[1]]
X <- checkScaleX(X, kind = scaleX.chk)
terms <- attr(fr, "terms")
varnames.random <- attr(terms, "varnames.random")
varnames.fixed <- attr(terms, "varnames.fixed")
varnames.bart <- attr(terms, "varnames.bart")
varnames.mixed <- varnames.bart %in% varnames.random | varnames.bart %in% varnames.fixed
if (any(varnames.mixed) && as.integer(verbose) > -1L) {
warning("variable(s) '", paste0(varnames.bart[varnames.mixed], collapse = "', '"),
"' that appear in both parametric and nonparametric are not identifiable; ",
"model will fit but some results may be uninterpretable")
}
if (length(varnames.random) == 0L && length(varnames.fixed) == 0L)
stop("no parametric component detected in formula; consider using dbarts package instead")
result <- list(fr = fr, X = X, bartData = bartData, reTrms = reTrms, formula = formula,
terms = terms)
if (length(reTrms) > 0L)
result$wmsgs <- c(Nlev = wmsgNlev, Zdims = wmsgZdims, Zrank = wmsgZrank)
result
}
tryResult <- tryCatch(lme4ns <- asNamespace("lme4"), error = function(e) e)
if (!inherits(tryResult, "error")) {
`%ORifNotInLme4%` <- function(a, b) {
mc <- match.call()
if (is.symbol(mc[["a"]])) a <- deparse(mc[["a"]])
if (!is.null(lme4ns[[a]])) lme4ns[[a]] else b
}
} else {
`%ORifNotInLme4%` <- function(a, b) b
}
glmerControl <- glmerControl %ORifNotInLme4% function (optimizer = c("bobyqa", "Nelder_Mead"), restart_edge = FALSE,
boundary.tol = 1e-05, calc.derivs = TRUE, use.last.params = FALSE,
sparseX = FALSE, standardize.X = FALSE, check.nobs.vs.rankZ = "ignore",
check.nobs.vs.nlev = "stop", check.nlev.gtreq.5 = "ignore",
check.nlev.gtr.1 = "stop", check.nobs.vs.nRE = "stop", check.rankX = c("message+drop.cols",
"silent.drop.cols", "warn+drop.cols", "stop.deficient",
"ignore"), check.scaleX = c("warning", "stop", "silent.rescale",
"message+rescale", "warn+rescale", "ignore"), check.formula.LHS = "stop",
check.conv.grad = .makeCC("warning", tol = 0.002, relTol = NULL),
check.conv.singular = .makeCC(action = "message", tol = formals(isSingular)$tol),
check.conv.hess = .makeCC(action = "warning", tol = 1e-06),
optCtrl = list(), mod.type = "glmer", tolPwrss = 1e-07, compDev = TRUE,
nAGQ0initStep = TRUE, check.response.not.const = "stop")
{
stopifnot(is.list(optCtrl))
if (mod.type == "glmer" && length(optimizer) == 1) {
optimizer <- replicate(2, optimizer)
}
c.opts <- paste0(mod.type, "Control")
merOpts <- getOption(c.opts)
if (!is.null(merOpts)) {
nn <- names(merOpts)
nn.ok <- .get.checkingOpts(names(merOpts))
if (length(nn.ignored <- setdiff(nn, nn.ok)) > 0) {
warning("some options in ", shQuote(sprintf("getOption('%s')",
c.opts)), " ignored : ", paste(nn.ignored, collapse = ", "))
}
for (arg in nn.ok) {
if (do.call(missing, list(arg)))
assign(arg, merOpts[[arg]])
}
}
check.rankX <- match.arg(check.rankX)
check.scaleX <- match.arg(check.scaleX)
me <- sys.function()
chk.cconv(check.conv.grad, me)
chk.cconv(check.conv.singular, me)
chk.cconv(check.conv.hess, me)
if (mod.type == "glmer" && use.last.params && calc.derivs) {
warning("using ", shQuote("use.last.params"), "=TRUE and ",
shQuote("calc.derivs"), "=TRUE with ", shQuote("glmer"),
" will not give backward-compatible results")
}
ret <- namedList(optimizer, restart_edge, boundary.tol, calc.derivs,
use.last.params, checkControl = namedList(check.nobs.vs.rankZ,
check.nobs.vs.nlev, check.nlev.gtreq.5, check.nlev.gtr.1,
check.nobs.vs.nRE, check.rankX, check.scaleX, check.formula.LHS),
checkConv = namedList(check.conv.grad, check.conv.singular,
check.conv.hess), optCtrl = optCtrl)
if (mod.type == "glmer") {
ret <- c(ret, namedList(tolPwrss, compDev, nAGQ0initStep))
ret$checkControl <- c(ret$checkControl, namedList(check.response.not.const))
}
class(ret) <- c(c.opts, "merControl")
ret
}
tolPwrss <- tolPwrss %ORifNotInLme4% NULL
compDev <- compDev %ORifNotInLme4% NULL
nAGQ0initStep <- nAGQ0initStep %ORifNotInLme4% NULL
check.response.not.const <- check.response.not.const %ORifNotInLme4% NULL
lmerControl <- lmerControl %ORifNotInLme4%
function (optimizer = "nloptwrap", restart_edge = TRUE, boundary.tol = 1e-05,
calc.derivs = TRUE, use.last.params = FALSE, sparseX = FALSE,
standardize.X = FALSE, check.nobs.vs.rankZ = "ignore", check.nobs.vs.nlev = "stop",
check.nlev.gtreq.5 = "ignore", check.nlev.gtr.1 = "stop",
check.nobs.vs.nRE = "stop", check.rankX = c("message+drop.cols",
"silent.drop.cols", "warn+drop.cols", "stop.deficient",
"ignore"), check.scaleX = c("warning", "stop", "silent.rescale",
"message+rescale", "warn+rescale", "ignore"), check.formula.LHS = "stop",
check.conv.grad = .makeCC("warning", tol = 0.002, relTol = NULL),
check.conv.singular = .makeCC(action = "message", tol = formals(isSingular)$tol),
check.conv.hess = .makeCC(action = "warning", tol = 1e-06),
optCtrl = list(), mod.type = "lmer")
{
stopifnot(is.list(optCtrl))
if (mod.type == "glmer" && length(optimizer) == 1) {
optimizer <- replicate(2, optimizer)
}
c.opts <- paste0(mod.type, "Control")
merOpts <- getOption(c.opts)
if (!is.null(merOpts)) {
nn <- names(merOpts)
nn.ok <- .get.checkingOpts(names(merOpts))
if (length(nn.ignored <- setdiff(nn, nn.ok)) > 0) {
warning("some options in ", shQuote(sprintf("getOption('%s')",
c.opts)), " ignored : ", paste(nn.ignored, collapse = ", "))
}
for (arg in nn.ok) {
if (do.call(missing, list(arg)))
assign(arg, merOpts[[arg]])
}
}
check.rankX <- match.arg(check.rankX)
check.scaleX <- match.arg(check.scaleX)
me <- sys.function()
chk.cconv(check.conv.grad, me)
chk.cconv(check.conv.singular, me)
chk.cconv(check.conv.hess, me)
if (mod.type == "glmer" && use.last.params && calc.derivs) {
warning("using ", shQuote("use.last.params"), "=TRUE and ",
shQuote("calc.derivs"), "=TRUE with ", shQuote("glmer"),
" will not give backward-compatible results")
}
ret <- namedList(optimizer, restart_edge, boundary.tol, calc.derivs,
use.last.params, checkControl = namedList(check.nobs.vs.rankZ,
check.nobs.vs.nlev, check.nlev.gtreq.5, check.nlev.gtr.1,
check.nobs.vs.nRE, check.rankX, check.scaleX, check.formula.LHS),
checkConv = namedList(check.conv.grad, check.conv.singular,
check.conv.hess), optCtrl = optCtrl)
if (mod.type == "glmer") {
ret <- c(ret, namedList(tolPwrss, compDev, nAGQ0initStep))
ret$checkControl <- c(ret$checkControl, namedList(check.response.not.const))
}
class(ret) <- c(c.opts, "merControl")
ret
}
checkArgs <- checkArgs %ORifNotInLme4% function (type, ...)
{
l... <- list(...)
if (isTRUE(l...[["sparseX"]]))
warning("sparseX = TRUE has no effect at present", call. = FALSE)
if (length(l... <- list(...))) {
if (!is.null(l...[["family"]])) {
warning("calling lmer with family() is deprecated: please use glmer() instead",
call. = FALSE)
type <- "glmer"
}
if (!is.null(l...[["method"]])) {
msg <- paste("Argument", sQuote("method"), "is deprecated.")
if (type == "lmer")
msg <- paste(msg, "Use the REML argument to specify ML or REML estimation.")
else if (type == "glmer")
msg <- paste(msg, "Use the nAGQ argument to specify Laplace (nAGQ=1) or adaptive",
"Gauss-Hermite quadrature (nAGQ>1). PQL is no longer available.")
warning(msg, call. = FALSE)
l... <- l...[names(l...) != "method"]
}
if (length(l...)) {
warning("extra argument(s) ", paste(sQuote(names(l...)),
collapse = ", "), " disregarded", call. = FALSE)
}
}
}
checkCtrlLevels <- checkCtrlLevels %ORifNotInLme4% function (cstr, val, smallOK = FALSE)
{
bvals <- c("message", "warning", "stop", "ignore")
if (smallOK)
bvals <- outer(bvals, c("", "Small"), paste0)
if (!is.null(val) && !val %in% bvals)
stop("invalid control level ", sQuote(val), " in ", cstr,
": valid options are {", paste(sapply(bvals, sQuote),
collapse = ","), "}")
invisible(NULL)
}
checkFormulaData <- checkFormulaData %ORifNotInLme4% function (formula, data, checkLHS = TRUE, checkData = TRUE, debug = FALSE)
{
nonexist.data <- missDataFun(data)
wd <- tryCatch(eval(data), error = identity)
if (wrong.data <- inherits(wd, "simpleError")) {
wrong.data.msg <- wd$message
}
bad.data <- nonexist.data || wrong.data
if (bad.data || debug) {
varex <- function(v, env) exists(v, envir = env, inherits = FALSE)
allvars <- all.vars(as.formula(formula))
allvarex <- function(env, vvec = allvars) all(vapply(vvec,
varex, NA, env))
}
if (bad.data) {
if (allvarex(environment(formula))) {
stop("bad 'data', but variables found in environment of formula: ",
"try specifying 'formula' as a formula rather ",
"than a string in the original model", call. = FALSE)
}
else {
if (nonexist.data) {
stop("'data' not found, and some variables missing from formula environment",
call. = FALSE)
}
else {
stop("bad 'data': ", wrong.data.msg)
}
}
}
else {
denv <- if (is.null(data)) {
if (!is.null(ee <- environment(formula))) {
ee
}
else {
parent.frame(2L)
}
}
else list2env(data)
}
if (debug) {
cat("Debugging parent frames in checkFormulaData:\n")
glEnv <- 1L
while (!identical(parent.frame(glEnv), .GlobalEnv)) {
glEnv <- glEnv + 1L
}
for (i in 1:glEnv) {
OK <- allvarex(parent.frame(i))
cat("vars exist in parent frame ", i)
if (i == glEnv)
cat(" (global)")
cat(" ", OK, "\n")
}
cat("vars exist in env of formula ", allvarex(denv),
"\n")
}
stopifnot(!checkLHS || length(as.formula(formula, env = denv)) ==
3)
return(denv)
}
mkReTrms <- mkReTrms %ORifNotInLme4% function (bars, fr, drop.unused.levels = TRUE, reorder.terms = TRUE,
reorder.vars = FALSE)
{
if (!length(bars))
stop("No random effects terms specified in formula",
call. = FALSE)
stopifnot(is.list(bars), vapply(bars, is.language, NA), inherits(fr,
"data.frame"))
names(bars) <- barnames(bars)
term.names <- vapply(bars, safeDeparse, "")
blist <- lapply(bars, mkBlist, fr, drop.unused.levels, reorder.vars = reorder.vars)
nl <- vapply(blist, `[[`, 0L, "nl")
if (reorder.terms) {
if (any(diff(nl) > 0)) {
ord <- rev(order(nl))
blist <- blist[ord]
nl <- nl[ord]
term.names <- term.names[ord]
}
}
Ztlist <- lapply(blist, `[[`, "sm")
Zt <- do.call(rbind, Ztlist)
names(Ztlist) <- term.names
q <- nrow(Zt)
cnms <- lapply(blist, `[[`, "cnms")
nc <- lengths(cnms)
nth <- as.integer((nc * (nc + 1))/2)
nb <- nc * nl
if (sum(nb) != q) {
stop(sprintf("total number of RE (%d) not equal to nrow(Zt) (%d)",
sum(nb), q))
}
boff <- cumsum(c(0L, nb))
thoff <- cumsum(c(0L, nth))
Lambdat <- Matrix::t(do.call(Matrix::sparseMatrix, do.call(rbind, lapply(seq_along(blist),
function(i) {
mm <- matrix(seq_len(nb[i]), ncol = nc[i], byrow = TRUE)
dd <- diag(nc[i])
ltri <- lower.tri(dd, diag = TRUE)
ii <- row(dd)[ltri]
jj <- col(dd)[ltri]
data.frame(i = as.vector(mm[, ii]) + boff[i], j = as.vector(mm[,
jj]) + boff[i], x = as.double(rep.int(seq_along(ii),
rep.int(nl[i], length(ii))) + thoff[i]))
}))))
thet <- numeric(sum(nth))
ll <- list(Zt = Matrix::drop0(Zt), theta = thet, Lind = as.integer(Lambdat@x),
Gp = unname(c(0L, cumsum(nb))))
ll$lower <- -Inf * (thet + 1)
ll$lower[unique(Matrix::diag(Lambdat))] <- 0
ll$theta[] <- is.finite(ll$lower)
Lambdat@x[] <- ll$theta[ll$Lind]
ll$Lambdat <- Lambdat
fl <- lapply(blist, `[[`, "ff")
fnms <- names(fl)
if (length(fnms) > length(ufn <- unique(fnms))) {
fl <- fl[match(ufn, fnms)]
asgn <- match(fnms, ufn)
}
else asgn <- seq_along(fl)
names(fl) <- ufn
attr(fl, "assign") <- asgn
ll$flist <- fl
ll$cnms <- cnms
ll$Ztlist <- Ztlist
ll
}
findbars <- findbars %ORifNotInLme4% function (term)
{
fb <- function(term) {
if (is.name(term) || !is.language(term))
return(NULL)
if (term[[1]] == as.name("("))
return(fb(term[[2]]))
stopifnot(is.call(term))
if (term[[1]] == as.name("|"))
return(term)
if (length(term) == 2)
return(fb(term[[2]]))
c(fb(term[[2]]), fb(term[[3]]))
}
expandSlash <- function(bb) {
makeInteraction <- function(x) {
if (length(x) < 2)
return(x)
trm1 <- makeInteraction(x[[1]])
trm11 <- if (is.list(trm1))
trm1[[1]]
else trm1
list(substitute(foo:bar, list(foo = x[[2]], bar = trm11)),
trm1)
}
slashTerms <- function(x) {
if (!("/" %in% all.names(x)))
return(x)
if (x[[1]] != as.name("/"))
stop("unparseable formula for grouping factor",
call. = FALSE)
list(slashTerms(x[[2]]), slashTerms(x[[3]]))
}
if (!is.list(bb))
expandSlash(list(bb))
else unlist(lapply(bb, function(x) {
if (length(x) > 2 && is.list(trms <- slashTerms(x[[3]])))
lapply(unlist(makeInteraction(trms)), function(trm) substitute(foo |
bar, list(foo = x[[2]], bar = trm)))
else x
}))
}
modterm <- expandDoubleVerts(if (inherits(term, "formula"))
term[[length(term)]]
else term)
expandSlash(fb(modterm))
}
RHSForm <- RHSForm %ORifNotInLme4% function (form, as.form = FALSE)
{
rhsf <- form[[length(form)]]
if (as.form)
reformulate(deparse(rhsf))
else rhsf
}
factorize <- factorize %ORifNotInLme4% function (x, frloc, char.only = FALSE)
{
for (i in all.vars(RHSForm(x))) {
if (!is.null(curf <- frloc[[i]]))
frloc[[i]] <- makeFac(curf, char.only)
}
return(frloc)
}
checkNlevels <- checkNLevels %ORifNotInLme4% function (flist, n, ctrl, allow.n = FALSE)
{
stopifnot(is.list(ctrl), is.numeric(n))
nlevelVec <- vapply(flist, function(x) nlevels(factor(x,
exclude = NA)), 1)
cstr <- "check.nlev.gtr.1"
checkCtrlLevels(cstr, cc <- ctrl[[cstr]])
if (doCheck(cc) && any(nlevelVec < 2)) {
wstr <- "grouping factors must have > 1 sampled level"
switch(cc, warning = warning(wstr, call. = FALSE), stop = stop(wstr,
call. = FALSE), stop(gettextf("unknown check level for '%s'",
cstr), domain = NA))
}
else wstr <- character()
cstr <- "check.nobs.vs.nlev"
checkCtrlLevels(cstr, cc <- ctrl[[cstr]])
if (doCheck(cc) && any(if (allow.n) nlevelVec > n else nlevelVec >=
n)) {
wst2 <- gettextf("number of levels of each grouping factor must be %s number of observations",
if (allow.n)
"<="
else "<")
switch(cc, warning = warning(wst2, call. = FALSE), stop = stop(wst2,
call. = FALSE))
}
else wst2 <- character()
cstr <- "check.nlev.gtreq.5"
checkCtrlLevels(cstr, cc <- ctrl[[cstr]])
if (doCheck(cc) && any(nlevelVec < 5)) {
wst3 <- "grouping factors with < 5 sampled levels may give unreliable estimates"
switch(cc, warning = warning(wst3, call. = FALSE), stop = stop(wst3,
call. = FALSE), stop(gettextf("unknown check level for '%s'",
cstr), domain = NA))
}
else wst3 <- character()
c(wstr, wst2, wst3)
}
checkZdims <- checkZdims %ORifNotInLme4% function (Ztlist, n, ctrl, allow.n = FALSE)
{
stopifnot(is.list(Ztlist), is.numeric(n))
cstr <- "check.nobs.vs.nRE"
checkCtrlLevels(cstr, cc <- ctrl[[cstr]])
term.names <- names(Ztlist)
rows <- vapply(Ztlist, nrow, 1L)
cols <- vapply(Ztlist, ncol, 1L)
stopifnot(all(cols == n))
if (doCheck(cc)) {
unique(unlist(lapply(seq_along(Ztlist), function(i) {
ww <- wmsg(cols[i], rows[i], allow.n, "number of observations",
"number of random effects", sprintf(" for term (%s)",
term.names[i]))
if (ww$unident) {
switch(cc, warning = warning(ww$wstr, call. = FALSE),
stop = stop(ww$wstr, call. = FALSE), stop(gettextf("unknown check level for '%s'",
cstr), domain = NA))
ww$wstr
}
else character()
})))
}
else character()
}
checkZrank <- checkZrank %ORifNotInLme4% function (Zt, n, ctrl, nonSmall = 1e+06, allow.n = FALSE)
{
stopifnot(is.list(ctrl), is.numeric(n), is.numeric(nonSmall))
cstr <- "check.nobs.vs.rankZ"
if (doCheck(cc <- ctrl[[cstr]])) {
checkCtrlLevels(cstr, cc, smallOK = TRUE)
d <- dim(Zt)
doTr <- d[1L] < d[2L]
if (!(grepl("Small", cc) && prod(d) > nonSmall)) {
rankZ <- Matrix::rankMatrix(if (doTr)
t(Zt)
else Zt, method = "qr", sval = numeric(min(d)))
ww <- wmsg(n, rankZ, allow.n, "number of observations",
"rank(Z)")
if (ww$unident) {
switch(cc, warningSmall = , warning = warning(ww$wstr,
call. = FALSE), stopSmall = , stop = stop(ww$wstr,
call. = FALSE), stop(gettextf("unknown check level for '%s'",
cstr), domain = NA))
ww$wstr
}
else character()
}
else character()
}
else character()
}
subbart <- function(term)
{
if (is.name(term) || !is.language(term)) return(term)
if (length(term) == 2) {
if (term[[1]] == as.name("bart"))
term[[1]] <- as.name("(")
term[[2]] <- subbart(term[[2]])
return(term)
}
for (j in 2:length(term)) term[[j]] <- subbart(term[[j]])
term
}
subbars <- subbars %ORifNotInLme4% function (term)
{
if (is.name(term) || !is.language(term))
return(term)
if (length(term) == 2) {
term[[2]] <- subbars(term[[2]])
return(term)
}
stopifnot(length(term) >= 3)
if (is.call(term) && term[[1]] == as.name("|"))
term[[1]] <- as.name("+")
if (is.call(term) && term[[1]] == as.name("||"))
term[[1]] <- as.name("+")
for (j in 2:length(term)) term[[j]] <- subbars(term[[j]])
term
}
chkRank.drop.cols <- chkRank.drop.cols %ORifNotInLme4% function (X, kind, tol = 1e-07, method = "qr.R")
{
stopifnot(is.matrix(X))
kinds <- eval(formals(lmerControl)[["check.rankX"]])
if (!kind %in% kinds)
stop(gettextf("undefined option for 'kind': %s", kind))
if (kind == "ignore")
return(X)
p <- ncol(X)
if (kind == "stop.deficient") {
if ((rX <- Matrix::rankMatrix(X, tol = tol, method = method)) <
p)
stop(gettextf(sub("\n +", "\n", "the fixed-effects model matrix is column rank deficient (rank(X) = %d < %d = p);\n the fixed effects will be jointly unidentifiable"),
rX, p), call. = FALSE)
}
else {
qr.X <- qr(X, tol = tol, LAPACK = FALSE)
rnkX <- qr.X$rank
if (rnkX == p)
return(X)
msg <- sprintf(ngettext(p - rnkX, "fixed-effect model matrix is rank deficient so dropping %d column / coefficient",
"fixed-effect model matrix is rank deficient so dropping %d columns / coefficients"),
p - rnkX)
if (kind != "silent.drop.cols")
(if (kind == "warn+drop.cols")
warning
else message)(msg, domain = NA)
contr <- attr(X, "contrasts")
asgn <- attr(X, "assign")
keep <- qr.X$pivot[seq_len(rnkX)]
dropped.names <- colnames(X[, -keep, drop = FALSE])
X <- X[, keep, drop = FALSE]
if (Matrix::rankMatrix(X, tol = tol, method = method) < ncol(X))
stop(gettextf("Dropping columns failed to produce full column rank design matrix"),
call. = FALSE)
if (!is.null(contr))
attr(X, "contrasts") <- contr
if (!is.null(asgn))
attr(X, "assign") <- asgn[keep]
attr(X, "msgRankdrop") <- msg
attr(X, "col.dropped") <- setNames(qr.X$pivot[(rnkX +
1L):p], dropped.names)
}
X
}
checkScaleX <- checkScaleX %ORifNotInLme4% function (X, kind = "warning", tol = 1000)
{
kinds <- eval(formals(lmerControl)[["check.scaleX"]])
if (!kind %in% kinds)
stop(gettextf("unknown check-scale option: %s", kind))
if (is.null(kind) || kind == "ignore")
return(X)
cont.cols <- apply(X, 2, function(z) !all(z %in% c(0, 1)))
col.sd <- apply(X[, cont.cols, drop = FALSE], 2L, sd)
sdcomp <- outer(col.sd, col.sd, "/")
logcomp <- abs(log(sdcomp[lower.tri(sdcomp)]))
logsd <- abs(log(col.sd))
if (any(c(logcomp, logsd) > log(tol))) {
wmsg <- "Some predictor variables are on very different scales:"
if (kind %in% c("warning", "stop")) {
wmsg <- paste(wmsg, "consider rescaling")
switch(kind, warning = warning(wmsg, call. = FALSE),
stop = stop(wmsg, call. = FALSE))
}
else {
wmsg <- paste(wmsg, "auto-rescaled (results NOT adjusted)")
X[, cont.cols] <- sweep(X[, cont.cols, drop = FALSE],
2, col.sd, "/")
attr(X, "scaled:scale") <- setNames(col.sd, colnames(X)[cont.cols])
if (kind == "warn+rescale")
warning(wmsg, call. = FALSE)
}
}
else wmsg <- character()
structure(X, msgScaleX = wmsg)
}
missDataFun <- missDataFun %ORifNotInLme4% function(d) {
ff <- sys.frames()
ex <- substitute(d)
ii <- rev(seq_along(ff))
foundAnon <- FALSE
for(i in ii) {
ex <- eval(substitute(substitute(x, env=sys.frames()[[n]]),
env = list(x = ex, n=i)))
if (is.symbol(ex) && grepl("^\\.\\.[0-9]+$",safeDeparse(ex))) {
foundAnon <- TRUE
break
}
}
return(!foundAnon && is.symbol(ex) && !exists(deparse(ex)))
}
safeDeparse <- safeDeparse %ORifNotInLme4% function(x, collapse=" ") paste(deparse(x, 500L), collapse=collapse)
makeFac <- makeFac %ORifNotInLme4% function (x, char.only = FALSE)
{
if (!is.factor(x) && (!char.only || is.character(x)))
factor(x)
else x
}
expandDoubleVerts <- expandDoubleVerts %ORifNotInLme4% function (term)
{
expandDoubleVert <- function(term) {
frml <- formula(substitute(~x, list(x = term[[2]])))
newtrms <- paste0("0+", attr(terms(frml), "term.labels"))
if (attr(terms(frml), "intercept") != 0)
newtrms <- c("1", newtrms)
as.formula(paste("~(", paste(vapply(newtrms, function(trm) paste0(trm,
"|", deparse(term[[3]])), ""), collapse = ")+("),
")"))[[2]]
}
if (!is.name(term) && is.language(term)) {
if (term[[1]] == as.name("(")) {
term[[2]] <- expandDoubleVerts(term[[2]])
}
stopifnot(is.call(term))
if (term[[1]] == as.name("||"))
return(expandDoubleVert(term))
term[[2]] <- expandDoubleVerts(term[[2]])
if (length(term) != 2) {
if (length(term) == 3)
term[[3]] <- expandDoubleVerts(term[[3]])
}
}
term
}
barnames <- barnames %ORifNotInLme4% function (bars)
vapply(bars, function(x) safeDeparse(x[[3]]), "")
mkBlist <- mkBlist %ORifNotInLme4% function (x, frloc, drop.unused.levels = TRUE, reorder.vars = FALSE)
{
frloc <- factorize(x, frloc)
if (is.null(ff <- tryCatch(eval(substitute(makeFac(fac),
list(fac = x[[3]])), frloc), error = function(e) NULL)))
stop("couldn't evaluate grouping factor ", deparse(x[[3]]),
" within model frame:", " try adding grouping factor to data ",
"frame explicitly if possible", call. = FALSE)
if (all(is.na(ff)))
stop("Invalid grouping factor specification, ", deparse(x[[3]]),
call. = FALSE)
if (drop.unused.levels)
ff <- factor(ff, exclude = NA)
nl <- length(levels(ff))
mm <- model.matrix(eval(substitute(~foo, list(foo = x[[2]]))),
frloc)
if (reorder.vars) {
mm <- mm[colSort(colnames(mm)), ]
}
sm <- Matrix::fac2sparse(ff, to = "d", drop.unused.levels = drop.unused.levels)
sm <- Matrix::KhatriRao(sm, t(mm))
dimnames(sm) <- list(rep(levels(ff), each = ncol(mm)), rownames(mm))
list(ff = ff, sm = sm, nl = nl, cnms = colnames(mm))
}
doCheck <- doCheck %ORifNotInLme4% function (x)
{
is.character(x) && !any(x == "ignore")
}
nobars_ <- nobars_ %ORifNotInLme4% function (term)
{
if (!anyBars(term))
return(term)
if (isBar(term))
return(NULL)
if (isAnyArgBar(term))
return(NULL)
if (length(term) == 2) {
nb <- nobars_(term[[2]])
if (is.null(nb))
return(NULL)
term[[2]] <- nb
return(term)
}
nb2 <- nobars_(term[[2]])
nb3 <- nobars_(term[[3]])
if (is.null(nb2))
return(nb3)
if (is.null(nb3))
return(nb2)
term[[2]] <- nb2
term[[3]] <- nb3
term
}
nobars <- nobars %ORifNotInLme4% function (term)
{
nb <- nobars_(term)
if (inherits(term, "formula") && length(term) == 3 && is.symbol(nb)) {
nb <- reformulate("1", response = deparse(nb))
}
if (is.null(nb)) {
nb <- if (inherits(term, "formula"))
~1
else 1
}
nb
}
anyBars <- anyBars %ORifNotInLme4% function (term)
{
any(c("|", "||") %in% all.names(term))
}
isBar <- isBar %ORifNotInLme4% function(term)
{
if (is.call(term)) {
if ((term[[1]] == as.name("|")) || (term[[1]] == as.name("||"))) {
return(TRUE)
}
}
FALSE
}
nobart_ <- function(term)
{
if (length(term) == 1L)
return(if (term == as.name("bart")) NULL else term)
if (length(term) == 2L) {
if (term[[1]] == as.name("bart")) return(NULL)
nb <- nobart_(term[[2]])
if (is.null(nb)) return(NULL)
term[[2]] <- nb
return(term)
}
nb2 <- nobart_(term[[2]])
nb3 <- nobart_(term[[3]])
if (is.null(nb2)) return(nb3)
if (is.null(nb3)) return(nb2)
term[[2]] <- nb2
term[[3]] <- nb3
term
}
nobart <- function(term)
{
nb <- nobart_(term)
if (inherits(term, "formula") && length(term) == 3 && is.symbol(nb)) {
nb <- reformulate("1", response = deparse(nb))
}
if (is.null(nb)) {
nb <- if (inherits(term, "formula"))
~1
else 1
}
nb
}
allbart_ <- function(term, inbart)
{
if (length(term) == 1L)
return(if (term == as.name("bart") || !inbart) NULL else term)
if (length(term) == 2L) {
if (term[[1]] == as.name("bart")) return(allbart_(term[[2L]], TRUE))
ab <- allbart_(term[[2]], inbart)
if (is.null(ab)) return(NULL)
term[[2L]] <- ab
return(term)
}
# should be a + b, or some other binary operator
ab2 <- allbart_(term[[2L]], inbart)
ab3 <- allbart_(term[[3L]], inbart)
if (is.null(ab2)) return(ab3)
if (is.null(ab3)) return(ab2)
term[[2L]] <- ab2
term[[3L]] <- ab3
term
}
allbart <- function(term)
{
ab <- allbart_(term, FALSE)
if (inherits(term, "formula") && length(term) == 3 && is.symbol(ab)) {
ab <- reformulate("1", response = deparse(ab))
}
if (is.null(ab)) {
ab <- if (inherits(term, "formula"))
~1
else 1
}
ab
}
`RHSForm<-` <- `RHSForm<-` %ORifNotInLme4% function (formula, value)
{
formula[[length(formula)]] <- value
formula
}
reOnly <- function (f, response = FALSE)
{
found_bars <- findbars(f)
if (is.null(found_bars)) return(~0)
reformulate(paste0("(", vapply(findbars(f), safeDeparse,
""), ")"), response = if (response && length(f) == 3L)
f[[2]])
}
mkVarCorr <- mkVarCorr %ORifNotInLme4% function (sc, cnms, nc, theta, nms)
{
ncseq <- seq_along(nc)
thl <- split(theta, rep.int(ncseq, (nc * (nc + 1))/2))
if (!all(nms == names(cnms)))
warning("nms != names(cnms) -- whereas lme4-authors thought they were --\n",
"Please report!", immediate. = TRUE)
ans <- lapply(ncseq, function(i) {
Li <- diag(nrow = nc[i])
Li[lower.tri(Li, diag = TRUE)] <- thl[[i]]
rownames(Li) <- cnms[[i]]
val <- tcrossprod(sc * Li)
stddev <- sqrt(diag(val))
corr <- t(val/stddev)/stddev
diag(corr) <- 1
structure(val, stddev = stddev, correlation = corr)
})
if (is.character(nms)) {
if (anyDuplicated(nms))
nms <- make.names(nms, unique = TRUE)
names(ans) <- nms
}
structure(ans, sc = sc)
}
wmsg <- wgmsg %ORifNotInLme4% function (n, cmp.val, allow.n, msg1 = "", msg2 = "", msg3 = "")
{
if (allow.n) {
unident <- n < cmp.val
cmp <- "<"
rstr <- ""
}
else {
unident <- n <= cmp.val
cmp <- "<="
rstr <- " and the residual variance (or scale parameter)"
}
wstr <- sprintf("%s (=%d) %s %s (=%d)%s; the random-effects parameters%s are probably unidentifiable",
msg1, n, cmp, msg2, cmp.val, msg3, rstr)
list(unident = unident, wstr = wstr)
}
.makeCC <- .makeCC %ORifNotInLme4% function (action, tol, relTol, ...)
{
stopifnot(is.character(action), length(action) == 1)
if (!is.numeric(tol))
stop("must have a numeric 'tol' component")
if (length(tol) != 1 || tol < 0)
stop("'tol' must be number >= 0")
mis.rt <- missing(relTol)
if (!mis.rt && !is.null(relTol))
stopifnot(is.numeric(relTol), length(relTol) == 1, relTol >=
0)
c(list(action = action, tol = tol), if (!mis.rt) list(relTol = relTol),
list(...))
}
isSingular <- isSingular %ORifNotInLme4% function (x, tol = 1e-04)
{
lwr <- getME(x, "lower")
theta <- getME(x, "theta")
any(theta[lwr == 0] < tol)
}
.get.checkingOpts <- .get.checkingOpts %ORifNotInLme4% function (nms)
nms[grepl("^check\\.(?!conv|rankX|scaleX)", nms, perl = TRUE)]
chk.cconv <- chk.cconv %ORifNotInLme4% function (copt, callingFn)
{
cnm <- deparse(substitute(copt))
if (is.character(copt)) {
def <- eval(formals(callingFn)[[cnm]])
def$action <- copt
assign(cnm, def, envir = sys.frame(sys.parent()))
}
else chk.convOpt(copt)
}
namedList <- namedList %ORifNotInLme4% function (...)
{
L <- list(...)
snm <- sapply(substitute(list(...)), deparse)[-1]
if (is.null(nm <- names(L)))
nm <- snm
if (any(nonames <- nm == ""))
nm[nonames] <- snm[nonames]
setNames(L, nm)
}
colSort <- colSort %ORifNotInLme4% function (x)
{
termlev <- vapply(strsplit(x, ":"), length, integer(1))
iterms <- split(x, termlev)
iterms <- sapply(iterms, sort, simplify = FALSE)
ilab <- "(Intercept)"
if (ilab %in% iterms[[1]]) {
iterms[[1]] <- c(ilab, setdiff(iterms[[1]], ilab))
}
unlist(iterms)
}
isAnyArgBar <- isAnyArgBar %ORifNotInLme4% function (term)
{
if ((term[[1]] != as.name("~")) && (term[[1]] != as.name("("))) {
for (i in seq_along(term)) {
if (isBar(term[[i]]))
return(TRUE)
}
}
FALSE
}
chk.convOpt <- chk.convOpt %ORifNotInLme4% function (opt)
{
cnm <- deparse(nm <- substitute(opt))[[1]]
if (!is.list(opt))
stop("check.conv* option ", cnm, " must be a list")
if (!is.character(opt$action))
stop("check.conv* option ", cnm, " has no 'action' string")
if (!is.numeric(tol <- opt$tol))
stop("check.conv* option ", cnm, " must have a numeric 'tol' component")
if (length(tol) != 1 || tol < 0)
stop("check.conv* option ", cnm, "$tol must be number >= 0")
if (!is.null(relTol <- opt$relTol))
stopifnot(is.numeric(relTol), length(relTol) == 1, relTol >=
0)
invisible()
}
getME <- getME %ORifNotInLme4% function (object, name, ...)
UseMethod("getME")
levelfun <- function (x, nl.n, sample_new_levels, Sigma)
{
old_levels <- dimnames(x)[["group"]]
if (!all(nl.n %in% old_levels)) {
nl.n.comb <- union(old_levels, nl.n)
d <- dim(x)
dn <- dimnames(x)
dnn <- names(dn)
# newx: num_predictors x num_levels x num_samples x num_chains
newx <- array(0, c(d[1L], length(nl.n.comb), d[3L:4L]),
dimnames = list(dn[[1L]], nl.n.comb, dn[[3L]], dn[[4L]]))
newx[,old_levels,,] <- x
if (sample_new_levels) {
new_levels <- nl.n.comb[!(nl.n.comb %in% old_levels)]
L <- apply(Sigma, c(3L, 4L), function(x) t(base::chol(x)))
n_predictors <- d[1L]
n_groups <- length(new_levels)
n_samples <- d[3L]
n_chains <- d[4L]
L <- array(Sigma, c(n_predictors, n_predictors, n_samples * n_chains))
L <- bdiag(lapply(seq_len(n_samples * n_chains), function(i) t(base::chol(L[,,i,drop = FALSE]))))
# L: block diagonal where each block is a sample of the covariance
# matrix for the random effects at that level
# (p x p) x (n_samp x n_chain)
u <- matrix(rnorm(n_predictors * n_predictors * n_groups * n_samples * n_chains),
n_predictors * n_predictors * n_samples * n_chains,
n_groups)
# L %*% u: (n_predictors x n_samp x n_chain) x n_groups
newx[,new_levels,,] <- aperm(array(as.vector(Matrix::crossprod(L, u)), c(n_predictors, n_samples, n_chains, n_groups)), c(1L, 4L, 2L, 3L))
}
x <- newx
names(dimnames(x)) <- dnn
old_levels <- dimnames(x)[["group"]]
}
if (!all(r.inn <- old_levels %in% nl.n)) {
x <- x[,r.inn,,,drop = FALSE]
}
x
}
# can't use lme4 version, as we need to use stan4bartFit terms and model.frame
get.orig.levs <- function (object, FUN = levels, newdata = NULL, sparse = FALSE, ...)
{
terms <- terms(object, ...)
frame <- model.frame(object, ...)
isFac <- vapply(frame, is.factor, FUN.VALUE = TRUE)
isFac[attr(terms, "response")] <- FALSE
frame <- frame[isFac]
hasSparse <- any(grepl("sparse", names(formals(FUN))))
orig_levs <- if (any(isFac) && hasSparse) lapply(frame, FUN, sparse = sparse) else if(any(isFac) && !hasSparse) lapply(frame, FUN)
if (!is.null(newdata)) {
for (n in names(frame)) {
orig_levs[[n]] <- c(orig_levs[[n]],
setdiff(unique(as.character(newdata[[n]])), orig_levs[[n]]))
}
}
if (!is.null(orig_levs)) attr(orig_levs, "isFac") <- isFac
orig_levs
}
formula.stan4bartFit <- function(x, ...)
{
dots_list <- list(...)
type <- match.arg(dots_list$type, c("all", "fixed", "random", "bart"))
if (is.null(formula <- x$formula)) {
if (!grepl("stan4bart$", deparse(x$call[[1]])))
stop("can't find formula stored in model frame or call")
form <- as.formula(formula(x$call))
}
if (type == "fixed") {
RHSForm(formula) <- nobart(nobars(RHSForm(formula)))
} else if (type == "random") {
formula <- reOnly(formula, response = TRUE)
} else if (type == "bart") {
RHSForm(formula) <- allbart(nobars(RHSForm(formula)))
}
formula
}
terms.stan4bartFit <- function(x, ...)
{
dots_list <- list(...)
type <- match.arg(dots_list$type, c("all", "fixed", "random", "bart"))
terms <- attr(x$frame, "terms")
if (type == "all")
return(terms)
if (type == "fixed") {
tt <- terms.formula(formula(x, type = "fixed"), data = x$frame)
attr(tt, "predvars") <- attr(terms, "predvars.fixed")
} else if (type == "random") {
tt <- terms.formula(subbars(formula(x, type = "random")), data = x$frame)
attr(tt, "predvars") <- attr(terms, "predvars.random")
} else if (type == "bart") {
tt <- terms.formula(formula(x, type = "bart"), data = x$frame)
pred_vars <- attr(terms, "predvars.bart")
attr(tt, "predvars") <- pred_vars
# bart, when used with a . variable expansion, can end up
# trying to build a frame with the offset variable. This
# doesn't work with the dbarts model frame functions, so we
# strip out anything that isn't in predvars
term_labels <- attr(tt, "term.labels")
extra_terms <- term_labels[term_labels %not_in% as.character(pred_vars)[-1L]]
if (length(extra_terms) > 0L)
tt <- strip_extra_terms(tt, extra_terms)
}
# Possibly re-order terms in case they apepared in a different order in
# the original, which happens if, for example, a variable appears in the
# fixed formula and again in the bart one (or in a . expression)
vnames <- as.character(attr(tt, "variables"))
pvnames <- as.character(attr(tt, "predvars"))
attr(tt, "variables")[vnames %in% pvnames] <- attr(tt, "variables")[match(pvnames, vnames)]
tt
}
model.frame.stan4bartFit <- function(formula, ...)
{
dots_list <- list(...)
type <- match.arg(dots_list$type, c("all", "fixed", "random", "bart"))
frame <- formula$frame
if (type == "fixed") {
# can return a lot of useless junk
frame <- frame[as.character(attr(terms(formula), "predvars.fixed"))[-1L]]
} else if (type == "random") {
frame <- frame[as.character(attr(terms(formula), "predvars.random"))[-1L]]
} else if (type == "bart") {
frame <- frame[attr(terms(formula), "varnames.bart")]
}
frame
}
rm(`%ORifNotInLme4%`)
if (exists("lme4ns")) rm("lme4ns")
make_glmerControl <- function(..., ignore_lhs = FALSE, ignore_x_scale = FALSE) {
glmerControl(check.nlev.gtreq.5 = "ignore",
check.nlev.gtr.1 = "stop",
check.nobs.vs.rankZ = "ignore",
check.nobs.vs.nlev = "ignore",
check.nobs.vs.nRE = "ignore",
check.formula.LHS = if (ignore_lhs) "ignore" else "stop",
check.scaleX = if (ignore_x_scale) "ignore" else "warning",
...)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.