Welcome to OStack Knowledge Sharing Community for programmer and developer-Open, Learning and Share
Welcome To Ask or Share your Answers For Others

Categories

0 votes
290 views
in Technique[技术] by (71.8m points)

lme4 - problem with bootMer CI: upper and lower limits are identical

I'm having the hardest time generating confidence intervals for my glmer poisson model. After following several very helpful tutorials (such as https://drewtyre.rbind.io/classes/nres803/week_12/lab_12/) as well as stackoverflow posts, I keep getting very strange results, i.e. the upper and lower limits of the CI are identical.

Here is a reproducible example containing a response variable called "production," a fixed effect called "Treatment_Num" and a random effect called "Genotype":

df1 <- data.frame(production=c(15,12,10,9,6,8,9,5,3,3,2,1,0,0,0,0), Treatment_Num=c(1,1,1,1,2,2,2,2,3,3,3,3,4,4,4,4), Genotype=c(1,1,2,2,1,1,2,2,1,1,2,2,1,1,2,2))

#run the glmer model

df1_glmer <- glmer(production ~ Treatment_Num +(1|Genotype),
                     data = df1, family = poisson(link = "log"))

#make an empty data set to predict from, that contains the explanatory variables but no response

require(magrittr)
df_empty <- df1 %>%
  tidyr::expand(Treatment_Num, Genotype)

#create new column containing predictions

df_empty$PopPred <- predict(df1_glmer, newdata = df_empty, type="response",re.form = ~0)

#function for bootMer

myFunc_df1_glmer <- function(mm) {
  predict(df1_glmer, newdata = df_empty, type="response",re.form=~0)
}

#run bootMer

require(lme4)
    merBoot_df1_glmer <- bootMer(df1_glmer, myFunc_df1_glmer, nsim = 10)

#get confidence intervals out of it

predCL <- t(apply(merBoot_df1_glmer$t, MARGIN = 2, FUN = quantile, probs = c(0.025, 0.975)))

#enter lower and upper limits of confidence interval into df_empty

df_empty$lci <- predCL[, 1]
df_empty$uci <- predCL[, 2]

#when viewing df_empty the problem becomes clear: the lci and uci are identical!

df_empty

df_empty

Any insights you can give me will be much appreciated!


与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome To Ask or Share your Answers For Others

1 Answer

0 votes
by (71.8m points)

Ignore my comment!

The issue is with the function you created to pass to bootMer(). You wrote:

myFunc_df1_glmer <- function(mm) {
  predict(df1_glmer, newdata = df_empty, type="response",re.form=~0)
}

The argument mm should be a fitted model object derived from the bootstrapped data. However, you don't pass this object to predict(), but rather the original model object. If you change the function to:

myFunc_df1_glmer <- function(mm) {
  predict(mm, newdata = df_empty, type="response",re.form=~0)
         #^^ pass in the object created by bootMer 
}

then it works:

> df_empty

# A tibble: 8 x 5

  Treatment_Num Genotype PopPred   lci   uci
          <dbl>    <dbl>   <dbl> <dbl> <dbl>
1             1        1  12.9   9.63  15.7 
2             1        2  12.9   9.63  15.7 
3             2        1   5.09  3.87   5.89
4             2        2   5.09  3.87   5.89
5             3        1   2.01  1.20   2.46
6             3        2   2.01  1.20   2.46
7             4        1   0.796 0.361  1.14
8             4        2   0.796 0.361  1.14

As an aside -- how many genotypes in your actual data? If less than 5-7 you might do better using a straight up glm() with genotype as a factor using sum-to-zero contrasts.


与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome to OStack Knowledge Sharing Community for programmer and developer-Open, Learning and Share
Click Here to Ask a Question

...