R-Conomics

Introducción

En la última entrada analizamos los gráficos de ranking, los cuales se basan principalmente en ordenar los datos de mayor a menor (o al revés) para poder observar cuál de ellos es el que ha tenido una mayor relevancia/valor en un periodo \(t\). En esta ocasión analizaré algo más estadístico, ya que voy a centrarme en los gráficos de distribución de frecuencias, los cuales tienen como característica principal el mostrarnos donde se encuentra la media de la información que tenemos, mostrándonos las zonas en las que tienden a agruparse los datos. Esto es esencialmente útil cuando trabajamos con análisis de información, ya que generalmente necesitamos distribuciones de frecuencia de tipo normal, por lo que verlas gráficamente nos puede ayudar a entender mejor la información con la que estamos trabajando.

Histograma

El gráfico más representativo en este caso será un histograma, ya que directamente mediante el uso de barras nos muestra la distribución que tienen nuestros datos, es decir, su frecuencia por cada valor, por lo que podremos observar con cierta claridad la zona en la que se agrupan mayormente los datos.

Vamos a ver un ejemplo de ello con los datos de gapminder, pero antes vamos a plantearlo de manera lógica. Supongamos que queremos comparar la expectativa de vida de 1957 con el de 2007 para saber si este se ha incrementado a lo largo del tiempo. Lo que podríamos hacer es graficar un histograma para cada año, por lo que sabríamos hacía que valor tiende a agruparse dicha expectativa para cada año respectivamente y, de esa manera, sabriamos si hubo un cambio notable en la calidad de vida de las personas. Voy a hacer el ejercicio a continuación:

# Cargamos los paquetes
library(echarts4r)
library(hrbrthemes)
library(dplyr)
library(gapminder)
library(ggplot2)
library(pyramid)
library(beeswarm)
library(viridis)
# Observamos los datos que usaremos
head(gapminder)
## # A tibble: 6 x 6
##   country     continent  year lifeExp      pop gdpPercap
##   <fct>       <fct>     <int>   <dbl>    <int>     <dbl>
## 1 Afghanistan Asia       1952    28.8  8425333      779.
## 2 Afghanistan Asia       1957    30.3  9240934      821.
## 3 Afghanistan Asia       1962    32.0 10267083      853.
## 4 Afghanistan Asia       1967    34.0 11537966      836.
## 5 Afghanistan Asia       1972    36.1 13079460      740.
## 6 Afghanistan Asia       1977    38.4 14880372      786.
# Creamos un dataset para cada año requerido
calidad_1957 <- gapminder %>% 
  filter(year == 1957)

calidad_2007 <- gapminder %>% 
  filter(year == 2007)


# Creamos un histograma para cada año
## Para 1957
calidad_1957 %>% 
  e_charts() %>% 
  e_histogram(lifeExp) %>% 
  e_title("Distribución de la expectativa de vida para 1957") %>% 
  e_axis_labels(x = "Edad", y = "Frecuencia") %>% 
  e_tooltip(trigger = "axis") 
## Para 2007
calidad_2007 %>% 
  e_charts() %>% 
  e_histogram(lifeExp) %>% 
  e_title("Distribución de la expectativa de vida para 2007") %>% 
  e_tooltip(trigger = "axis") 

Muy bien, ahora ¿qué podemos inferir de los gráficos obtenidos? Podemos observar que para 1957 el mayor pico en los datos fue en 42.5, lo que implica que la expectativa de vida promedio tendía a un valor cercano a ello (aunque eso ya depende de la cantidad de valor atípicos que tuviesemos), mientras que para 2007 el pico mayor en los datos (y por mucho) fue de 72.5 años, lo que implica que la expectativa de vida se incrementó bastante en un periodo de 50 años, aumentando su frecuencia en 30 años.

Diagrama de densidad

Un diagrama de densidad o curva de densidad es básicamente lo mismo que un histograma, solo que en una versión mucho más suavizada, ya que no se utilizan barras, si no una simple linea que de igual manera nos muestra la distribución de los datos observados, pero esto en porcentajes que deben dar un total de 1. De hecho, mediante un diagrama de densidad o una función de densidad de probabilidad se puede hacer uso de pruebas de hipótesis, aunque eso lo veremos después.

Los diagramas de densidad tienen un valor total de 1 ya que son, en esencia, una forma de distribuir las observaciones de una población o muestra igual a 1, y como son un área bajo la curva que puede ser calculada entonces tendríamos en este caso una integral definida para una variable aleatoria \(X\) de la siguiente manera:

\[ \alpha < X < \beta = \int_{\alpha}^{\beta} f(x) \, dx \]

Donde \(\alpha\) es el valor más pequeño observado y \(\beta\) el más grande. Lo que estoy diciendo con esto es que el diagrama de densidad nos dice el valor representativo en un ratio que va de 0 a 1 de la frecuencia de los datos observados, por lo que aquel valor más cercano a 1 será el que tenga una mayor frecuencia. Vamos a hacer el diagrama de densidad respectivo para los datos anteriormente graficados:

# Creamos un diagrama de densidad para cada año
## Para 1957
calidad_1957 %>% 
  e_charts() %>% 
  e_density(lifeExp) %>% 
  e_title("Distribución de la expectativa de vida para 1957") %>% 
  e_axis_labels(x = "Edad", y = "Frecuencia") %>% 
  e_tooltip(trigger = "axis")
## Para 2007
calidad_2007 %>% 
  e_charts() %>% 
  e_density(lifeExp) %>% 
  e_title("Distribución de la expectativa de vida para 2007") %>% 
  e_tooltip(trigger = "axis") 

Como podemos observar, la distribución de los datos es la misma, pero ahora tenemos una curva mucho más suavizada que nos muestra de mejor manera los picos de frecuencias.

Es posible combinar los diagramas de densidad con los histogramas, esto lo podemos hacer de la siguiente manera:

# Gráficos combinados

## Para 1957
calidad_1957 %>% 
  e_charts() %>% 
  e_histogram(lifeExp) %>% 
  e_title("Distribución de la expectativa de vida para 1957") %>% 
  e_axis_labels(x = "Edad", y = "Frecuencia") %>% 
  e_density(lifeExp, name = "density", areaStyle = list(opacity = .4), 
            smooth = TRUE, y_index = 1) %>% 
  e_tooltip(trigger = "axis")
## Para 2007
calidad_2007 %>% 
  e_charts() %>% 
  e_histogram(lifeExp) %>% 
  e_title("Distribución de la expectativa de vida para 2007") %>% 
  e_density(lifeExp, name = "density", areaStyle = list(opacity = .4), 
            smooth = TRUE, y_index = 1) %>% 
  e_tooltip(trigger = "axis") 

Y como podemos observar, los datos son similares tanto para el diagrama de densidad como para el histograma. Es importante no olvidar considerar el argumento y_index = 1 dentro de la función de densidad que hicimos para el gráfico, ya que esto nos creará un doble eje Y y podremos comparar de manera simétrica la información.

Puntos de Cleveland

Estos gráficos nos dicen la diferencia que puede existir entre el valor de una variable, por ejemplo, a lo largo de dos periodos de tiempo o para dos observaciones. Lo que hacen es graficar el valor inicial y el valor final, para luego unir dichos valores con una linea entre ambos y, de esta forma, poder visualizar mejor la distancia que existe entre ambos. La página de R Graph Gallery muestra un ejemplo interesante de este tipo de gráficos utilizando, pero utiliza datos ficticios sin sentido aparente, lo que a mi parecer no termiaría de dar a entender el uso que se le puede dar a estos gráficos.

Para poder dejar más en claro el uso que se le puede dar a este tipo de gráficos, se me ocurre utilizar la base gapminder para ver como ha evolucionado la expectativa de vida. Esto lo podríamos hacer creando dos bases filtradas con los datos para cada año respectivo, para luego unirlos en una base final de la siguiente manera:

# Base 1
gap_1 <- gapminder %>% 
  filter(year == 1957) %>% 
  rename(lifeExp1 = lifeExp) %>% 
  top_n(10)

# Base 2
gap_2 <- gapminder %>% 
  filter(year == 2007) %>% 
  rename(lifeExp2 = lifeExp) %>% 
  top_n(10)

# Leemos ambas
attach(gap_1)
attach(gap_2)

# Creamos un data frame con todos los datos
gap_3 <- data.frame(country, lifeExp1, lifeExp2)

# Agregamos una media, reducimos observaciones y reordenamos los datos del data frame
gap_3 <- gap_3 %>% 
  rowwise() %>% 
  mutate(media = mean(c(lifeExp1, lifeExp2))) %>% 
  arrange(media)

Ahora que ya tenemos los datos que necesitamos, vamos a crear el gráfico:

# Gráfico de puntos de Cleveland
ggplot(gap_3) +
  # Para crear la barra del centro
  geom_segment( aes(x=country, xend=country, y=lifeExp1, yend=lifeExp2), color="grey") +
  # Para crear el punto verde inicial
  geom_point( aes(x=country, y=lifeExp1), color=rgb(0.2,0.7,0.1,0.9), size=3 ) +
  # Para crear el punto rojo final
  geom_point( aes(x=country, y=lifeExp2), color=rgb(0.7,0.2,0.1,0.9), size=3 ) +
  # Para invertir las coordenadas
  coord_flip()+
  # Fondo y tipo de fuente
  theme_ipsum() +
  # Quitamos la leyenda de la derecha
  theme(
    legend.position = "none",
  ) +
  # Agregamos nombre a los ejes
  xlab("Países") +
  ylab("Expectativa de vida")+
  labs(
    title = "Ejercicio de gráfico para ver la expectativa de vida.",
    subtitle = "Y su cambio en dos periodos de tiempo."
  )

Notese que usamos la media para reordenar los datos, ya que esto nos da una mejor homogeneidad en la información.

gráficos de caja (Boxplots)

Unos de los más representativos que podemos usar cuando trabajamos cone estadística y, sobre todo, análisis de la varianza (ANOVA, MANOVA, ANCOVA), ya que éstos gráficos nos dicen la distribución que tienen varias observaciones clasificadas por diferentes características.

Podríamos utilizar el ejemplo de la base de datos iris para ver como se distribuyen los pétalos de las flores por especie, lo que haremos de la siguiente manera:

# Observamos los datos que usaremos
head(iris)
##   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1          5.1         3.5          1.4         0.2  setosa
## 2          4.9         3.0          1.4         0.2  setosa
## 3          4.7         3.2          1.3         0.2  setosa
## 4          4.6         3.1          1.5         0.2  setosa
## 5          5.0         3.6          1.4         0.2  setosa
## 6          5.4         3.9          1.7         0.4  setosa
# Gráfico de caja
ggplot(iris, aes(x=Species, y=Petal.Length, fill=Species)) +
  geom_boxplot() +
  scale_fill_viridis(discrete = TRUE, alpha=0.6, option="A")  + 
  theme_ipsum() + 
  xlab("Especie") + 
  ylab("Distribución")+ 
  theme(plot.title = element_text(size = 11, face = "bold", hjust = 0)) + 
  theme(plot.subtitle = element_text(size = 10, hjust = 0)) + 
  theme(plot.caption = element_text(size = 10, hjust = 0)) +
  theme(plot.margin = unit(c(1,1,1,1), "cm")) +
  labs(
    title = "Ejemplo de un gráfico de caja"
  )

Bien, ahora vamos a interpretarlo. Como podemos observar, el gráfico nos dice que la especie setosa es la que tiene en promedio un pétalo de tamaño más pequeño, mientras que las de la especie virgínica tienen un tamaño mucho mayor, llegando a más de 6. La lina dentro de la caja nos indica donde se encuentra la media, mientras que el resto de la caja nos muestra cuantiles de los datos, además de que las lineas que sobresalen de las cajas nos dicen los valores mínimos y máximos observados y, por último, los puntos son valores atípicos.

Podemos hacer lo mismo con un gráfico del paquete echarts4r de la siguiente manera:

# Boxplot con echarts4r
iris %>% 
  group_by(Species) %>% 
  e_charts() %>% 
  e_boxplot(Petal.Length)  %>% 
  e_title("Gráfico de cajas") %>% 
  e_axis_labels(x = "Especie", y = "Valor") %>% 
  e_tooltip(trigger = "axis")

Personalmente recomiendo hacerlo con éste último paquete, ya que al pasar el puntero sobre cada caja nos dice información estadística de ella, como la media, cuartiles, outliers (valores atípicos), etc.

Gráficos de violín

Estos gráficos nos dicen básicamente lo mismo que uno de caja, pero la diferencia radica en que nos muestran de mejor manera la distribución que tienen los datos, por lo que si encontramos una variable con datos muy homogeneos el gráfico de violín tenderá a ser más delgado, mientras que si encontramos otra con una forma más leptocúrtica el gráfico será más gordo. No me detendré demasiado a explicarlo ya que la teoría es la misma que con los gráficos de cajas. La forma de graficarlo es la misma en ggplot2, solo que en lugar de utilizar geom_boxplot() utilizaremos geom_violin(), de la siguiente manera:

# Gráfico de violín
ggplot(iris, aes(x=Species, y=Petal.Length, fill=Species)) +
  geom_violin() +
  scale_fill_viridis(discrete = TRUE, alpha=0.6, option="A")  + 
  theme_ipsum() + 
  xlab("Especie") + 
  ylab("Distribución")+ 
  theme(plot.title = element_text(size = 11, face = "bold", hjust = 0)) + 
  theme(plot.subtitle = element_text(size = 10, hjust = 0)) + 
  theme(plot.caption = element_text(size = 10, hjust = 0)) +
  theme(plot.margin = unit(c(1,1,1,1), "cm")) +
  labs(
    title = "Ejemplo de un gráfico de violín"
  )

Piramide poblacional

Estos gráficos son usados casi exclusivamente para ver como ha variado la distribución de la población para ambos sexos. En este caso contamos con un paquete bastante útil para graficar lo que necesitamos, llamado pyramid. El paquete en cuestión nos graficará automáticamente una piramide de población para una base de datos que nos diga la cantidad de de personas que tenemos por edad, divididas por sexo. Por suerte el mismo paquete nos incluye bases de datos para practicar, como la que usaremos a continuación:

# Cargamos los datos del paquete
datos <- GunmaPop2010

# Leemos la base
head(datos)
##   Males Females Ages
## 1  8135    7612    0
## 2  8555    7952    1
## 3  8869    8176    2
## 4  8607    8397    3
## 5  8688    8438    4
## 6  9001    8458    5

Ahora, podemos observar que tenemos una columna para edad, una para cantidad de hombres y otra para cantidad de mujeres por edad, por lo que ya sabemos cómo están estructurados los datos que vamos a usar. Vamos con el gráfico:

# Para cambiar tamaño del gráfico
par( mar=c(2,2,2,2))

# El gráfico
pyramid(data.frame(M=datos[,1], F=datos[,2], A=datos[,3]), 
        Clab="Edad", Llab="Hombres", Rlab="Mujeres", Cstep=10, AxisFM="d",
        main="Distribución de población por edad")

Nótese que la función pyramid() nos pide primero un data frame que incluya, en primer lugar, el eje izquierdo, luego el eje derecho y, finalmente, las edades del centro, mientras que resto de funciones son solo para títulos de los ejes y título general. La función Cstep sirve para saber en cuantas partes se va a dividir el rango que tenemos de edades.

Curva cumulativa

Estos gráficos son utilizados generalmente para saber que tan desigual puede ser una distribución de datos, siendo que el eje Y nos muestra la frecuencia acumulada de los datos (la cual podemos calcular con la función cumsum() si tenemos una base cualquiera con frecuencias) y el eje X nos muestra un ratio que va de 0 a 1, como es el caso de la función de densidad que vimos anteriormente.

En el portal de Tidyverse nos muestran un ejemplo de una curva acumulativa, haciendo uso de una distribución normal creada de manera aleatoria con la función rnorm(), usando la función stat_ecdf() para crear el gráfico respectivo. Se hizo de la siguiente manera:

# Datos de distribución ficticios
ac <- c(rnorm(100, 0, 3))
g <- gl(2, 100)
datos_fic <- data.frame(ac, g)

# Vemos los datos creados
head(datos_fic)
##          ac g
## 1  0.415885 1
## 2  3.194220 1
## 3  2.639286 1
## 4 -1.133012 1
## 5 -0.486785 1
## 6 -4.519169 1
# Creamos el gráfico
datos_fic %>% 
  ggplot(aes(ac)) + stat_ecdf(geom = "step")

Aunque me parece algo simple, sería bueno agregarle formato, como un mejor fondo, título, nombre a ejes, etc.

# Gráfico con mejor formato
ggplot(datos_fic, aes(ac)) + stat_ecdf(geom = "step", col = "darkblue") + 
  theme_bw() + 
  theme(legend.position = "none") +
  theme(legend.title = element_blank()) +
  guides(col = guide_legend(nrow = 1, byrow = TRUE)) + 
  xlab("Frecuencia acumulada") + 
  ylab("Ratio")+ 
  theme(plot.title = element_text(size = 11, face = "bold", hjust = 0)) + 
  theme(plot.subtitle = element_text(size = 10, hjust = 0)) + 
  theme(plot.caption = element_text(size = 10, hjust = 0)) +
  theme(plot.margin = unit(c(1,1,1,1), "cm")) +
  labs(
    title = "Gráfico de curva cumulativa",
    subtitle = "Ejemplo de distribución normal",
    caption = "Fuente: Elaboración propia con ayuda de Tidyverse. \nNotas: Espero que te ayude este ejercicio."
  )

Gráfico beeswarm

Estos gráficos son muy parecidos a los histogramas de frecuencias usados al principio, solo que en lugar de utilizar barras lo que hacemos es graficar todos los puntos en una distribución de frecuencias. No suelen ser tan comunes ni utilizados para el análisis de datos, ya que son poco agradables a la vista, pero se pueden hacer con el paquete beeswarm. Usaré los datos de expectativa de vida que filtré al principio para este ejercicio a continuación:

beeswarm(calidad_1957$lifeExp, col = "darkblue", main = "Gráfico de frecuencia de puntos",
         cex=1.5, xlab = "Frecuencia", ylab = "Expectativa de vida")

Es un gráfico bastante simple, pero de igual manera puede ser ilustrativo cuando queremos ver más específicamente la distribución de nuestros datos.

¡Muchas gracias por leerme! Sería de gran ayuda que compartieras ésta entrada con quien creas que pueda serle útil. Me ayudaría mucho que ésto llegue a más gente.

Referencias

[1] Urias, H. Q., & Salvador, B. R. P. (2014). Estadística para ingeniería y ciencias. Grupo Editorial Patria.

[2] Wickham, H., Chang, W., & Wickham, M. H. (2016). Package ‘ggplot2’. Create Elegant Data Visualisations Using the Grammar of Graphics. Version, 2(1), 1-189.

[3] Kuan, J. (2015). Learning Highcharts 4. Packt Publishing Ltd.

[4] Bivand, R. S., Pebesma, E. J., Gomez-Rubio, V., & Pebesma, E. J. (2013). Applied spatial data analysis with R (Vol. 2). New York: Springer.

[5] Smith, A. (2018). Visual Vocabulary. 4 de marzo de 2020, de Financial Times. Sitio web: https://github.com/ft-interactive/chart-doctor/tree/master/visual-vocabulary

[6] Coene, J. (2017). Graphs. 26 de febrero de 2021, de echarts4r Sitio web: https://echarts4r.john-coene.com/articles/map.html

R-Conomics
Todos los derechos reservados