25-resumenes.Rmd

```{r echo = FALSE}
library(knitr)

# Color text
colorize <- function(x, color) {
  
  if (knitr::is_latex_output()) {
    
    sprintf("\\textcolor{%s}{%s}", color, x)
    
  } else if (knitr::is_html_output()) {
    
    sprintf("<span style='color: %s;'>%s</span>", color, x)
    
  } else { x }
}
```


# Resumiendo y relacionado datos {#resumenes}

:::: {.blackbox data-latex=""}

Scripts usados:

* [**script25.R**](https://github.com/dadosdelaplace/courses-ECI-2022/blob/main/scripts/script25.R): resumiendo y relacionando datos.  Ver en <https://github.com/dadosdelaplace/courses-ECI-2022/blob/main/scripts/script25.R>
::::


Ahora que **ya sabemos depurar y transformar los datos**, y algunas de las características de los tipos de variables, en esta sección vamos a aprender dos cosas básicas a la hora de trabajar con datos:

* Realizar <mark>**resúmenes numéricos**</mark> de los datos, de forma general pero también cuando queremos calcular estadísticas desagregadas por grupos.
* <mark>**Medidas de centralización, dispersión y posición**</mark>
* <mark>**Relacionar tablas**</mark> entre sí (los famosos _join_).


## Resúmenes numéricos (summarise y skimr) y por grupos (group_by)

Antes de pasar a ver **cómo generar nuestras propias estadísticas de los datos**, veamos la funcionalidad de un maravilloso paquete llamado `{skimr}`, que nos permite tener en un vistazo un <mark>**resumen numérico**</mark> muy completo de nuestros datos, con **histograma/diagrama de barras incluido**.

```{r}
library(skimr)
starwars %>% skim()
```

Dicho resumen nos proporciona:

- Variables de tipo **caracter**: `n_missing` (número de ausentes), `complete_rate` (proporción de datos sin ausentes), `min/max` y `n_unique` (número de valores únicos).

- Variables de tipo **lista**: `n_missing` (número de ausentes), `complete_rate` (proporción de datos sin ausentes), `n_unique` (número de valores únicos) y `min_length/max_length` (longitud mínimo/máxima de las listas).

- Variables de tipo **numérico**: `n_missing` (número de ausentes), `complete_rate` (proporción de datos sin ausentes), `mean/sd` (media y cuasidesviación típica), `p0/p25/p50/p75/p100` (cuartiles, percentiles 0%-25%-50%-75%-100%, valores que nos dividen nuestro conjunto en 4 trozos) e `hist` (una especie de histograma/barras sencillo).

Este resumen no solo podemos visualizarlo sino que podemos guardarlo para exportarlo por ejemplo en un `.csv`. Sin embargo, aunque el resumen es bastante completo, muchas veces querremos <mark>**generar nuestras propias estadísticas o resúmenes numéricos**</mark>, y para eso vamos a aplicar `summarise()`, que nos calculará estadísticas de nuestros datos.


### Medidas de centralización y dispersión

Por ejemplo, vamos a <mark>**calcular las medidas de centralización (media-mediana-moda)**</mark>. Para la media y la mediana basta con usar las funciones correspondientes `mean()` y `median()`, dentro de `summarise()`

```{r}
starwars %>% summarise(media = mean(mass),
                       mediana = median(mass))
```

Al contrario que las demás órdenes de `{tidyverse}`, la función `summarise()` no nos devuelve la tabla original modificada, sino un resumen de los datos, con las funciones que le hayamos indicdo. Como ves nos **devuelve un dato ausente** ya que al existir datos ausentes en la variables, **la media también lo es**. Para evitar ese problema **podemos hacer dos cosas**: eliminar antes los ausentes, o indicarle en la propia media y mediana que haga el cálculo ignorando los valores `NA`.

```{r}
# Primero eliminamos NA
starwars %>% drop_na(mass) %>%
  summarise(media = mean(mass), mediana = median(mass))

# Al realizar el cálculo los ignora
starwars %>%
  summarise(media = mean(mass, na.rm = TRUE),
            mediana = median(mass, na.rm = TRUE))
```

Las <mark>**medidas de centralización**</mark> son aquellos parámetros o valores que nos informe en torno a que <mark>**valores**</mark> se concentran los datos. La <mark>**media**</mark> $\overline{x}$ es una medida de centralización basada en el valor que nos **minimiza el promedio de desviaciones al cuadrado**: la media es el valor que está más cerca de todos los puntos a la vez (con la distancia Euclídea, definida como la diferencia al cuadrado)

$$\overline{x} = \frac{1}{N} \sum_{i=1}^{N} x_i, \qquad \overline{x} = \arg \min_{x} \frac{1}{N} \sum_{i=1}^{N} \left(x_i - x \right)^2$$ 

La media <mark>**solo se puede calcular para variables cuantitativas**</mark>: solo podemos calcular medias de números.

```{r}
# Media de todas las cuantitativas
starwars_nueva %>%
  summarise(media = across(where(is.numeric), mean, na.rm = TRUE))
```

&nbsp;

La <mark>**mediana**</mark> se define como el valor que ocupa el <mark>**centro de los datos cuando los ordenamos de menor a mayor**</mark>, un valor que nos deja por debajo al menos el 50% y por encima al menos el 50%.

```{r}
# Media y mediana de todas las cuantitativas
starwars_nueva %>%
  summarise(media = across(where(is.numeric), mean, na.rm = TRUE),
            mediana = across(where(is.numeric), median, na.rm = TRUE))
```

La mediana también se puede calcular para <mark>**variables cualitativas ordinales (categorías que se puedan ordenar)**</mark>, además de para variables cuantitativas.. Vamos a construir una variable que sea «muy bajo-bajo-medio-alto-muy alto» según la estatura, y calcular la mediana para dichas categorías.

```{r}
categorias <- c("muy bajo", "bajo", "medio", "alto", "muy alto")
starwars_talla <- 
  starwars %>%
  mutate(talla =
           cut(height, breaks = c(-Inf, 80, 120, 160, 190, Inf),
               labels = categorias))

# Mediana (hay que pasárselo como número)
starwars_talla %>%
  summarise(mediana =
              categorias[median(as.numeric(talla), na.rm = TRUE)])
```

Para la **moda tenemos algún problema mayor** ya que no hay un función en los paquetes `{base}` para su cálculo directo. La <mark>**moda**</mark> se define como el <mark>**valor más repetido**</mark>: si tenemos algo discreto o cualitativo es la barra más alta.

```{r}
library(tidyverse)

# Carga
netflix <-
  read_csv('https://raw.githubusercontent.com/elartedeldato/datasets/main/netflix_titles.csv')

# Películas y series de insti
netflix_hs <- netflix %>%
  filter(str_detect(toupper(description), "HIGH SCHOOL"))
netflix_hs 

# Manipulamos fechas
library(lubridate)
netflix_final <- 
  netflix_hs %>%
  mutate(year = year(mdy(date_added))) %>%
  filter(!is.na(year))

netflix_resumen <- 
  netflix_final %>%
  group_by(year) %>%
  count() %>%
  ungroup()

library(sysfonts)
library(showtext)
font_add_google(family = "Bebas Neue", name = "Bebas Neue")
font_add_google(family = "Permanent Marker", name = "Permanent Marker")
showtext_auto()

ggplot(netflix_resumen, aes(x = year, y = n)) +
  geom_col(fill = "red") +
  scale_x_continuous(breaks = netflix_resumen$year) +
  theme_void() +
  theme(legend.position = "none",
        plot.title = element_text(family = "Bebas Neue",
                                  color = "red", size = 80)) +
  labs(title = "NETFLIX",
       subtitle = "Películas y series de instituto",
       caption = "Basada en El Arte del Dato (https://elartedeldato.com) | Datos: Kaggle") +
  theme(panel.background = element_rect(fill = "black"),
        plot.background = element_rect(fill = "black", 
                                       color = "black"),
        panel.grid.major.y =
          element_line(size = 0.1, color = "white"),
        plot.subtitle = element_text(family = "Permanent Marker",
                                     size = 21, color = "white"),
        plot.caption =  element_text(family = "Permanent Marker",
                                     color = "white", size = 19),
        axis.text = 
          element_text(size = 15, family = "Permanent Marker",
                       color = "white"),
        plot.margin = margin(t = 4, r = 4, b = 4, l = 8, "pt")) +
  annotate("text", label = "Moda", 
           x = 2020.7, y = 47, hjust = 0.2, vjust = 0, family = "Permanent Marker", size = 6, color = 'white', angle = 10) +
  annotate("curve", x = 2020.3, y = 48, xend = 2020, yend = 44,
           color = "white")
```


Para calcularla podemos usar el paquete `{modeest}`: la función `mfv()` nos calcula la **moda exacta** de una variable numérica discreta o cuantitativa (busca los valores más repetidos), la función `mlv()` nos calcula la **moda estimada** de una variable numérica continua.

```{r, warning = FALSE}
library(modeest)

# Media y mediana y moda de mass y n_films
resumen <- 
  starwars_nueva %>%
  summarise(media = across(c(mass, n_films), mean, na.rm = TRUE),
            mediana = across(c(mass, n_films), median, na.rm = TRUE),
            moda_n_films = mfv(n_films, na_rm = TRUE),
            moda_mass = mlv(mass, na.rm = TRUE))
resumen
```

Fíjate que el resumen ha agrupado todas las modas y medianas que hemos calculado a la vez: **la salida de nuestro resumen es una lista, con tablas a su vez dentro**.

&nbsp;

Veamos un ejemplo sintético de <mark>**distribuciones multimodales**</mark>, generando con `rnorm()` (genera datos aleatorias provenientes de una distribución normal) un conjunto <mark>**unimodal**</mark>, <mark>**bimodal**</mark> y <mark>**trimodal**</mark>, los 3 con la misma media.

```{r}
n <- 900

# Datos
datos <-
  tibble("unimod" = rnorm(n, 5, 1),
         "bimod" = c(rnorm(n/2, 3, 1), rnorm(n/2, 7, 1)),
         "trimod" = c(rnorm(n/3, 2, 1), rnorm(n/3, 5, 1),
                      rnorm(n/3, 8, 1)))

# Tidy data
datos <-
  pivot_longer(datos, cols = everything(),
               names_to = "tipo", values_to = "values")

# Estadisticas
datos %>% 
  group_by(tipo) %>% 
  summarise(media = mean(values),
            mediana = median(values))
```

  
Sólo utilizando la media estos casos son indistiguibles, ya que en todos los casos la media y mediana coinciden, pero si podemos hacerlo visualizando sus densidades.

```{r}
library(ggridges)
ggplot(datos,
       aes(y = tipo, x = values, fill = tipo)) +
  geom_density_ridges(alpha = 0.7, height = 1) +
  scale_fill_manual(values = met.brewer("Klimt")) +
  labs(fill = "Distribución",
       y = "Distribución",
       title = "DISTRIBUCIONES MULTIMODALES")
```

Con `stat_density_ridges()` podemos además marcar algunos valores, como por ejemplo los cuartiles.

```{r}
ggplot(datos,
       aes(y = tipo, x = values, fill = tipo)) +
  geom_density_ridges(alpha = 0.7, height = 1) +
  stat_density_ridges(quantile_lines = TRUE,
                      quantiles = c(0.25, 0.5, 0.75),
                      color = "black", alpha = .8,
                      size = 0.7) + 
  scale_fill_manual(values = met.brewer("Klimt")) +
  labs(fill = "Distribución",
       y = "Distribución",
       title = "DISTRIBUCIONES MULTIMODALES")
```



De la misma manera podemos pedirle que nos calcule <mark>**medidas de dispersión (varianza, desv. típica y cv)**</mark>, que nos indican cómo de dispersos están los datos respecto a un centro, normalmente la media, y <mark>**medidas de localización (percentiles)**</mark> (valores que nos parten los datos en trozos iguales). Dada una proporción de los datos $0 < p < 1$, el <mark>**cuantil de orden $p$**</mark> de una variable cuantitativa, denotado como $Q_p$, es el valor más pequeño tal que su frecuencia relativa acumulada es mayor o igual que $p$: es el número más pequeño que deja a su izquierda (incluyéndolo a él) como mínimo la fracción $p$ de los datos (la propia mediana es el cuantil $Q_{0.5}$.



```{r, warning = FALSE}
starwars %>%
  summarise(media = mean(mass, na.rm = TRUE),
            mediana = median(mass, na.rm = TRUE),
            moda = mlv(mass, na.rm = TRUE),
            var = var(mass, na.rm = TRUE),
            sd = sd(mass, na.rm = TRUE),
            cv = sd / abs(media),
            p13 = quantile(mass, probs = c(0.13), na.rm = TRUE),
            p87 = quantile(mass, probs = c(0.87), na.rm = TRUE))
``` 

Haciendo uso de cosas ya vistas podemos, por ejemplo, <mark>**calcular el coeficiente de variación (CV)**</mark> de todas las variables numéricas para poder decidir cual es más o menos dispersa.

```{r, warning = FALSE}
resumen <- 
  starwars %>%
  summarise(media = across(where(is.numeric), mean, na.rm = TRUE),
            var = across(where(is.numeric), var, na.rm = TRUE),
            cv = sqrt(var) / abs(media))
resumen$media
resumen$var
resumen$cv # más homogénea la altura por tener menos CV.
``` 

### Agrupando datos: group_by

Una de las funcionalidades más potentes es la opción de <mark>**añadir antes una agrupación**</mark> con `group_by()`. Esta función per se no cambia los datos sino que cambia la forma en la que se aplicarán las funciones posteriores, realizándose por **desagregadas grupos**.

Imagina que queremos calcular la **media de altura y peso** de cada personaje PERO desagregada por cada una de las clases que tenemos en la variable `sex`.

```{r}
starwars %>% group_by(sex) %>%
  summarise(media_altura = mean(height, na.rm = TRUE),
            media_peso = mean(mass, na.rm = TRUE)) %>%
  ungroup()
```

Lo que obtenemos no es la media de todos los personajes sino la **media desagregada por grupo**, grupo marcado por la variable `sex`. Las **agrupaciones pueden estar en función de varias variables** a la vez.

```{r}
starwars %>% group_by(sex, gender) %>%
  summarise(media_altura = mean(height, na.rm = TRUE),
            media_peso = mean(mass, na.rm = TRUE)) %>%
  ungroup()
```

Dicha función también es muy útil cuando queremos **realizar un filtro de registros en base al número de cada clase**: por ejemplo, vamos a filtrar los registros que pertenezcan a una clase de `sex` que tenga al menos 10 individuos dentro de dicho grupo.

```{r}
starwars %>%
  group_by(sex) %>%
  count() %>%
  ungroup()
```

Si te fijas solo deberíamos filtrar los registros que sean `female` y `male`. Vamos a hacerlo de forma automática en base a dicho **umbral**.

```{r}
# Podemos filtrar por grupos solo aquellos que superen un
# un umbral mínimo
starwars %>% 
  group_by(sex) %>% 
  filter(n() > 10) %>%
  ungroup()
```

Para comprobar rápidamente que efectivamente solo ha filtrado aquellos grupos con más de 10 elementos en ellos podemos **añadir `count()`**.

```{r}
starwars %>% 
  group_by(sex) %>% 
  filter(n() > 10) %>%
  count() %>%
  ungroup()
```

Los resúmenes y las agrupaciones podemos **combinarlas de todas las formas que nos imaginemos**, por ejemplo, calculando la media desagregada por sexo y género pero solo de las variables numéricas.

```{r}
starwars %>% # doble agrupación
  group_by(sex, gender) %>%
  summarise(across(where(is.numeric), mean, na.rm = TRUE)) %>%
  ungroup()
```

### Rowwise: operaciones por filas

Por último veamos un ejemplo de una **opción muy útil usada antes de una operación que es** `rowwise()`: toda operación que venga después se aplicará <mark>**en cada fila por separado**</mark>. Para el ejemplo vamos a definir un fichero _dummy_ de notas en mates, lengua y dibujo de 50 alumnos. Para ello usamos la función `sample()`: le decimos que seleccione aleatoriamente `size = 50` notas de un conjunto de notas posibles (`0:10`), y que lo haga con reemplazamiento (`replace = TRUE`, es decir, que dos alumnos puedan tener la misma nota)

```{r}
notas <- tibble("mates" = sample(0:10, size = 50, replace = TRUE),
                "lengua" = sample(0:10, size = 50, replace = TRUE),
                "dibujo" = sample(0:10, size = 50, replace = TRUE))
notas
```

Una vez que tenemos 3 notas aleatorias por cada uno de los 50 alumnos, **¿qué sucede si yo quiero calcular la nota media del curso?**

```{r}
notas %>% mutate(media_curso = mean(c(mates, lengua, dibujo)))
```

Como ves si aplicamos la media de las tres variables, en cada fila el valor de `media_curso` es idéntico ya que nos ha hecho la **media global**: ha tomado las 50 filas, las 3 columnas, y ha hecho la media de 150 datos. Cuando en realidad a nosotros nos gustaría **sacar una media por registro**, que para cada alumno tengamos la media de las asignaturas.

```{r}
# Aplicado por fila
notas %>% rowwise() %>%
  mutate(media_curso = mean(c(mates, lengua, dibujo)))
```

Solos nos falta poner la guinda a lo aprendido a esta introducción: vamos a ver como podemos <mark>**relacionar dos conjuntos de datos distintos entre sí**</mark>.



## Relacionar datos (joins)

Una de las **opciones más comunes para trabajar con datos** es tener a nuestra disposición diversas tablas, con alguno o varios de los campos en común, y nos interesa a veces <mark>**cruzar datos**</mark> de ambos para tener una información más completa con el conjunto de las tablas. Es lo que se conoce en ciencia de datos e informática como **hacer un _join_ de tablas**. Para ese cruce será **indispensable que haya uno o varios campos clave**, campos que sirvan para identificar unívocamente cada registro (por ejemplo, DNI).

```{r joins, echo = FALSE, fig.cap = "Esquema con los principales tipos de join, extraído de <https://estradawebgroup.com/Post/-Que-es-y-para-que-sirve-SQL-Joins-/4278>"}
knitr::include_graphics("./img/sql-joins.jpg")
```

Existen principalmente **cuatro tipos de cruces** si pensamos en cruzar un conjunto `A` con otro conjunto `B`:

- **Inner join**: solo nos quedamos con las **filas que tengan correspondencia en ambas tablas** (personas cuyo DNI aparezca en ambas tablas, por ejemplo).

- **Left (outer) join**: nos quedamos con **todas las filas de `A`**, buscando que registros de dicha tabla están también en `B`, completando los datos de esta tabla para esos registros.

- **Right (outer) join**: nos quedamos con **todas las filas de `B`**, buscando que registros de dicha tabla están también en `A`, completando los datos de esta tabla para esos registros.

- **Full join**: nos quedamos con **todas las filas de `A` y `B`**, tengan o no correspondencia en la otra tabla (si no está en una de ellas, las columnas correspondientes quedarán como campo ausente).


Esos **campos clave (keys)** serán las columnas que usaremos para definir los cruces. Para los ejemplos usaremos las tablas del paquete `{nycflights13}`.

```{r}
library(nycflights13)
```

Dicho paquete cuenta con las siguientes **tablas**:

* `airlines`:	nombre de la aerolínea (con su abreviatura).
* `airports`: datos de aeropuertos (nombres, longitud, latitud, altitud, etc).
* `flights`: datos de vuelos (con `tailnum` como marca de idenfiticación).
* `planes`:	datos de los aviones.
* `weather`: datos meteorológicos horarios de las estaciones LGA, JFK y  EWR.

Veamos un ejemplo: imagina que queremos **completar en la tabla de vuelos los datos de cada una de las aerolíneas que operan dichos vuelos**. Vamos a seleccionar unas pocas columnas para que sea más fácil de visualizar.

```{r}
# Seleccionamos antes columnas para que sea más corto
flights_filtrada <- flights %>%
  select(year:day, arr_time, carrier:dest)
flights_filtrada
```

Queremos **TODAS las filas de los vuelos**, todos sus registros, pero **añadiendo la información** que tenemos de la aerolínea que opere los vuelos, así que haremos un **LEFT JOIN** de `flights` vs `airlines`. El campo común que nos permite cruzarla, la **clave (key)** es el código abreviado de las aerolíneas (variable `carrier`).

```{r}
# Mismas filas pero con una nueva columna: siempre que sea
# posible el cruce tendrá la info de la aerolínea
l_join_flights_airlines <-
  flights_filtrada %>% left_join(airlines, by = "carrier")
l_join_flights_airlines
```

Si te fijas ahora tenemos en `l_join_flights_airlines` las mismas 336776 filas pero con una columna más: la tabla `airlines` tenía 2 columnas, una la común con `flights` y en otra la nueva que se ha incorporado a la tabla. **¿Y si en `flights` había algún vuelo operado por alguna aerolínea que no estuviese en `airlines`?**

```{r}
# Filtramos
l_join_flights_airlines %>% filter(is.na(name))
```

En este caso todos los registros de la primera tabla tenían su correspondencia en la segunda tabla. Veamos un ejemplo donde suceda eso, quitándole algunas filas a `airlines`, quitando las aerolíneas con códigos `"B6"`.

```{r}
airlines_filtrada <- airlines %>% 
  filter(!(carrier %in% c("B6")))
l_join <- flights_filtrada %>% left_join(airlines_filtrada, by = "carrier")
dim(l_join %>% filter(is.na(name)))
```

En el anterior ejemplo tenemos 54 635 filas que cuyo nombre está ausente, es decir, 54 635 filas de `flights` que no tienen correspondencia en la tabla filtrada `airlines_filtrada`. Ahora las **filas que no han encontrado su _match_ en la segunda tabla están como campo ausente**.

```{r}
l_join
```

```{r left-join, echo = FALSE, fig.cap = "Esquema del left join, extraído de <https://r4ds.had.co.nz/relational-data.html#mutating-joins>"}
knitr::include_graphics("./img/left_join.jpg")
```

La misma lógica para los **RIGHT JOIN** y los **FULL JOIN**. En este último tendremos bastantes datos ausentes, ya que todos los registros que no estén en ambas tablas tendrán una parte de sus variables como `NA`


```{r}
tabla1 <- tibble("key_1" = 1:7, "valor_1" = paste0("x", 1:7))
tabla2 <- tibble("key_2" = c(0, 1, 2, 5, 6, 9, 10),
                 "valor_2" = paste0("x", c(0, 1, 2, 5, 6, 9, 10)))

# Left
tabla1 %>% left_join(tabla2, by = c("key_1" = "key_2"))

# Right
tabla1 %>% right_join(tabla2, by = c("key_1" = "key_2"))

# Full
tabla1 %>% full_join(tabla2, by = c("key_1" = "key_2"))
```

```{r right-join, echo = FALSE, fig.cap = "Esquema del right join, extraído de <https://r4ds.had.co.nz/relational-data.html#mutating-joins>"}
knitr::include_graphics("./img/right_join.jpg")
```

```{r full-join, echo = FALSE, fig.cap = "Esquema del full join, extraído de <https://r4ds.had.co.nz/relational-data.html#mutating-joins>"}
knitr::include_graphics("./img/full_join.jpg")
```

Un tipo de join especial son los **INNER JOIN**. Hasta ahora la tabla resultante tenía al menos tantas filas como tuviese la tabla más pequeña en el cruce. En el caso de los **INNER JOIN** vamos a **reducir el tamaño** ya que solo nos quedaremos con **aquellos registros que podamos encontrar en ambas**, de forma que **el cruce nunca generará datos ausentes**.

```{r}
# Inner
tabla1 %>% inner_join(tabla2, by = c("key_1" = "key_2"))
```



```{r inner-join, echo = FALSE, fig.cap = "Esquema del inner join, extraído de <https://r4ds.had.co.nz/relational-data.html#mutating-joins>"}
knitr::include_graphics("./img/inner_join.png")
```

## 📝 Ejercicios

en construcción (cuarteto de anscombe)