Este es el primer de 2 juegos de notas introductorias a R con un enfoque en análisis de datos. A diferencia de otros recursos, no se pretende dar una introducción general a R sino mostrar las herramientas más importantes para comenzar a utilizar R en análisis de datos. En este primer juego se cubre la sección de visualización. Más adelante aprenderemos a usar R para manipulación de datos y modelación. Estas notas siguen material e ideas de Hadley Wickham y en particular el libro R for Data Science. Las notas están ordenadas como sigue:
Un estándar mínimo para el análisis de datos es la reproducibilidad: debe ser posible reproducir el análisis en todos sus pasos y en cualquier momento. Para ello los pasos del análisis deben estar documentados apropiadamente, de manera que las decisiones importantes puedan ser entendidas claramente.
Estos dos principios generalmente implican que debemos trabajar escribiendo código, más que usando interfaces gráficas de point and click. Esto permite crear programas reproducibles que son fácilmente documentados, y tiene otras consecuencias positivas como la facilidad de comunicación (compartir código), la posibilidad de trabajar con versiones que documenten la historia del desarrollo, respaldos fáciles del trabajo, e incluso el uso de lenguajes de programación más flexibles que integren nuestro trabajo en procesos de producción de reportes o monitoreo.
Los scripts de R son oportunos para llevar a cabo análisis reproducbiles, pero
hay más herramientas que nos ayudan a documentar y compartir nuestro trabajo:
Los paquetes rmarkdown y knitr se utilizan para generar documentos en formato pdf, html o word que integran texto, código de R y resultados producidos por el código.
Packrat: Sistema para administrar dependencias de paquetes en R.
Organizar los análisis para ser reproducibles no es trivial pero es una buena práctica que te agradecerán los que consulten o utilicen tu trabajo (incluído tu yo del futuro), puedes leer más recomendaciones para lograr análisis reproducibles en initial steps toward reproducible research. También es conveniente usar un controlador de versiones este es un buen tutorial para Git y Github con R.
Para comenzar se debe descargar R, esta descarga incluye R básico y un editor de textos para escribir código. Después de descargar R se recomienda descargar RStudio (gratis y libre).
Rstudio es un ambiente de desarrollo integrado para R: incluye una consola, un editor de texto y un conjunto de herramientas para administrar el espacio de trabajo cuando se utiliza R.
Algunos shortcuts útiles:
En el editor
En la consola
La mejor manera de usar R para análisis de datos es aprovechando la gran cantidad de paquetes que aportan funcionalidad adicional. Desde Rstudio podemos instalar paquetes (Tools - > Install packages o usar la función install.packages("nombre_paquete")
). Una vez instalados, podemos cargarlos a nuestra sesión de R mediante library
. Por ejemplo, para cargar el paquete readr
hacemos:
# print(read_csv)
# Error in print(read_csv) : object 'read_csv' not found
library(tidyverse)
print(read_csv)
## function (file, col_names = TRUE, col_types = NULL, locale = default_locale(),
## na = c("", "NA"), quoted_na = TRUE, quote = "\\"", comment = "",
## trim_ws = TRUE, skip = 0, n_max = Inf, guess_max = min(1000,
## n_max), progress = show_progress())
## {
## tokenizer <- tokenizer_csv(na = na, quoted_na = TRUE, quote = quote,
## comment = comment, trim_ws = trim_ws)
## read_delimited(file, tokenizer, col_names = col_names, col_types = col_types,
## locale = locale, skip = skip, comment = comment, n_max = n_max,
## guess_max = guess_max, progress = progress)
## }
## <environment: namespace:readr>
read_csv
es una función que aporta el paquete readr
, que a su vez está incluido en el tidyverse.
Los paquetes se instalan una sola vez, sin embargo, se deben cargar (ejecutar library(tidyverse)
) en cada sesión de R que los ocupemos.
En estas notas utilizaremos la colección de paquetes incluídos en el tidyverse. Estos paquetes de R están diseñados para ciencia de datos, y para funcionar juntos como parte de un flujo de trabajo.
La siguiente imagen tomada de Why the tidyverse (Joseph Rickert) indica que paquetes del tidyverse se utilizan para cada etapa del análisis de datos.
knitr::include_graphics("imagenes/tidyverse.png")
En R se puede trabajar con distintas estructuras de datos, algunas son de una sola dimensión y otras permiten más, como indica el diagrama de abajo:
nosotros trabajaremos principalmente con vectores y data frames.
Comenzamos viendo algunas operaciones básicas con vectores.
a <- c(5, 2, 4.1, 7, 9.2)
a
## [1] 5.0 2.0 4.1 7.0 9.2
a[1]
## [1] 5
a[2]
## [1] 2
a[2:4]
## [1] 2.0 4.1 7.0
Las operaciones básicas con vectores son componente a componente:
b <- a + 10
b
## [1] 15.0 12.0 14.1 17.0 19.2
d <- sqrt(a)
d
## [1] 2.236068 1.414214 2.024846 2.645751 3.033150
a + d
## [1] 7.236068 3.414214 6.124846 9.645751 12.233150
10 * a
## [1] 50 20 41 70 92
a * d
## [1] 11.180340 2.828427 8.301867 18.520259 27.904982
Y podemos crear secuencias como sigue:
ejemplo_1 <- 1:10
ejemplo_1
## [1] 1 2 3 4 5 6 7 8 9 10
ejemplo_2 <- seq(0, 1, 0.25)
ejemplo_2
## [1] 0.00 0.25 0.50 0.75 1.00
Para calcular características de vectores usamos funciones:
# media del vector
mean(a)
## [1] 5.46
# suma de sus componentes
sum(a)
## [1] 27.3
# longitud del vector
length(a)
## [1] 5
También podemos construir vectores de caracteres:
frutas <- c('manzana', 'manzana', 'pera', 'plátano', 'fresa')
frutas
## [1] "manzana" "manzana" "pera" "plátano" "fresa"
Podemos juntar vectores del mismo tamaño en tablas, que se llaman data.frame
. Por ejemplo:
tabla <- data_frame(n = 1:5, valor = a, fruta = frutas) # la función data_frame de tibble es más conveniente que data.frame de R base.
tabla
## # A tibble: 5 x 3
## n valor fruta
## <int> <dbl> <chr>
## 1 1 5.0 manzana
## 2 2 2.0 manzana
## 3 3 4.1 pera
## 4 4 7.0 plátano
## 5 5 9.2 fresa
Los data frames son estructuras rectangulares donde cada columna es del mismo tipo (e.g. categórica o factor, numérica, caracter) pero columnas distintas pueden tener diferentes tipos.
library(ggplot2)
head(diamonds)
## # A tibble: 6 x 10
## carat cut color clarity depth table price x y z
## <dbl> <ord> <ord> <ord> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1 0.23 Ideal E SI2 61.5 55 326 3.95 3.98 2.43
## 2 0.21 Premium E SI1 59.8 61 326 3.89 3.84 2.31
## 3 0.23 Good E VS1 56.9 65 327 4.05 4.07 2.31
## 4 0.29 Premium I VS2 62.4 58 334 4.20 4.23 2.63
## 5 0.31 Good J SI2 63.3 58 335 4.34 4.35 2.75
## 6 0.24 Very Good J VVS2 62.8 57 336 3.94 3.96 2.48
La instrucción str
nos describe el tipo de variables en el data.frame:
str(diamonds)
## Classes 'tbl_df', 'tbl' and 'data.frame': 53940 obs. of 10 variables:
## $ carat : num 0.23 0.21 0.23 0.29 0.31 0.24 0.24 0.26 0.22 0.23 ...
## $ cut : Ord.factor w/ 5 levels "Fair"<"Good"<..: 5 4 2 4 2 3 3 3 1 3 ...
## $ color : Ord.factor w/ 7 levels "D"<"E"<"F"<"G"<..: 2 2 2 6 7 7 6 5 2 5 ...
## $ clarity: Ord.factor w/ 8 levels "I1"<"SI2"<"SI1"<..: 2 3 5 4 2 6 7 3 4 5 ...
## $ depth : num 61.5 59.8 56.9 62.4 63.3 62.8 62.3 61.9 65.1 59.4 ...
## $ table : num 55 61 65 58 58 57 57 55 61 61 ...
## $ price : int 326 326 327 334 335 336 336 337 337 338 ...
## $ x : num 3.95 3.89 4.05 4.2 4.34 3.94 3.95 4.07 3.87 4 ...
## $ y : num 3.98 3.84 4.07 4.23 4.35 3.96 3.98 4.11 3.78 4.05 ...
## $ z : num 2.43 2.31 2.31 2.63 2.75 2.48 2.47 2.53 2.49 2.39 ...
Para lograr una programación eficiente en R es importante conocer las técnicas de indexar data frames:
# extraemos los primeros cinco renglones
diamonds[1:5, ]
## # A tibble: 5 x 10
## carat cut color clarity depth table price x y z
## <dbl> <ord> <ord> <ord> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1 0.23 Ideal E SI2 61.5 55 326 3.95 3.98 2.43
## 2 0.21 Premium E SI1 59.8 61 326 3.89 3.84 2.31
## 3 0.23 Good E VS1 56.9 65 327 4.05 4.07 2.31
## 4 0.29 Premium I VS2 62.4 58 334 4.20 4.23 2.63
## 5 0.31 Good J SI2 63.3 58 335 4.34 4.35 2.75
# extraemos los primeros cinco renglones y las columnas 2,4,6
diamonds[1:5, c(2, 4, 6)]
## # A tibble: 5 x 3
## cut clarity table
## <ord> <ord> <dbl>
## 1 Ideal SI2 55
## 2 Premium SI1 61
## 3 Good VS1 65
## 4 Premium VS2 58
## 5 Good SI2 58
# también podemos extraer columnase usando $: extraemos la columna x
head(diamonds$x)
## [1] 3.95 3.89 4.05 4.20 4.34 3.94
# ¿Que extraemos con las siguientes 2 instrucciones?
diamonds[diamonds$x == diamonds$y, ]
## # A tibble: 17 x 10
## carat cut color clarity depth table price x y z
## <dbl> <ord> <ord> <ord> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1 0.30 Ideal H VVS2 62.5 54 567 4.30 4.30 2.70
## 2 0.27 Very Good F VVS1 62.0 55 591 4.16 4.16 2.59
## 3 1.00 Very Good H VS2 63.3 53 5139 0.00 0.00 0.00
## 4 1.14 Fair G VS1 57.5 67 6381 0.00 0.00 0.00
## 5 1.00 Premium E VS2 60.0 60 6600 6.43 6.43 3.89
## 6 1.00 Premium E VS2 60.0 60 6720 6.43 6.43 3.89
## 7 1.22 Premium G SI2 62.4 61 6969 6.79 6.79 4.23
## 8 1.56 Ideal G VS2 62.2 54 12800 0.00 0.00 0.00
## 9 1.20 Premium D VVS1 62.1 59 15686 0.00 0.00 0.00
## 10 2.25 Premium H SI2 62.8 59 18034 0.00 0.00 0.00
## 11 0.32 Ideal D VVS2 62.1 54 858 4.40 4.40 2.74
## 12 0.42 Ideal H VVS1 62.8 57 1108 4.79 4.79 3.01
## 13 0.61 Premium G SI1 60.8 60 1255 5.42 5.42 3.31
## 14 0.48 Ideal F VS2 62.4 54 1279 5.03 5.03 3.15
## 15 0.51 Premium F SI1 61.4 59 1421 5.13 5.13 3.16
## 16 0.71 Good F SI2 64.1 60 2130 0.00 0.00 0.00
## 17 0.71 Good F SI2 64.1 60 2130 0.00 0.00 0.00
diamonds[-(1:53929), c("carat", "price")]
## # A tibble: 11 x 2
## carat price
## <dbl> <int>
## 1 0.71 2756
## 2 0.71 2756
## 3 0.71 2756
## 4 0.70 2757
## 5 0.70 2757
## 6 0.72 2757
## 7 0.72 2757
## 8 0.72 2757
## 9 0.70 2757
## 10 0.86 2757
## 11 0.75 2757
como vemos arriba para indexar los data frames tenemos que indicar filas y columnas, en el lado izquierdo de los corchetes se indica (con un vector) que filas queremos extraer, y en el lado derecho se indican las columnas: diamonds[filas, columnas]
. También vale la pena notar que diamonds$x
regresa la columna x como vector, es decir, diamonds$x
es de una sola dimensión.
En R los datos faltantes se expresan como NA
, ¿qué regresan las siguientes expresiones?
5 + NA
NA / 2
sum(c(5, 4, NA))
mean(c(5, 4, NA))
NA < 3
NA == 3
NA == NA
Las expresiones anteriores regresan NA
, el hecho que la media de un vector que incluye NAs o su suma regrese NAs se debe a que el default en R es propagar los valores faltantes, esto es, si deconozco el valor de una de las componentes de un vector, también desconozco la suma del mismo; sin embargo, muchas funciones tienen un argumento na.rm para removerlos,
sum(c(5, 4, NA), na.rm = TRUE)
## [1] 9
mean(c(5, 4, NA), na.rm = TRUE)
## [1] 4.5
El manejo de datos faltantes en R utiliza una lógica ternaria (como SQL):
NA == NA
## [1] NA
La expresión anterior puede resultar confusa, una manera de pensar en esto es considerar los NA como no sé, por ejemplo si no se la edad de Juan y no se la edad de Esteban, la pregunta a ¿Juan tiene la misma edad que Esteban? es no sé (NA).
edad_Juan <- NA
edad_Esteban <- NA
edad_Juan == edad_Esteban
## [1] NA
edad_Jose <- 32
# Juan es menor que José?
edad_Juan < edad_Jose
## [1] NA
?mean
o help()
a <- 10
a
## [1] 10
(a <- 15)
## [1] 15
qplot(carat, price, data = diamonds, colour = clarity)
Las funciones son una base importante de la programación en R. Veamos un ejemplo:
wtd_mean <- function(x, wt = rep(1, length(x))) {
sum(x * wt) / sum(wt)
}
wtd_mean(1:10)
## [1] 5.5
wtd_mean(1:10, 10:1)
## [1] 4
Las funciones de R tienen tres partes:
body(wtd_mean)
## {
## sum(x * wt)/sum(wt)
## }
formals(wtd_mean)
## $x
##
##
## $wt
## rep(1, length(x))
environment(wtd_mean)
## <environment: R_GlobalEnv>
environment(qplot)
## <environment: namespace:ggplot2>
Veamos mas ejemplos, ¿qué regresan las siguientes funciones?
# 1
x <- 5
f <- function(){
y <- 10
c(x = x, y = y)
}
rm(x, f)
# 2
x <- 5
g <- function(){
x <- 20
y <- 10
c(x = x, y = y)
}
# 3
x <- 5
h <- function(){
y <- 10
i <- function(){
z <- 20
c(x = x, y = y, z = z)
}
i()
}
# 4 ¿qué ocurre si la corremos por segunda vez?
j <- function(){
if (!exists("a")){
a <- 5
} else{
a <- a + 1
}
print(a)
}
x <- 0
y <- 10
# 5 ¿qué regresa k()? ¿y k()()?
k <- function(){
x <- 1
function(){
y <- 2
x + y
}
}
Las reglas de búsqueda determinan como se busca el valor de una variable libre en una función. A nivel lenguaje R usa lexical scoping, una alternativa es dynamic scoping. En R (lexical scoping) los valores de los símbolos se basan en como se anidan las funciones cuando fueron creadas y no en como son llamadas. Esto es, en R no importa como son las llamadas a una función para saber como se va a buscar el valor de una variable. Una consecuencia de las reglas de búsqueda es que todos los objetos deben ser guardados en memoria.
Utilizaremos el paquete ggplot2 y cubriremos los siguientes puntos:
# install.packages("ggplot2") # sólo se hace una vez
library(tidyverse) # Cargamos el paquete en nuestra sesión
Usaremos el conjunto de datos mpg que se incluye en R, puedes encontrar información de esta base de datos tecleando ?mpg
.
?mpg
glimpse(mpg)
## Observations: 234
## Variables: 11
## $ manufacturer <chr> "audi", "audi", "audi", "audi", "audi", "audi", "...
## $ model <chr> "a4", "a4", "a4", "a4", "a4", "a4", "a4", "a4 qua...
## $ displ <dbl> 1.8, 1.8, 2.0, 2.0, 2.8, 2.8, 3.1, 1.8, 1.8, 2.0,...
## $ year <int> 1999, 1999, 2008, 2008, 1999, 1999, 2008, 1999, 1...
## $ cyl <int> 4, 4, 4, 4, 6, 6, 6, 4, 4, 4, 4, 6, 6, 6, 6, 6, 6...
## $ trans <chr> "auto(l5)", "manual(m5)", "manual(m6)", "auto(av)...
## $ drv <chr> "f", "f", "f", "f", "f", "f", "f", "4", "4", "4",...
## $ cty <int> 18, 21, 20, 21, 16, 18, 18, 18, 16, 20, 19, 15, 1...
## $ hwy <int> 29, 29, 31, 30, 26, 26, 27, 26, 25, 28, 27, 25, 2...
## $ fl <chr> "p", "p", "p", "p", "p", "p", "p", "p", "p", "p",...
## $ class <chr> "compact", "compact", "compact", "compact", "comp...
Y comencemos con nuestra primera gráfica:
ggplot(data = mpg) +
geom_point(mapping = aes(x = displ, y = hwy))
En ggplot2 se inicia una gráfica con la instrucción ggplot()
, debemos especificar explicitamente que base de datos usamos, este es el primer argumento en la función ggplot. Una vez que creamos la base añadimos capas, y dentro de aes() escribimos las variables que queremos graficar y el atributo de la gráfica al que queremos mapearlas.
La función geom_point()
añade una capa de puntos, hay muchas funciones geometrías incluídas en ggplot2: geom_line()
, geom_boxplot()
, geom_histogram
,… Cada una acepta distintos argumentos para mapear las variables en los datos a características estéticas de la gráfica. En el ejemplo de arriba mapeamos displ
al eje x, hwy
al eje y, pero geom_point()
nos permite representar más variables usando la forma, color y/o tamaño del punto. Esta flexibilidad nos permite entender o descubrir patrones más interesantes en los datos.
ggplot(mpg) +
geom_point(aes(x = displ, y = hwy, color = class))
Experimenta con los aesthetics color (color), tamaño (size) y forma (shape).
¿Qué diferencia hay entre las variables categóricas y las continuas?
¿Qué ocurre cuando combinas varios aesthetics?
El mapeo de las propiedades estéticas se denomina escalamiento y depende del tipo de variable, las variables discretas (por ejemplo, genero, escolaridad, país) se mapean a distintas escalas que las variables continuas (variables numéricas como edad, estatura, etc.), los defaults para algunos atributos son (estos se pueden modificar):
Discreta | Continua | |
---|---|---|
Color (color) |Arcoiris de colores |Gradiente de colores Tamaño ( size) |Escala discreta de tamaños |Mapeo lineal entre el área y el valor Forma ( shape) |Distintas formas |No aplica Transparencia ( alpha`) |
No aplica | Mapeo lineal a la transparencia |
Los geoms controlan el tipo de gráfica
p <- ggplot(mpg, aes(x = displ, y = hwy))
p + geom_line() # en este caso no es una buena gráfica
¿Qué problema tiene la siguiente gráfica?
p <- ggplot(mpg, aes(x = cty, y = hwy))
p + geom_point()
p + geom_jitter()
¿Cómo podemos mejorar la siguiente gráfica?
ggplot(mpg, aes(x = class, y = hwy)) +
geom_point()
Intentemos reodenar los niveles de la variable clase
ggplot(mpg, aes(x = reorder(class, hwy), y = hwy)) +
geom_point()
Podemos probar otros geoms.
ggplot(mpg, aes(x = reorder(class, hwy), y = hwy)) +
geom_jitter()
ggplot(mpg, aes(x = reorder(class, hwy), y = hwy)) +
geom_boxplot()
También podemos usar más de un geom!
ggplot(mpg, aes(x = reorder(class, hwy), y = hwy)) +
geom_jitter() +
geom_boxplot()
Lee la ayuda de reorder y repite las gráficas anteriores ordenando por la mediana de hwy.
¿Cómo harías para graficar los puntos encima de las cajas de boxplot?
Veamos ahora como hacer páneles de gráficas, la idea es hacer varios múltiplos de una gráfica donde cada múltiplo representa un subconjunto de los datos, es una práctica muy útil para explorar relaciones condicionales.
En ggplot podemos usar facet_wrap() para hacer paneles dividiendo los datos de acuerdo a las categorías de una sola variable
ggplot(mpg, aes(x = displ, y = hwy)) +
geom_jitter() +
facet_wrap(~ cyl)
También podemos hacer una cuadrícula de 2 dimensiones usando facet_grid(filas~columnas)
ggplot(mpg, aes(x = displ, y = hwy)) +
geom_jitter() +
facet_grid(.~ class)
ggplot(mpg, aes(x = displ, y = hwy)) +
geom_jitter() +
facet_grid(drv ~ class)
Los páneles pueden ser muy útiles para entender relaciones en nuestros datos. En la siguiente gráfica es difícil entender si existe una relación entre radiación solar y ozono
data(airquality)
ggplot(airquality, aes(x = Solar.R, y = Ozone)) +
geom_point()
## Warning: Removed 42 rows containing missing values (geom_point).
Veamos que ocurre si realizamos páneles separando por velocidad del viento
library(Hmisc)
airquality$Wind.cat <- cut2(airquality$Wind, g = 3)
ggplot(airquality, aes(x = Solar.R, y = Ozone)) +
geom_point() +
facet_wrap(~ Wind.cat)
Podemos agregar un suavizador (loess) para ver mejor la relación de las variables en cada panel.
ggplot(airquality, aes(x = Solar.R, y = Ozone)) +
geom_point() +
facet_wrap(~ Wind.cat) +
geom_smooth(span = 3)
## `geom_smooth()` using method = 'loess'
Escribe algunas preguntas que puedan contestar con estos datos.
En ocasiones es necesario realizar transformaciones u obtener subconjuntos de los datos para poder responder preguntas de nuestro interés.
library(babynames)
glimpse(babynames)
## Observations: 1,858,689
## Variables: 5
## $ year <dbl> 1880, 1880, 1880, 1880, 1880, 1880, 1880, 1880, 1880, 188...
## $ sex <chr> "F", "F", "F", "F", "F", "F", "F", "F", "F", "F", "F", "F...
## $ name <chr> "Mary", "Anna", "Emma", "Elizabeth", "Minnie", "Margaret"...
## $ n <int> 7065, 2604, 2003, 1939, 1746, 1578, 1472, 1414, 1320, 128...
## $ prop <dbl> 0.072384329, 0.026679234, 0.020521700, 0.019865989, 0.017...
Supongamos que queremos ver la tendencia del nombre “John”, para ello debemos generar un subconjunto de la base de datos.
babynames_John <- babynames[babynames$name == "John", ]
ggplot(babynames_John, aes(x = year, y = prop)) +
geom_point()
ggplot(babynames_John, aes(x = year, y = prop, color = sex)) +
geom_line()
La preparación de los datos es un aspecto muy importante del análisis y suele ser la fase que lleva más tiempo. Es por ello que el siguiente tema se enfocará en herramientas para hacer transformaciones de manera eficiente.
Tarea. Explora la base de datos gapminder, estos datos están incluidos en el paquete del mismo nombre, para acceder a ellos basta con cargar el paquete:
# install.packages("gapminder")
library(gapminder)
gapminder
## # A tibble: 1,704 x 6
## country continent year lifeExp pop gdpPercap
## <fctr> <fctr> <int> <dbl> <int> <dbl>
## 1 Afghanistan Asia 1952 28.801 8425333 779.4453
## 2 Afghanistan Asia 1957 30.332 9240934 820.8530
## 3 Afghanistan Asia 1962 31.997 10267083 853.1007
## 4 Afghanistan Asia 1967 34.020 11537966 836.1971
## 5 Afghanistan Asia 1972 36.088 13079460 739.9811
## 6 Afghanistan Asia 1977 38.438 14880372 786.1134
## 7 Afghanistan Asia 1982 39.854 12881816 978.0114
## 8 Afghanistan Asia 1987 40.822 13867957 852.3959
## 9 Afghanistan Asia 1992 41.674 16317921 649.3414
## 10 Afghanistan Asia 1997 41.763 22227415 635.3414
## # ... with 1,694 more rows
realiza al menos 3 gráficas y explica las relaciones que encuentres. Debes usar lo que revisamos en estas notas: al menos una de las gráficas debe ser de páneles, realiza una gráfica con datos de México, y (opcional)si lo consideras interesante, puedes crear una variable categórica utilizando la función cut2 del paquete Hmisc.