library(ggplot2)
class(diamonds)

## [1] "tbl_df"     "tbl"        "data.frame"

diamonds

## # A tibble: 53,940 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
##  7  0.24 Very Good     I    VVS1  62.3    57   336  3.95  3.98  2.47
##  8  0.26 Very Good     H     SI1  61.9    55   337  4.07  4.11  2.53
##  9  0.22      Fair     E     VS2  65.1    61   337  3.87  3.78  2.49
## 10  0.23 Very Good     H     VS1  59.4    61   338  4.00  4.05  2.39
## # ... with 53,930 more rows

Summarizing data

Las funciones que cubriremos son:

• filter(),
• select(),
• arrange(),
• mutate(),
• summarise(), y
• group_by(),

estas serán comparadas con las versiones equivalentes en R estándar.

Obteniendo subconjunto de datos usando filter()

Así como las opciones básicas de R para obtener subconjuntos de datos a través de expresiones condicionales como el uso de corchetees doblesdf[,], filter() crear un subconjunto de datos basado en algún criterio definido por el usuario.

Supongamos que se desea crear un subconjunto de diamantes diamonds manteniendo todos aquellos de color D,E o F con un corte de calidad (cut quality) Premium’ or ‘Ideal’ y un peso de más de 3 carat. Una forma estándar muy básica sería:

Usando filter() notemos que no se requiere repetir el nombre del data frame:

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

filter(diamonds, carat > 3 &
         cut %in% c("Premium", "Ideal") &
         color %in% c("D", "E", "F"))

## # A tibble: 2 x 10
##   carat     cut color clarity depth table price     x     y     z
##   <dbl>   <ord> <ord>   <ord> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1  3.01 Premium     F      I1  62.2    56  9925  9.24  9.13  5.73
## 2  3.05 Premium     E      I1  60.9    58 10453  9.26  9.25  5.66

Seleccinando columnas usando select()

dplyr::select(diamonds, carat, cut, color, clarity)

## # A tibble: 53,940 x 4
##    carat       cut color clarity
##    <dbl>     <ord> <ord>   <ord>
##  1  0.23     Ideal     E     SI2
##  2  0.21   Premium     E     SI1
##  3  0.23      Good     E     VS1
##  4  0.29   Premium     I     VS2
##  5  0.31      Good     J     SI2
##  6  0.24 Very Good     J    VVS2
##  7  0.24 Very Good     I    VVS1
##  8  0.26 Very Good     H     SI1
##  9  0.22      Fair     E     VS2
## 10  0.23 Very Good     H     VS1
## # ... with 53,930 more rows

Nótese que tanto el uso de c() para combinar los nombres de varias columnas como el uso de corchetes dobles se vuelve obsoleto. Además es posible omitir el criterio de conteo de columnas (2:5) y expresar la selección de columnas contiguas a través de los nombres.

dplyr::select(diamonds, carat:clarity, price)

## # A tibble: 53,940 x 5
##    carat       cut color clarity price
##    <dbl>     <ord> <ord>   <ord> <int>
##  1  0.23     Ideal     E     SI2   326
##  2  0.21   Premium     E     SI1   326
##  3  0.23      Good     E     VS1   327
##  4  0.29   Premium     I     VS2   334
##  5  0.31      Good     J     SI2   335
##  6  0.24 Very Good     J    VVS2   336
##  7  0.24 Very Good     I    VVS1   336
##  8  0.26 Very Good     H     SI1   337
##  9  0.22      Fair     E     VS2   337
## 10  0.23 Very Good     H     VS1   338
## # ... with 53,930 more rows

Note: otras funciones útiles para definir el criterio de selección: starts_with(), ends_with(), matches() and contains().

Pipelining

Ahora, supongamos que se desea seleccionar algunas columnas a partir del subconjunto de datos creado anteriormente. Normalmente usaríamos dos pasos con resultados intermedios innecesarios:

## first, create subset
diamonds_sub <- filter(diamonds, carat > 3 &
                         cut %in% c("Premium", "Ideal") &
                         color %in% c("D", "E", "F"))

## second, select columns
select(diamonds_sub, carat:clarity, price)

o quizas funciones anidadas que usualmnete son difíciles de interpretar:

## all-in-one nested solution
select(
  filter(diamonds, carat > 3 &
           cut %in% c("Premium", "Ideal") &
           color %in% c("D", "E", "F")), 
  carat:clarity, price
)

dplyr introduce el operador %>% el cual creado para interconectar los pasos de procesamiento en los datos eliminando el uso de variables intermedias o funciones anidadas como se ilustró anteriormente. Podemos pensar en %>% como una traducción de “entonces” o “luego” que conecta dos partes de una oración:

diamonds %>%
  filter(carat > 3 &
           cut %in% c("Premium", "Ideal") &
           color %in% c("D", "E", "F")) %>%
  dplyr::select(carat:clarity, price)

## # A tibble: 2 x 5
##   carat     cut color clarity price
##   <dbl>   <ord> <ord>   <ord> <int>
## 1  3.01 Premium     F      I1  9925
## 2  3.05 Premium     E      I1 10453

Este encadenamiento se vuelve aún más útil caundo se efectúan múltiples operaciones a la vez, reduciendo el código y haciéndolo fácil de leer e interpretar. La analogía es la lectura de la página de un libro, que se efectúa de arriba-abajo e izquierda-derecha.

Reordenar filas usando arrange()

arrange() facilita la rapidez en el reordenamiento de filas basada en criterios multi-columna.

Asumamos un ejemplo con nuestro data set diamonds. Se quiere obtener un listado del subconjunto de datos referentes a los atributos de color, precio y carat (color, price, carat) ordenados en base a los siguientes criterios.

1. De acuerdo al color con el mejor al inicio (color D)
2. Sobre el orden anterior, es decir sobre diamantes de igual color ordenarlos con el menor precio al inicio
3. Finalmente de cada conjunto de color y precio iguales, se desea ordenar sus pesos de manera descendente.

La forma trivial:

Usando la librería dplyr:

diamonds %>%
  dplyr::select(color, price, carat) %>%
  arrange(color, price, desc(carat))

## # A tibble: 53,940 x 3
##    color price carat
##    <ord> <int> <dbl>
##  1     D   357  0.23
##  2     D   357  0.23
##  3     D   361  0.32
##  4     D   362  0.23
##  5     D   367  0.24
##  6     D   367  0.20
##  7     D   367  0.20
##  8     D   367  0.20
##  9     D   373  0.24
## 10     D   373  0.23
## # ... with 53,930 more rows

Añadir nuevas columnas usando mutate()

mutate() permite añadir nuevas variables a un data frame existente. Además tiene la ventaja de poder utilizar las columnas (variables) recién creadas, así:

diamonds %>%
  dplyr::select(color, carat, price) %>%
  mutate(ppc = price / carat, 
         ppc_rnd = round(ppc, 2))

## # A tibble: 53,940 x 5
##    color carat price      ppc ppc_rnd
##    <ord> <dbl> <int>    <dbl>   <dbl>
##  1     E  0.23   326 1417.391 1417.39
##  2     E  0.21   326 1552.381 1552.38
##  3     E  0.23   327 1421.739 1421.74
##  4     I  0.29   334 1151.724 1151.72
##  5     J  0.31   335 1080.645 1080.65
##  6     J  0.24   336 1400.000 1400.00
##  7     I  0.24   336 1400.000 1400.00
##  8     H  0.26   337 1296.154 1296.15
##  9     E  0.22   337 1531.818 1531.82
## 10     H  0.23   338 1469.565 1469.57
## # ... with 53,930 more rows

Resumir valores usando summarise()

Asumamos que se desea calcular el mínimo, media y máximo precio por color de diamante. Usando la versión muy conocida de R estándar:

Ahora, usando la solución de la función group_by, la cual es sumamente útil para trabajar sobre subconjuntos de datos:

diamonds %>% 
  group_by(color) %>%
  summarise(MIN = min(price), MEAN = mean(price), MAX = max(price))

## # A tibble: 7 x 4
##   color   MIN     MEAN   MAX
##   <ord> <dbl>    <dbl> <dbl>
## 1     D   357 3169.954 18693
## 2     E   326 3076.752 18731
## 3     F   342 3724.886 18791
## 4     G   354 3999.136 18818
## 5     H   337 4486.669 18803
## 6     I   334 5091.875 18823
## 7     J   335 5323.818 18710

microbenchmark

El paquete microbenchmark es útil para ejecutar secciones pequeñas de código y evaluar su rendimiento, así como comparar la velocidad de varias funciones que tienen el mismo objetivo. La función del mismo nombre de este paquete ejecuta el código múltiples veces (100 por defecto) y provee un resumen de estadísticas del tiempo utilizado para ejecutar el código. La ventaja de este paquete es que toma en cuenta el tiempo invertido en el procesamiento de la función de microbenchmark por sí misma y de esta manera utiliza iteraciones de “warm-up” antes de ejecutar las iteraciones para el procesamiento estadístico.

You can use the times argument in microbenchmark to customize how many iterations are used. For example, if you are working with a function that is a bit slow, you might want to run the code fewer times when benchmarking (although with slower or more complex code, it likely will make more sense to use a different tool for profiling, like profvis).

Se puede incluir múltiples líneas de cpodigo dentro de una única llamada de microbenchmark. Sin embargo, para obtener benchmarks independientes de cada línea de código se debe separarlas con una coma:

library(microbenchmark)
microbenchmark(a <- rnorm(1000), 
               b <- mean(rnorm(1000)))

Ahora trabajemos en un ejemplo donde vamos a comparar dos funciones que pretenden obtener los mismo resultados. Supoongamos que se necesita una función que puede identificar los días que cumplen con dos condiciones: (a) la temperatura iguala o excede un umbral de temperatura (e.g., 25 grados centígrados) y (b) la temperatura iguala o excede la temperatura más alta registrada antes de ese día.

Los datos de entrada lucen como el data frame abajo que incluye una columna denominada “temp” que indica la temperatura media diaria registrada:

date	temp
2015-07-01	26.5
2015-07-02	27.2
2015-07-03	28.0
2015-07-04	26.9
2015-07-05	27.5
2015-07-06	25.9
2015-07-07	28.0
2015-07-08	28.2

y el data frame de salida se espera tenga una columna adicional de datos binarios denominda “record_temp” indicando si las dos condiciones se han cumplido:

date	temp	record_temp
2015-07-01	26.5	FALSE
2015-07-02	27.2	TRUE
2015-07-03	28.0	TRUE
2015-07-04	26.9	FALSE
2015-07-05	27.5	FALSE
2015-07-06	25.9	FALSE
2015-07-07	28.0	TRUE
2015-07-08	28.2	TRUE

Abajo se plantean dos ejemplos de funciones para este problema.

# Function that uses a loop 
find_records_1 <- function(datafr, threshold){
  highest_temp <- c()
  record_temp <- c()
  for(i in 1:nrow(datafr)){
    highest_temp <- max(highest_temp, datafr$temp[i])
    record_temp[i] <- datafr$temp[i] >= threshold & 
      datafr$temp[i] >= highest_temp
  }
  datafr <- cbind(datafr, record_temp)
  return(datafr)
}

# Function that uses tidyverse functions
find_records_2 <- function(datafr, threshold){
  datafr <- datafr %>%
    mutate_(over_threshold = ~ temp >= threshold,
            cummax_temp = ~ temp == cummax(temp),
            record_temp = ~ over_threshold & cummax_temp) %>%
    select_(.dots = c("-over_threshold", "-cummax_temp"))
  return(as.data.frame(datafr))
}

If you apply the two functions to the small example data set, you can see that they both create the desired output:

example_data <- data_frame(date = c("2015-07-01", "2015-07-02",
                                    "2015-07-03", "2015-07-04",
                                    "2015-07-05", "2015-07-06",
                                    "2015-07-07", "2015-07-08"),
                           temp = c(26.5, 27.2, 28.0, 26.9, 
                                    27.5, 25.9, 28.0, 28.2))

(test_1 <- find_records_1(example_data, 27))

##         date temp record_temp
## 1 2015-07-01 26.5       FALSE
## 2 2015-07-02 27.2        TRUE
## 3 2015-07-03 28.0        TRUE
## 4 2015-07-04 26.9       FALSE
## 5 2015-07-05 27.5       FALSE
## 6 2015-07-06 25.9       FALSE
## 7 2015-07-07 28.0        TRUE
## 8 2015-07-08 28.2        TRUE

(test_2 <- find_records_2(example_data, 27))

##         date temp record_temp
## 1 2015-07-01 26.5       FALSE
## 2 2015-07-02 27.2        TRUE
## 3 2015-07-03 28.0        TRUE
## 4 2015-07-04 26.9       FALSE
## 5 2015-07-05 27.5       FALSE
## 6 2015-07-06 25.9       FALSE
## 7 2015-07-07 28.0        TRUE
## 8 2015-07-08 28.2        TRUE

all.equal(test_1, test_2)

## [1] TRUE

The performance of these two functions can be compared using microbenchmark:

library(microbenchmark)
record_temp_perf <- microbenchmark(find_records_1(example_data, 27), 
                                   find_records_2(example_data, 27))
record_temp_perf

## Unit: milliseconds
##                              expr       min        lq      mean    median
##  find_records_1(example_data, 27)  1.630206  1.808161  2.137233  1.971466
##  find_records_2(example_data, 27) 11.470211 12.793909 13.836272 13.790387
##         uq       max neval
##   2.225151  4.407326   100
##  14.566731 18.791127   100

It’s useful to check next to see if the relative performance of the two functions is similar for a bigger data set. The chicagoNMMAPS data set from the dlnm package includes temperature data over 15 years in Chicago, IL. Here are the results when we benchmark the two functions with that data (note, this code takes a minute or two to run):

library(dlnm)

## This is dlnm 2.3.4. For details: help(dlnm) and vignette('dlnmOverview').

data("chicagoNMMAPS")

record_temp_perf_2 <- microbenchmark(find_records_1(chicagoNMMAPS, 27), 
                                     find_records_2(chicagoNMMAPS, 27))
record_temp_perf_2

## Unit: milliseconds
##                               expr       min        lq      mean    median
##  find_records_1(chicagoNMMAPS, 27) 196.71548 203.66758 211.27056 208.27330
##  find_records_2(chicagoNMMAPS, 27)  13.27619  14.24635  15.35658  15.14507
##         uq       max neval
##  212.71677 322.16791   100
##   16.07052  22.96398   100

While the function with the loop (find_records_1) performed better with the very small sample data, the function that uses tidyverse functions (find_records_2) performs much, much better with a larger data set.

The microbenchmark function returns an object of the “microbenchmark” class. This class has two methods for plotting results, autoplot.microbenchmark and boxplot.microbenchmark. To use the autoplot method, you will need to have ggplot2 loaded in your R session.

library(ggplot2)
# For small example data
autoplot(record_temp_perf)

# For larger data set
autoplot(record_temp_perf_2)

profvis

Una vez que se ha identificado el código cuya ejecución es más lenta, el siguiente paso es comprender que partes del código están causando los cuellos de botella. La función profvis del paquete profvis es la herramienta adecuada para esto la cual usa la función RProf base de R y visualiza un entorno interactivo de navegación. El profiling se efectúa escribiendo los resultados de la pila de llamadas (sentencias) cada 10 milisegundos.

Para usar esta herramienta, únicamente incluya el código en llaves en el caso de tener varias sentencias. Probemos a identificar cual fue el problema con la función find_records_1, para ejecutar el profiling:

library(profvis)
library(dlnm)
datafr <- chicagoNMMAPS
threshold <- 25

profvis({
  highest_temp <- c()
  record_temp <- c()
  for(i in 1:nrow(datafr)){
    highest_temp <- max(highest_temp, datafr$temp[i])
    record_temp[i] <- datafr$temp[i] >= threshold & 
      datafr$temp[i] >= highest_temp
  }
  datafr <- cbind(datafr, record_temp)
})

La salida de profvis proporciona dos opciones de visualización, “Flame Graph” o “Data”.

Es posible usar el agumento interval para personalizar el intervalo de muestreo (interval = 0.01 -> 10 milisegundos). Se sugiere eviar el uso de periodos de intervalo muy cortos (e.g., menores a 5 milisegundos). Esto puede provocar resultados no precisos. Incluso si el profiling corresponde a un código muy rápido es mejor utilizar microbenchmark el cual produce resultados más precisos a intervalos más cortos.

###References

Advance Programming with R Florian https://oer.uni-marburg.de/data/mriliasmooc/lm_data/lm_2050/speeding-up-iteration-procedures.html

Link from the Data Sciente Toolbox course Coursera https://www.youtube.com/watch?v=jWjqLW-u3hc

Enlace a Cheatsheet https://www.rstudio.com/wp-content/uploads/2015/02/data-wrangling-cheatsheet.pdf

Advance Programming with R (from https://bookdown.org/rdpeng/RProgDA/functions.html)