ANÁLISIS DE SERIES DE TIEMPO

Analizaremos una serie de tiempo de Google Trends, la cual nos dice la frecuencia con la que se buscó la palabra “pan” en México los últimos cinco años.

library(ggplot2)
library(ggthemes)
library(dplyr)
library(lubridate)
library(forecast)
library(janitor) 
library(readxl)
library(dynlm)
library(stats)

setwd(file.path("C:","Users","Nath","Documents","Semestre 2022-1","series_de_tiempo","ejercicios", "data"))

pan_data<-read.csv("C://Users//Nath//Documents//Semestre 2022-1//series_de_tiempo//ejercicios//data//pandata.csv", header = TRUE)%>%
  janitor::clean_names() 

pan.ts <- ts(pan_data$freq_busqueda, start = c(2016,10), frequency = 48)
plot(pan.ts)

Se observa que la serie presenta tendencia aparentemente lineal y estacionalidad. Ahora analicemos su descomposición para ver si la serie puede representarse en su forma simple \(X_t = T_t + E_t + I_t +\epsilon_t\)

pan.ts.des <-  decompose(pan.ts)
plot(pan.ts.des)

La tendencia puede verse creciente y decreciente entre el 2018 y pasado 2019.

shapiro.test(x=as.integer(pan.ts.des$random))

## 
##  Shapiro-Wilk normality test
## 
## data:  as.integer(pan.ts.des$random)
## W = 0.92507, p-value = 6.041e-09

Hecha la prueba de Shapiro, se rechaza que los errores se distribuyen normales. Por lo tanto la serie no puede modelarse como \(X_t = T_t + E_t + I_t +\epsilon_t\).

MODELO AR(p)

Para hacer un modelo autorregresivo es necesario que la serie sea estacionaria, es decir, que su media y varianza sean constantes. Para saber si nuestra serie es estacionaria, nos ayudaremos de las funciones de autocorrelación:

acf(pan.ts)

pacf(pan.ts)

El primer gráfico no luce como si se tratase de una serie estacionaria, por lo tanto tendremos que usar transformaciones que le obliguen a serlo. Nos ayudaremos de una función que nos dirá cuántos rezagos debemos hacer para que, al aplicar diferencias a la serie, tenga media constante.

ndiffs(pan.ts)

## [1] 1

La función nos dice que basta con un rezago:

pan.ts.dif <- diff(pan.ts, 1)
plot(pan.ts.dif)

Podemos observar que la serie ya está centrada en cero, es decir, su media ya parece ser constante. Por lo tanto, ya podemos aplicar el modelo autorregresivo.

MODELO 1: AR(1)

ar_1 <- ar(pan.ts.dif, order.max = 1)
ar_1

## 
## Call:
## ar(x = pan.ts.dif, order.max = 1)
## 
## Coefficients:
##       1  
## -0.1701  
## 
## Order selected 1  sigma^2 estimated as  39.15

plot(forecast(ar_1, 48))

Desde aquí podemos ver que el modelo AR(1) no ajusta nada a la serie de diferencias, pero veamos cómo se ve en la serie original:

res <- forecast(ar_1, 48) #Predicción de un año
forecasted <- as.vector(res$mean)
base <- pan_data$freq_busqueda[261]
pan_data$freq_busqueda[261]+forecasted[1]

## [1] 51.04307

undiff <- c()

for(i in 1:48){
  val <- base+forecasted[i]
  undiff[i] <- val
  base <- val
}

serie_completa <- c(pan_data$freq_busqueda, undiff)

plot(serie_completa, type = "l")

Hagamos uso del gráfico de autocorrelación para intentar determinar el parámetro p que ajustará mejor nuestro modelo:

pacf(pan.ts.dif)

Intentemos con \(p=3\)

MODELO 2: AR(3)

ar_3 <- ar(pan.ts.dif, order.max = 3)
ar_3

## 
## Call:
## ar(x = pan.ts.dif, order.max = 3)
## 
## Coefficients:
##       1        2        3  
## -0.2143  -0.1713  -0.1346  
## 
## Order selected 3  sigma^2 estimated as  37.92

plot(forecast(ar_3, 48))

De igual manera se ve que no ajusta nada de nada…

MODELO 3: AR(9)

ar_9 <- ar(pan.ts.dif, order.max = 9)
ar_9

## 
## Call:
## ar(x = pan.ts.dif, order.max = 9)
## 
## Coefficients:
##       1        2        3        4        5        6        7        8  
## -0.2549  -0.2209  -0.2045  -0.0956  -0.1069  -0.2001  -0.1492  -0.1063  
##       9  
## -0.1532  
## 
## Order selected 9  sigma^2 estimated as  36.58

plot(forecast(ar_9, 48))

A pesar de que la predicción sigue siendo malísima, vemos que hay un poquito más de movimiento que en las series anteriores. Por curiosidad veamo cómo se ve en su serie original:

res <- forecast(ar_9, 48) #Predicción de un año
forecasted <- as.vector(res$mean)
base <- pan_data$freq_busqueda[261]
pan_data$freq_busqueda[261]+forecasted[1]

## [1] 49.80237

undiff <- c()

for(i in 1:48){
  val <- base+forecasted[i]
  undiff[i] <- val
  base <- val
}

serie_completa <- c(pan_data$freq_busqueda, undiff)

plot(serie_completa, type = "l")

Este mal ajuste del modelo se debe a que, a pesar de haber centrado la serie, ésta no es estrictamente estacionaria como lo pide el proceso AR().

HANDS-ON 1

ANÁLISIS DE SERIES DE TIEMPO

MODELO AR(p)