I have a dataset with 72 individuals from 5 separate groups, with repeated measures - each individual was sampled 4 times. All data is binomial, and in most cases there are more 0s than 1s.
I have a lot of potential predictor variables BacteriaA
,BacteriaB
etc, (9+ depending on how I group the data), and a number of outcome variables I want to test DiseaseA
, Disease B
, etc.
data <- read.csv("https://pastebin.com/raw/gwEFqh79 ")
cols <- colnames(data)
data[cols] <- lapply(data[cols], factor)
str(data)
'data.frame': 500 obs. of 18 variables:
$ Group : Factor w/ 5 levels "A","B","C","D",..: 1 1 1 1 1 1 1 1 1 1 ...
$ ID : Factor w/ 73 levels "E1","E10","E11",..: 1 1 1 1 2 2 2 2 2 2 ...
$ Time : Factor w/ 4 levels "1","2","3","4": 2 4 2 4 1 2 3 1 2 3 ...
$ DiseaseA : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
$ DiseaseB : Factor w/ 2 levels "0","1": 2 1 2 1 1 1 1 1 1 1 ...
$ DiseaseC : Factor w/ 2 levels "0","1": 1 2 1 2 1 1 1 1 1 1 ...
$ DiseaseD : Factor w/ 2 levels "0","1": 2 2 2 2 1 1 1 1 1 1 ...
$ DiseaseE : Factor w/ 2 levels "0","1": 2 2 2 2 1 1 1 1 1 1 ...
$ DiseaseF : Factor w/ 2 levels "0","1": 2 2 2 2 1 1 1 1 1 1 ...
$ BacteriaA: Factor w/ 2 levels "0","1": 2 1 1 1 1 2 1 2 1 1 ...
$ BacteriaB: Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
$ BacteriaC: Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 2 1 1 ...
$ BacteriaD: Factor w/ 2 levels "0","1": 2 1 2 1 2 2 1 1 2 1 ...
$ BacteriaE: Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
$ BacteriaF: Factor w/ 2 levels "0","1": 1 1 1 2 1 1 1 1 1 1 ...
$ BacteriaG: Factor w/ 2 levels "0","1": 1 2 2 2 1 2 2 1 2 2 ...
$ BacteriaH: Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
$ BacteriaI: Factor w/ 2 levels "0","1": 2 1 1 1 1 1 1 1 1 1 ...
I have run a model using glmer
because I believe I need a mixed effects model, and I have binomial data.
model1 <- glmer(DiseaseA ~ (1|Group) + (1|ID) + Time + BacteriaA + BacteriaB + BacteriaC + BacteriaD + BacteriaE + BacteriaF + BacteriaG + BacteriaH + BacteriaI, data = data, family = 'binomial'(link = "logit"))
summary(model1)
Model failed to converge with max|grad| = 0.0337228 (tol = 0.001, component 1)Model failed to converge with max|grad| = 0.0337228 (tol = 0.001, component 1)Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
Family: binomial ( logit )
Formula: DiseaseA ~ (1 | Group) + (1 | ID) + Time + BacteriaA + BacteriaB +
BacteriaC + BacteriaD + BacteriaE + BacteriaF + BacteriaG + BacteriaH + BacteriaI
Data: data
AIC BIC logLik deviance df.resid
426.2 489.4 -198.1 396.2 485
Scaled residuals:
Min 1Q Median 3Q Max
-1.2452 -0.3965 -0.2829 -0.1610 4.8949
Random effects:
Groups Name Variance Std.Dev.
ID (Intercept) 1.298e+00 1.139325
Group (Intercept) 3.841e-05 0.006198
Number of obs: 500, groups: ID, 73; Group, 5
Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -1.5414 0.4423 -3.485 0.000493 ***
Time2 0.6009 0.4518 1.330 0.183535
Time3 0.4617 0.4782 0.966 0.334287
Time4 0.4136 0.4577 0.904 0.366223
BacteriaA1 -0.5733 0.4593 -1.248 0.211939
BacteriaB1 -0.6688 0.5057 -1.323 0.185963
BacteriaC1 -0.6867 0.4371 -1.571 0.116228
BacteriaD1 -0.7494 0.3458 -2.167 0.030236 *
BacteriaE1 -0.1796 0.7715 -0.233 0.815964
BacteriaF1 -1.3041 0.6412 -2.034 0.041972 *
BacteriaG1 -0.7576 0.3260 -2.324 0.020133 *
BacteriaH1 -10.1035 124.1785 -0.081 0.935154
BacteriaI1 -0.1340 0.7987 -0.168 0.866712
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation matrix not shown by default, as p = 13 > 12.
Use print(x, correlation=TRUE) or
vcov(x) if you need it
convergence code: 0
Model failed to converge with max|grad| = 0.0337228 (tol = 0.001, component 1)
I am super new to modelling, never mind mixed models, and I'm just really stuck on what I should do next.
Is this the right model? How do I work out the odds ratios from this? What does Model failed to converge with max|grad|
mean?
And most of all, how do I choose a model from this? Is there a stepwise way of dealing with so many variables?