Similar questions have been posed here https://stats.stackexchange.com/questions/13000/random-effect-nested-under-fixed-effect-model-in-r and here https://stats.stackexchange.com/questions/190120/nesting-random-effect-within-fixed-effect-using-lmer-of-lme4-in-r however, I have trouble applying those answers to my particular problem.
The Data
Plots were established in 23 different forest locations. Each plot has 4 sub plots. The response I'm interested in is the number of sub plots that contain a successful tree. To avoid pseudo-replication, this is modeled at the plot level as a binomial with x successes (x can be an integer [0:4]) in 4 trials. I have two fixed effects that I am interested in: Fencing and average seedling size. Fencing is a stand-level variable, and avg. seedling size is measured at the plot level. Because multiple plots were nested within the same location, my comittee members want me to include location as a random effect to account for lack of complete independence between plots.
The following code simulates my data:
location<-as.character(LETTERS[rep(seq(from=1,to=23),30)])
plot<-rep(seq(from=1,to=30),23)
success<-sample(0:4,690,replace=T)
n_subplots<-rep(4,690)
df<-data.frame(location,plot,success,n_subplots)
df$FENCED<-ifelse(df$location %in% c("A","B","C","D","E","F","G"),
1,0)
df$seedling_size<-df$success*.05+rnorm(1,0,2.5)+.1*df$FENCED
df$success<-ifelse(df$success>2 & df$FENCED ==1, df$success-2,df$success)
df$success<-ifelse(df$location %in% c("A","J","L"),df$success-1,df$success)
df$success<-ifelse(df$success<0,0,df$success)
df<-df[order(df$location,df$plot),]
The model
I want to test the fencing effect by itself as well as a possible interaction between fencing and seedling size (fence protects small seedlings from deer). This is the model I am using:
model<-glmer(cbind(success,n_subplots)~FENCED*seedling_size +
(FENCED|location),
family=binomial(),data=df,nAGQ=1)
summary(model)
Summary Results
Random effects:
Groups Name Variance Std.Dev. Corr
location (Intercept) 6.412e-08 0.0002532
FENCED 2.714e-02 0.1647464 0.30
Number of obs: 690, groups: location, 23
Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 13.9404 0.8173 17.057 < 2e-16 ***
FENCED -7.1809 1.4386 -4.991 5.99e-07 ***
seedling_size 12.9083 0.7257 17.787 < 2e-16 ***
FENCED:seedling_size -5.2464 1.3464 -3.897 9.76e-05 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) FENCED sdlng_
FENCED -0.568
seedling_sz 0.999 -0.567
FENCED:sdl_ -0.538 0.996 -0.539
convergence code: 0
Model failed to converge with max|grad| = 0.00261073 (tol = 0.001, component 1)
The Questions
Is this the correct syntax to specify the model I want? How to interpret the results? Fenced appears in both the random summary and the fixed effects. Finally, What might be causing the convergence warning, and how to address it?
At least for convergence issues, I did run the model on all the different optimizers as suggested here:https://www.rdocumentation.org/packages/lme4/versions/1.1-17/topics/convergence The effect estimates and significances were similar between optimizers, and all but 1 were "OK" or "worked". Whatever that means specifically to me is unclear, other than it is a not a problem with the method of approximation. However, doing something like the "Richardson Extrapolation" is a bit over my head.
