I have created an example data set here
My data represent the presence/absence of a particular animal species (data$outcome) and measurements of trees. I would like to investigate how tree species (data$var1), tree diameter (data$var2) and tree height (data$var3) influence the presence/absence of my animal (data$outcome). Due to the study design I also have two nested grouping variables/nested random effect, habitat (data$id1) and individual animal nested within habitat (data$id2), and a weighting variable which represents the relative abundance of a particular tree species within a particular habitat. This weighting variable is important to be able to determine that a particular tree species was more or less important than another even after taking their relative abundance in the landscape into account. I have used binomial GLMMs in lem4 to try and investigate my questions.
I can run a model which includes tree diameter and hight as explanatory variables. However, when I include tree species as an explanatory variable in my model I get model convergence problems. This is not surprising as tree species is a factor with 32 levels and some levels are only associated with absence records or only associated with presence records. Hence the inclusion of tree species in my model is likely introducing two problems; 1) it is making my models very complex and difficult to estimate; and 2) it introduces separation into my models.
If I remove all levels of tree species where there are <=1 occurrence of a particular tree species amongst the presence/absence records and then re-run my model (model3) I still have convergence issues, however these can be solved by using a different optimiser (model4).
My question is: After running model2 and finding convergence problems, likely due to model complexity and separation issues, what are my options going forward/how can I proceed? Is it possible to not filter out levels of tree species, like I have done, but instead do something else. I am very interested in the effect of tree species on the presence/absence of my animal so would prefer not to remove levels of tree species if possible.
There is a good troubleshooting in lme4 page here
# Load library
library(lme4); library(AICcmodavg)
# Load data
data <- read.csv("data.csv", header = T)
# Classify variable class
data$id1 <- factor(data$id1)
data$id2 <- factor(data$id2)
data$var1 <- factor(data$var1)
data$var2 <- as.numeric(data$var2)
data$var3 <- as.numeric(data$var3)
data$var4 <- as.numeric(data$var4)
data$weights <- as.numeric(data$weights)
data$outcome <- as.numeric(data$outcome)
# Scale variables
data$s.var2 <- scale(data$var2, scale = T, center = F)
data$s.var3 <- scale(data$var3, scale = T, center = F)
# Model 1
model1 <- glmer(outcome ~ s.var3 + s.var2 + (1 | id2/id1), data = data, family = binomial, weights = weights)
Warning message:
In eval(family$initialize, rho) : non-integer #successes in a binomial glm!
AICc(model1)
[1] 23588.41
# Model 2
model2 <- glmer(outcome ~ var1 + s.var3 + s.var2 + (1 | id2/id1), data = data, family = binomial, weights = weights)
Warning messages:
1: In eval(family$initialize, rho) :
non-integer #successes in a binomial glm!
2: In (function (fn, par, lower = rep.int(-Inf, n), upper = rep.int(Inf, :
failure to converge in 10000 evaluations
3: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
unable to evaluate scaled gradient
4: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge: degenerate Hessian with 8 negative eigenvalues
# Observe where separation might be occurring
table(data$var1, data$outcome)
0 1
A. binervata 6 0
A. costata 4 1
A. decurrens 144 7
A. floribunda 7 2
A. implexa 14 0
A. littoralis 36 0
A. mearnsii 18 0
A. parramattensis 4 0
B. serrata 22 0
C. gummifera 243 1
E. agglomerata 3 0
E. amplifolia 10 1
E. beyeriana 23 7
E. bosistoana 1 0
E. crebra 30 3
E. eugenioides 8 1
E. fibrosa 69 5
E. globoidea 116 18
E. longifolia 33 18
E. moluccana 4 3
E. oblonga 33 3
E. paniculata 19 1
E. pilularis 37 12
E. piperita 53 1
E. punctata 92 44
E. resinifera 54 3
E. sclerophylla 36 0
E. sieberi 31 0
E. sparsifolia 1 0
E. tereticornis 49 5
M. decora 4 0
M. styphelioides 56 0
# Filter out and drop levels from var1 where there are <= 1 occurance of a particular species amongst the presence or absence records
data.2 <- filter(data, var1 != "A. binervata" &
var1 != "A. implexa" &
var1 != "A. littoralis" &
var1 != "A. mearnsii" &
var1 != "A. parramattensis" &
var1 != "B. serrata" &
var1 != "E. agglomerata" &
var1 != "E. bosistoana" &
var1 != "E. sclerophylla" &
var1 != "E. sieberi" &
var1 != "E. sparsifolia" &
var1 != "M. decora" &
var1 != "M. styphelioides" &
var1 != "A. costata" &
var1 != "C. gummifera" &
var1 != "E. amplifolia" &
var1 != "E. eugenioides" &
var1 != "E. paniculata" &
var1 != "E. piperita")
data.2 <- droplevels(data.2)
# Model 3
model3 <- glmer(outcome ~ var1 + s.var3 + s.var2 + (1 | id2/id1), data = data.2, family = binomial, weights = weights)
Warning messages:
1: In eval(family$initialize, rho) :
non-integer #successes in a binomial glm!
2: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.0544203 (tol = 0.002, component 1)
# Model 4
model4 <- glmer(outcome ~ var1 + s.var3 + s.var2 + (1 | id2/id1), data = data.2, family = binomial, weights = weights, control = glmerControl(optimizer ="bobyqa"))
Warning message:
In eval(family$initialize, rho) : non-integer #successes in a binomial glm!
