Cuándo tratar un número como un recuento

Aug 24 2020

Estoy probando las diferencias de grupo en el número de días que los participantes usaron un medicamento en los 28 días anteriores. Estos son los datos, pero tengo problemas para decidir qué enfoque usar: regresión gaussiana estándar o regresión binomial agregada. He hecho preguntas similares antes en CV (por ejemplo, aquí ) pero todavía estoy un poco inseguro.

He proporcionado el código R para la replicabilidad, pero, por supuesto, cualquiera que quiera opinar, usuario de R o no, es más que bienvenido.

df <- data.frame(group = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
                 baseline = as.integer(c(28, 28, 28, 28, 28, 28, 28, 28, 28, 28, 12, 28, 28, 28, 28, 28, 28, 24, 28, 28, 28, 28, 28, 28, 28, 28, 20, 28, 28, 24, 24, 28, 28, 28, 28, 28, 28, 28, 24, 28, 28, 28, 28, 28, 16, 28)),
                 outcome = as.integer(c(28, 0, 28, 0, 0, NA, NA, 16, 28, 10, 12, 0, 28, 12, 0, 0, 28, 8, 0, 28, 28, 0, 4, NA, NA, 0, NA, 28, NA, 20, 1, 3, 28, 26, NA, 0, 20, 16, 16, 0, NA, 3, 0, 1, 20, 0)),
                 coverage = 28)

groupes el tratamiento que recibieron los participantes; baselineel número de días que utilizaron en los 28 días previos a comenzar el estudio; outcomeel número de días que utilizaron durante el estudio de 28 días ( coveragees el número de días del ensayo).

Aquí están las estadísticas resumidas:

library(tidyverse)

df %>%
  group_by(group) %>%
    drop_na(outcome) %>%
      summarise(mean = mean(outcome, na.rm = T),
                sd = sd(outcome, na.rm = T),
                median = median(outcome, na.rm = T),
                firstQuartile = quantile(outcome, probs = 0.25, na.rm = T),
                thirdQuartile = quantile(outcome, probs = 0.75, na.rm = T),
                tally = n()) 

# output
# group  mean    sd median firstQuartile thirdQuartile tally
# <dbl> <dbl> <dbl>  <int>         <dbl>         <dbl> <int>
#     0  10.7  11.3      3             0            20    17
#     1  12.3  12.3     10             0            28    21

Y la distribución de los resultados en cada grupo

ggplot(df, aes(x = outcome, group = group)) + geom_histogram() + facet_wrap(~group) + scale_x_continuous(breaks = seq(0,28,7))

Como es habitual en los datos sobre el consumo de sustancias, los resultados se distribuyen de forma bastante bimodal.

Cuando analizo el resultado, regresando los días utilizados, tratados como variable continua, sobre el tratamiento groupy los baselinedías utilizados ...

summary(contMod <- lm(formula = outcome ~ group + baseline, 
                      data = df, 
                      na.action = na.exclude))

# output
# Coefficients:
#             Estimate Std. Error t value Pr(>|t|)
# (Intercept)  17.7783    16.0047   1.111    0.274
# group         1.7020     3.9248   0.434    0.667
# baseline     -0.2660     0.5919  -0.449    0.656

El modelo indica que no hay diferencias significativas entre los grupos en la media de días utilizados al controlar los días de referencia utilizados. Sin embargo, los residuos del modelo son muy anormales:

hist(resid(contMod))

La gráfica de cuantiles-cuantiles cuenta la misma historia

plot(contMod,2)

Entonces, para mí, parece que la regresión gaussiana estándar puede no ser apropiada para modelar estos datos.

Dado que los datos son básicamente recuentos enteros de ocurrencias de un evento binario (usado el día x versus no usado el día x ) dentro de un número conocido de 'ensayos' (28 días). Me preguntaba si una regresión binomial agregada podría ser una forma más apropiada de modelar los datos.

summary(contMod <- glm(formula = cbind(outcome, coverage-outcome) ~ group + baseline, 
                       data = df, 
                       family = binomial,
                       na.action = na.exclude))

# output
# Coefficients:
#             Estimate Std. Error z value Pr(>|z|)  
# (Intercept)  0.54711    0.50908   1.075   0.2825  
# group        0.25221    0.12634   1.996   0.0459 *
# baseline    -0.03866    0.01886  -2.050   0.0403 *

Ahora la diferencia de grupo es significativa cuando se controla la línea de base.

Una diferencia tan dramática en los resultados de dos modelos diferentes del mismo me sorprende bastante. Por supuesto, sabía que era posible, pero nunca antes lo había encontrado en la naturaleza.

Así que tengo varias preguntas para los usuarios inteligentes de CV.

1. ¿Es la regresión binomial agregada una mejor manera de modelar estos datos dada la extrema no normalidad tanto del resultado como de los residuos del modelo? ¿Y si es apropiado lo hice correctamente? E, incluso si lo hice correctamente, ¿hay otra forma aún mejor (no paramétrica, por ejemplo? La prueba de Kruskal-Wallis kruskal.test(outcome ~ group, data = df)arrojó resultados similares a los de Gauss,$\chi^2 = 0.07, p = 0.80$, pero no controla la línea de base).

2. ¿Cómo interpreto el resultado de la regresión logística agregada? Si el resultado fuera un proceso de Bernoulli, usaría una regresión logística binaria simple e interpretar los resultados sería sencillo, exponencializaría el coeficiente del grupo y eso representaría cuánto mayores son las probabilidades de usar en el día en cuestión en el 1grupo que en el 0grupo. Pero estoy menos seguro de cómo se informaría el resultado del binomio agregado.

Ayuda muy apreciada como siempre.

Respuestas

2 NickCox Aug 25 2020 at 16:14

Estás haciendo una pregunta sobre los métodos aquí, pero prefiero comenzar una respuesta a partir de tus datos y lo que quieres saber.

Aquí hay una versión de sus datos que puede ser útil para las personas que no usan R habitualmente; las líneas de apertura y cierre son específicamente para Stata, pero los usuarios de otro software deberían poder adaptarse según sus necesidades. Los puntos son el código genérico de Stata para las faltas numéricas y corresponden a NA en R.

clear
input byte(id group baseline outcome coverage)
 1 1 28 28 28
 2 1 28  0 28
 3 1 28 28 28
 4 1 28  0 28
 5 1 28  0 28
 6 1 28  . 28
 7 1 28  . 28
 8 1 28 16 28
 9 1 28 28 28
10 1 28 10 28
11 1 12 12 28
12 1 28  0 28
13 1 28 28 28
14 1 28 12 28
15 1 28  0 28
16 1 28  0 28
17 1 28 28 28
18 1 24  8 28
19 1 28  0 28
20 1 28 28 28
21 1 28 28 28
22 1 28  0 28
23 1 28  4 28
24 1 28  . 28
25 0 28  . 28
26 0 28  0 28
27 0 20  . 28
28 0 28 28 28
29 0 28  . 28
30 0 24 20 28
31 0 24  1 28
32 0 28  3 28
33 0 28 28 28
34 0 28 26 28
35 0 28  . 28
36 0 28  0 28
37 0 28 20 28
38 0 28 16 28
39 0 24 16 28
40 0 28  0 28
41 0 28  . 28
42 0 28  3 28
43 0 28  0 28
44 0 28  1 28
45 0 16 20 28
46 0 28  0 28
end

El núcleo del problema es comparar la media outcomede dos valores de group. Una distracción es que baselinevaría y parece que lo más sencillo al menos al principio es ignorar los casos para los que no son 28 días baseline. No es obvio para mí que agregar baselinecomo predictor sea la mejor manera de ajustar las variaciones baseline; una alternativa es escalar outcomea una fracción de baseline, pero es probable que también sea confuso.

La comparación de medias puede, naturalmente, volver a plantearse como un problema de regresión. Aquí es un gráfico con la línea de regresión superpuesta para la regresión de outcomesobre grouppara baseline28 días. Con esta simplificación, la línea simplemente conecta las dos medias del grupo.

Estoy rotando sus histogramas y tratando los datos como lo que son, puntos de datos en un problema de regresión comparando medias. La acumulación de resultados idénticos es solo una convención gráfica y no afecta los resultados de la regresión.

Aunque se refiere a la "regresión gaussiana", la condición ideal de los errores gaussianos o normales es el aspecto menos importante de la regresión lineal. El texto reciente de Gelman y amigos

https://www.cambridge.org/core/books/regression-and-other-stories

incluso desaconseja los gráficos cuantílicos normales de residuos como una pérdida de tiempo. Yo no iría tan lejos, pero es un punto de vista serio.

Un vistazo al gráfico y los resultados de la regresión apuntan a una diferencia de 2,9 días; Mi suposición laica es que una diferencia de esa magnitud sería clínica o científicamente interesante, pero los resultados de la regresión muestran que la muestra es demasiado pequeña para confirmarla como significativa a niveles convencionales. Si le preocupa que tal indicación sea demasiado dependiente de la suposición implícita de errores normales, un poco de arranque de esos resultados de regresión implica que una diferencia de cero está dentro de casi cualquier intervalo de confianza para la diferencia de medias. Retirarme a Kruskal-Wallis me parece un consejo de desesperación; ¿Por qué utilizar la tecnología de la década de 1950 cuando la tecnología de la década de 1970 (bootstrap) está disponible y le permite concentrarse en la diferencia de medios que es de interés principal?

En general, es una muy buena idea ser sensible a si sus datos se cuentan o se miden; pensar en sus distribuciones condicionales; y observar si un resultado está necesariamente acotado. En este caso particular, estos resultados de regresión simple implican que poco importa lo que asume o es lo que se asume o es ideal para los métodos utilizados. La diferencia entre las medias parece interesante pero no es convencionalmente significativa y esa indicación es robusta para cualquier cosa que haga a través del análisis.

Sin embargo, si trato de hacer coincidir su regresión binomial, pero concentrándome en baselineigual a 28, encuentro de manera similar que es suficiente para cambiar la diferencia a convencionalmente significativa. Al principio no entendí por qué hay una diferencia tan grande en la indicación.

Pero deberíamos preocuparnos por lo que se supone sobre las distribuciones. Observo que los binomios no pueden tener forma de U. Primero dudé de que ese fuera el problema, pero ese pensamiento era visceral, no lógico. Si repite el análisis con errores estándar robustos (Eicker-Huber-White), la importancia se evapora.

En resumen, al aplicar la regresión binomial en lugar de la regresión simple, está reemplazando un supuesto distributivo que no muerde, aunque parece bastante incorrecto, por un supuesto distributivo que sí muerde. Ese es mi diagnóstico.

FWIW, el uso de días aquí como un recuento de números enteros es en parte natural (las personas tienen ritmos diarios que siguen, a veces de forma rígida y a veces de forma flexible) y en parte una convención (en principio, los datos también pueden estar basados en horas del día, lo que produce días fraccionarios) .

Finalmente, la comparación de medias no es el único juego en la ciudad. Observo que en el grupo 0 solo 2 de 13, pero en el grupo 1 7 de 19 personas informaron los 28 días completos. Esas diferencias afectaron naturalmente los medios, pero los detalles también pueden ser importantes.

Nitty-gritty sigue como salida de Stata. La gente R espera que seamos lo suficientemente inteligentes como para decodificar la salida R si somos lo suficientemente perversos como para no usarla (no para usarla rutinariamente, en mi caso) y les devuelvo el cumplido. El minimalismo de la salida de R es admirable, excepto que no mostrar el tamaño de muestra utilizado ni siquiera en el resumen predeterminado me desconcierta.