0

I have a simple experimental design. Each experiment consists of 6 treatments organized in 6 blocks with each treatment replicated once per block. The whole experiment was repeated 6 times. The outcome variable is an insect count for each subject. During each trial a certain amount of insects was released per block and then the number of insects caught in each treatment was counted. There where trials with 50 and 100 insects per block, so this number was used as offset variable. In my dataset the blocks have the same coding for each trial, however their not the same. The blocks are nested in trial. I'm aiming to fit a poisson model with offset, but for simplification I fitted a gaussian LMM without offset and intercept only.

My problem: My model has singularity, i.e. the variance of "block" is estimated to be zero. However, only if I fit the model with the lme4 package, not with the nlme package!

I can't see, why this is. The model is fairly simple. The data does not contain many zeroes.

Here is my data:

example <- structure(list(count = c(19L, 23L, 10L, 25L, 12L, 1L, 12L, 9L, 
31L, 7L, 9L, 6L, 10L, 48L, 4L, 5L, 6L, 7L, 6L, 19L, 44L, 2L, 
12L, 2L, 13L, 6L, 35L, 16L, 5L, 4L, 9L, 18L, 34L, 15L, 3L, 3L, 
5L, 15L, 40L, 13L, 14L, 1L, 9L, 15L, 40L, 3L, 6L, 9L, 12L, 47L, 
6L, 6L, 3L, 4L, 7L, 12L, 46L, 8L, 5L, 4L, 6L, 4L, 42L, 13L, 13L, 
7L, 6L, 22L, 41L, 13L, 1L, 1L, 3L, 23L, 17L, 26L, 6L, 0L, 20L, 
30L, 17L, 1L, 9L, 4L, 9L, 30L, 1L, 4L, 20L, 1L, 4L, 13L, 41L, 
4L, 12L, 4L, 23L, 1L, 32L, 9L, 7L, 11L, 34L, 14L, 13L, 2L, 2L, 
10L, 1L, 7L, 8L, 12L, 4L, 0L, 2L, 6L, 9L, 0L, 4L, 4L, 8L, 12L, 
0L, 2L, 15L, 8L, 1L, 2L, 15L, 1L, 12L, 4L, 4L, 1L, 5L, 18L, 11L, 
2L, 8L, 3L, 15L, 6L, 2L, 3L, 12L, 23L, 26L, 18L, 12L, 1L, 25L, 
25L, 31L, 3L, 8L, 8L, 11L, 38L, 7L, 9L, 9L, 12L, 6L, 10L, 40L, 
5L, 14L, 1L, 8L, 3L, 38L, 28L, 6L, 10L, 20L, 15L, 22L, 14L, 7L, 
4L, 4L, 18L, 26L, 15L, 22L, 0L, 11L, 20L, 36L, 8L, 5L, 6L, 13L, 
49L, 12L, 3L, 10L, 11L, 5L, 7L, 67L, 6L, 2L, 4L, 17L, 5L, 40L, 
27L, 5L, 6L, 17L, 20L, 25L, 17L, 6L, 4L), trial = c(1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 
5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 
5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 
5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 
6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 
6L, 6L, 6L, 6L, 6L), block = structure(c(1L, 1L, 1L, 1L, 1L, 
1L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 
4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 1L, 
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 
3L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 
6L, 6L, 6L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 
3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 
5L, 6L, 6L, 6L, 6L, 6L, 6L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 
2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 
5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 1L, 1L, 1L, 1L, 1L, 
1L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 
4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 1L, 
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 
3L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 
6L, 6L, 6L), levels = c("A", "B", "C", "D", "E", "F"), class = "factor"), 
    treat = structure(c(1L, 2L, 5L, 3L, 4L, 6L, 2L, 3L, 1L, 6L, 
    5L, 4L, 4L, 5L, 6L, 1L, 3L, 2L, 5L, 1L, 4L, 6L, 3L, 2L, 2L, 
    6L, 3L, 5L, 4L, 1L, 2L, 4L, 3L, 6L, 5L, 1L, 1L, 2L, 5L, 3L, 
    4L, 6L, 2L, 3L, 1L, 6L, 5L, 4L, 4L, 5L, 6L, 1L, 3L, 2L, 5L, 
    1L, 4L, 6L, 3L, 2L, 2L, 6L, 3L, 5L, 4L, 1L, 2L, 4L, 3L, 6L, 
    5L, 1L, 1L, 2L, 5L, 3L, 4L, 6L, 2L, 3L, 1L, 6L, 5L, 4L, 4L, 
    5L, 6L, 1L, 3L, 2L, 5L, 1L, 4L, 6L, 3L, 2L, 2L, 6L, 3L, 5L, 
    4L, 1L, 2L, 4L, 3L, 6L, 5L, 1L, 1L, 2L, 5L, 3L, 4L, 6L, 2L, 
    3L, 1L, 6L, 5L, 4L, 4L, 5L, 6L, 1L, 3L, 2L, 5L, 1L, 4L, 6L, 
    3L, 2L, 2L, 6L, 3L, 5L, 4L, 1L, 2L, 4L, 3L, 6L, 5L, 1L, 1L, 
    2L, 5L, 3L, 4L, 6L, 2L, 3L, 1L, 6L, 5L, 4L, 4L, 5L, 6L, 1L, 
    3L, 2L, 5L, 1L, 4L, 6L, 3L, 2L, 2L, 6L, 3L, 5L, 4L, 1L, 2L, 
    4L, 3L, 6L, 5L, 1L, 1L, 2L, 5L, 3L, 4L, 6L, 2L, 3L, 1L, 6L, 
    5L, 4L, 4L, 5L, 6L, 1L, 3L, 2L, 5L, 1L, 4L, 6L, 3L, 2L, 2L, 
    6L, 3L, 5L, 4L, 1L, 2L, 4L, 3L, 6L, 5L, 1L), levels = c("P", 
    "BH", "GWC", "CC", "WF_u", "Control"), class = "factor"), 
    offset = c(100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 
    50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 
    50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 50L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
    100L, 100L)), row.names = c(NA, -216L), class = c("tbl_df", 
"tbl", "data.frame"))

Here is my code:

# The final model I have in mind
glmer1 <- glmer(count ~ treat+ (1|trial/block) + offset(log(offset)), family = poisson, data = example); summary(glmer1) # nested ranef

# Simplified models to show my issues
lmer1 <- lme(count ~ 1, random = ~1|trial/block, data = example); summary(lmer1) # no singularity
lmer2 <- lmer(count ~ 1 + (1|trial/block),  data = example); summary(lmer2) # singularity

Any advice why I've got singularity and how I should deal with it?

Pharcyde
  • 397
  • 1
  • 3
  • 14

0 Answers0