How to do a generalized linear model with multiple dependent variables in R?

Question

I have six dependent variables (count data) and several independent variables, I see that in a MMR the script goes like this:

my.model <- lm(cbind(DV1,DV2,DV3,DV4,DV5,DV6) ~ IV1 + IV2 + ... + IVn)

But, since my data are counts, I want to use a generalized linear model and I tried this:

my.model <- glm(cbind(DV1,DV2,DV3,DV4,DV5,DV6) ~ IV1 + IV2 + ... + IVn, family="poisson")

And appears this error message:

Error in glm.fit(x = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,  : 
  (subscript) logical subscript too long`

Can anyone explain me this error message or a way to solve my problem?

Answer 1

The short answer is that glm doesn't work like that. The lm will create mlm objects if you give it a matrix, but this is not widely supported in the generics and anyway couldn't easily generalize to glm because users need to be able to specify dual column dependent variables for logistic regression models.

The solution is to fit the models separately. Assume your IVs and DVs live in a data.frame called dd and are labelled the way they are in your question. The following code makes a list of fitted models indexed by the name of the dependent variable they use:

models <- list()
dvnames <- paste("DV", 1:6, sep='')
ivnames <- paste("IV", 1:n, sep='') ## for some value of n

for (y in dvnames){
  form <- formula(paste(y, "~", ivnames))
  models[[y]] <- glm(form, data=dd, family='poisson') 
}

To examine the results, just wrap your usual functions in a lapply, like this:

lapply(models, summary) ## summarize each model

There are no doubt more elegant ways to do this in R, but that should work.

Answer 2

I was told Multivariate Generalized Linear (Mixed) Models exists that address your problem. I'm not an expert about it, but I would have a look to SABRE documentation and this book on multivariate GLMs. Maybe they help...