R: Tabulación cruzada y prueba chi-cuadrado con gmodels

Problema

En R deseamos crear una tabla de referencias cruzadas similar a la producida en SAS (SAS Proc Freq) o crosstabs SPSS.

Solución

Utilizamos la función CrossTable del paquete gmodels.

SAS

CrossTable(mtcars$am, mtcars$gear, dnn = c("am", "gear"), expected = TRUE, format = "SAS")

   Cell Contents
|-------------------------|
|                       N |
|              Expected N |
| Chi-square contribution |
|           N / Row Total |
|           N / Col Total |
|         N / Table Total |
|-------------------------|

 
Total Observations in Table:  32 

 
             | gear 
          am |         3 |         4 |         5 | Row Total | 
-------------|-----------|-----------|-----------|-----------|
           0 |        15 |         4 |         0 |        19 | 
             |     8.906 |     7.125 |     2.969 |           | 
             |     4.169 |     1.371 |     2.969 |           | 
             |     0.789 |     0.211 |     0.000 |     0.594 | 
             |     1.000 |     0.333 |     0.000 |           | 
             |     0.469 |     0.125 |     0.000 |           | 
-------------|-----------|-----------|-----------|-----------|
           1 |         0 |         8 |         5 |        13 | 
             |     6.094 |     4.875 |     2.031 |           | 
             |     6.094 |     2.003 |     4.339 |           | 
             |     0.000 |     0.615 |     0.385 |     0.406 | 
             |     0.000 |     0.667 |     1.000 |           | 
             |     0.000 |     0.250 |     0.156 |           | 
-------------|-----------|-----------|-----------|-----------|
Column Total |        15 |        12 |         5 |        32 | 
             |     0.469 |     0.375 |     0.156 |           | 
-------------|-----------|-----------|-----------|-----------|

 
Statistics for All Table Factors


Pearson's Chi-squared test 
------------------------------------------------------------
Chi^2 =  20.94467     d.f. =  2     p =  2.830889e-05 

SPSS

CrossTable(mtcars$am, mtcars$gear, dnn = c("am", "gear"), expected = TRUE, format = "SPSS")

   Cell Contents
|-------------------------|
|                   Count |
|         Expected Values |
| Chi-square contribution |
|             Row Percent |
|          Column Percent |
|           Total Percent |
|-------------------------|

Total Observations in Table:  32 

             | gear 
          am |        3  |        4  |        5  | Row Total | 
-------------|-----------|-----------|-----------|-----------|
           0 |       15  |        4  |        0  |       19  | 
             |    8.906  |    7.125  |    2.969  |           | 
             |    4.169  |    1.371  |    2.969  |           | 
             |   78.947% |   21.053% |    0.000% |   59.375% | 
             |  100.000% |   33.333% |    0.000% |           | 
             |   46.875% |   12.500% |    0.000% |           | 
-------------|-----------|-----------|-----------|-----------|
           1 |        0  |        8  |        5  |       13  | 
             |    6.094  |    4.875  |    2.031  |           | 
             |    6.094  |    2.003  |    4.339  |           | 
             |    0.000% |   61.538% |   38.462% |   40.625% | 
             |    0.000% |   66.667% |  100.000% |           | 
             |    0.000% |   25.000% |   15.625% |           | 
-------------|-----------|-----------|-----------|-----------|
Column Total |       15  |       12  |        5  |       32  | 
             |   46.875% |   37.500% |   15.625% |           | 
-------------|-----------|-----------|-----------|-----------|

 
Statistics for All Table Factors


Pearson's Chi-squared test 
------------------------------------------------------------
Chi^2 =  20.94467     d.f. =  2     p =  2.830889e-05 


 
       Minimum expected frequency: 2.03125 
Cells with Expected Frequency < 5: 3 of 6 (50%)

Notas

Argumentos:

expected - si es TRUE para incluir la prueba chi-cuadrado de Pearson.
format - indicamos el formato de la salida de datos, SAS —por defecto— o SPSS.
dnn - nombres de las dimensiones en el resultado.

Referencias

• CrossTable
• gmodels
• Tablas de contingencia en R
• Tablas de proporciones en R
• Transformar tabla de contingencia en tabla de de frecuencias

Project Euler - Problema 5 en R

Continuamos con los problemas planteados en Project Euler.

Problema - Smallest multiple

2520 is the smallest number that can be divided by each of the numbers from 1 to 10 without any remainder.
What is the smallest positive number that is evenly divisible by all of the numbers from 1 to 20?

2520 es el número más pequeño que puede ser dividido por cada uno de los números del 1 al 10 sin ningún resto.
¿Cuál es el número positivo más pequeño que es divisible por todos los números del 1 al 20?

Solución

require(numbers)
mLCM(1:20)

[1] 232792560

Notas

Usamos la función mLCM para calcular el mínimo común múltiplo de los números del 1 al 20.

Otra alternativa manual sería calcular los factores primos de cada uno de los números. Después inspeccionamos visualmente los factores comunes de todos ellos.

x <- matrix(1:20)
apply(x, 1, primeFactors)

[[1]]
[1] 1

[[2]]
[1] 2

[[3]]
[1] 3

[[4]]
[1] 2 2

[[5]]
[1] 5

[[6]]
[1] 2 3

[[7]]
[1] 7

[[8]]
[1] 2 2 2

[[9]]
[1] 3 3

[[10]]
[1] 2 5

[[11]]
[1] 11

[[12]]
[1] 2 2 3

[[13]]
[1] 13

[[14]]
[1] 2 7

[[15]]
[1] 3 5

[[16]]
[1] 2 2 2 2

[[17]]
[1] 17

[[18]]
[1] 2 3 3

[[19]]
[1] 19

[[20]]
[1] 2 2 5

Los factores serían: (2^4)*(3^2)*5*7*11*13*17*19 = 232792560

Para calcular los factores primos únicos

unique(unlist(apply(x, 1, primeFactors)))

[1]  1  2  3  5  7 11 13 17 19

Referencias

• Problem 5

Transformar tabla de contingencia en tabla de de frecuencias

Problema

Deseamos transformar una tabla de contingencia en una tabla de frecuencias.

# Tabla de contigencia
tabla <- table(mtcars[, c("am", "gear")])

   gear
am   3  4  5
  0 15  4  0
  1  0  8  5

Tabla de frecuencias

tabla <- table(mtcars[, c("am", "gear")])
cuenta <- as.data.frame(tabla)
cuenta

  am gear Freq
1  0    3   15
2  1    3    0
3  0    4    4
4  1    4    8
5  0    5    0
6  1    5    5

Transformar tabla de frecuencias en tabla de contingencia

ftable(xtabs(Freq ~ am + gear, data = cuenta)) 

   gear  3  4  5
am              
0       15  4  0
1        0  8  5

Es el equivalente a:

ftable(mtcars[, c("am", "gear")])

Referencias

• Tablas de contingencia en R
• Tablas de proporciones en R

Tablas de proporciones en R

Problema

Deseamos crear una tabla de proporciones.

Solución

• Una dimensión

tabla <- table(mtcars$am)
prop.table(tabla)

      0       1 
0.59375 0.40625

• Dos dimensiones

tabla <- table(mtcars[, c("am", "gear")])
prop.table(tabla)

   gear
am        3       4       5
  0 0.46875 0.12500 0.00000
  1 0.00000 0.25000 0.15625

La función prop.table tiene dos argumentos:
• x, tabla creada anteriormente con la función table
• margin, con tres posibles valores:
  Nulo, calcula x/suma(x), como en el ejemplo anterior.
  1, calcula la proporción por filas.
  2, calcula la proporción por columnas.


# Por filas
prop.table(tabla, 1)

   gear
am          3         4         5
  0 0.7894737 0.2105263 0.0000000
  1 0.0000000 0.6153846 0.3846154

# Por columnas
prop.table(tabla, 2)

   gear
am          3         4         5
  0 1.0000000 0.3333333 0.0000000
  1 0.0000000 0.6666667 1.0000000

• Tres dimensiones

tabla <- table(mtcars[, c("am", "gear", "cyl")])
prop.table(tabla)

, , cyl = 4

   gear
am        3       4       5
  0 0.03125 0.06250 0.00000
  1 0.00000 0.18750 0.06250

, , cyl = 6

   gear
am        3       4       5
  0 0.06250 0.06250 0.00000
  1 0.00000 0.06250 0.03125

, , cyl = 8

   gear
am        3       4       5
  0 0.37500 0.00000 0.00000
  1 0.00000 0.00000 0.06250

• Tabla de contingencia plana

En caso anterior sería más claro crear una tabla de contingencia plana (Flat Contingency Table).

tabla <- ftable(mtcars[, c("am", "gear", "cyl")])
prop.table(tabla)

        cyl       4       6       8
am gear                            
0  3        0.03125 0.06250 0.37500
   4        0.06250 0.06250 0.00000
   5        0.00000 0.00000 0.00000
1  3        0.00000 0.00000 0.00000
   4        0.18750 0.06250 0.00000
   5        0.06250 0.03125 0.06250

• Tabla en porcentajes

Empleamos la función round.

round(prop.table(tabla)*100, 2)

         cyl     4     6     8
am gear                      
0  3         3.12  6.25 37.50
   4         6.25  6.25  0.00
   5         0.00  0.00  0.00
1  3         0.00  0.00  0.00
   4        18.75  6.25  0.00
   5         6.25  3.12  6.25

round(prop.table(tabla, 1)*100, 2) # Por fila, am y gear.

        cyl     4     6     8
am gear                      
0  3         6.67 13.33 80.00
   4        50.00 50.00  0.00
   5          NaN   NaN   NaN
1  3          NaN   NaN   NaN
   4        75.00 25.00  0.00
   5        40.00 20.00 40.00

round(prop.table(tabla, 2)*100, 2) # Por columna, cilindro

        cyl     4     6     8
am gear                      
0  3         9.09 28.57 85.71
   4        18.18 28.57  0.00
   5         0.00  0.00  0.00
1  3         0.00  0.00  0.00
   4        54.55 28.57  0.00
   5        18.18 14.29 14.29

Referencias

• Tablas de contingencia en R
• Transformar tabla de contigencia en tabla de de frecuencias

Tablas de contingencia en R

Problema

Deseamos crear una tabla de contingencia para uno o varios factores.

Solución

• Una dimensión

table(mtcars$am)

 0  1 
19 13 

• Dos dimensiones

table(mtcars$am, mtcars$gear)

Como hemos especificado dos vectores, no aparece el nombre de las dimensiones: am y gear.

     3  4  5
  0 15  4  0
  1  0  8  5

table(mtcars[, c("am", "gear")]) 
tabla <- table(mtcars[, 9:10])
# O usando el argumento dnn:
table(mtcars$am, mtcars$gear, dnn = c("am", "gear"))

   gear
am   3  4  5
  0 15  4  0
  1  0  8  5

• Tres dimensiones

table(mtcars[, c("am", "gear", "cyl")])

, , cyl = 4

   gear
am   3  4  5
  0  1  2  0
  1  0  6  2

, , cyl = 6

   gear
am   3  4  5
  0  2  2  0
  1  0  2  1

, , cyl = 8

   gear
am   3  4  5
  0 12  0  0
  1  0  0  2

• Tabla de contingencia plana

En el caso anterior sería más claro crear una tabla de contingencia plana (Flat Contingency Table).

ftable(mtcars[, c("am", "gear", "cyl")])

        cyl  4  6  8
am gear             
0  3         1  2 12
   4         2  2  0
   5         0  0  0
1  3         0  0  0
   4         6  2  0
   5         2  1  2

Podemos usar los argumentos row.vars y col.vars para especificar las filas y columnas que deseamos en la tabla de contingencia. Si ninguno es especificado, la última variable es usada para las columnas. Como en nuestro caso cyl.

ftable(mtcars[, c("am", "gear", "cyl")], col.vars = c(1, 2))

     am    0        1      
    gear  3  4  5  3  4  5
cyl                       
4         1  2  0  0  6  2
6         2  2  0  0  2  1
8        12  0  0  0  0  2

Alternativa

La función xtabs crea tablas de contingencia mediante el formato de fórmula, con cada variable separada por un +.

# Una dimensión
xtabs(~ am, mtcars)
# Dos dimensiones
xtabs(~ am + gear, mtcars)
# Tres dimensiones
xtabs(~ am + gear + cyl, mtcars)
# Tabla de contingencia plana
ftable(xtabs(~ am + gear + cyl, mtcars))

Referencias

• Tablas de proporciones en R
• Transformar tabla de contigencia en tabla de de frecuencias

Modo predeterminado de apertura de Ms Access: compartido o exclusivo

Problema

Deseamos cambiar el modo predeterminado de apertura de nuestros ficheros de Ms Access.

Solución

1. En la ficha Archivo, o en clic el botón de Microsoft Office (dependiendo de la versión) y, a continuación, clic en Opciones.

2. En el panel izquierdo del cuadro de diálogo Opciones de Access, clic en Configuración de cliente. En el panel derecho, en la sección Avanzadas, en Modo predeterminado de apertura, seleccionamos Exclusivo.

3. Clic en Aceptar y salimos de Access. A partir de ese momento, cualquier fichero de Access se abrirá en modo exclusivo.

Notas

Por defecto, el modo predeterminado de apertura es el compartido. En el ejemplo anterior lo hemos cambiado a exclusivo. Para deshacer el cambio, repetimos el procedimiento y elegimos Compartido.

Referencias

• Abrir fichero de Ms Access en modo exclusivo

Abrir fichero de Ms Access en modo exclusivo

Si queremos abrir un fichero de Access en modo exclusivo seguimos los siguientes pasos:

1. En la ficha Archivo, clic en Abrir.

2. En el cuadro de diálogo Abrir buscamos y seleccionamos el fichero de Access, clic en la flecha situada junto al botón Abrir y, a continuación, clic en Abrir en modo exclusivo.

3. Tendremos acceso a la base de datos en modo exclusivo y ningún usuario podrá abrirla mientras la mantengamos abierta.

Referencias

Etiquetar selectivamente los puntos de un diagrama de dispersión con ggplot2

Problema

Tenemos un diagrama de dispersión.

qplot(dpi, sr, data = LifeCycleSavings)

Si etiquetáramos todos los puntos obtendríamos un gráfico con muchas de las etiquetas superpuestas.

qplot(dpi, sr, data = LifeCycleSavings, label = rownames(data))+
  geom_text(vjust = 1.5)

Solución

Opción 1

LifeCycleSavings$lab <- rownames(LifeCycleSavings)
newlab <- subset(LifeCycleSavings, lab %in% c("Canada", "Sweden", "United Kingdom", "United States", "Iceland", "Japan", "Spain", "Netherlands", "Switzerland", "Denmark"))

1. lab, creamos una nueva columna con los nombres de las filas.
2. newlab, subconjunto de las filas que deseamos etiquetar.
3. Utilizamos lab dentro de aes y las nuevas, newlab, en geom_text.

• qplot

qplot(dpi, sr, data = LifeCycleSavings, label = lab)+  
  geom_text(data = newlab , vjust = 1.5)

• ggplot

ggplot(LifeCycleSavings, aes(x = dpi, y = sr, label = lab)) + 
  geom_point() + 
  geom_text(data = newlab , vjust = 1.5)

Opción 2

Name <- rownames(LifeCycleSavings)
idx <- Name %in% c("Canada", "Sweden", "United Kingdom", "United States", "Iceland", "Japan", "Spain", "Netherlands", "Switzerland", "Denmark") 
Name[!idx] <- NA

1. Name, vector con los nombres de las etiquetas: nombres de las filas.
2. idx, vector lógico (verdadero o falso) con TRUE para las etiquetas que deseamos conservar.
3. Asignamos al vector Name el valor NA para aquellos nombres no presentes en el vector anterior.

• qplot

qplot(dpi, sr, data = LifeCycleSavings, label = Name)+  
  geom_text(vjust = 1)

• ggplot

p <- ggplot(LifeCycleSavings, aes(x = dpi, y = sr, label = Name))
p +  geom_point() + geom_text(vjust = 1.5)

Opción3

Con el paquete base graphics

plot(sr ~ dpi, 
     LifeCycleSavings, 
     xlab = 'Real Per-Capita Disposable Income', 
     ylab = 'Aggregate Personal Savings', 
     main = 'Intercountry Life-Cycle Savings Data')
text(LifeCycleSavings$dpi,
     LifeCycleSavings$sr,
     Name, pos = 1)

Referencias

• R Graphics Cookbook de Winston Chang
• geom_text
• Set or Query Graphical Parameters
• Add Text to a Plot
• Identify Points in a Scatter Plot
• Regresión lineal simple para principiantes en R
• Regresión lineal simple en R con ggplot2
• Regresión lineal simple para principiantes en Excel