Hola mundo en esta entrada vamos a hacer una introducción a series de tiempo en R.
si recuerdas en otra entrada estuvimos hablando de modelos estadístico o problemáticas de machine learning.
Pronóstico de series de tiempo
Bien, hablar de series de tiempo es hablar de pronosticar. En este caso, o en la mayoría de los casos, el objetivo siempre es saber qué va a pasar en el futuro: estimar o predecir, pronosticar cantidades a futuro.
Ejemplo de esto son temperaturas, precios de acciones y, cómo vamos a ver en este ejemplo, la cantidad de bicicletas que estén en uso en una determinada hora.
Las series de tiempo entrarían dentro de lo que son un modelo estadístico que toma como primicia que el futuro se puede predecir a través del pasado. Un pronóstico de serie de tiempo lo que nos da es una idea de lo que pudiera pasar después si sigue pasando lo que ha venido ocurriendo, valga la redundancia.
Datos Ecobici
En este caso, vamos a hablar sobre el sistema de bicicletas públicas de la Ciudad de México. Esto, digamos que en varios puntos de la ciudad hay diferentes estaciones, en donde hay aparcadas o estacionadas bicicletas que los turistas y habitantes pueden usar luego de que se den de alta.
Lo que vamos a estudiar en este vídeo, es cuantos usuarios por hora existen. Entonces, para esto vamos a utilizar datos de febrero. Ahora, vamos a ver cómo obtener los datos.
Por ejemplo vamos a utilizar datos del lunes 24 de febrero, martes 25 y miércoles 26. En donde, vamos a agrupar por hora la cantidad de personas que están usando bicicletas.
Lo que vamos a hacer es pronosticar o predecir lo que pasaría el jueves y el viernes. Entonces, lo que está acá en negro serían los datos que tenemos, es decir el pasado. Por otra parte, lo azul sería el futuro.
Esto es, algo que todavía no hemos visto y queremos predecir en base al pasado.
Obtención de datos de Ecobici
Vamos a ir a la web oficial de ecobici e información de servicio -> https://www.ecobici.cdmx.gob.mx/es/informacion-del-servicio/open-data . Bajamos a descarga de archivos histórico.
Para este ejercicio vamos a seleccionar el año 2020 y el mes de febrero. Luego, descargar y guardamos el archivo. NOTA: el archivo puede ser pesado, si tarda un poco no te preocupes. Luego abrimos el CSV para examinarlo.
Donde nos va a interesar la Fecha y la Hora de retiro. Para esto vamos a hacer un análisis en dplyr, como hemos visto en entradas anteriores.
Procesamiento de datos
En este caso, nos interesan los viajes por hora. Vamos a trabajar con las fechas y luego vamos a tomar la última semana de febrero para hacer un conteo por hora. Vamos a ver cómo hacer esto en R.
Muy bien, lo primero acá es que vamos a utilizar los paquetes
library(tidyverse)
library(lubridate)
En este caso, lubridate que es un paquete especial para todo lo que tiene que ver con fechas en R viene muy bien para trabajar series de tiempo de forma conjunta.
Lectura del csv
Para esto, usamos la función file.choose() y buscamos la ruta donde está el archivo. Recuerda que en esta otra entrada vimos como cargar datos CSV en R.
# leer csv
file.choose()
raw_data <- read_csv("C:\\ruta_carpeta\\2020-02.csv")
Luego de esto deberíamos tener en R un CSV con la misma estructura que el excel de ecobici.
Procesamiento en dplyr
La idea ahora es poder pasar de datos que están a nivel de usuario. A datos que estén a nivel de horas. Para eso podemos utilziar este código de manipulación en dplyr
viajes_diarios <- raw_data %>%
mutate(fecha_hora = dmy_hms(paste(Fecha_Retiro, Hora_Retiro))) %>%
filter(fecha_hora >= as.Date('2020-02-24'),
fecha_hora <= as.Date('2020-02-27')) %>%
group_by(horas = floor_date(fecha_hora, unit = 'hour')) %>%
summarise(conteo = n())
Donde fíjate el uso de la función dmy_hms de lubridate para transformar la fecha y hora de retiro. El resto son filtros para la última semana de febrero y una agrupación por hora usando la función floor_date también de lubridate. Por último, nos quedamos con los conteos para modelarlos con la serie temporal.
Si todo ha salido bien, deberíamos poder ver este dataframe
Si observas, hay un salto en las horas de 00, 01, 02, a 05! Esto es porque no hay datos para las 03 y 04 horas. Para el análisis de serie de tiempo tenemos que corregir este error. Esto podemos hacerlo utilizando el siguiente código
# rellenando los ceros
horas_completas <- data.frame(
horas = seq(
floor_date(min(viajes_diarios$horas), unit = 'hour'),
floor_date(max(viajes_diarios$horas), unit = 'hour'),
by = 'hour'
))
Luego tenemso que hacer un left_join para cruzar las tablas. Esto es posible utilizando el código:
# left join con horas
viajes_hora <- horas_completas %>%
group_by(horas_readondeadas = floor_date(horas, unit = 'hour')) %>%
left_join(viajes_diarios) %>%
mutate(conteo = ifelse(is.na(conteo), 0, conteo))
Si todo ha ido bien, el dataframe debería verse ahora como
Ahora sí incluye los ceros para las 03 y 04 horas. Perfecto!
Gráfica inicial
Podemos hacer una gráfica con los datos que hemos obtenidos en ggplot2 utilizando el siguiente código.
# grafica de lineas inicial
ggplot(data = viajes_diarios,
aes(x = horas,
y = conteo)) +
geom_line()
Eso arroja la gráfica
Modelo Arima Estacional en R
Antes de continuar, vamos a ver algunos aspectos teorícos de este modelo de serie temporal.
En este caso, vamos a trabajar el modelo auto regresivo integrado de media móvil estacional o arima estacional o sarima. También se le conoce este es un modelo que tiene varias partes digamos unidas y así puede hacer una mejor predicción o pronóstico. Estas partes son
- AR: Autoregresivo – cómo una regresión sobre sí mismo
- I: Integrado – diferencias con respecto a valores posteriores
- MA – Media movil (Moving Average en inglés) – tiene que ver con los «errores» de la predicción
- S: Estacional – para fenómenos que se repiten.
Parámetros: (p,d,q) (p_s, d_s, q_s) [s]
Con esto en mente, volvamos a R a ver cómo se hace esto
Implementación en R
Para hacer series de tiempo en R. Necesitamos un tipo de objeto que se llama ts (por time series). Podemos crearlo utilizando el siguiente código.
# creando objeto ts para modelo
conteo_ts <- ts(viajes_hora$conteo,
start = 1,
frequency = 24)
Ahora para impementar el modelo utilizamos el paquete forecast. Lo cargamos con el código (previo instalación del paquete)
# haciendo modelo ARIMA
# install.packages('forecast')
library(forecast)
Y ahora vendría el ajuste del modelo. Esto podemos lograrlo con este comando
ajuste <- auto.arima(y = conteo_ts)
Y eso básicamente sería todo… Podemos hacer un summary del modelo con:
summary(ajuste)
Una explicación de este summary la puedes ver en el video asociado a esta entrada.
Ahora, para hacer las predicciones podemos utilizar este código
predicciones <- forecast(ajuste)
Y para finalizar podemos obtener la gráfica con el comando
autoplot(predicciones)
Resultando en esta gráfica
Espero esta entrada te ayude a realizar pronósticos de serie de tiempo en R!
Te mando un saludo donde quiera que estés