original dummyA dummyB dummyC
1 A 1 0 0
2 B 0 1 0
3 C 0 0 1
4 A 1 0 0
5 B 0 1 0
6 C 0 0 1
7 A 1 0 0
8 B 0 1 0
9 C 0 0 1
10 A 1 0 0
11 B 0 1 0
12 C 0 0 1
13 A 1 0 0
14 B 0 1 0
15 C 0 0 1
16 A 1 0 0
17 B 0 1 0
18 C 0 0 1
Lo que aquí se realizó fue crear variables indicadoras (dummy variables) a partir de nuestra variable politómica. En el data.frame podemos ver que la variable A es 1 en dummyA y 0 en dummyB y dummy C. La B es 1 en dummyB y 0 en dummyA y dummy C. Y la C es 1 en dummyC y 0 en dummyA y dummyB. Esta relación se puede observar de forma más clara en la siguiente tabla:
En este ejemplo utilizaremos la base de datos big5. Queremos comparar los niveles de neuroticismo en distintos países. Es decir, queremos examinar si las personas en distintos países tienen distintos niveles de neuroticismo.
La hipótesis a considerar son:
\[ H_0: \mu_1 = \mu_2 = \mu_3 \]
\[ H_a: \text{No todos los promedios poblacionales son iguales} \]
Entonces, la variable cuantitativa que usaremos en nuestro análisis será creada en base a preguntas de una escala de neuroticismo.
# Cargamos la libraría 'psych'library(psych)# Cargamos los datosbig5 <-read.csv("https://david-ti.github.io/introstats/data/big_five.csv")# Visualizamos los datoshead(big5)
vars n mean sd median trimmed mad min max range skew kurtosis se
X1 1 500 3.04 0.84 3.1 3.04 0.89 1 5 4 -0.01 -0.52 0.04
Veamos los grupos que compararemos en nuestro análisis.
# Cargamos la libraría 'dplyr'library(dplyr)
Attaching package: 'dplyr'
The following objects are masked from 'package:stats':
filter, lag
The following objects are masked from 'package:base':
intersect, setdiff, setequal, union
# Variable Categórica (país) table(big5$country)
AE AO AT AU BD BE BG BH BR CA CH CM CN DE DK DZ EG ES FI FR GB GH HK HR HU ID
1 1 2 79 2 2 1 1 1 10 3 1 2 5 4 1 15 3 4 3 40 4 1 1 2 5
IE IL IN IT JP KE LK LS MA MM MU MY NG NL NO NZ PH PK PL PT QA RO RS RW SA SE
4 1 31 7 3 9 2 1 4 1 4 9 9 4 2 21 22 9 1 2 1 4 1 1 3 4
SG SI TH TN TR TZ UG US ZA ZM ZW
3 1 2 2 1 1 2 89 42 1 2
# Seleccionamos tres paísesbig5 <-filter(big5, country =='PH'|country =='IN'|big5$country =='ZA')# Graficamosbarplot(table(big5$country), ylab ='Frecuencia', xlab ='Países')
El gráfico de barras anterior muestra la distribución de los países seleccionados en el análisis (PH, IN, ZA). Este paso permite observar el tamaño de cada grupo antes de hacer la comparación de los niveles de neuroticismo.
Comenzamos usando describeBy(y,x).
# Obtenemos los descriptivos por paísdescribeBy(big5$media_neuro, big5$country)
Descriptive statistics by group
group: IN
vars n mean sd median trimmed mad min max range skew kurtosis se
X1 1 31 3.15 0.76 3.2 3.16 0.59 1.5 4.8 3.3 -0.04 -0.32 0.14
------------------------------------------------------------
group: PH
vars n mean sd median trimmed mad min max range skew kurtosis se
X1 1 22 3.5 0.77 3.4 3.53 0.67 1.8 4.9 3.1 -0.2 -0.52 0.16
------------------------------------------------------------
group: ZA
vars n mean sd median trimmed mad min max range skew kurtosis se
X1 1 42 2.78 0.79 2.7 2.72 0.74 1.2 4.7 3.5 0.65 -0.06 0.12
Podemos también generar un boxplot con los mismos datos y visualizar las distribuciones para cada país.
# Generamos un boxplotboxplot(big5$media_neuro~big5$country, xlab ='País', ylab ='Neuroticismo')
Con los datos anteriores entonces podemos ilustrar el uso del modelo de regresión, evaluando si existen diferencias entre los grupos.
Tomando en cuenta que tenemos 3 grupos que debemos comparar, las comparaciones serían:
PH-ZA
PH-IN
ZA-IN
En primer lugar, crearemos las variables Dummy como aprendimos anteriormente.
Ahora que tenemos las variables Dummy, realizamos una primera regresión en la que fijaremos como intercepto a PH. Al fijarlo como intercepto, la media de PH en media_neuro será el valor de referencia para comparar con IN y ZA.
regresion1<-lm(media_neuro ~1+ dummyZA + dummyIN, data = big5)summary(regresion1)
Call:
lm(formula = media_neuro ~ 1 + dummyZA + dummyIN, data = big5)
Residuals:
Min 1Q Median 3Q Max
-1.70000 -0.53333 -0.04516 0.38575 1.91667
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.5000 0.1657 21.118 < 2e-16 ***
dummyZA -0.7167 0.2046 -3.503 0.000712 ***
dummyIN -0.3548 0.2167 -1.637 0.104958
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.7774 on 92 degrees of freedom
Multiple R-squared: 0.1216, Adjusted R-squared: 0.1025
F-statistic: 6.366 on 2 and 92 DF, p-value: 0.002574
A través de esta regresión podemos observar la comparación entre el intercepto (PH) con ZA y el intercepto (PH) con IN, por lo que nos falta la comparación entre ZA e IN.
Para evaluar esta relación realizamos una segunda regresión, pero esta vez fijaremos ZA como intercepto.
regresion2<-lm(media_neuro ~1+ dummyPH + dummyIN, data = big5)summary(regresion2)
Call:
lm(formula = media_neuro ~ 1 + dummyPH + dummyIN, data = big5)
Residuals:
Min 1Q Median 3Q Max
-1.70000 -0.53333 -0.04516 0.38575 1.91667
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.7833 0.1199 23.204 < 2e-16 ***
dummyPH 0.7167 0.2046 3.503 0.000712 ***
dummyIN 0.3618 0.1841 1.966 0.052346 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.7774 on 92 degrees of freedom
Multiple R-squared: 0.1216, Adjusted R-squared: 0.1025
F-statistic: 6.366 on 2 and 92 DF, p-value: 0.002574
Ahora que tenemos ambas regresiones, observemos que:
El valor F en las regresiones:
F-statistic: 6.366
El tamaño del efecto en las regresiones:
Multiple R-squared: 0.1216
Respecto a las diferencias de grupos:
La unica diferencia significativa que vemos en relación a la variable neuroticismo es entre los grupos PH y ZA, en los que se observa un p-value menor a 0.05.
Uso de ANOVA y Post Hoc: Ejemplo
También podemos ilustrar el uso de ANOVA. A través de esta prueba estadística veremos si existen diferencias estadísticamente significativas entre uno o más de estos grupos.
# Se genera un objeto, y se utiliza aov(y~x, data=basededatos)anova1 =aov(media_neuro~country, data=big5)# Ahora para saber el reporte del anova usamos summary del objeto (anova1)summary(anova1)
Df Sum Sq Mean Sq F value Pr(>F)
country 2 7.69 3.847 6.366 0.00257 **
Residuals 92 55.60 0.604
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Contrastaremos el p-value (pr(>F)) con nuestro alfa determinado por el grado de confianza de nuestra prueba estadística.
En este caso, el valor p es menor que nuestro alfa, entonces rechazamos la hipótesis nula.
El rechazo de la hipótesis nula indica que hay diferencias significativas en los niveles de neuroticismo entre al menos dos de los países analizados, lo que sugiere que no todos comparten la misma distribución para dicha escala.
ANOVA es solo un análisis global. Por ende, no nos dice entre cuáles grupos hay diferencias.
¿Cuál es el tamaño del efecto?
Podemos calcular el tamaño del efecto global en base a los resultados obtenidos anteriormente.
# Calculamos el tamaño del efecto dividiendo la suma de cuadrados intergrupal por la suma de cuadrados total (inter + intra)t_efecto =7.69/ (7.69+55.60)t_efecto
[1] 0.1215042
El valor del tamaño del efecto en este caso corresponde a 0.12 aproximadamente, lo cual nos indica que un 12% de la varianza de las personas en la variable neuroticismo puede ser explicada por la pertenencia a los distintos países.
Post Hoc
A través de análisis adicionales, conocidos como pruebas Post Hoc, podemos determinar entre qué grupos existen diferencias. Hay muchas pruebas distintas que pueden ser usadas para identificar qué diferencias específicas existen entre grupos. En este curso vamos a usar la prueba de Tukey.
# Realizamos una prueba Tukey considerando un 95% de confianzaTukeyHSD(anova1, conf.level = .95)
Tukey multiple comparisons of means
95% family-wise confidence level
Fit: aov(formula = media_neuro ~ country, data = big5)
$country
diff lwr upr p adj
PH-IN 0.3548387 -0.1614045 0.87108194 0.2352320
ZA-IN -0.3618280 -0.8003229 0.07666696 0.1265576
ZA-PH -0.7166667 -1.2040410 -0.22929237 0.0020403
Acá podemos visualizar los contrastes pareados entre todos los grupos. Podemos identificar aquellos contrastes que son significativos a través de un p-value < 0.05 (p adj).
En este caso, podemos visualizar que existen diferencias estadísticamente significativas entre los países ZA y PH.
Ahora que hemos realizado el ANOVA y las pruebas post hoc, podemos destacar lo siguiente:
El valor F es igual en ANOVA que en las regresiones: F-statistic: 6.366
El tamaño del efecto que calculamos en ANOVA es el mismo que encontramos en las regresiones: Multiple R-squared: 0.1216
Respecto a las diferencias de grupos:
La única diferencia significativa en los niveles de neuroticismo se encuentra entre los países PH y ZA, con un p-value menor a 0.05. Este resultado coincide tanto con los valores obtenidos en la regresión como con los hallazgos de las pruebas post hoc.
