3/12/2021

Análisis de varianza de un factor

El análisis de varianza de un factor sirve para contrastar varios grupos de una variable cuantitativa.

La variable nominal u ordinal que define los grupos se denomina variable independiente o factor. La variable continua con los datos se considera la variable dependiente.

El estadístico utilizado para poner a prueba la igualdad de medias se llama \(F\) y refleja el parecido entre las medias que se están comparando tomando en cuenta la variabilidad de las medias de cada grupo.

\[F={\sigma_b^2 \over \sigma_w^2}\]

El numerador indica la varianza poblacional entre las medias de cada grupo mientras que el denominador es la varianza poblacional dentro de cada grupo.

\[F={{n\Sigma(\bar{x}-\bar{x_G})^2 \over k-1} \over {SC_1+SC_2+ SC_k \over N-k}}={{SC_b \over gl_b} \over {SC_w \over gl_w}}\]

Si las medias poblacionales son iguales la estimación \(\sigma_b^2\) reflejará la misma variación que la estimación \(\sigma_w^2\) y el cociente \(F\) se aproximará al valor de uno. Cuanto más diferentes sean las medias mayor será el valor de \(F\).

El estadístico \(F\) se interpreta de forma similar a la prueba \(t\), si el valor p asociado a \(F\) es menor a .05 rechazamos \(H_0\) y asumimos que al menos una de las medias grupales es diferente.

Existen ciertos supuestos que deben cumplirse para utilizar el análisis de varianza:

  • Muestro aleatorio o mínimamente asignación aleatoria a los grupos.
  • Independencia de los puntajes en la variable respuesta.
  • Muestras con distribución normal para cada uno de los grupos en la variable explicativa.
  • Homogeneidad de varianza en todos los grupos de la variable explicativa.

Si tomamos en cuenta un grupo de 15 estudiantes con diferente nivel de habilidad en el uso de la computadora (baja, media, alta) y sus calificaciones en el examen de estadística. ¿Existen diferencias en las calificaciones del examen de estadística de acuerdo al nivel de habilidad en el uso de la computadora?

grupo <- gl(3, 5, labels = c("hb", "hm", "ha"), ordered = T)
cal <- c(7, 8, 7.5, 9, 8, 7, 8, 8.5, 9, 8.5, 8, 9, 8.5, 9, 10)
datos <- cbind.data.frame(grupo, cal)

Exploremos los datos

attach(datos)
tapply(cal, grupo, mean)
 hb  hm  ha 
7.9 8.2 8.9 
tapply(cal, grupo, var)
   hb    hm    ha 
0.550 0.575 0.550 
tapply(cal, grupo, length)
hb hm ha 
 5  5  5 

boxplot(cal ~ grupo)

Para probar la homogeneidad de varianza contamos con la prueba de Bartlett, si el valor p es mayor al nivel de significancia asumimos que existe homogeneidad de varianzas.

bartlett.test(cal ~ grupo, data = datos)
    Bartlett test of homogeneity of variances

data:  cal by grupo
Bartlett's K-squared = 0.0023827, df = 2, p-value = 0.9988

Calculamos la prueba anova

oneway.test(cal ~ grupo, var.equal = T)
    One-way analysis of means

data:  cal and grupo
F = 2.3582, num df = 2, denom df = 12, p-value = 0.1368

Corrección de Welch para desigualdad de var.

oneway.test(cal ~ grupo)
    One-way analysis of means (not assuming equal variances)

data:  cal and grupo
F = 2.2065, num df = 2.0000, denom df = 7.9991, p-value = 0.1725

Otra opción para el análisis

aov.r <- aov(cal ~ grupo)
summary(aov.r)
            Df Sum Sq Mean Sq F value Pr(>F)
grupo        2  2.633  1.3167   2.358  0.137
Residuals   12  6.700  0.5583               

Prueba post hoc

TukeyHSD(aov.r)
  Tukey multiple comparisons of means
    95% family-wise confidence level

Fit: aov(formula = cal ~ grupo)

$grupo
      diff        lwr      upr     p adj
hm-hb  0.3 -0.9607832 1.560783 0.8041587
ha-hb  1.0 -0.2607832 2.260783 0.1281361
ha-hm  0.7 -0.5607832 1.960783 0.3336871

Prueba de Kruskal-Wallis

Cuando no se cumplen los supuestos del anova de un factor se puede utilizar la prueba no paramétrica de Kruskal-Wallis.

kruskal.test(cal ~ grupo)
    Kruskal-Wallis rank sum test

data:  cal by grupo
Kruskal-Wallis chi-squared = 3.7735, df = 2, p-value = 0.1516

Análisis de regresión simple

La regresión lineal es otro método estadístico utilizado para conocer la asociación entre variables. A diferencia de la correlación, en este análisis se valora la contribución de la variable independiente (x) sobre la variable dependiente (y).

Supuestos de la regresión

La regresión lineal al igual que otros modelos estadísticos requiere que se cumplan ciertas condiciones para su aplicación.

Linealidad.- La ecuación de la regresión adopta la forma de una recta conjuntando los valores observados en la variable independiente el origen de la recta y los residuos. Si los datos no se comportan linealmente podemos hablar de un problema de especificación.

Independencia.- Los residuos deben ser independientes, es decir, no debe existir un patrón en los residuos.

Homocedasticidad.- Para cada valor de la variable independiente la varianza de los residuos es constante.

Normalidad.- Para cada valor de la variable independiente los residuos se distribuyen de forma normal.

La recta de regresión

La regresión simple brinda una medida de qué tanto se parecen los datos observados a una línea recta representada por la ecuación:

\[y_i= \beta_0+ \beta_1 x_i+ \epsilon\]

Donde

\(Y_i\) es el valor de la variable dependiente para la observación i

\(\beta_0\) es la constante al origen

\(\beta_1\) es el coeficiente de aportación la variable independiente sobre la variable dependiente

\(x_i\) es el valor de la variable independiente para la observación iésima

\(\epsilon\)es el término de error o residuo

El parámetro \(\beta_0\), conocido como la “ordenada en el origen”, indica cuál es el valor de \(y\) cuando \(x\) = 0. El parámetro \(\beta_1\), conocido como la “pendiente”, indica cuánto aumenta \(y\) por cada aumento de una unidad en \(x\). El problema consiste en obtener estimaciones de estos coeficientes a partir de una muestra de observaciones sobre las variables \(y\) y \(x\).

Existen diferentes métodos para determinar la recta de la regresión, el más utilizado es el método de mínimos cuadrados en el cual se elevan al cuadrado las diferencias entre los valores observados \(y\) y los valores predichos por la regresión \(y\)’ como se muestra en la siguiente ecuación.

\[\Sigma(y_i-y')^2\]

Coeficiente de determinación

Se simboliza como \(R^2\) e indica la cantidad de variación en la variable dependiente que se explica por los cambios en la variable independiente. Mientras mayor es la proporción de la varianza de \(y\) explicada por el modelo, mayor será el coeficiente de determinación (\(R^2\)). El rango de \(R^2\) puede variar entre 0, que indica un nulo ajuste, hasta 1, que señala un ajuste perfecto.

Bondad de ajuste

Siempre existe la posibilidad de encontrar una línea de ajuste para un conjunto de datos, no obstante, esto no quiere decir que sea la mejor opción para explicar esos datos.

En primer lugar el obtener un coeficiente de determinación de .01 nos indica que el ajuste de los datos al modelo de regresión sólo explica una pequeña cantidad de la varianza en los datos.

En segundo lugar, el estadístico \(F\) del ANOVA permite contrastar la \(H_0\) en que el valor de \(R^2\) es cero, cuando el nivel crítico es menor al nivel de significancia establecido, se rechaza \(H_0\).

Como ejemplo tomemos la edad y pulso máximo de 15 personas. Se tiene como estándar que el pulso máximo para una persona está dado por \(Max=220-edad\).

edad <- c(18, 23, 25, 35, 65, 54, 34, 56, 72, 19, 23, 42, 18, 39, 37)
pulso <- c(202, 186, 187, 180, 156, 169, 174, 172, 153, 199, 193, 174, 198, 183,
    178)
plot(edad, pulso)  # crear la gráfica
abline(lm(pulso ~ edad))  # añadir la línea de regresión

library(UsingR)
res.lm <- simple.lm(edad, pulso)

# Resumen del modelo
summary(res.lm)
Call:
lm(formula = y ~ x)

Residuals:
    Min      1Q  Median      3Q     Max 
-8.9258 -2.5383  0.3879  3.1867  6.6242 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 210.04846    2.86694   73.27  < 2e-16 ***
x            -0.79773    0.06996  -11.40 3.85e-08 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 4.578 on 13 degrees of freedom
Multiple R-squared:  0.9091,    Adjusted R-squared:  0.9021 
F-statistic:   130 on 1 and 13 DF,  p-value: 3.848e-08

# Coeficientes de intercepto y pendiente
coef(res.lm)
(Intercept)           x 
210.0484584  -0.7977266 
# residuos
resid.lm <- resid(res.lm)
summary(resid.lm)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
-8.9258 -2.5383  0.3879  0.0000  3.1867  6.6242 

# probando los supuestos del modelo
plot(res.lm)

# Modelo con intervalos de confianza
simple.lm(edad, pulso, show.ci = TRUE, conf.level = 0.9)

Call:
lm(formula = y ~ x)

Coefficients:
(Intercept)            x  
   210.0485      -0.7977  

# Predicción de datos no contenidos en la muestra
predict(res.lm, data.frame(x = c(50, 60)))
       1        2 
170.1621 162.1849 

Regresión Seigel

Es una modificación del método Kendall-Theil, este método es robusto ante la presencia de datos atípicos en la variable dependiente. Se basa en el cálculo de todas las líneas de regresión entre cada par de puntos y utiliza la mediana de las pendientes para generar la pendiente del modelo.

library(mblm)
reg_seigel <- mblm(pulso ~ edad)
summary(reg_seigel)
Call:
mblm(formula = pulso ~ edad)

Residuals:
   Min     1Q Median     3Q    Max 
-7.128 -1.012  1.051  5.149  8.587 

Coefficients:
            Estimate      MAD V value Pr(>|V|)    
(Intercept) 207.2348   5.0352     120 0.000725 ***
edad         -0.7679   0.1293       0 0.000725 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 4.958 on 13 degrees of freedom

Regresión de cuantiles

Este método de regresión utiliza la mediana u otro cuantil en lugar de la media para generar el modelo y no hace supuestos sobre la distribución de los datos. Este método es robusto ante la presencia de datos atípicos.

library(quantreg)
reg_cuant <- rq(pulso ~ edad, tau = 0.5)
summary(reg_cuant)
Call: rq(formula = pulso ~ edad, tau = 0.5)

tau: [1] 0.5

Coefficients:
            coefficients lower bd  upper bd 
(Intercept) 211.77551    200.13810 215.03918
edad         -0.81633     -1.03508  -0.60778

mod_nulo <- rq(pulso ~ 1, tau = 0.5)
anova(reg_cuant, mod_nulo)
Quantile Regression Analysis of Deviance Table

Model 1: pulso ~ edad
Model 2: pulso ~ 1
  Df Resid Df F value    Pr(>F)    
1  1       13  69.142 1.461e-06 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

library(rcompanion)
nagelkerke(reg_cuant)
$Models
                              
Model: "rq, pulso ~ edad, 0.5"
Null:  "rq, pulso ~ 1, 0.5"   

$Pseudo.R.squared.for.model.vs.null
                             Pseudo.R.squared
McFadden                             0.283743
Cox and Snell (ML)                   0.904171
Nagelkerke (Cragg and Uhler)         0.904403

$Likelihood.ratio.test
 Df.diff LogLik.diff  Chisq    p.value
      -1     -17.589 35.178 3.0094e-09

$Number.of.observations
         
Model: 15
Null:  15

$Messages
[1] "Note: For models fit with REML, these statistics are based on refitting with ML"

$Warnings
[1] "None"