Importando datos de precios alrededor del mundo

Cada vez que voy al supermercado, mi cartera se queja. Pero, ¿qué tan costosa es la comida alrededor del mundo? En este proyecto exploraremos series de tiempo de los precios de la comida en Ruanda, de United Nations Humanitarian Data Exchange Global Food Price Database. La agricultura representa un 30% de la economía de Rwanda, y alrededor de 60% de sus ingresos por exportaciones (CIA World Factbook), por lo que el precio de la comida es muy importante para el estilo de vida de muchos ruandeses. En este proyecto manipularemos, visualizaremos, y haremos proyecciones del costo de la papa en Ruanda. Igualmente crearemos funciones para optimizar este proceso y hacer análisis de otros rubros.

library(readr)
library(dplyr)
potato_prices <- read_csv2('potatoes.csv')
glimpse(potato_prices)
## Observations: 4,739
## Variables: 18
## $ adm0_id            <int> 205, 205, 205, 205, 205, 205, 205, 205, 205...
## $ adm0_name          <chr> "Rwanda", "Rwanda", "Rwanda", "Rwanda", "Rw...
## $ adm1_id            <int> 2587, 2587, 2587, 2587, 2587, 2587, 2587, 2...
## $ adm1_name          <chr> "Kigali City/Umujyi wa Kigali", "Kigali Cit...
## $ mkt_id             <int> 1054, 1054, 1054, 1054, 1054, 1054, 1054, 1...
## $ mkt_name           <chr> "Gahanga", "Gahanga", "Gahanga", "Gahanga",...
## $ cm_id              <int> 148, 148, 148, 148, 148, 148, 148, 148, 148...
## $ cm_name            <chr> "Potatoes (Irish) - Retail", "Potatoes (Iri...
## $ cur_id             <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0...
## $ cur_name           <chr> "RWF", "RWF", "RWF", "RWF", "RWF", "RWF", "...
## $ pt_id              <int> 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15,...
## $ pt_name            <chr> "Retail", "Retail", "Retail", "Retail", "Re...
## $ um_id              <int> 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5...
## $ um_name            <chr> "KG", "KG", "KG", "KG", "KG", "KG", "KG", "...
## $ mp_month           <int> 7, 8, 9, 10, 11, 12, 1, 2, 3, 4, 5, 6, 7, 8...
## $ mp_year            <int> 2010, 2010, 2010, 2010, 2010, 2010, 2011, 2...
## $ mp_price           <dbl> 136.2500, 178.3333, 213.7500, 125.0000, 135...
## $ mp_commoditysource <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,...
potato_prices <- read_csv2("potatoes.csv", 
                          col_types=cols_only(adm1_name='c', 
                                              mkt_name='c', 
                                              cm_name='c', 
                                              mp_month='i', 
                                              mp_year='i', 
                                              mp_price='d'))

Limpiando nuestra data

Muchas de las columnas en la data de “potato” no son útiles para nuestro análisis. Por ejemplo, la columna adm1_name siempre es “Rwanda”, y cur_name siempre es “RWF” (abreviatura para Franco Ruandés, para entrar en contexto, 1000 RFW son alrededor de 1 USD). Similarmente, no necesitamos ninguna de las columnas de ID.

Asimismo, las columnas que necesitamos tienen nombres extraños, por ejemplo, adm1_id no parece significar “Region”, y mkt_name no parece estar tan claro como “Market”. Uno de los errores en análisis de datos es no conocer el significado de la variable, así que nombrar una variable claramente es una buena solución para evitar esto.

potato_prices_renamed <- plyr::rename(potato_prices, c('adm1_name'='region',
                                                       'mkt_name'='market',
                                                       'cm_name'='commodity_kg',
                                                       'mp_month'='month',
                                                       'mp_year'='year',
                                                       'mp_price'='price_rwf'))
head(potato_prices_renamed)
## # A tibble: 6 x 6
##   region                market  commodity_kg         month  year price_rwf
##   <chr>                 <chr>   <chr>                <int> <int>     <dbl>
## 1 Kigali City/Umujyi w~ Gahanga Potatoes (Irish) - ~     7  2010      136.
## 2 Kigali City/Umujyi w~ Gahanga Potatoes (Irish) - ~     8  2010      178.
## 3 Kigali City/Umujyi w~ Gahanga Potatoes (Irish) - ~     9  2010      214.
## 4 Kigali City/Umujyi w~ Gahanga Potatoes (Irish) - ~    10  2010      125 
## 5 Kigali City/Umujyi w~ Gahanga Potatoes (Irish) - ~    11  2010      135 
## 6 Kigali City/Umujyi w~ Gahanga Potatoes (Irish) - ~    12  2010      111.

Tratando con fechas

Es muy común en el análisis de datos, que la data no esté en la forma que nosotros queremos. Por ejemplo, pudimos notar que las variables de mes y año están dadas como integers. Como realizaremos análisis de series de tiempo, debemos crear una variable que contenga esas fechas en tal formato.

library(lubridate)

potato_prices_cleaned <- potato_prices_renamed %>% 
  mutate(date=ymd(paste(year, month, '01'))) %>% 
  select(-year, -month)

head(potato_prices_cleaned)
## # A tibble: 6 x 5
##   region                 market commodity_kg          price_rwf date      
##   <chr>                  <chr>  <chr>                     <dbl> <date>    
## 1 Kigali City/Umujyi wa~ Gahan~ Potatoes (Irish) - R~      136. 2010-07-01
## 2 Kigali City/Umujyi wa~ Gahan~ Potatoes (Irish) - R~      178. 2010-08-01
## 3 Kigali City/Umujyi wa~ Gahan~ Potatoes (Irish) - R~      214. 2010-09-01
## 4 Kigali City/Umujyi wa~ Gahan~ Potatoes (Irish) - R~      125  2010-10-01
## 5 Kigali City/Umujyi wa~ Gahan~ Potatoes (Irish) - R~      135  2010-11-01
## 6 Kigali City/Umujyi wa~ Gahan~ Potatoes (Irish) - R~      111. 2010-12-01

¿Repetir el proceso?

Las papas no son el único alimento que consumen los ruandeses, por más versátiles que sean. Los ruandeses tienen más variedad culinaria, esto significa que necesitaremos importar datos de diferentes rubros. Si queremos hacer una misma tarea en R, podríamos copiar y pegar nuestro código y cambiar ciertas cosas. Esto es una idea terrible, ya que podría afectarnos en algún momento por distintas razones. Una mejor idea es escribir una función, así evitaríamos copiar y pegar tantos códigos y tener un código más leíble.

read_price_data <- function(commodity){
  file_name = paste0(commodity, '.csv')
  
  prices <- read_csv2(file_name,  col_types = cols_only(adm1_name = col_character(),
                                                       mkt_name = col_character(),
                                                       cm_name = col_character(),
                                                       mp_month = col_integer(),
                                                       mp_year = col_integer(),
                                                       mp_price = col_double()))
  prices_renamed <- prices %>% rename(region = adm1_name, 
                                      market = mkt_name,
                                      commodity_kg = cm_name,
                                      month = mp_month,
                                      year = mp_year,
                                      price_rwf = mp_price)
  prices_cleaned <- prices_renamed %>% 
    mutate(date = ymd(paste(year, month, "01"))) %>% 
    select(-month, -year)
}

pea_prices <- read_price_data("peas")
glimpse(pea_prices)
## Observations: 1,958
## Variables: 5
## $ region       <chr> "Kigali City/Umujyi wa Kigali", "Kigali City/Umuj...
## $ market       <chr> "Gahanga", "Gahanga", "Gahanga", "Gahanga", "Gaha...
## $ commodity_kg <chr> "Peas (fresh) - Retail", "Peas (fresh) - Retail",...
## $ price_rwf    <dbl> 600.0000, 800.0000, 600.0000, 450.0000, 700.0000,...
## $ date         <date> 2011-01-01, 2011-04-01, 2011-05-01, 2011-12-01, ...

Visualizando el precio de las papas

El primer paso para analizar unos datos es visualizándola. En este caso, tenemos algunos precios y algunas fechas, así que grafiquemos cómo han variado esos precios en el tiempo.

library(ggplot2)
ggplot(potato_prices_cleaned, aes(x=date, y=price_rwf, group=market)) +
    geom_line(alpha=0.2) +
    labs(title = "Precio de la papa en el tiempo")

Parece haber una tendencia en los precios de las papas, incrementándose hasta 2013, y después disminuye la variación. Observamos también la estacionalidad: los precios son más bajos en diciembre y enero, y tienen un pico en agosto. En algunos años muestra picos en abril y mayo.

Más funciones…

Al igual que importando y limpiando la data, si queremos hacer muchos gráficos similares, necesitamos crear una función que lo haga por nosotros.

plot_price_vs_time <- function(prices, commodity){
  prices %>% 
    ggplot(aes(date, price_rwf, group = market)) +
    geom_line(alpha = 0.2) +
    ggtitle(paste0("Precio de ", commodity, " en el tiempo"))
}

plot_price_vs_time(pea_prices, 'arvejas')

Preparándonos para predecir el futuro (parte 1)

Aunque es útil observar cómo los precios han variado en el pasado, es más excitante predecir cómo cambiarán en el futuro. Antes de hacer eso, debemos preparar nuestra data.

Las bases de datos de cada rubro son variadas: en vez de ser una sola serie de tiempo, consisten en una serie de tiempo para cada mercado. El modo “fancy” de analizar esto es tratarla como una sola serie de tiempo jerárquica. La forma fácil que intentaremos aca, es tomar el costo medio alrededor de los mercados durante cada tiempo y analizar la serie de tiempo resultante.

Observando las gráficas de las papas y arvejas, podemos ver que ocasionalmente hay un gran pico en los precios. Probablemente se deba a algún problema logistico donde la comida no estaba disponible fácilmente en un mercado en particular, o el comprador parecía un turista y fue “timado”. La consecuencia de estas dispersiones es que es una mala idea usar el costo promedio de cada punto en el tiempo, en cambio, utilizar la mediana tiene más sentido debido a que es más robusta contra los puntos de dispersión.

potato_prices_summarized <- potato_prices_cleaned %>% 
  group_by(date) %>% 
  summarize(median_price_rwf=median(price_rwf))
potato_prices_summarized
## # A tibble: 122 x 2
##    date       median_price_rwf
##    <date>                <dbl>
##  1 2008-01-01             98.8
##  2 2008-02-01            100  
##  3 2008-03-01             95  
##  4 2008-04-01             98.1
##  5 2008-05-01             95.6
##  6 2008-06-01            109. 
##  7 2008-07-01            116. 
##  8 2008-08-01            125  
##  9 2008-09-01            138. 
## 10 2008-10-01            131. 
## # ... with 112 more rows

Preparándonos para predecir el futuro (parte 2)

Ya tenemos nuestro reusmen de precios resumidos a través de la mediana de cada mercado en el tiempo. Ahora debemos convertirlo en un objeto de series de tiempo.

library(magrittr)

max_date <- max(potato_prices_summarized$date)
min_date <- min(potato_prices_summarized$date)

potato_time_series <- ts(potato_prices_summarized$median_price_rwf, 
                         end=c(year(max_date), month(max_date)),
                         start=c(year(min_date), month(min_date)),
                         frequency=12)

potato_time_series
##           Jan      Feb      Mar      Apr      May      Jun      Jul
## 2008  98.7500 100.0000  95.0000  98.1250  95.6250 109.3750 115.8333
## 2009 120.0000 122.5000 130.0000 131.2500 135.0000 124.3125 125.8333
## 2010 109.6875 113.5000 131.2500 132.0833 140.4167 147.3750 141.2500
## 2011 105.7000 108.1750 118.8750 145.0143 148.6667 148.0500 137.4048
## 2012 150.5625 175.8750 186.7778 188.1250 182.9643 162.7500 180.0000
## 2013 154.3333 157.0000 171.2500 187.5000 177.0000 202.2500 210.0000
## 2014 138.3333 158.7500 186.2500 198.2500 191.0000 189.3333 182.5000
## 2015 136.2500 157.6071 178.0000 190.2778 180.0000 168.3333 180.0000
## 2016 300.0000 300.0000 285.0000 260.0000 207.5000 200.0000 220.0000
## 2017 250.0000 270.0000 250.0000 200.0000 150.0000 160.0000 190.0000
## 2018 280.0000 250.0000  98.7500 100.0000  95.0000  98.1250  95.6250
##           Aug      Sep      Oct      Nov      Dec
## 2008 125.0000 137.5000 131.2500 127.5000 113.7500
## 2009 144.2500 181.2500 170.0000 150.2500 112.0000
## 2010 159.7500 181.8000 162.5000 151.5000 122.5000
## 2011 137.2619 141.6667 144.2000 133.1750 141.5000
## 2012 200.0000 230.0000 213.5000 169.2500 144.0000
## 2013 233.1875 241.3333 237.5000 176.7083 140.0000
## 2014 187.6191 200.0000 183.1309 150.0000 133.9286
## 2015 202.1250 223.5000 217.5000 216.1250 190.0000
## 2016 250.0000 265.0000 250.0000 240.0000 250.0000
## 2017 200.0000 230.0000 230.0000 230.0000 250.0000
## 2018 109.3750 115.8333

Y más funciones…

Hicimos una preparación de data extensa. Sería una buena idea no tener que volverlo a hacer, ¿verdad?. Creemos una función que haga la preparación por nosotros.

create_price_time_series <- function(prices){
  prices_summarized <- prices %>%
    group_by(date) %>% 
    summarize(median_price_rwf = median(price_rwf))
  
  time_series <- prices_summarized %$% 
    ts(median_price_rwf, 
       start = c(year(min(date)), month(min(date))), 
       end   = c(year(max(date)), month(max(date))), 
       frequency = 12)
}

## Utilizandola en la data de arvejas
pea_time_series <- create_price_time_series(pea_prices)
pea_time_series
##            Jan       Feb       Mar       Apr       May       Jun       Jul
## 2011  536.1667  700.0000  950.0000  720.0000  591.5000  597.8572  666.3572
## 2012  655.0000  950.0000 1272.1667 1132.0000  920.8750  819.4167  714.2857
## 2013  633.3333  781.6334  829.9875  975.0000  904.1250  789.9444  806.8000
## 2014  695.5000 1025.0000 1166.6250 1066.5000  825.0000  808.3333  809.5714
## 2015  800.0000 1035.7619 1100.0000 1051.8889  950.0000  866.6667  800.0000
##            Aug       Sep       Oct       Nov       Dec
## 2011  758.5000  938.8333 1506.2500  775.0000  546.1042
## 2012  788.1250  990.7222 1413.7500  952.8333  650.0000
## 2013 1000.0000 1158.2500 1316.7500  914.0000  623.8571
## 2014 1000.0000 1000.0000 1666.6667  700.0000  633.3333
## 2015  900.0000 1166.6667 1550.0000 1062.2500  804.2500

El precio futuro de las papas

Toda la preparación está lista para empezar a proyectar. La pregunta más frecuente a la hora de hacer proyecciones es “¿cómo se que puedo confiar en una proyección?”. Recordemos que la data de las papas y de las arvejas tenían una fuerte estacionalidad. Para datos agrícolas, una buena proyección debería mostrar una forma similar durante las estaciones.

library(forecast)
potato_price_forecast <- forecast(potato_time_series)
potato_price_forecast
##          Point Forecast    Lo 80    Hi 80     Lo 95    Hi 95
## Oct 2018      110.36168 94.21006 126.5133 85.659913 135.0634
## Nov 2018      101.48131 80.63033 122.3323 69.592491 133.3701
## Dec 2018       93.84696 70.38658 117.3073 57.967408 129.7265
## Jan 2019      102.04683 72.77320 131.3205 57.276680 146.8170
## Feb 2019      106.16736 72.31819 140.0165 54.399514 157.9352
## Mar 2019      109.60651 71.54139 147.6716 51.390939 167.8221
## Apr 2019      110.86834 69.50418 152.2325 47.607318 174.1294
## May 2019      106.34779 64.14936 148.5462 41.810856 170.8847
## Jun 2019      113.44907 65.93769 160.9604 40.786688 186.1115
## Jul 2019      121.22081 67.96123 174.4804 39.767309 202.6743
## Aug 2019      131.78452 71.33196 192.2371 39.330308 224.2387
## Sep 2019      144.64439 75.64204 213.6467 39.114409 250.1744
## Oct 2019      137.25523 69.38006 205.1304 33.449113 241.0614
## Nov 2019      125.71871 61.45766 189.9798 27.439912 223.9975
## Dec 2019      115.82366 54.77522 176.8721 22.458123 209.1892
## Jan 2020      125.48639 57.42399 193.5488 21.393930 229.5788
## Feb 2020      130.09542 57.61439 202.5765 19.245243 240.9456
## Mar 2020      133.85433 57.37229 210.3364 16.885145 250.8235
## Apr 2020      134.95139 55.98163 213.9212 14.177558 255.7252
## May 2020      129.03820 51.80279 206.2736 10.916829 247.1596
## Jun 2020      137.23181 53.30853 221.1551  8.882234 265.5814
## Jul 2020      146.19652 54.94138 237.4517  6.633831 285.7592
## Aug 2020      158.47847 57.60240 259.3545  4.201829 312.7551
## Sep 2020      173.45692 60.95799 285.9558  1.404658 345.5092
autoplot(potato_price_forecast, main='Proyección del precio de la papa')

La función final

La proyección parece mostrar la misma estacionalidad que poseía el precio de la papa en el pasado, por lo que podemos confiar en ella. Si queremos hacer más proyecciones de los rubros, nos vendría bien escribir una última función.

plot_price_forecast <- function(time_series, commodity){
  
  price_forecast <- forecast(time_series)
  autoplot(price_forecast, 
           main = paste0("Proyeccion del precio de ", commodity))
}

##Probamos la función en la data de arvejas
plot_price_forecast(pea_time_series, 'arvejas')

¿Hacerlo todo de nuevo?

Se realizó un esfuerzo en escribir todos estos códigos para analizar el precio de la papa en Ruanda. Al tener nuestras funciones creadas, podemos fácilmente hacer el análisis en otro rubro.

commodity <- "granos"
bean_prices <- read_price_data('beans')
plot_price_vs_time(bean_prices, 'granos')

bean_time_series <- create_price_time_series(bean_prices)
plot_price_forecast(bean_time_series, 'granos')