Econometría II: SARIMA

Departamento de Economía

Carlos A. Yanes G.

2024-04-02

Paquetes con que se trabaja la sesión

Los paquetes que vamos a utilizar en la sesión de hoy son:

Note

No olvide que debe desde luego debe tener presente en su documento de R Markdown los siguientes paquetes:

library(pacman)
p_load(readxl, TSstudio, tidyverse, stats, urca, forecast, ggfortify, ggplot2, tseries, fpp2))

Preambulo

SARIMA

  • Nos acercamos a modelos un poco mas superiores que los anteriormente vistos.

  • Llega la hora de trabajar con modelos estacionales, estos desde luego intentaran mostrar la referencia estacional en las series de tiempo.

  • Se presenta la opción de Pronostico Rápido. Se denomina auto.arima

Datos

Cartera comercial de bancos

Producto de industria Manufacturera Perú (op)

Proceso 1

Nuevamente Cartera comercial

  • La vez pasada supimos que no era una serie Ruido Blanco

  • Eso nos afecta el Pronostico de esa serie. Tal vez podamos mirar otro modelo. Puede ser estacional, como de otro orden. Eso vamos hacer hoy!!.

Cargamos datos

Code
bd <- read_excel("cartera.xlsx")
bd <- bd |> select(Cartera)
cartera <- ts(bd, frequency=12, start=c(2016,1))
Code
bd2 <- read_excel("Pind01.xlsx")
bd2 <- bd2 |> select(pind)
pind<- ts(bd2, frequency=12, start=c(1992,1))

Code
cartera %>% diff(lag=12) %>% diff() %>%
  ggtsdisplay()

Code
cartera %>% diff() %>%
  ggtsdisplay()

  • Este modelo ya diferenciado estacionalmente y -resumido (vamos paso adelante)- nos dice que el ACF cae lentamente y tenemos dos picos uno de AR y uno de MA en el lag 12 de la función de PACF. Esta parte nos indica en conformidad un ARIMA (1,1,1) para ambos correlogramas. Sin embargo los candidatos pueden ser los siguientes:

  • Los modelos que pueden ser candidatos son:

\[\begin{equation} (1,1,0)(0,1,1)_{12} \end{equation} \]

\[ (1,1,0)(1,1,1)_{12} \]

\[ (1,1,1)(1,1,1)_{12} \] (opcional)

Test de Raices Unitarias

Code
dif_cartera_12<- diff(cartera, lag = 12)
dcartera12<-diff(dif_cartera_12, 1)
dfuller<-ur.df(dcartera12, lags=0, type='trend')
summary(dfuller)

############################################### 
# Augmented Dickey-Fuller Test Unit Root Test # 
############################################### 

Test regression trend 


Call:
lm(formula = z.diff ~ z.lag.1 + 1 + tt)

Residuals:
     Min       1Q   Median       3Q      Max 
-14218.7  -2125.5    252.6   2011.5  13345.7 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 499.78587  870.03142   0.574    0.567    
z.lag.1      -0.35689    0.08565  -4.167 7.76e-05 ***
tt          -12.49149   18.08548  -0.691    0.492    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 3891 on 80 degrees of freedom
Multiple R-squared:  0.1783,    Adjusted R-squared:  0.1578 
F-statistic: 8.681 on 2 and 80 DF,  p-value: 0.0003872


Value of test-statistic is: -4.1668 5.7879 8.6812 

Critical values for test statistics: 
      1pct  5pct 10pct
tau3 -4.04 -3.45 -3.15
phi2  6.50  4.88  4.16
phi3  8.73  6.49  5.47
Code
testkp<-ur.kpss(dcartera12, type = "tau", lags = "short")
summary(testkp)

####################### 
# KPSS Unit Root Test # 
####################### 

Test is of type: tau with 3 lags. 

Value of test-statistic is: 0.1413 

Critical value for a significance level of: 
                10pct  5pct 2.5pct  1pct
critical values 0.119 0.146  0.176 0.216

Code
sarima_1<-Arima(cartera, order=c(1,1,0), seasonal=c(0,1,1))
checkresiduals(sarima_1)


    Ljung-Box test

data:  Residuals from ARIMA(1,1,0)(0,1,1)[12]
Q* = 33.345, df = 17, p-value = 0.01019

Model df: 2.   Total lags used: 19
Code
sarima_2<-Arima(cartera, order=c(1,1,0), seasonal=c(1,1,1))
checkresiduals(sarima_2)


    Ljung-Box test

data:  Residuals from ARIMA(1,1,0)(1,1,1)[12]
Q* = 33.492, df = 16, p-value = 0.006356

Model df: 3.   Total lags used: 19
Code
sarima_3<-Arima(cartera, order=c(1,1,1), seasonal=c(1,1,1))
checkresiduals(sarima_3)


    Ljung-Box test

data:  Residuals from ARIMA(1,1,1)(1,1,1)[12]
Q* = 20.836, df = 15, p-value = 0.1422

Model df: 4.   Total lags used: 19
  • El mejor modelo por lo pronto es el (opcional) o sarima de orden (1,1,1) en ambas partes.

  • Tenemos que igual revisar cada uno en test de ruido blanco.

Test de Ruido para cada modelo

Code
checkresiduals(sarima_1, plot=FALSE)

    Ljung-Box test

data:  Residuals from ARIMA(1,1,0)(0,1,1)[12]
Q* = 33.345, df = 17, p-value = 0.01019

Model df: 2.   Total lags used: 19

    Ljung-Box test

data:  Residuals from ARIMA(1,1,0)(1,1,1)[12]
Q* = 33.492, df = 16, p-value = 0.006356

Model df: 3.   Total lags used: 19

    Ljung-Box test

data:  Residuals from ARIMA(1,1,1)(1,1,1)[12]
Q* = 20.836, df = 15, p-value = 0.1422

Model df: 4.   Total lags used: 19

Selección final

Usamos el criterio AKAIKE

Code
aicc <- c(
  Arima(cartera, order=c(1,1,0), seasonal=c(0,1,1))$aicc,
  Arima(cartera, order=c(1,1,0), seasonal=c(1,1,1))$aicc,
  Arima(cartera, order=c(1,1,0), seasonal=c(1,1,0))$aicc,
  Arima(cartera, order=c(1,1,1), seasonal=c(1,1,1))$aicc
  )
  • Primer modelo SARIMA (1,1,0)(0,1,1): 1594.73.
  • Segundo modelo SARIMA (1,1,0)(1,1,1): 1596.92.
  • Tercer modelo SARIMA (1,1,0)(1,1,0): 1608.75.
  • Cuarto modelo SARIMA (1,1,1)(1,1,1): 1592.52.

Tenemos modelo 🟢

Resultado

Code
autoplot(forecast(sarima_3, h=12))

Code
summary(sarima_3)
Series: cartera 
ARIMA(1,1,1)(1,1,1)[12] 

Coefficients:
         ar1      ma1    sar1     sma1
      0.8476  -0.3788  0.0035  -0.9999
s.e.  0.0821   0.1295  0.1222   0.1677

sigma^2 = 6876625:  log likelihood = -790.87
AIC=1591.75   AICc=1592.52   BIC=1603.9

Training set error measures:
                   ME     RMSE      MAE        MPE      MAPE       MASE
Training set 162.8221 2381.479 1591.688 0.03376412 0.3177587 0.04273688
                    ACF1
Training set -0.03705195
Code
pb<-forecast(sarima_3, h=12)
pb
         Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
Feb 2024       650712.2 647139.0 654285.4 645247.5 656177.0
Mar 2024       654306.3 647957.6 660654.9 644596.9 664015.6
Apr 2024       657170.1 647968.6 666371.6 643097.7 671242.6
May 2024       660444.0 648345.5 672542.5 641941.0 678947.1
Jun 2024       663814.4 648815.1 678813.7 640874.9 686753.8
Jul 2024       665011.7 647138.1 682885.4 637676.3 692347.1
Aug 2024       666769.8 646069.3 687470.3 635111.1 698428.5
Sep 2024       669826.0 646360.3 693291.7 633938.2 705713.7
Oct 2024       672219.8 646059.7 698379.9 632211.4 712228.3
Nov 2024       676258.5 647480.6 705036.3 632246.6 720270.4
Dec 2024       677967.3 646651.9 709282.7 630074.5 725860.1
Jan 2025       675372.5 641601.7 709143.3 623724.5 727020.5

Auto arima

Miremos como trabaja la IA

Code
auto.arima(cartera)
Series: cartera 
ARIMA(2,1,0)(2,0,0)[12] 

Coefficients:
         ar1     ar2    sar1    sar2
      0.4441  0.3050  0.3046  0.2771
s.e.  0.0969  0.0972  0.0927  0.1010

sigma^2 = 8655529:  log likelihood = -903.43
AIC=1816.85   AICc=1817.52   BIC=1829.67

Sin pasos

Code
auto.arima(cartera, 
  stepwise=FALSE, approximation=FALSE)
Series: cartera 
ARIMA(2,1,0)(2,0,0)[12] 

Coefficients:
         ar1     ar2    sar1    sar2
      0.4441  0.3050  0.3046  0.2771
s.e.  0.0969  0.0972  0.0927  0.1010

sigma^2 = 8655529:  log likelihood = -903.43
AIC=1816.85   AICc=1817.52   BIC=1829.67

Gracias por su atención!!