Seleção entre um modelo binomial inflado de zero, OLRE e beta-binomial
Preciso de ajuda para decidir qual dos modelos a seguir se encaixa melhor nos dados que tenho. Esta foi uma pesquisa em que os participantes relataram proporções de sucessos (definidos como n / m) nas condições A e B. O modelo prevê as proporções pela conditionvariável binária e contínua xe zvariáveis (variando de 1 a 7), bem como efeitos aleatórios para cada um subjecte 13 tipos de task. Esta é a distribuição das proporções
Portanto, o modelo é definido como
mod_b0 <- glmmTMB(n/m ~ x*condition + z*condition + (1|subject) + (1|task), weights = m, family = binomial)
summary(mod_b0)
AIC BIC logLik deviance df.resid
22830.4 22883.7 -11407.2 22814.4 5781
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
task (Intercept) 0.2094 0.4576
subject (Intercept) 1.5546 1.2468
Number of obs: 5789, groups: task, 13; subject, 225
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -3.44713 0.25706 -13.410 < 2e-16 ***
x 0.38560 0.03690 10.449 < 2e-16 ***
conditionB -1.36826 0.20133 -6.796 1.08e-11 ***
z -0.07328 0.02276 -3.220 0.00128 **
x:conditionB 0.17682 0.03807 4.644 3.41e-06 ***
conditionB:z 0.12544 0.02512 4.994 5.91e-07 ***
O teste de resíduos por DHARMa(N = 1000 simulações) sugere que não há superdispersão, que há inflação zero e que o modelo não se ajusta bem aos dados.
Tentei três soluções:
- Modelo binomial com inflação zero
- Modelo binomial OLRE
- Modelo beta-binomial
Aqui estão os resultados de todos os três.
Modelo binomial com inflação zero
mod_bzi <- glmmTMB(n/m ~ x*condition + z*condition + (1|task) + (1|subject),
data = dx, family = binomial, weights = m, ziformula = ~ 1 + condition*z)
summary(mod_bzi)
AIC BIC logLik deviance df.resid
17949.0 18029.0 -8962.5 17925.0 5777
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
task (Intercept) 0.09208 0.3034
subject (Intercept) 1.95087 1.3967
Number of obs: 5789, groups: task, 13; subject, 225
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.65838 0.29974 -8.869 < 2e-16 ***
x 0.40498 0.04874 8.309 < 2e-16 ***
conditionB -1.31011 0.26986 -4.855 1.21e-06 ***
z -0.01559 0.02852 -0.547 0.5847
x:conditionB 0.14559 0.05150 2.827 0.0047 **
conditionB:z 0.19289 0.03291 5.861 4.59e-09 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Zero-inflation model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.393898 0.084827 -4.644 3.42e-06 ***
conditionB 0.307062 0.126750 2.423 0.0154 *
z 0.034095 0.034146 0.999 0.3180
conditionB:z -0.003092 0.046014 -0.067 0.9464
Observe que as linhas de regressão no gráfico correto não são significativamente diferentes das linhas de quantis se o número de simulações for 250!
Agora vemos uma ligeira subdispersão.
Modelo OLRE
mod_OLRE <- glmmTMB(n/m ~ x*condition + z*condition + (1|task) + (1|subject) + (1|obs_id),
data = dx, family = binomial, weights = m)
AIC BIC logLik deviance df.resid
15588.2 15648.1 -7785.1 15570.2 5780
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
task (Intercept) 0.4361 0.6604
subject (Intercept) 3.0721 1.7527
obs_id (Intercept) 4.8962 2.2127
Number of obs: 5789, groups: task, 13; subject, 225; obs_id, 5789
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -4.46870 0.55951 -7.987 1.38e-15 ***
x 0.43727 0.09152 4.778 1.77e-06 ***
conditionB -2.65037 0.53953 -4.912 9.00e-07 ***
z -0.17483 0.06014 -2.907 0.003650 **
x:conditionB 0.35813 0.10186 3.516 0.000438 ***
conditionB:z 0.21831 0.06827 3.198 0.001384 **
Novamente, não há mais inflação zero, mas há alguma subdispersão.
Modelo beta-binomial
mod_bb <- glmmTMB(n/m ~ x*condition + z*condition + (1|task) + (1|subject),
data = dx, family = betabinomial(link = "logit"), weights = m)
AIC BIC logLik deviance df.resid
15305.4 15365.4 -7643.7 15287.4 5780
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
task (Intercept) 0.2267 0.4761
subject (Intercept) 0.9929 0.9965
Number of obs: 5789, groups: task, 13; subject, 225
Overdispersion parameter for betabinomial family (): 1.54
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.51074 0.33909 -7.404 1.32e-13 ***
x 0.24238 0.05426 4.467 7.94e-06 ***
conditionB -1.31799 0.32146 -4.100 4.13e-05 ***
z -0.08722 0.03508 -2.486 0.01291 *
x:conditionB 0.17975 0.06081 2.956 0.00312 **
conditionB:z 0.09051 0.04010 2.257 0.02400 *
Aqui, há mais subdispersão nos modelos anteriores.
Minhas conclusões e perguntas
- Pela aparência da distribuição residual, parece-me que o modelo beta-binomial faz o melhor trabalho para contabilizar os dados. Todos os modelos têm alguns problemas com níveis mais altos de preditores, pois há menos casos para esses valores. Assim, não é de admirar que os ajustes sejam um pouco mais pobres nesse segmento da trama.
- Os valores de AIC são mais baixos para o modelo beta-binomial. No entanto, não tenho certeza se posso comparar o AIC para modelos com diferentes distribuições do critério. Se sim, então esse seria outro argumento para escolher o modelo beta-binomial.
- Os coeficientes são um tanto semelhantes nos modelos beta-binomial e binomial inflado a zero. O modelo OLRE tem alguns coeficientes bastante diferentes. De acordo com Harrison (2014) , os modelos beta-binomiais tendem a produzir estimativas mais confiáveis do que OLRE. Portanto, eu ficaria com aquele.
- Você concorda com minhas conclusões de que o modelo beta-binomial é o melhor de todos os propostos?
- Existe alguma outra maneira de melhorar o ajuste dos modelos que não pensei?
- Posso tentar ajustar o parâmetro de inflação zero no modelo beta-binomial para obter um melhor ajuste, embora nenhuma inflação zero tenha sido diagnosticada pelo DHARMa?
- Existe alguma outra maneira de testar o ajuste dos modelos?
- A subdispersão é "problemática" para o modelo beta-binomial? De acordo com o FAQ do GLMM , a dispersão é um problema apenas para modelos com variância fixa, como binomiais ou de Poisson.
Respostas
Você concorda com minhas conclusões de que o modelo beta-binomial é o melhor de todos os propostos?
Sim, você parece ter feito um trabalho completo nessa análise. Sua opinião sobre se é normal comparar esses modelos com a AIC é boa. Lembro-me de ter lido informações conflitantes sobre esse ponto, mas rapidamente encontrei uma referência que apóia a ideia de que está tudo bem:
Hardin, JW e Hilbe, JM, 2014. Estimativa e teste de modelos de regressão binomial e beta-binomial com e sem inflação zero. The Stata Journal, 14 (2), páginas 292-303.https://journals.sagepub.com/doi/pdf/10.1177/1536867X1401400204
Existe alguma outra maneira de melhorar o ajuste dos modelos que não pensei?
Você pode observar a precisão preditiva usando uma abordagem de treinar / validar / testar.
Posso tentar ajustar o parâmetro de inflação zero no modelo beta-binomial para obter um melhor ajuste, embora nenhuma inflação zero tenha sido diagnosticada pelo DHARMa?
Valeria a pena tentar, mas dada a saída do DHARMa provavelmente não vai melhorar as coisas.
Existe alguma outra maneira de testar o ajuste dos modelos?
Novamente, eu sugeriria olhar para as previsões.
A subdispersão é "problemática" para o modelo beta-binomial? De acordo com o FAQ do GLMM, a dispersão é um problema apenas para modelos com variância fixa, como binomiais ou de Poisson.
A dispersão insuficiente e excessiva é "tratada" por modelos beta-binomiais, portanto, não deve ser um problema.