lmer giving worse performance than lm

Question

I am having trouble training a model for nested data about house prices. Lets say my data looks like following:

  logPrice bedCount bathCount                city
 0.6517920        4       2-3        Redwood City
 0.4402192        1       1-2 South San Francisco
 0.5922396        2       1-2           San Mateo
 0.4606918        3       1-2 South San Francisco
 0.7592523    5plus       3-4           San Mateo
 0.4710397        1       1-2        Redwood City

bedCount, bathCount and city are factors.

As a baseline, I trained a simple linear model ignoreing nested structure of the data (houses are nested within cities).

lm.model = lm(formula = logPrice ~ 1 + bedCount + bathCount + city)

which corresponds to following assumption:

logPrice$_i = \beta_0 + \beta_1\cdot$ bedCount$_i + \beta_2\cdot$ bathCount$_i + \beta_{3,j[i]}\cdot I$(city$_{j[i]}) + \epsilon_i$

where

$\epsilon_i \sim N(0, \sigma^2_{logPrice})$ and $I$(city$_{j[i]}$) is the indicator function for city of the $i^{th}$ house (which is 1).

Now, I trained a 2-level hierarchical model:

lmer.model = lmer(formula = logPrice ~ 1 + bedCount + bathCount + (1 | city))

which corresponds to the following assumption:

logPrice$_i = \beta_0 + \beta_1\cdot$ bedCount$_i + \beta_2\cdot$ bathCount$_i + \beta_{3,j[i]}\cdot I$(city$_{j[i]}) + \epsilon_i$

where $\epsilon_i \sim N(0, \sigma^2_{logPrice})$ and $\beta_{3,j} \sim N(0, \sigma^2_{\beta_3})$

Now, on the training data, lm.model gives me lesser average RMSE than lmer.model which shouldn't happen because linear regression is a special case of multilevel linear regression (I didn't care to check average RMSE on test data because that on training data itself should be lower for 2nd model than that for 1st model). In fact, my data has multiple levels (houses nested within subdivisions, which are nested within zipcodes, which in turn are nested within cities) and the performance gets worse as I add more and more levels to the model (i.e. model with random effect (1 | subdivision) does worse than that with random effect (1 | zipcode) + (1 | zipcode:subdivision), which in turn does worse than a model with random effect (1 | city) + (1 | city:zipcode) + (1 | city:zipcode:subdivision)).

What am I missing?

Answer 1

I just saw this while looking for the answer to another question, and thought an answer to this might be informative for others, even if it is a little old... I'm not as expert on this as some people on here, so please pull me up if I got this wrong.

What the random effect does, compared with the fixed effect, is to impose further assumptions on how the means of the factor levels will be distributed (specifically, that they are normally distributed). So, whereas the fixed effect is free to just go straight for the data, the random effect has certain constraints on it. What this means in practice, is that the fixed effect will (I think?) always be closer to your training data than the random effect (because it can just go and sit right on the mean of the values for that level). But, the random effect should (IF our assumptions are correct) capture the underlying mechanics better, and so (as if by magic) should be closer to your test data than a fixed effect would have been.

See my answer here for an illustration of an extreme case, where it is clear that the random effect comes out miles away from the training data, but should be expected (on average) to be closer to future examples.

This is really very much like any other overfitting example. The more assumptions we can make, the more constrained the model is. Which means it can't get as close to the training data, but (IF our assumptions are good) should get closer to future examples. Much like with a linear vs quadratic regression: the quadratic model can get closer to the data, but IF the true relationship is linear, the linear model will be a better predictor of future cases.