Tutorial de tablas de contingencia en R

Cómo hacer una mesa

En primer lugar, obtengamos algunos datos. El paquete MASS contiene datos sobre 93 coches a la venta en EE.UU. en 1993. Se almacenan en el objeto Cars93 e incluyen 27 características para cada coche, algunas de las cuales son categóricas. Así que vamos a cargar el paquete MASS y ver el tipo de vehículos incluidos en cars93:

library(MASS)
Cars93$Type

##  [1] Small   Midsize Compact Midsize Midsize Midsize Large   Large  
##  [9] Midsize Large   Midsize Compact Compact Sporty  Midsize Van    
## [17] Van     Large   Sporty  Large   Compact Large   Small   Small  
## [25] Compact Van     Midsize Sporty  Small   Large   Small   Small  
## [33] Compact Sporty  Sporty  Van     Midsize Large   Small   Sporty
## [41] Sporty  Small   Compact Small   Small   Sporty  Midsize Midsize
## [49] Midsize Midsize Midsize Large   Small   Small   Compact Van    
## [57] Sporty  Compact Midsize Sporty  Midsize Small   Midsize Small  
## [65] Compact Van     Midsize Compact Midsize Van     Large   Sporty
## [73] Small   Compact Sporty  Midsize Large   Compact Small   Small  
## [81] Small   Compact Small   Small   Sporty  Midsize Van     Small  
## [89] Van     Compact Sporty  Compact Midsize
## Levels: Compact Large Midsize Small Sporty Van

Ahí tenemos 6 tipos de coches. La función table nos dice cuántos tenemos de cada tipo:

table(Cars93$Type)

##
## Compact   Large Midsize   Small  Sporty     Van
##      16      11      22      21      14       9

prop.table lo convierte en fracciones:

prop.table(table(Cars93$Type))

##
##    Compact      Large    Midsize      Small     Sporty        Van
## 0.17204301 0.11827957 0.23655914 0.22580645 0.15053763 0.09677419

Lo mismo ocurre con el origen de los coches:

table(Cars93$Origin)

##
##     USA non-USA
##      48      45

prop.table(table(Cars93$Origin))

##
##      USA  non-USA
## 0.516129 0.483871

Cómo hacer una tabla de contingencias

Genial, hemos visto que nuestro conjunto de datos contiene un número similar de coches estadounidenses y no estadounidenses y que los tipos más frecuentes son los medianos y los pequeños. Sin embargo, ¿quizás difieran en el tipo los estadounidenses y los no estadounidenses?

Veamos los tipos de coches en función de su origen. Podemos volver a utilizar table, pero ahora con dos argumentos. La primera se convertirá en variable de fila y la segunda en variable de columna:

table(Cars93$Type, Cars93$Origin)

##          
##           USA non-USA
##   Compact   7       9
##   Large    11       0
##   Midsize  10      12
##   Small     7      14
##   Sporty    8       6
##   Van       5       4

Ahora, vimos lo que todo el mundo sabe: A los estadounidenses les encantan los vehículos grandes. La tabla anterior muestra la distribución conjunta de dos variables categóricas (Type y Origin). Estas tablas se denominan tablas de contingencia.

Cómo obtener márgenes de la tabla de contingencia

rowSums y colSums se explican por sí solas.

(tab1<-table(Cars93$Type, Cars93$Origin))

##          
##           USA non-USA
##   Compact   7       9
##   Large    11       0
##   Midsize  10      12
##   Small     7      14
##   Sporty    8       6
##   Van       5       4

rowSums(tab1)

## Compact   Large Midsize   Small  Sporty     Van
##      16      11      22      21      14       9

colSums(tab1)

##     USA non-USA
##      48      45

Cómo obtener porcentajes de una tabla de contingencia

prop.table anidado con table da frecuencias

prop.table(table(Cars93$Type, Cars93$Origin))

##          
##                  USA    non-USA
##   Compact 0.07526882 0.09677419
##   Large   0.11827957 0.00000000
##   Midsize 0.10752688 0.12903226
##   Small   0.07526882 0.15053763
##   Sporty  0.08602151 0.06451613
##   Van     0.05376344 0.04301075

La conversión a porcentajes es simplemente multiplicar por 100:

prop.table(table(Cars93$Type, Cars93$Origin))*100

##          
##                 USA   non-USA
##   Compact  7.526882  9.677419
##   Large   11.827957  0.000000
##   Midsize 10.752688 12.903226
##   Small    7.526882 15.053763
##   Sporty   8.602151  6.451613
##   Van      5.376344  4.301075

Obsérvese que se trata de una distribución de probabilidad conjunta, a partir de la cual podemos ver que, por ejemplo, alrededor del 7,5% de los coches son pequeños y de origen estadounidense.

Más a menudo, nos interesa la distribución de una variable dentro de grupos creados por otra. Aquí parece interesante la distribución de los tipos de coches entre los estadounidenses y (por separado) los no estadounidenses. Para obtenerlo, utilizamos el argumento margin de la función prop.table. Indica en qué lugar de las filas (margin=1) o de las columnas (margin=2) se encuentra la variable de agrupación:

prop.table(table(Cars93$Type, Cars93$Origin), margin=2)*100

##          
##                 USA   non-USA
##   Compact 14.583333 20.000000
##   Large   22.916667  0.000000
##   Midsize 20.833333 26.666667
##   Small   14.583333 31.111111
##   Sporty  16.666667 13.333333
##   Van     10.416667  8.888889

Ahora podemos ver fácilmente que los coches pequeños son dos veces más frecuentes en la parte no estadounidense que en la estadounidense de nuestro conjunto de datos.

Observe también que los porcentajes suman 100 en las columnas, mientras que en la tabla de distribución conjunta (la que no tiene el argumento margin ), 100 era la suma de toda una tabla.

(tab2<-prop.table(table(Cars93$Type, Cars93$Origin), margin=2)*100)

##          
##                 USA   non-USA
##   Compact 14.583333 20.000000
##   Large   22.916667  0.000000
##   Midsize 20.833333 26.666667
##   Small   14.583333 31.111111
##   Sporty  16.666667 13.333333
##   Van     10.416667  8.888889

colSums(tab2)

##     USA non-USA
##     100     100

Prueba Chi-cuadrado

La pregunta más común que surge de las tablas de contingencia es si las variables de fila y columna son independientes. La forma más básica de responder es realizar una prueba chi-cuadrado. En este tutorial se explica con todo detalle. Comprobemos si Type y Origin son independientes:

chisq.test(Cars93$Type, Cars93$Origin)

## Warning in chisq.test(Cars93$Type, Cars93$Origin): Chi-squared
## approximation may be incorrect

##
##  Pearson's Chi-squared test
##
## data:  Cars93$Type and Cars93$Origin
## X-squared = 14.08, df = 5, p-value = 0.01511

Aparentemente, no lo son, pero también recibimos el aviso de Chi-squared approximation may be incorrect. Esto se debe a que la estadística chi-cuadrado sigue la distribución chi-cuadrado sólo aproximadamente. Cuantas más observaciones tengamos, mejor será la aproximación. La función chisq.test lanza la advertencia anterior siempre que uno de los recuentos esperados sea inferior a 5 (para saber qué es un "recuento esperado", véase el tutorial enlazado más arriba).

Prueba exacta de Fisher

La prueba exacta de Fisher es una alternativa a la prueba chi-cuadrado utilizada principalmente cuando una aproximación chi-cuadrado no es satisfactoria. Vamos a ejecutarlo:

fisher.test(Cars93$Type, Cars93$Origin)

##
##  Fisher's Exact Test for Count Data
##
## data:  Cars93$Type and Cars93$Origin
## p-value = 0.007248
## alternative hypothesis: two.sided

Los resultados son bastante similares a los de chi-cuadrado, pero no siempre es así. Uno de los principales inconvenientes de la prueba exacta de Fisher es que para tablas grandes (o muestras grandes) resulta ineficaz desde el punto de vista computacional.

Prueba G

Otra alternativa es la denominada prueba G. Su estadístico también tiene una distribución aproximada chi-cuadrado, pero para muestras pequeñas, esta aproximación es más cercana que la que utiliza la prueba chi-cuadrado. Para la prueba G podemos utilizar la función GTest del paquete DescTools. Los resultados vuelven a ser bastante similares a los de las dos pruebas anteriores: Type y Origin no son independientes.

library(DescTools)
GTest(Cars93$Type, Cars93$Origin)

##
##  Log likelihood ratio (G-test) test of independence without
##  correction
##
## data:  Cars93$Type and Cars93$Origin
## G = 18.362, X-squared df = 5, p-value = 0.002526

Corrección de Yates

En las tablas de contingencia 2x2, la prueba chi-cuadrado puede mejorarse con la corrección de continuidad de Yates. Simplemente resta 0,5 de cada término | Observado - Esperado | en una estadística chi-cuadrado. Si te sientes perdido, consulta de nuevo este tutorial. Además, R aplica automáticamente la corrección de Yates siempre que es necesario. Veamos la disponibilidad de versiones con transmisión manual de coches estadounidenses y no estadounidenses:

(tab3<-table(Cars93$Man.trans.avail, Cars93$Origin))

##      
##       USA non-USA
##   No   26       6
##   Yes  22      39

chisq.test(tab3)

##
##  Pearson's Chi-squared test with Yates' continuity correction
##
## data:  tab3
## X-squared = 15.397, df = 1, p-value = 8.712e-05

La corrección de Yates se aplicó automáticamente, y R nos lo comunicó.

Tabla de 3 (o más) dimensiones

Veamos Man.trans.avail, Origin, y Type a la vez. table divide el conjunto de datos según la variable proporcionada como tercera:

table(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)

## , ,  = Compact
##
##      
##       USA non-USA
##   No    2       0
##   Yes   5       9
##
## , ,  = Large
##
##      
##       USA non-USA
##   No   11       0
##   Yes   0       0
##
## , ,  = Midsize
##
##      
##       USA non-USA
##   No    9       4
##   Yes   1       8
##
## , ,  = Small
##
##      
##       USA non-USA
##   No    0       0
##   Yes   7      14
##
## , ,  = Sporty
##
##      
##       USA non-USA
##   No    0       0
##   Yes   8       6
##
## , ,  = Van
##
##      
##       USA non-USA
##   No    4       2
##   Yes   1       2

ftable ofrece una visión más compacta:

ftable(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)

##              Compact Large Midsize Small Sporty Van
##                                                    
## No  USA            2    11       9     0      0   4
##     non-USA        0     0       4     0      0   2
## Yes USA            5     0       1     7      8   1
##     non-USA        9     0       8    14      6   2

Cochran-Mantel-Haenszel Test

A partir de los resultados anteriores de table, debe quedar claro que la relación entre Origin y Man.trans.avail difiere de Type. En el caso de los coches pequeños y deportivos, por ejemplo, no existe ninguna asociación: todos los coches pequeños o deportivos de EE.UU., así como de otros países, vienen en versión manual. Por otro lado, la mayoría de los coches estadounidenses de tamaño medio no tienen versión manual, mientras que la mayoría de los coches de este tipo que no son estadounidenses sí la tienen. Para tener en cuenta posibles relaciones diferentes en los estratos, puede utilizarse la prueba de Cochran-Mantel-Haenszel:

mantelhaen.test(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)

##
##  Mantel-Haenszel chi-squared test with continuity correction
##
## data:  Cars93$Man.trans.avail and Cars93$Origin and Cars93$Type
## Mantel-Haenszel X-squared = 8.0153, df = 1, p-value = 0.004638
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
##   2.226531 76.891307
## sample estimates:
## common odds ratio
##          13.08438

El tercer argumento que se pasa a mantelhaen.test identifica los estratos. Compare los resultados anteriores con los obtenidos sin estratificación (ejemplo anterior de corrección de Yates). La asociación sigue presente, pero las pruebas son mucho más débiles.

Medidas de asociación

Una vez descubiertas algunas asociaciones entre variables, es hora de medir su fuerza. Existe una gran variedad de medidas posibles. Muchas de ellas se describen aquí. Ahora me centraré en dos de los más utilizados.

V de Cramer

V se basa en el estadístico chi-cuadrado:

$$ V = \sqrt{\frac{\chi^2/N}{\min(C-1, R-1)}}, $$ where:

• N es el total general de la tabla de contingencia (suma de todas sus celdas),
• C es el número de columnas
• R es el número de filas.

V ∈ [0; 1]. Cuanto mayor es V, más fuerte es la relación entre las variables. V= 0 puede interpretarse como independencia (ya que V= 0 si y sólo si ^χ2= 0). El principal inconveniente de V es la falta de una interpretación precisa. ¿V= 0,6 es una asociación fuerte, moderada o débil?

CramerV de DescTools puede calcularlo por nosotros:

CramerV(Cars93$Type, Cars93$Origin)

## [1] 0.3890967

Una vez más: ¿se trata de una asociación fuerte, moderada o débil?

Goodman y Kruskal lambda

La lambda de Goodman y Kruskal es un ejemplo de medida basada en la reducción proporcional de la variación. Este tipo de medidas pretende imitar el^R2 (coeficiente de determinación de la regresión lineal):

• toman valores de [0, 1]
• expresan una fracción de variación explicada por la variable independiente.

Utilicemos una variable-columna como independiente. La fórmula lambda de Goodman y Kruskal es entonces: $$\lambda = \frac{L - \sum_j L_j}{L},$$ donde:

• L_j es la suma de las frecuencias no modales de la columna j-ésima,
• L es una suma de frecuencias no modales en la columna "total".

Esto se explica mejor con un ejemplo. Veamos Type y Origin, este último independiente:

table(Cars93$Type, Cars93$Origin)

##          
##           USA non-USA
##   Compact   7       9
##   Large    11       0
##   Midsize  10      12
##   Small     7      14
##   Sporty    8       6
##   Van       5       4

El modo de la columna US es 11, por lo que _L1= 7 + 10 + 7 + 8 + 5 = 37. El modo de la columna no estadounidense es 14, por lo que _L2= 9 + 0 + 12 + 6 + 4 = 31. La columna "Total" es

rowSums(table(Cars93$Type, Cars93$Origin))

## Compact   Large Midsize   Small  Sporty     Van
##      16      11      22      21      14       9

El modo es 22, así que L= 16 + 11 + 21 + 14 + 9 = 71 y $$\lambda = \frac{71 - (37+31)}{71}=0.042$$

En su lugar, podemos utilizar la función Lambda del paquete DescTools:

Lambda(Cars93$Type, Cars93$Origin, direction='row')

## [1] 0.04225352

direction determina dónde (en row o en column) se encuentra la variable dependiente.

Origin sólo explica alrededor del 4% de la variación en Type. Obsérvese en la fórmula que lambda define la variación como una dicotomía entre pertenecer o no pertenecer al grupo mayor.

Cabe mencionar que lambda es cero siempre que la categoría modal de cada columna sea la misma. Echa un vistazo a esta tabla:

(lambdaTab<-cbind(c(0,100), c(49,51)))

##      [,1] [,2]
## [1,]    0   49
## [2,]  100   51

chisq.test(lambdaTab)

##
##  Pearson's Chi-squared test with Yates' continuity correction
##
## data:  lambdaTab
## X-squared = 62.279, df = 1, p-value = 2.981e-15

La categoría modal es la misma para cada columna (es la segunda fila). Así que lambda debe ser cero a pesar de la asociación significativa y visible:

Lambda(lambdaTab, direction='row')

## [1] 0

Si desea aprender más sobre estadística en R, siga el curso Statistical Modeling in R (Part 1) de DataCamp.