De acuerdo con lo expuesto en el apartado anterior, para conocer la relación entre dos variables cualitativas (nominales u ordinales) se pueden utilizar:
Tablas de contingencia.
Gráficos.
Pruebas de asociación, por ejemplo, la prueba de \(\chi^2\).
Supongamos que queremos analizar las preferencias políticas en una región en la que existen 3 opciones. Un analísta político afirma que pueden existir diferencias en la preferencia de voto entre el ámbito urbano y el ámbito rural, por lo que el objetivo es determinar si esa hipótesis es verdad.
Los datos de partida son el número total de votos de cada una de las 3 opciones políticas (“P1”, “P2”, “P3”) en municipios rurales y urbanos.
datos <- matrix(c(1200, 900, 400, 680, 250, 650),
ncol=3,
byrow=TRUE)
colnames(datos) <- c("Conservador","Liberal","Socialdemócrata")
rownames(datos) <- c("Rural","Urbano")
El análisis de dos variables cualitativas se realiza sobre una tabla
de frecuencias, que es creada transformando la matriz de datos en una
tabla con as.table()
. En caso de trabajar con datos
procedentes de un dataframe, se utilizaría la función
table()
, como se mostró en el apartado 4.
tabla1 <- as.table(datos)
Como no volveremos a utlizar los datos originales, los borramos del Global Environment.
rm(datos)
Las frecuencias relativas, en términos de proporciones (tantos por uno), son:
prop.table(tabla1)
## Conservador Liberal Socialdemócrata
## Rural 0.29411765 0.22058824 0.09803922
## Urbano 0.16666667 0.06127451 0.15931373
En la imagen anterior también aparecen una columna a la derecha y un fila debajo con la denominación de “Total”; son las Frecuencias marginales. Además, la figura incluya una cifra situada en la esquina inferior derecha, denominada Gran total.
Para añadir las frecuencias marginales a la tabla de contingencia, se
debe recurrir a la función addmargins()
cuya sintaxis
básica es:
addmargins(tabla, margin=…, FUN=…)
donde:
tabla es una tabla de contingencia.
margin es un argumento puede tomar diferentes valores:
1 si queremos una nueva fila con las marginales de cada columna.
2 si queremos una nueva columna con las marginales de cada fila.
c(1,2), que es el valor por defecto para tablas de contingencia bidimensionales (y por lo tanto no hace falta especificarlo), si queremos las marginales por filas y por columnas.
FUN es la función que se aplica a las filas o
columnas para obtener el valor marginal. Por defecto es la suma
(sum()
), que es la función que nos interesa, y por tanto
tampoco hace falta especificarlo.
El resultado es otro objeto de la clase table
al que se
le han añadido una o varias filas o columnas, que contienen las
frecuencias marginales, tanto absolutas como relativas
addmargins(tabla1) # Frecuencias marginales absolutas
## Conservador Liberal Socialdemócrata Sum
## Rural 1200 900 400 2500
## Urbano 680 250 650 1580
## Sum 1880 1150 1050 4080
addmargins(prop.table(tabla1)) # Frecuencias marginales relativas
## Conservador Liberal Socialdemócrata Sum
## Rural 0.29411765 0.22058824 0.09803922 0.61274510
## Urbano 0.16666667 0.06127451 0.15931373 0.38725490
## Sum 0.46078431 0.28186275 0.25735294 1.00000000
En caso de que necesitásemos únicamente la fila o la columna de
frecuencias marginales, las funciones colSums()
y
rowSums()
crean tablas sólo con la suma por columnas o por
filas, respectivamente.
colSums(tabla1)
## Conservador Liberal Socialdemócrata
## 1880 1150 1050
rowSums(tabla1)
## Rural Urbano
## 2500 1580
📝 ACTIVIDAD DE EVALUACIÓN CONTINUA: :
Crea una tabla de contingencia relacionando las variables
género
y uso_movil_let
del dataframe
alumnos_UC
.
La primera actividad consistirá en importar el citado dataframe.
A continuación, se creará una tabla de contingencia que contendrá además las frecuencias marginales absolutas y otra que contendrá las frecuencias marginales relativas.
tabla_uc <- table(alumnos_uc$genero, alumnos_uc$uso_movil_let) # Creación de la tabla
addmargins(tabla_uc) # Frecuencias marginales absolutas
addmargins(prop.table(tabla_uc)) # Frecuencias marginales relativas
Para la visualización de las tablas de contingencia es posible utilizar diferentes tipos de gráficos.
La función balloonplot()
(paquete
gplots) dibuja una matriz gráfica donde cada celda
contiene un punto cuyo tamaño refleja la magnitud relativa del
componente correspondiente.
library("gplots")
balloonplot(tabla1,
main ="Opciones políticas según ámbito geográfico", # Título del gráfico
xlab ="Ámbito geográfico", # Etiqueta del eje X
ylab="Opción política", # Etiqueta del eje Y
label = TRUE, # Dibuja el valor numérico de los componentes de la tabla.
show.margins = TRUE) # Dibuja las frecuencias marginales.
Como ya sabemos, se obtiene cuando se aplica la función
barplot()
a una tabla con dos variables (bidimensional),
con dos posibilidades dependiendo del argumento beside
:
si es “FALSE” barras apiladas.
si “TRUE” barras yuxtapuestas.
OPCION 1: barras apiladas (por defecto)
par(mar = c(5, 5, 4, 11))
barplot(tabla1,
main = "Gráfico de barras apiladas",
sub = "Subtítulo",
xlab = "Etiqueta eje X",
ylab = "Etiqueta eje Y",
col = c("red", "blue"),
axes = TRUE,
ylim = c(0, 2000),
legend.text = rownames(tabla1), # Texto de la leyenda
args.legend = list(x = "topright", inset = c(-0.30, 0))) # Posición de la leyenda
par(mar = c(5, 4, 4, 2))
OPCION 2: barras yuxtapuestas
barplot(tabla1,
main = "Gráfico de barras yuxtapuestas",
ylim = c(0,1500),
xlab = "Etiqueta eje X",
ylab = "Etiqueta eje Y",
col = c("red", "blue"),
axes = TRUE,
beside = TRUE)
Es un diagrama que facilita la representación visual de la asociación
entre dos o más variables categóricas. En este caso, cada entrada de la
tabla de frecuencias corresponde a una región rectangular de área
proporcional a su valor, utilizando para ello la función
mosaicplot()
del paquete graphics.
Como argumentos adicionales se incluye:
shade
usado para colorear el gráfico.
Color azul: indica que los valores observados son más elevados que los esperados si estos fueran aleatorios
Color rojo: indica que los valores observados son más bajos que los esperados.
las = 2
dispone verticalmente las
etiquetas.
library("graphics")
mosaicplot(tabla1,
shade = TRUE,
las = 1)
La leyenda representada en el gráfico de mosaico muestra cuáles son las celdas más relevantes para los resultados de la prueba de chi-cuadrado. Contiene las mismas gamas de colores que en el caso del paquete graphics.
Existe otro paquete, denominado vcd, que también
puede ser usado para elaborar un gráfico de mosaico con la función
mosaic()
library("vcd")
Para la elaboración del diagrama de mosaico hay que señalar que la
primera variable mencionada en la sintaxis corresponde a las columnas
(vertical), mientras que la segunda variable corresponde a las filas
(horizontal). Esta función trabaja tanto con tables
como
con dataframes.
OPCION A: a partir de una tabla.
mosaic(tabla1, shade = TRUE,
labeling= labeling_border(rot_labels = c(90,0,0,0),
just_labels = c("left",
"center",
"center",
"center")))
OPCIÓN B: a partir de un dataframe.
#create data
df <- data.frame(producto = rep(c('TV', 'Radio', 'PC'), times=c(30, 14, 26)),
region = rep(c('A', 'B', 'C', 'D'), times=c(20, 15, 10, 25)))
mosaic(~producto + region, gp = shading_max,
data = df)
Los spineplots son una generalización de los gráficos de barras apilados. Un gráfico de mosaico es un spineplot para dos o más variables.
ATENCIÓN: se requiere una tabla de frecuencias
spineplot(tabla1)
Este gráfico puede ser personalizado, por ejemplo, trasponiendo las variables
spineplot(t(tabla1)) # Trasposición de variables
También puede modificarse el color de los cuadros y el borde de éstos
spineplot(tabla1,
col = c("red", "yellow", "blue"),
border = "darkblue") # Color del borde de las barras
Por defecto, la función spineplot crea un spineplot, pero el
argumento off = 0
crea un espinograma, eliminando la
distancia entre las barras.
spineplot(tabla1, off = 0)
📝 ACTIVIDAD DE EVALUACIÓN CONTINUA: :
A partir de la tabla creada a partir de las variables
genero
y `uso_movil_let”, elabora
ballplot
.balloonplot(tabla_uc,
main ="Relación entre género y uso de móvil",
xlab ="Género",
ylab="Uso de móvil",
label = FALSE,
show.margins = FALSE)
barplot(tabla_uc,
main = "Relación entre género y uso de móvil",
sub = "Subtítulo",
ylim = c(0,25),
xlab = "Categorías de uso",
ylab = "Número de alumnos",
axes = TRUE,
beside = TRUE,
col = c("red", "blue"),
legend.text = rownames(tabla_uc),
args.legend = list(x = "topright"))
mosaicplot(tabla_uc,
shade = TRUE,
las = 1)
Una vez tabulados los datos y visualizadas las posibles diferencias entre variables, el tercer paso es verificar estadísticamente esas posibles diferentes. Para ello recurrimos a varias pruebas.
La prueba \(\chi^2\) de Pearson tiene varios usos:
Opción 1: Evaluar si hay una asociación significativa entre dos o más variables cualitativas en una tabla de frecuencias.
Opción 2: determinar si un variable categórica sigue una distribución hipotética de valores esperados obtenidos de acuerdo a una hipótesis de proporcionalidad (bondad de ajuste).
En el caso que nos ocupa se utiliza la primera opción, siendo denominada prueba de “contraste de independencia” para tablas de contingencia bidimensionales. Consiste en decidir si las dos variables de la tabla tienen distribuciones independientes (es decir, si la distribución de probabilidades conjunta es igual al producto de las probabilidades marginales).
Si la proporción de individuos en cada columna varía entre las diversas filas y viceversa, se dice que existe asociación entre las dos variables.
Si no existe asociación se dice que ambas variables son independientes.
La prueba de \(\chi^2\) presenta algunas características singulares:
Test no paramétrico: no requiere el conocimiento de ciertas características de la población de la que se han extraído las muestras.
Las categorías deben ser por lo menos dos y excluyentes (ningún individuo puede estar en más de una categoría a la vez.
Además, requiere varios pasos.
⚠️ IMPORTANTE:
PRIMER PASO: definir las hipótesis nula y alternativa
Hipótesis nula (H0): las variables son independientes (no hay relación entre ellas).
Hipótesis alternativa (H1): las variables no son independientes (están relacionadas).
SEGUNDO PASO: decidir el valor de α
. Este valor
indica la probabilidad de obtener un estadístico de chi-cuadrado tan
extremo como el observado, asumiendo que la hipótesis nula es verdadera,
es decir, decidir qué riesgo se asume para llegar a una conclusión
errónea. Normalmente α=0,05 para las pruebas de independencia, es decir,
existe un riesgo del 5 % de concluir que las dos variables son
independientes, cuando en realidad no lo son. En este caso:
Si el valor-p α>0,05 se acepta la Hipótesis nula (H0): las dos variables son independientes (no hay relación entre ellas).
Si el valor-p α<0,05 se rechaza la hipótesis nula y se acepta Hipótesis alternativa (H1): las dos variables están relacionadas (no son independientes).
Para obtener de forma manual el valor de \(\chi^2\) comenzamos calculando cálculo de los valores esperados de cada celda de la tabla bajo la hipótesis nula, de acuerdo con la fórmula
El test de independencia usa las frecuencias absolutas esperadas bajo la hipótesis nula de independencia, que se obtienen, para cada celda (i,j), multiplicando la frecuencia marginal de la fila i por la de la columna j y dividiendo por el tamaño de la muestra. En nuestro ejemplo estas frecuencias esperadas son
Por último, se calcula el estadístico comparando los valores observados (O) y esperados (E), según la siguiente fórmula:
Un valor elevado de χ2 indica que los valores observados son muy distintos de los valores esperados, y por lo tanto mayor será la evidencia en contra de la hipótesis nula.
Cuanto más pequeño sea el valor del estadístico, más se parecen los datos observados a los datos esperados.
Si el resultado de la prueba chi-cuadrado es 0, los valores observados y los esperados son exactamente iguales.
Si el resultado es significativo, hay una asociación real entre ambas variables (es decir confirmamos la Hipótesis Alternativa). Si no es significativo, podríamos decir que no hay evidencia suficiente para sugerir una asociación (confirmaría la Hipótesis Nula)
chiq.test()
para el cálculo de \(\chi^2\)La sintaxis de la prueba chi-cuadrado en R es la siguiente
chisq.test(tabla1)
##
## Pearson's Chi-squared test
##
## data: tabla1
## X-squared = 382.76, df = 2, p-value < 2.2e-16
Como el p-valor
es inferior a 0.05, se rechaza la
hipótesis nula y se acepta la alternativa: las variables objeto de
estudio están relacionadas.
La función chisq.test()
ofrece una serie de outputs
adicionales, como son los siguientes
chisq.test(tabla1)$observed # Valores observados
## Conservador Liberal Socialdemócrata
## Rural 1200 900 400
## Urbano 680 250 650
chisq.test(tabla1)$expected # Valores esperados bajo la hipótesis nula
## Conservador Liberal Socialdemócrata
## Rural 1151.9608 704.6569 643.3824
## Urbano 728.0392 445.3431 406.6176
chisq.test(tabla1)$parameter # Grados de libertad
## df
## 2
chisq.test(tabla1)$p.value # Probabilidad
## [1] 7.686217e-84
qchisq(0.95, # Valor critico de chi-cuadrado, al 95% con 1 grado de libertad
1) # Grados de libertad, calculados a partir de n-1 (numero de categorias menos uno).
## [1] 3.841459
En el caso de tablas formadas por 2 filas y 2 columnas, se debe aplicar una corrección de continuidad.
chisq.test(tabla1,
correct=FALSE) # Indicador lógico para aplicar una corrección de continuidad
##
## Pearson's Chi-squared test
##
## data: tabla1
## X-squared = 382.76, df = 2, p-value < 2.2e-16
Si algunas frecuencias absolutas esperadas fueran inferiores a 5, la
aproximación del p-valor por una distribución \(\chi^2\) podría no ser adecuada. En este
caso, al ser las variables cualitativas, no podemos recurrir al
agrupamiento de valores consecutivos, puesto que no tienen orden. Si se
da esta situación, lo mejor es recurrir a simular el p-valor usando el
argumento simulate.p.value=TRUE
.
set.seed(123)
chisq.test(tabla1,
correct=FALSE,
simulate.p.value=TRUE, # Calcula los valores de probabilidad mediante una simulación Monte Carlo
B=5000)$p.value # Número de réplicas simuladas y probabilidad
## [1] 0.00019996
Si se desea conocer qué celdas contribuyen más a los valores de la prueba, se calculan los resíduos de Pearson correspondientes a cada celda (residuos estandarizados)
round(chisq.test(tabla1)$residuals, 3)
## Conservador Liberal Socialdemócrata
## Rural 1.415 7.359 -9.595
## Urbano -1.780 -9.257 12.070
Finalmente, es posible visualizar estos resíduos con la función
corrplot()
del paquete homónimo: el tamaño del círculo es
proporcional a la contribución de cada celda.
library(corrplot)
corrplot(chisq.test(tabla1)$residuals,is.cor = FALSE)
Otra opción en caso de frecuencias menores a 5 en tablas de contingencia es el test de Fisher.
#?fisher.test
La prueba exacta de Fisher utiliza las siguientes hipótesis nula y alternativa:
Hipótesis nula H0: las dos variables son independientes.
Hipótesis alternativa H1: las dos variables no son independientes (existe relación entre ellas).
fisher.test(tabla1)
##
## Fisher's Exact Test for Count Data
##
## data: tabla1
## p-value < 2.2e-16
## alternative hypothesis: two.sided
En la prueba exacta de Fisher, la hipótesis nula es que las dos columnas son independientes (o, de manera equivalente, que la odds ratio -razón de posibilidades- es igual a 1).
Para determinar si las dos columnas son independientes, se utiliza el
p-value
de la prueba. Si el valor es inferior a 0.05, no
rechazamos la hipótesis nula: no tenemos evidencia suficiente para decir
que existe una asociación estadísticamente significativa entre ambas
variables.
fisher.test(tabla1,
simulate.p.value=TRUE,
B=2^3)
##
## Fisher's Exact Test for Count Data with simulated p-value (based on 8
## replicates)
##
## data: tabla1
## p-value = 0.1111
## alternative hypothesis: two.sided
El grado de asociación entre dos variables se puede cuantificar adicionalmente empleando otros coeficientes:
El coeficiente phi es el valor de la chi cuadrado entre el número de observaciones. En general, cuanto más lejos de 0, más fuerte será la relación entre las dos variables.
Un valor próximo a 0 indica independencia entre los factores.
Valores próximos a -1 indica una relación perfectamente negativa entre las dos variables.
Valores próximos a 1 indica una relación perfectamente positiva entre las dos variables.
El coeficiente de contingencia también es una medida de la intensidad de la relación basado en la chi cuadrado y toma valores entre 0 (independencia) y 1 (dependencia).
La V de Cramer es menos susceptible a valores muestrales.
0 implica independencia
1 una relación perfecta entre los factores. Habitualmente valores superiores a 0,30 ya nos están indicando que hay una posible relación entre las variables; un valor de 0.5 o superior ya indicaría una gran relación entre ambas variables.. En este caso tenemos valores muy próximos a 0 en todos los estadísticos que nos ofrece, ambos factores son independientes. Saludos.
library(vcd)
assocstats(tabla1)
## X^2 df P(> X^2)
## Likelihood Ratio 386.57 2 0
## Pearson 382.76 2 0
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.293
## Cramer's V : 0.306
📝 ACTIVIDAD DE EVALUACIÓN CONTINUA: :
PRIMER EJERCICIO. En un análisis sobre prevalencia de enfermedades cardiacas, se pregunta a 100 personas su grado de actividada física. Los resultados son los siguientes:
Activo | Sedentario | |
---|---|---|
Con enfermedad | 20 | 30 |
Sin enfermedad | 40 | 10 |
ejercicio1 <- matrix(c(20, 30, 40, 10),
ncol=2,
byrow=TRUE)
colnames(ejercicio1) <- c("Activo","Sedentario")
rownames(ejercicio1) <- c("Con_enfermedad","Sin_enfermedad")
tabla_ejercicio1 <- as.table(ejercicio1)
addmargins(tabla_ejercicio1) # Frecuencias marginales absolutas
tabla_ejercicio1
.mosaic(tabla_ejercicio1, shade = TRUE)
# Hipótesis Nula (H0): El nivel de actividad física y la incidencia de enfermedades cardíacas son independientes entre sí.
# Hipótesis Alternativa (H1): El nivel de actividad física y la incidencia de enfermedades cardíacas están asociados.
chisq.test(tabla_ejercicio1,
correct=FALSE) # Indicador lógico para aplicar una corrección de continuidad
fisher.test(tabla_ejercicio1,
simulate.p.value=TRUE,
B=2^3)
# Dado que el valor p es menor que el nivel de significancia comúnmente utilizado de 0.05, se rechaza la hipótesis nula: hay una asociación
# significativa entre el nivel de actividad física y la incidencia de enfermedades cardíacas.
SEGUNDO EJERCICIO
zonas_verdes.RData
. A partir de ese
dataframe convierte la variable obras
en un factor y
transforma las etiquetas originales (0 = “sin obras”, 1, = “en
obras”).load("D:/G2040/TEMA_3_ESDA_Graficos/datos/zonas_verdes.Rdata")
zonas_verdes$obras <- as.factor(zonas_verdes$obras)
tabla2 <- table(zonas_verdes$obras, zonas_verdes$barrio)
rownames(tabla2) <- c("sin obras","con obras")
tabla2
prop.table(tabla2) * 100
addmargins(tabla2, margin = 1)
addmargins(tabla2, margin = 2)
par(mar = c(5, 5, 4, 11))
barplot(tabla2,
main = "Gráfico de barras apiladas",
sub = "Subtítulo",
xlab = "Etiqueta eje X",
ylab = "Etiqueta eje Y",
col = c("red", "blue"),
axes = TRUE,
ylim = c(0, 20),
legend.text = rownames(tabla2), # Texto de la leyenda
args.legend = list(x = "topright", inset = c(-0.30, 0))) # Posición de la leyenda
par(mar = c(5, 4, 4, 2))
chisq.test(tabla2)
TERCER EJERCICIO
periodo
deberá ser
convertida en factor para su uso posterior.zonas_verdes <- within(zonas_verdes,{
periodo <- NA # Crea la variable grupo, pero vacía
periodo[fecha_inauguracion > "1990-01-01"] <- "Nuevos" # Los restantes comandos rellenan esa variable de manera ordenada
periodo[fecha_inauguracion <= "1990-01-01"] <- "Antiguos"})
zonas_verdes$periodo <- as.ordered(zonas_verdes$periodo)
zonas_verdes$obras <- factor(zonas_verdes$obras,
levels = c(0,1),
labels = c("Sin obras", "Con obras"))
tabla3 <- addmargins(table(zonas_verdes$obras, zonas_verdes$periodo), margin = c(1,2))
tabla3
mosaicplot(tabla3,
shade = TRUE,
las = 1)
chisq.test(tabla3)
assocstats(tabla3)
CUARTO EJERCICIO
El dataset “mtcars” contiene los detalles técnicos de los conches disponibles en un concesionario de automóviles. De ese dataset nos interesan las siguientes variables:
vs: tipo de motor (0 = en V, 1 = lineal)
ame: tipo de transmisión (0 = automática, 1 = manual).
gear: número de marchas (excluyendo la marcha atrás).
carb: número de carburadores.
Carga este dataset en la memoria
data(mtcars)
str(mtcars)
## 'data.frame': 32 obs. of 11 variables:
## $ mpg : num 21 21 22.8 21.4 18.7 18.1 14.3 24.4 22.8 19.2 ...
## $ cyl : num 6 6 4 6 8 6 8 4 4 6 ...
## $ disp: num 160 160 108 258 360 ...
## $ hp : num 110 110 93 110 175 105 245 62 95 123 ...
## $ drat: num 3.9 3.9 3.85 3.08 3.15 2.76 3.21 3.69 3.92 3.92 ...
## $ wt : num 2.62 2.88 2.32 3.21 3.44 ...
## $ qsec: num 16.5 17 18.6 19.4 17 ...
## $ vs : num 0 0 1 1 0 1 0 1 1 1 ...
## $ am : num 1 1 1 0 0 0 0 0 0 0 ...
## $ gear: num 4 4 4 3 3 3 3 4 4 4 ...
## $ carb: num 4 4 1 1 2 1 4 2 2 4 ...
vsAm <- table(mtcars$vs, mtcars$am)
vsAm
gearCarb <- table(mtcars$gear,mtcars$carb)
gearCarb
cylGear <- table(mtcars$cyl,mtcars$gear)
cylGear
chisq.test(vsAm)
chisq.test(gearCarb)
chisq.test(cylGear)
mosaicplot(cylGear,
shade = TRUE,
las = 1)
corrplot(chisq.test(cylGear)$residuals,is.cor = FALSE)