Is it possible to specify a lmer model without any fixed effects?

Question

In R, how do I specify lmer model without global fixed effect? For example, if I say something like

lmer(y ~ (1 | group) + (0 + x | group), data = my_df)

the fitted model will be

$y_{ij} = a + \alpha_{i} + \beta_{i} x_{ij}$

How do I fit model

$y_{ij} = \alpha_{i} + \beta_{i} x_{ij}$?

Answer 1

As @Mike Lawrence mentioned the obvious thing to do when defining a model without fixed effects is something in the form of:

lmer(y ~ -1 + (1|GroupIndicator))

which is actually quite straightforward; one defines no intercept or an X matrix. The basic reason which this doesn't work out is that as @maxTC pointed out "lme4 package is dedicated to mixed models only".

In particular what lmer() fitting does is calculate the profiled deviance by solving the penalized least square regression between the $\hat{y}$ and ${y}$ as well as the spherical random effects $u$ and $0$ (Eq. (11), Ref.(2)). Computationally this optimization procedure computes the Cholesky decomposition of the corresponding system exploiting the system's block structure (Eq. (5), Ref.(1)). Setting no global fixed effects practically distorts that block structure in a way that the code of lmer() can't cope. Among other things the conditional expected value of $u$ is based on $\hat{\beta}$'s but solving for $\hat{\beta}$ asks the solution of a matrix system that never existed (the matrix $R_{XX}$ in Ref.(1), or $L_X$ in Ref.(2)). So you get an error like:

Error in mer_finalize(ans) : 
  Cholmod error 'invalid xtype' at file:../Cholesky/cholmod_solve.c, line 970

cause after all there was nothing to solve for in the first place.

Assuming you don't want to re-write lmer() profiled deviance cost function the easiest solution is based on the CS-101 axiom: garbage in, garbage out.

 N = length(y); Garbage <- rnorm(N); 
 lmer(y ~ -1 + Garbage + (1|GroupIndicator));

So what we do is define a variable $Garbage$ that is just noise; as before lmer() is instructed to use no fixed intercept but only the X matrix defined us (in this case the single column matrix Garbage). This extra Gaussian noise variable will be in expectation uncorrelated to our sample measurement errors as well as with your random effects variance. Needless to say the more structure your model has the smaller the probability of getting unwanted but statistically significant random correlations.

So lmer() has a placebo $X$ variable (matrix) to play with, you in expectation will get the associated $\beta$ to be zero and you didn't have to normalize your data in any way (centring them, whitening them etc.). Probably trying a couple random initialization of the placebo $X$ matrix won't hurt either. A final note for the "Garbage": using Gaussian noise wasn't "accidental"; it has the largest entropy amongst all random variables of equal variance so the least chance of providing an information gain.

Clearly, this is more a computational trick than a solution, but it allows the user to effectively specify an lmer model without global fixed effect. Apologies for hoping around the two references. In general I think Ref.(1) is the best bet for anyone for realize what lmer() is doing, but Ref.(2) is closer to the spirit of the actual code.

Here's a bit of code show-casing the idea above:

library(lme4)
N= 500;                 #Number of Samples
nlevA = 25;             #Number of levels in the random effect 
set.seed(0)             #Set the seed
e = rnorm(N); e = 1*(e - mean(e) )/sd(e); #Some errors

GroupIndicator =  sample(nlevA, N, replace=T)   #Random Nvel Classes 

Q = lmer( rnorm(N) ~ (1| GroupIndicator ));      #Dummy regression to get the matrix Zt easily
Z = t(Q@Zt);            #Z matrix

RA <-  rnorm(nlevA )                        #Random Normal Matrix
gammas =c(3*RA/sd(RA))                      #Colour this a bit

y = as.vector(  Z %*% gammas +  e )         #Our measurements are the sum of measurement error (e) and Group specific variance

lmer_native <- lmer(y ~ -1 +(1| GroupIndicator)) #No luck here.
Garbage <- rnorm(N)                              #Prepare the garbage

lmer_fooled <- lmer(y ~ -1 + Garbage+(1| GroupIndicator)) #OK...
summary(lmer_fooled)                             #Hey, it sort of works!

References:

1. Linear mixed models and penalized least squares by D.M. Bates and S. DebRoy, Journal of Multivariate Analysis, Volume 91 Issue 1, October 2004 (Link to preprint)
2. Computational methods for mixed models by Douglas Bates, June 2012. (Link to source )