mumm (version 0.2.1)

mumm: Fit Multiplicative Mixed Models with TMB

Description

Fit a multiplicative mixed-effects model to data with use of the Template Model Builder.

Usage

mumm(formula, data, cor = TRUE, start = c(), control = list())

Arguments

formula

a two-sided formula object describing the linear fixed-effects and random-effects part together with the multiplicative part. The response is on the left of a ~ operator and the terms which are separated by + operators are on the right. The random-effect terms are recognized by vertical bars "|", separating an expression for a model matrix and a grouping factor. The syntax for the multiplicative term is 'mp("random effect","fixed effect")'.

data

a data frame containing the variables in the formula.

cor

logical. If FALSE the random effect in the multiplicative term is assumed to be independent of the corresponding random main effect.

start

a numeric vector of starting values for the parameters in the model.

control

a list of control parameters passed on to the nlminb function used for the optimization.

Value

An object of class mumm.

Details

Fit a multiplicative mixed model via maximum likelihood with use of the Template Model Builder. A multiplicative mixed model is here considered as a model with a linear mixed model part and one multiplicative term. A multiplicative term is here defined as a product of a random effect and a fixed effect, i.e. a term that models a part of the interaction as a random coefficient model based on linear regression on a fixed main effect.

Examples

# NOT RUN {
set.seed(100)
sigma_e <- 1.5
sigma_a <- 0.8
sigma_b <- 0.5
sigma_d <- 0.7
nu <- c(8.2, 6.2, 2.3, 10.4, 7.5, 1.9)
nA <- 15
nP <- 6
nR <- 5
a <- rnorm(nA, mean = 0, sd = sigma_a)
b <- rnorm(nA, mean = 0, sd = sigma_b)
d <- rnorm(nA*nP, mean = 0, sd = sigma_d)
e <- rnorm(nA*nP*nR, mean = 0, sd = sigma_e)
Assessor <- factor(rep(seq(1,nA),each = (nP*nR)))
Product <- factor(rep(rep(seq(1,nP),each = nR), nA))
AssessorProduct <- (Assessor:Product)
y <- nu[Product] + a[Assessor] + b[Assessor]*(nu[Product]-mean(nu)) + d[AssessorProduct] + e
sim_data <- data.frame(y, Assessor, Product)
fit <- mumm(y ~ 1 + Product + (1|Assessor) + (1|Assessor:Product) +
             mp(Assessor,Product) ,data = sim_data)

# }