fit_with() is a pipe-friendly tool that applies a list of formulas to a fitting function such as stats::lm(). The list of formulas is typically created with formulas().

fit_with(data, .f, .formulas, ...)

Arguments

data

A dataset used to fit the models.

.f

A fitting function such as stats::lm(), lme4::lmer() or rstanarm::stan_glmer().

.formulas

A list of formulas specifying a model.

...

Additional arguments passed on to .f

Details

Assumes that .f takes the formula either as first argument or as second argument if the first argument is data. Most fitting functions should fit these requirements.

See also

Examples

# fit_with() is typically used with formulas(). disp_fits <- mtcars %>% fit_with(lm, formulas(~disp, additive = ~drat + cyl, interaction = ~drat * cyl, full = add_predictors(interaction, ~am, ~vs) )) # The list of fitted models is named after the names of the list of # formulas: disp_fits$full
#> #> Call: #> .f(formula = disp ~ drat * cyl + am + vs, data = data) #> #> Coefficients: #> (Intercept) drat cyl am vs drat:cyl #> -164.83 27.58 75.84 -38.22 -15.53 -6.97 #>
# Additional arguments are passed on to .f mtcars %>% fit_with(glm, list(am ~ disp), family = binomial)
#> [[1]] #> #> Call: .f(formula = am ~ disp, family = function (link = "logit") #> { #> linktemp <- substitute(link) #> if (!is.character(linktemp)) #> linktemp <- deparse(linktemp) #> okLinks <- c("logit", "probit", "cloglog", "cauchit", "log") #> if (linktemp %in% okLinks) #> stats <- make.link(linktemp) #> else if (is.character(link)) { #> stats <- make.link(link) #> linktemp <- link #> } #> else { #> if (inherits(link, "link-glm")) { #> stats <- link #> if (!is.null(stats$name)) #> linktemp <- stats$name #> } #> else { #> stop(gettextf("link \"%s\" not available for binomial family; available links are %s", #> linktemp, paste(sQuote(okLinks), collapse = ", ")), #> domain = NA) #> } #> } #> variance <- function(mu) mu * (1 - mu) #> validmu <- function(mu) all(is.finite(mu)) && all(mu > 0 & #> mu < 1) #> dev.resids <- function(y, mu, wt) .Call(C_binomial_dev_resids, #> y, mu, wt) #> aic <- function(y, n, mu, wt, dev) { #> m <- if (any(n > 1)) #> n #> else wt #> -2 * sum(ifelse(m > 0, (wt/m), 0) * dbinom(round(m * #> y), round(m), mu, log = TRUE)) #> } #> initialize <- expression({ #> if (NCOL(y) == 1) { #> if (is.factor(y)) y <- y != levels(y)[1L] #> n <- rep.int(1, nobs) #> y[weights == 0] <- 0 #> if (any(y < 0 | y > 1)) stop("y values must be 0 <= y <= 1") #> mustart <- (weights * y + 0.5)/(weights + 1) #> m <- weights * y #> if (any(abs(m - round(m)) > 0.001)) warning("non-integer #successes in a binomial glm!") #> } else if (NCOL(y) == 2) { #> if (any(abs(y - round(y)) > 0.001)) warning("non-integer counts in a binomial glm!") #> n <- y[, 1] + y[, 2] #> y <- ifelse(n == 0, 0, y[, 1]/n) #> weights <- weights * n #> mustart <- (n * y + 0.5)/(n + 1) #> } else stop("for the 'binomial' family, y must be a vector of 0 and 1's\nor a 2 column matrix where col 1 is no. successes and col 2 is no. failures") #> }) #> simfun <- function(object, nsim) { #> ftd <- fitted(object) #> n <- length(ftd) #> ntot <- n * nsim #> wts <- object$prior.weights #> if (any(wts%%1 != 0)) #> stop("cannot simulate from non-integer prior.weights") #> if (!is.null(m <- object$model)) { #> y <- model.response(m) #> if (is.factor(y)) { #> yy <- factor(1 + rbinom(ntot, size = 1, prob = ftd), #> labels = levels(y)) #> split(yy, rep(seq_len(nsim), each = n)) #> } #> else if (is.matrix(y) && ncol(y) == 2) { #> yy <- vector("list", nsim) #> for (i in seq_len(nsim)) { #> Y <- rbinom(n, size = wts, prob = ftd) #> YY <- cbind(Y, wts - Y) #> colnames(YY) <- colnames(y) #> yy[[i]] <- YY #> } #> yy #> } #> else rbinom(ntot, size = wts, prob = ftd)/wts #> } #> else rbinom(ntot, size = wts, prob = ftd)/wts #> } #> structure(list(family = "binomial", link = linktemp, linkfun = stats$linkfun, #> linkinv = stats$linkinv, variance = variance, dev.resids = dev.resids, #> aic = aic, mu.eta = stats$mu.eta, initialize = initialize, #> validmu = validmu, valideta = stats$valideta, simulate = simfun), #> class = "family") #> }, data = data) #> #> Coefficients: #> (Intercept) disp #> 2.6308 -0.0146 #> #> Degrees of Freedom: 31 Total (i.e. Null); 30 Residual #> Null Deviance: 43.23 #> Residual Deviance: 29.73 AIC: 33.73 #>