spaMM (version 4.6.1)

aliases: Variable aliases for multivariate-response fits

Description

Consider the toy_iris dataset defined in the Examples, noting in particular

u <- rnorm(10,sd=1)
id_y <- gl(10,5)
id_b <- rep(seq(10),5)

an the subsequent use of u[id_y] and u[id_b] in the expectations of Poisson draws. Suppose that we try to a multivariate-response model with two Poisson responses to it (yes, it is an artificial example, as a multivariate fit may not be necessary here, but if you know what pois4mlogit does, then a less artificial example can be defined from a similar dataset with binomial samples instead of Poisson ones). A fit by

fitmv(submodels = list(
    list(yellow ~ 1+(1|id_y), family = poisson()),
    list(blue ~ 1+(1|id_b), family = poisson())),
    data = toy_iris)

does not match the data-generating algorithm, because the fit of the (1|id_y) and (1|id_b) random effects does not take into account that the latent values of these random effects are sampled from a single vector u of 10 latent values. A matching fit should fit a single random-effect variance and produce a single predicted vector for u, rather than two distinct vectors of predicted values.

Conversely, a fit such as

fitmv(submodels = list(
    list(yellow ~ 1+(1|id_y), family = poisson()),
    list(blue ~ 1+(1|id_y), family = poisson())),
    data = toy_iris)

would fit a single random-effect variance and produce a single predicted vector, but the factor indices are incorrect in the second submodel.

We need a syntax such that a single variance and a single vector are fitted (as when the random-effect terms are identical in the two formulas), but where the factor indices are effectively id_y and id_b in the two submodels. The syntax of the proper fit in the Examples achieves this effect. The random effect is specified as (1|id) where id is *not* an actual factor in the data but is an alias for the variable id_y in the first formula and id_b in the second one. This interpretation of the (1|id) term is specified by the argument aliases=list(id=c("id_y","id_b")).

Examples

# Toy data
set.seed(123)
ssize <- 50L
u <- rnorm(10,sd=1)
id_y <- gl(10,5)
id_b <- rep(seq(10),5)
yellow <- rpois(ssize, lambda=exp(1+u[id_y]))
blue <- rpois(ssize, lambda=exp(1+u[id_b]))
toy_iris <- data.frame(
  id_y=id_y, id_b=id_b, yellow=yellow, blue=blue
)

# Fit
proper <- fitmv(submodels = list(
    list(yellow ~ 1+(1|id), family = poisson()),
    list(blue ~ 1+(1|id), family = poisson())),
    data = toy_iris, aliases=list(id=c("id_y","id_b")))
    
ranef(proper) # single vector of 10 values

# Rows 1-50 and 51-100 of Z matrix show 
# the distinct design of the two submodels:
get_matrix(proper,"ZA")