Comparación Modelos

• Modelo con datos individuales

reg<- lm(mathach~ses+female+sector, data=mlm)

• Modelo con datos agregados

reg_agg<- lm(mathach~ses+female+sector, data=agg_mlm)

Comparación Modelos

pacman::p_load(sjPlot,sjmisc,sjlabelled)
tab_model(reg, reg_agg, show.ci=F, show.se = T, dv.labels = c("Individual", "Agregado"))

Implicancias estimación individual/agregada

Componentes de la varianza

Componentes de la varianza

Correlación intra clase: ICC

• La correlación intra-clase ( ρ$\rho$ ) indica qué porcentaje de la varianza de la variable dependiente se debe a pertenencia a unidades de nivel 2

• Descomposición de la varianza en modelo nulo= Var y=τ00+σ2$Vary={\tau }_{00}+{\sigma }^2$

• Es decir, parte de la varianza se debe a los individuos ( σ2${\sigma }^2$ ) y parte al grupo ( τ00${\tau }_{00}$ )

Correlación intra clase: ICC

• Correlación intra-clase = ICC =

ρ=τ00τ00+σ2$$\rho =\frac{{\tau }_{00}}{{\tau }_{00}+{\sigma }^2}$$

• Una ICC baja indica baja variabilidad de la variable dependiente entre unidades de nivel 2, y por lo tanto, menores posibilidades de dar cuenta (explicar) de esa varianza con predictores de nivel 2.

librería lme4

• función lmer (linear mixed effects)

• forma general:

  • objeto <- lmer (depvar ~ predictor_1 + predictor_2 + predictor_n + (1 | cluster), data=data)

  • el objeto contiene la información de la estimación; para ver un resumen, summary(objeto), y de manera más presentable,screenreg(objeto)

Tipos de medidas de ajuste

1. Medidas relativas a la varianza de efectos aleatorios (tipo R2$R^2$)

2. Medidas de fit comparativo (deviance)

Ajuste por proporción de varianzas

Estimación de modelos para el ejemplo:

results_0: nulo results_1: agrega predictores individuales (nivel 1) results_2: agrega predictores grupales (nivel 2)

results_0 <-lmer(mathach ~ 1 + (1 | schoolid), data = mlm)
results_1 <-lmer(mathach ~ 1 + ses + female + (1 | schoolid), data = mlm)
results_2 <-lmer(mathach ~ 1 + sector + mnses + (1 | schoolid), data = mlm)

Ajuste por proporción de varianzas

• Se relacionan con el grado de varianza “explicada” (disminución en la(s) varianza(s) atribuida a la inclusión de predictores en el modelo de regresión)

• Controversia en la literatura de multinivel, no existe una medida única

• Las propuestas sugieren usualmente cálculos de R2$R^2$ para cada nivel

Bryck & Raudenbush R2 multinivel (1992)

• lógica general: calcular la diferencia entre componentes de la varianza entre los modelos estimados

• modelo base para la comparación: modelo nulo

• luego, al ir estimando nuevos models con más predictores, se compara (con el modelo nulo) en que medida los componentes de la varianza van disminuyendo

Componentes de la varianza del ejemplo modelos con HSB data

• Los componentes de la varianza van disminuyendo a medida que se ingresan predictores a los modelos

• Por ejemplo, en el modelo 1 el componente de la varianza individual σ2${\sigma }^2$ disminuye en comparación al modelo nulo: 39.148-36.813=2.335 - ¿Cómo interpretar esto?

Bryk & Raudenbush - R2 Nivel 1

Bryk & Raudenbush - R2 Nivel 1

Ej: R2L1=(39.148−36.813)/39.148=2.335/39.148=0.059

Bryck & Raudenbush - R2 Nivel 2

Bryck & Raudenbush - R2 Nivel 2

Ej: R2L2=(8.553−2.314)/8.553=6.239/8.553=0.73

• Recordar interpretación en relación a correlación intra-clase (para el caso de HSB data= 0.18): para el caso del R2 nivel 2 se está dando cuenta del 73% del 18%

Cálculo directo con R - modelo 1

función multilevel.r2, librería misty

pacman::p_load(misty)# for multilevel R2
misty::multilevel.r2(results_1, print = "RB")

##  R-Squared Measures for Multilevel and Linear Mixed Effects Models
## 
##   Reduction in Residual Variance (Raudenbush and Bryk, 2002)
## 
##     Within-Cluster R2: 0.060 
##    Between-Cluster R2: 0.479

Cálculo directo Bryk & Raudenbush con R - modelo 2

función multilevel.r2, librería misty

pacman::p_load(misty)# for multilevel R2
misty::multilevel.r2(results_2, print = "RB")

##  R-Squared Measures for Multilevel and Linear Mixed Effects Models
## 
##   Reduction in Residual Variance (Raudenbush and Bryk, 2002)
## 
##     Within-Cluster R2: -0.000 
##    Between-Cluster R2: 0.731

Nakagawa and Schielzeth R2 (2013)

• dos R2:

  • marginal: para los efectos fijos del modelo

  • condicional: para el modelo con efectos fijos y aleatorios

• utiliza la varianza de los valores predichos de la variable dependiente

Nakagawa and Schielzeth R2 (2013)

• marginal (efectos fijos): R2m=var(ˆYij)var(ˆYij)+τ00+σ2
• condicional (efectos fijos + aleatorios) R2c=var(ˆYij)+τ00var(ˆYij)+τ00+σ2

Cálculo directo Nakagawa con R - modelo 2

misty::multilevel.r2(results_2, print = "NS")

##  R-Squared Measures for Multilevel and Linear Mixed Effects Models
## 
##   Variance Partitioning (Nakagawa and Schielzeth, 2013; Johnson, 2014)
## 
##       Marginal R2: 0.130 
##    Conditional R2: 0.179

2. Ajuste comparativo

Deviance test

• El test o estadístico de deviance compara el ajuste de dos modelos basado en la log verosimilitud de cada modelo

• La hipótesis a contrastar es si predictores adicionales del modelo mejoran o no el ajuste

• Asume que los modelos son anidados, es decir, que un modelo con menos predictores puede ser derivado del modelo mayor mediante la fijación de ciertos coeficientes como 0.

• Deviance= −2∗LL (LL=Log Likelihood)

• Deviance test= deviance(anidado)−deviance(mayor)

2. Ajuste comparativo

Deviance test

• La distribución del estadístico de devianza es χ2, y los grados de libertad para calcular el valor crítico equivalen al número de parámetros extra en el modelo mayor

• Es decir, parámetros modelo mayor - parámetros modelo inicial (o anterior)

• Se utiliza con estimación ML en lugar de REML (restricted maximum likelihood).

Ejemplo deviance

results_0ml = lmer(mathach ~ 1 + (1 | schoolid), REML=FALSE)
results_1ml = lmer(mathach ~ 1 + minority + ses + (1 | schoolid), REML=FALSE)

• Deviance test= 47115.8-46443.64= 672.16, 2 DF (2 parámetros adicionales)

Ejemplo deviance

• Valor crítico χ2DF=2 para p<0.95=5.99

• Por lo tanto, se rechaza H0, es decir, las diferencias entre los modelos son distintas de 0 ( p<0.05 ). En otras palabras, el modelo con más parámetros presenta un mejor ajuste.

Ejemplo deviance II - de manera directa en R:

anova(results_0ml,results_1ml)

Y para reportar:

print(xtable::xtable(anova(results_0ml,results_1ml)), type="html")

Ajuste y modelos con pendiente aleatoria

• Los test de ajuste por proporción de varianza (R2) no son recomendados al momento de estimar modelos con pendiente aleatoria

• El test de deviance se utiliza al momento de reportar la significancia de los efectos aleatorios (ej: con pendiente aleatorio,comparando el mismo modelo con y sin aleatorización)

• Considerar que los componentes de la varianza son parámetros del modelo, es decir, se cuentan para la diferencia de grados de libertad

Tablas y descriptivos

• Tabla de descriptivos de variables L1 y L2

  • para esto conviene generar base de datos agregados para tabla descriptiva

  • se recomienda utilizar librería summarytools, función dfSummary (información aquí)

• Descripción del número de casos por nivel y también de otras variables relevantes (ej, dependiente e independientes comprometidas en hipótesis)

Descriptivos

Ejemplo tabla descriptiva (Ferron et al 2007, en O’Connell Cap.11)

Alternativas en R: Stargazer, summarytools

Reporte modelos

(O’Connell, cap.11)

Reporte modelos

Reporte modelos

• Efectos fijos, con valores t o se

• Asteriscos para niveles de significación

• Efectos aleatorios

• Fit: log likelihood

• Fit adicionales: AIC, BIC, deviance, R2

Reportes de relaciones bivariadas

• tablas cruzadas
• correlaciones / corrplot
• nubes de puntos

Revisar:

• https://multivariada.netlify.app/assignment/02-code/

• https://correlacional.netlify.app/assignment/03-practico

Sobre reporte reproducible

RMarkdown / Quarto

Tres cosas básicas a aprender:

• Markdown

• Insertar trozos (chunks) de código

• Opciones generales en encabezado YAML

Varios tutoriales ... revisar

• Guías de documentos dinámicos (con video) de clases de correlacional 2023:

  • Reportes dinámicos 1: RMarkdown
  • Reportes dinámicos 2: Tablas, gráficos y referencias
  • Reportes dinámicos 3: Quarto

• Tutorial Hola Quarto - quarto.org

• Video introductorio de Mine Centinkaya