I have a model in the form of:

y ~ A + B + random=~1|C + random= ~1|D, family=poisson, data=data)

I recently saw a thread on grokbase that suggested;

lme and, by extension, glmmPQL do not handle crossed random effects

easily.

You must create a factor of the same length as y, A, B, C, and D with

a single level

const = factor(rep(1, length(y)))

then use the non-obvious formulation

glmmPQL(y ~ A + B, random = list(const = pdBlocked(pdIdent(~ C - 1),

pdIdent(~ D - 1))))

When I run this script however – I get the output :

Error in pdConstruct.pdBlocked(object, form = form, nam = nam, data = data, :

'form' must be a list

Does anyone have suggestions on how to get this script to work for MASS (glmmPQL) or lme4?

Or any other approaches?