Kodiologist - 1 month ago 30

R Question

I'm using R's

`lme4`

I tried the following code:

`library(lme4)`

ilogit = function(v)

1 / (1 + exp(-v))

set.seed(20)

n = 100

x1 = rnorm(n)

x2 = rnorm(n)

x3 = rnorm(n)

x4 = rnorm(n)

x5 = rnorm(n)

y.p = ilogit(.5 + x1 - x2)

y = rbinom(n = n, size = 1, prob = y.p)

m1 = glm(

y ~ x1 + x2 + x3 + x4 + x5,

family = binomial)

print(round(d = 2, unname(coef(m1))))

m2 = glmer(

y ~ ((x1 + x2 + x3 + x4 + x5)|1),

family = binomial)

print(round(d = 2, unname(coef(m2))))

This yields:

`Loading required package: Matrix`

[1] 0.66 1.14 -0.78 -0.01 -0.16 0.25

Error: (p <- ncol(X)) == ncol(Y) is not TRUE

Execution halted

What did I do wrong? What's the right way to do this?

Answer

Looks like lme4 can't do this as-is. Here's what @amoeba said in stats.SE chat:

What Kodi wants to do is definitely a mixed model, in the sense of Bates et al. see e.g. eq (2) here https://cran.r-project.org/web/packages/lme4/vignettes/lmer.pdf As far as I can see, X and Z design matrices are equal in this case. However, there is no way one can use lme4 to fit this (without hacking into the code): it allows only particular Z matrices that arise from the model formulas of the type (formula|factor).

See https://stat.ethz.ch/pipermail/r-sig-mixed-models/2011q1/015581.html "We intend to allow lmer to be able to use more flexible model matrices for the random effects although, at present, that requires a certain amount of tweaking on the part of the user"

And https://stat.ethz.ch/pipermail/r-sig-mixed-models/2009q2/002351.html "I view the variance-covariance structures available in the lme4 package as being related to random-effects terms in the model matrix. A random-effects term is of the form (LMexpr | GrpFac). The expression on the right of the vertical bar is evaluated as a factor, which I call the grouping factor. The expression on the left is evaluated as a linear model expression."

That's all quotes from Bates. He does say "In future versions of lme4 I plan to allow for extensions of the unconditional variance-covariance structures." (in 2009) but I don't this was implemented.