I am trying to build a mixed model with random effect nested in fixed effect. Below is an example dataset.

example_data <- as.data.frame(list(  quality = c(4,4,2,3,3,4,1,3,2,3,2,1),cell.count = c(24,27,23,23,23,24,24, 24, 25,24, 22, 23),species = c("Pl","Pl" ,"Pl","Pl","Pl","Pl", "Pl","Pl","Pl","Pm","Pm","Pm"),mouse = c( "mouse1","mouse1","mouse1",
         "mouse15","mouse15","mouse15"), location = c("barn","barn","barn",  "basement","basement","basement","lab","lab", "lab", "lab", "lab", "lab"), month = c("aug", "aug","aug", "jun","jun","jun", "sep","sep","sep",  "sep","sep","sep"),  status = c("wild","wild","wild", "wild","wild","wild", "lab","lab","lab", "lab","lab","lab")
      )  )

I am modeling the dependant variable ‘cell.count’, there are multiple cell observations from a single mouse, which is coded with a unique ID. ‘quality’ is a score assigned to each cell. For the original data set

  • cell count, 517 observations
  • mouse, 20 levels
  • Species, 2 levels
  • location, 6 levels
  • month, 5 levels
  • status 2 levels

I’d like to build several mixed models from different subsets of the full dataset, that would be used to estimate particular effects separately for example.

m1 <- lmer(cell.count ~ species + quality + (1 | species/mouse), data = example_data[example_data$status == "lab",])
m2 <- lmer(cell.count ~ status + quality + (1 | status/mouse), data = example_data[example_data$species == "Pl",])

In the summary() output for m2, status (a fixed effect) is listed within the random effects section, The variance is estimated for mouse (within status), but variance for status is also estimated. The 'status' effect is estimated as a fixed effect as well, but this output indicates that there are multiple random effects. This complicates how we plan to use a LRT for the models (which require 1 random effect per model).

Random effects:
Groups       Name        Variance Std.Dev.
mouse:status (Intercept) 0.4220   0.650   
status       (Intercept) 0.0315   0.177   
Residual                 1.4001   1.183   
Number of obs: 9, groups:  mouse:status, 3; status, 2

My question is, Did I use the wrong syntax for coding a nested random effect? Or am I misunderstanding what a nested effect is?

I would like to use exactRLRT() to test the random effects for these models, but which these models, this function doesn't work since there are more than 1 random effect.

  • 1
    $\begingroup$ A factor should either be fixed, or random; not both. This should explain nested factors. $\endgroup$ – Robert Long Apr 17 at 9:52
  • $\begingroup$ Also, with only 2 levels each, status and Species should not be random. $\endgroup$ – Robert Long Apr 17 at 9:59
  • 1
    $\begingroup$ Thanks for linking to that answer. From my understanding the correct model in my case would be m0 <- lmer(cell.count ~ status + quality + (1 | mouse), data = ). Since the nesting is encoded in the data (unique mouse ID's), I do not need to explicitly tell lmer to nest the random factor 'mouse' within the fixed factors. $\endgroup$ – user1757654 Apr 17 at 14:38
  • $\begingroup$ Yes, that seems more appropriate $\endgroup$ – Robert Long Apr 17 at 19:03

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.