---
title: "Clase 4. Explorando datos cuantitativos"
subtitle: "Contenidos conceptuales"
output: 
  html_document:
    toc: TRUE
    toc_float: TRUE
author: "Germán Rosati y Laia Domenech Burin"
date: "`r format(Sys.time(), '%d %B, %Y')`"
theme: Boadilla
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo=FALSE)
#install.packages("gapminder")
```

```{r message=FALSE, warning=FALSE, include=FALSE, paged.print=FALSE}
knitr::opts_chunk$set(echo=TRUE, message=FALSE, warning=FALSE, tidy=TRUE)
```

```{r warning=TRUE, include=FALSE}
library(gapminder)
library(tidyverse)
library(eph)
library(gt)
df <- get_microdata(year = 2016, trimester = 3, type = "individual")

df <- df %>% calculate_poverty(canastas_reg_example)
df <- organize_labels(df=df, type='individual') # Vamos a etiquetar la base
df <- df %>% mutate_at(vars(REGION, AGLOMERADO, CH03, CH04, CH07:CAT_INAC), ~as.character(.))

df <- df %>%
        mutate(REGION = if_else(REGION == '44', 'Patagonia', REGION))
```

Pensemos, ahora, en variables como el ingreso laboral `P21` o la edad `CH06`. Ambas variables son cuantitativas en tanto podemos calcular la diferencia entre dos ingresos laborales o entre dos valores de edad. Pero variables como el código postal no refieren a magnitudes sino a zonas o áreas. Por lo tanto, no tiene sentido considerar la diferencia entre un código postal correspondiente a Villa Devoto (1417) y otro perteneciente a Villa Luro (1407).

Vamos a ver, entonces, ahora algunos métodos para realizar un análisis descriptivo de variables cuantitativas.

## Gráficos de dispersión (o scatterplots)
Un diagrama de dispersión proporciona una vista caso por caso de los datos de dos variables numéricas. En la Figura siguiente, se usó un diagrama de dispersión para examinar el promedio de ingresos laborales y no laborales de cada aglomerado 

```{r echo=FALSE}
df %>%
  group_by(AGLOMERADO) %>%
  summarise(media_ing_lab = mean(P21),
            media_ing_no_lab = mean(T_VI)) %>%
  ggplot() +
    geom_point(aes(x=media_ing_lab, y=media_ing_no_lab),
               size=2, color='lightblue') +
    labs(
        title='Medias de ingresos laborales y no laborales por aglomerado',
        x='Media de ingresos laborales',
         y='Media de ingresos no laborales') +
    scale_x_continuous(labels=scales::dollar_format(), limits=c(0,10000)) +
    scale_y_continuous(labels=scales::dollar_format(), limits=c(0,2500)) +
    theme_minimal()

```

En cualquier diagrama de dispersión, cada punto representa un solo caso. Como hay 32 aglomerados urbanos relevados en la EPH hay 32 puntos en la Figura 5.1.

Se ve que parece existir una cierta relación entre ambas variables. A medida que se incrementan los ingresos laborales lo mismo parece suceder con los no laborales. No obstante, existen algunos puntos en los que esta relación no se cumple. ¿Cuáles son?

A continuación hemos calculado el porcentaje de personas pobres y la media del ingreso total del hogar para cada aglomerado.

```{r echo=FALSE}
df %>%
  filter(ITF > 0) %>%
  group_by(AGLOMERADO) %>%
  summarise(media_ing_lab = mean(ITF)) %>%
  left_join(df %>%
  mutate(situacion = if_else(situacion %in% c('indigente', 'pobre'), 'pobre', 'no pobre')) %>%
  group_by(AGLOMERADO, situacion) %>%
  summarise(n=n()) %>%
  mutate(pobreza = n/sum(n)) %>%
  filter(situacion == 'pobre') %>%
  select(-n, -situacion)) %>%
  ggplot() +
    geom_point(aes(x=pobreza, y=media_ing_lab),
               size=2, color='lightblue') +
    geom_smooth(aes(x=pobreza, y=media_ing_lab),
                size=0.5,
               linetype='dashed', color='red', se=FALSE) +
    labs(
        title='Aglomerados según tasa de pobreza y media de ITF',
        x='Tasa de pobreza',
         y='Media de ingresos total del hogar') +
    scale_x_continuous(labels=scales::percent_format(), limits=c(0,0.6)) +
    scale_y_continuous(labels=scales::dollar_format(), limits=c(0,40000)) +
    theme_minimal()

```

¿Qué relación pueden ver acá? La relación es evidentemente no lineal, como lo destaca la línea discontinua. Esto es diferente de los diagramas de dispersión anteriores que hemos visto, que indican muy poca o ninguna curvatura en la tendencia.

## Medidas de tendencia central: media

```{r echo=FALSE}
df %>%
  #filter(P21 > 100) %>%
  mutate(CH06 = as.numeric(as.character(CH06))) %>%
  ggplot() +
    geom_histogram(aes(x=CH06), bins=100, color='white') +
    geom_vline(xintercept = mean(df$CH06), color='red') +
    theme_minimal() +
    labs(
            title = 'Histograma de edad',
            x='Edad en años simples',
            y='Frecuencia')

```

La media, a menudo llamada promedio, es una forma común de medir el centro de una distribución de datos. Para calcular la tasa de interés media, sumamos todas las tasas de interés y las dividimos por el número de observaciones.

La media de la muestra a menudo se etiqueta $\overline x$. Es útil pensar en la media como el punto de equilibrio de la distribución, y se muestra como una línea roja en la figura anterior.

La media se calcula de la siguiente forma: 
$$\overline x = \frac{x_{1} + x_{2} + ... x_{n}}{n} = \frac{\sum_{i=1}^n x_{i}}{n}$$

El conjunto de datos de la EPH representa una muestra de una población más grande (¿cuál?). Podríamos calcular una media para toda la población de la misma manera que la media de la muestra. Sin embargo, la media poblacional tiene una etiqueta especial: $\mu$. Es la letra griega mu y representa el promedio de todas las observaciones en la población. A veces un subíndice, como $_{x}$  se utiliza para representar a qué variable se refiere la media de la población ($\mu_{x}$).  A menudo es demasiado costoso medir la media de la población con precisión, por lo que a menudo estimamos $\mu$ utilizando $\overline x$.

Aquí estamos trabajando con un muestra, con lo cual no podemos calcular la media de edad exacta de la población... pero podemos hacer una _estimacion_ de ese valor usando los datos de la muestra. Esos 34.2 años que obtuvimos, son una estimación un tanto rústica de la media de edad de la población urbana ($\mu_{edad}$). Si bien no es perfecta, esta es nuestra mejor _estimación puntual_ de la edad media de toda la población urbana en Argentina. Más adelante, desarrollaremos herramientas para caracterizar la precisión de las estimaciones puntuales, como la media muestral. Como podrán intuir, las estimaciones puntuales basadas en muestras más grandes tienden a ser más precisas que las basadas en muestras más pequeñas.

### Medias ponderadas
Supongamos que queremos calcular [la esperanza de vida](https://es.wikipedia.org/wiki/Esperanza_de_vida) promedio en América durante 2007 de la siguiente tabla:

```{r echo=FALSE}
gapminder %>%
  filter(continent=="Americas" & year==2007) %>%
  select(country, continent, lifeExp, pop) %>%
  gt()
```

Podríamos tomar las esperanzas de vida de cada uno de los países en esta región y promediarla. ¿Qué problemas se les ocurre que tiene este enfoque? ¿Cómo podríamos hacerlo mejor?

El conjunto de datos es especial porque cada país en realidad representa a muchas personas individuales. Si tuviéramos que promediar simplemente la variable `lifeExp`, estaríamos tratando condados con 1.056.608 y 301.139.947 de residentes por igual en los cálculos. En su lugar, debemos calcular la expectativa de vida total de cada país, sumar los totales de todos los países y luego dividir por la cantidad de personas en todos los países. 

```{r echo=FALSE}
gapminder %>%
  filter(continent=="Americas" & year==2007) %>%
  select(country, continent, lifeExp, pop) %>%
  summarise(max = max(pop),
            min = min(pop),
            mean = mean(lifeExp),
            mean_weigted = weighted.mean(lifeExp, w=pop)) %>%
  gt()
```

Si completamos estos pasos encontraríamos que la expectativa de vida de la región es de 75.5 años. Si hubiéramos calculado la media simple el resultado hubiera sido solo 73.6. Dos años son una gran diferencia en términos de estimaciones poblacionales.

Este ejemplo usó lo que se llama una media ponderada, la formula es simple:

$$\overline x = \frac{w_{1}x_{1} + w_{1}x_{1} + ... + w_{n}x_{n}}{w_{1} + w_{2} + ... w_{n}} = \frac{\sum_{i=1}^n w_{i}x_{i}}{\sum_{i=1}^n w_{i}}$$ 


## Histogramas
Cuando tenemos muestras grandes y varaibles cuantitativas suele ser útil graficar cada caso como perteneciendo a un contenedor o a un "agrupamiento". Por ejemplo, en la EPH podríamos crear una tabla de conteos para los ingresos laborales (personas de 0 a $1000 años, de $1001 a $5000, etc.). Las observaciones que caen en el límite de un agrupamiento se asignan al contenedor inferior. Luego, a partir de este podemos generar un histograma:

```{r echo=FALSE}
df %>%
  mutate(P21 = as.numeric(as.character(P21))) %>%
  mutate(P21 = case_when(
    P21 >= 75000 ~ 75000, 
    P21 < 75000 ~ P21)) %>%
  filter(P21 > 0) %>%
  ggplot() +
    geom_histogram(aes(x=P21), color='white', fill='lightblue', bins = 50) +
    scale_x_continuous(breaks=seq(0,75000, 8000)) +
    theme_minimal() +
    labs(
            title = 'Histograma ingresos laborales',
            x='Ingresos laborales',
         y='Frec.')

``` 
Los histogramas son especialmente convenientes para comprender la forma de la distribución de datos. La figura anterior sugiere que la mayoría de las personas tenían en 2016 ingresos laborales menores a \$16000, mientras que solo unas pocas personas tienen ingresos superiores al \$32000. Cuando la distribución de una variable se desvía hacia la derecha de esta manera y tiene una cola derecha más larga, se dice que la forma está sesgada a la derecha.

La figura siguiente muestra un gráfico de densidad que es un histograma suavizado. Los detalles técnicos sobre cómo dibujar gráficos de densidad (precisamente cómo suavizar el histograma) están más allá del alcance de esta materia. Pero se nota que la forma, la escala y la dispersión de las observaciones se muestran de manera similar en un histograma que en un gráfico de densidad.

```{r}
df %>%
  mutate(P21 = as.numeric(as.character(P21))) %>%
  mutate(P21 = case_when(
    P21 >= 75000 ~ 75000, 
    P21 < 75000 ~ P21)) %>%
  filter(P21 > 0) %>%
  ggplot() +
    geom_density(aes(x=P21), color='white', fill='lightblue', bins = 50) +
    scale_x_continuous(breaks=seq(0,75000, 8000)) +
    theme_minimal() +
    labs(
            title = 'Density plot de ingresos laborales',
            x='Ingresos laborales',
         y='Frec.')
```

Se dice que las variables con la característica inversa (una cola larga y más delgada hacia la izquierda) están sesgadas (o son asimétricas) a la izquierda. También decimos que tal distribución tiene una larga cola izquierda. Las variables que muestran un final aproximadamente igual en ambas direcciones se denominan simétricas.

Además de observar si una distribución es asimétrica o simétrica, los histogramas se pueden usar para identificar modas. Una moda está representada por un pico prominente en la distribución. ¿Cuántos picos hay en la distribución de ingresos?

Se observa un pico prominente. Este tipo de distribuciones se denominan unimodal. Cualquier distribución con más de dos picos prominentes se denomina multimodal. 

Buscar modas no se trata de encontrar una respuesta clara y correcta sobre la cantidad de modas en una distribución, razón por la cual la noción de "prominente" no se define rigurosamente. La parte más importante de este examen es comprender mejor los datos.

La media se introdujo como un método para describir el centro de una variable. Ahora bien, la variabilidad en los datos también es importante. Aquí presentamos dos medidas de variabilidad: la varianza y la desviación estándar. Ambos son muy útiles en el análisis de datos, aunque sus fórmulas son un poco tediosas para calcularlas a mano. La desviación estándar es la más fácil de comprender de las dos, ya que describe aproximadamente qué tan lejos está la observación típica de la media.

## Varianza y desviación estándar
Llamamos a la distancia de una observación de su media su desviación. A continuación se muestran las desviaciones de los ingresos laborales para el cinco casos muestreados aleatoriamiente de la base de la EPH:

$$ edad_{1} - \overline x_{edad} = 17 - 34.22 =  -17.22$$
$$ edad_{2} - \overline x_{edad} = 50 - 34.22 =  15.78$$
$$ ... $$
$$ edad_{59550} - \overline x_{edad} = 22 - 34.22 =  -12.22$$

Si elevamos al cuadrado estos valores y los promediamos el resultado es igual a la varianza muestral, denotada como $s^2$:

$$ s^2 = \frac{(-12.22)^2 + 15.78^2 + ... + (-12.22)^2 +}{59550-1} = 496.29$$

Dividimos por $n-1$ en lugar de hacerlo por $n$ porque se trata de la varianza muestral (un estimador de la varianza poblacional). Es importante notar que elevar al cuadrado hace dos cosas: 

- elimina los valores negativos
- hace los valores grandes más grandes aún en términos relativos

La desviación estándar de la muestra se puede calcular como la raíz cuadrada de la suma de la distancia al cuadrado de cada valor desde la media dividida por el número de observaciones menos uno:

$$ s = \sqrt \frac{\sum_{i=1}(x_{i}^n - \overline x)^2}{n-1}$$
En nuestro caso, el desvío estándar es igual a 

$$s_{edad} = \sqrt 496.29  = 22.28$$
Aunque a menudo se omite, un subíndice de la variable en cuestión puede agregarse a la notación de la varianza y el desvío estándar como recordatorio de que estas son la varianza y la desviación estándar de las observaciones representadas por cierta variable (en nuestro caso, la edad)

La varianza es la distancia promedio al cuadrado de la media. La desviación estándar es la raíz cuadrada de la varianza. La desviación estándar es útil cuando se considera qué tan lejos se distribuyen los datos de la media.

La desviación estándar representa la desviación típica de las observaciones de la media. A menudo, alrededor del 68 % de los datos estarán dentro de una desviación estándar de la media y alrededor del 95 % estarán dentro de dos desviaciones estándar. Sin embargo, estos porcentajes no son reglas estrictas.

Al igual que la media, los valores de población para la varianza y la desviación estándar tienen símbolos especiales: $\sigma^2$  por la varianza y $\sigma$ para la desviación estándar.


## Boxplots, cuantiles y mediana

Un boxplot resume un dataset usando cinco números, cinco estadísticos.
```{r echo=FALSE}
df %>%
  mutate(P21 = as.numeric(as.character(P21))) %>%
  mutate(P21 = case_when(
    P21 >= 75000 ~ 75000, 
    P21 < 75000 ~ P21)) %>%
  filter(P21 > 0) %>%
  ggplot() +
    geom_boxplot(aes(x=P21)) +
    theme_minimal() +
    labs(x='Ingresos laborales',
         title="Boxplot ingresos laborales")
```

La línea oscura dentro del cuadro representa la _mediana_, que divide los datos por la mitad. El 50% de los datos caen por debajo de este valor y el 50% por encima. Dado que en el conjunto de datos de préstamo50 hay 50 observaciones (un número par), la mediana se define como el promedio de las dos observaciones centrales. Cuando hay un número impar de observaciones, habrá exactamente una observación que divide los datos en dos mitades y, en tal caso, esa observación es la mediana (no se necesita un promedio).

Si los datos están ordenados de menor a mayor, la mediana es la observación justo en el medio. Si hay un número par de observaciones, habrá dos valores en el medio y la mediana se tomará como su promedio.

El segundo paso en la construcción de un boxplot es dibujar un rectángulo para representar el 50% medio de los datos. La longitud de la caja se llama rango intercuartílico (o IQR por sus siglas en inglés para abreviar). Al igual que la desviación estándar, es una medida de la variabilidad de los datos. Cuanto más variables sean los datos, mayor será la desviación estándar y el IQR. Los dos límites de la caja se denominan primer cuartil (es decir, el valor que deja por debajo el 25% de los datos) y el tercer cuartil (es decir, el 75% de los datos caen por debajo de este valor), y estos a menudo se etiquetan $Q1$ $Q3$ respectivamente.

El rango intercuartílico $IQR$ es la longitud de la caja en un diagrama de caja. se calcula como $Q3-Q1$. Un $\alpha$ percentil es un número con $\alpha$\% de observaciones por encima y $100 - \alpha$\% de observaciones por debajo. Por ejemplo, el 90 avo. percentil de puntajes SAT es el valor del puntaje SAT con el 90% de los estudiantes por debajo de ese valor y el 10% de los estudiantes por encima de ese valor.

Extendiéndose desde la caja, los bigotes intentan capturar los datos fuera de la caja. Los bigotes de un diagrama de caja alcanzan los valores mínimo y máximo de los datos, a menos que haya puntos que se consideren inusualmente altos o inusualmente bajos, que el diagrama de caja identifica como posibles valores atípicos. Estos están etiquetados con un punto en el diagrama de caja. El propósito de etiquetar los puntos periféricos, en lugar de extender los bigotes a los valores mínimos y máximos observados, es ayudar a identificar cualquier observación que parezca inusualmente distante del resto de los datos. Hay una variedad de fórmulas para determinar si un punto de datos en particular se considera un valor atípico, y diferentes software estadísticos usan fórmulas diferentes. Una fórmula comúnmente utilizada es que cualquier observación más allá $1.5 \times IQR$ lejos del primer o tercer cuartil se considera un valor atípico. En cierto sentido, la caja es como el cuerpo del diagrama de caja y los bigotes son como sus brazos que intentan alcanzar el resto de los datos, hasta los valores atípicos.

Un _valor atípico_ es una observación que parece extrema en relación con el resto de los datos. Examinar los datos en busca de valores atípicos tiene muchos propósitos útiles, que incluyen

- identificando un fuerte sesgo en la distribución,
- identificar posibles errores de recopilación o ingreso de datos, y
- proporcionando información sobre propiedades interesantes de los datos.

Tenga en cuenta, sin embargo, que algunos conjuntos de datos tienen un sesgo naturalmente largo y los puntos periféricos no representan ningún tipo de problema en el conjunto de datos.