rm(list=ls())
cat("\014")
library(readxl)
library(ggplot2)
library(lubridate)
##
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, union
library(readr)
library(fitdistrplus)
## Loading required package: MASS
## Loading required package: survival
library(forecast)
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --
## v tibble 3.1.4 v dplyr 1.0.7
## v tidyr 1.1.4 v stringr 1.4.0
## v purrr 0.3.4 v forcats 0.5.1
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x lubridate::as.difftime() masks base::as.difftime()
## x lubridate::date() masks base::date()
## x dplyr::filter() masks stats::filter()
## x lubridate::intersect() masks base::intersect()
## x dplyr::lag() masks stats::lag()
## x dplyr::select() masks MASS::select()
## x lubridate::setdiff() masks base::setdiff()
## x lubridate::union() masks base::union()
library(tseries)
## Warning: package 'tseries' was built under R version 4.1.2
library(zoo)
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(smooth)
## Loading required package: greybox
## Package "greybox", v1.0.1 loaded.
##
## Attaching package: 'greybox'
## The following object is masked from 'package:tidyr':
##
## spread
## The following object is masked from 'package:forecast':
##
## forecast
## The following object is masked from 'package:lubridate':
##
## hm
## This is package "smooth", v3.1.3
Caso3 <- "C:\\Users\\user\\OneDrive - Pontificia Universidad Javeriana\\Documents\\Cata\\la javeriana\\ingenieria industrial\\9. Ingeneiria industrial noveno semestre\\Estocasticos\\Corte 3\\Caso 3\\Caso 3 - Estocásticos 2130 - Datos.xlsx"
D_CraftBeer <- read_excel(Caso3, sheet = "CraftBeer")
D_Whiskey <- read_excel(Caso3, sheet = "Whiskey")
S_CraftBeer <- ts(D_CraftBeer$`Relative Interest`, frequency = 12, start = c(2010,1))
S_Whiskey <- ts(D_Whiskey$`Relative Interest`, frequency = 12, start = c(2010,1))
#A=ajuste
#P=Prueba
SA_CraftBeer <- window(S_CraftBeer, start = c(2010, 1), end = c(2014, 12))
SP_CraftBeer <- window(S_CraftBeer, start = c(2015, 1))
SA_Whiskey <- window(S_Whiskey, start = c(2010, 1), end = c(2014, 12))
SP_Whiskey <- window(S_Whiskey, start = c(2015, 1))
autoplot(SA_CraftBeer)+ theme_minimal()+ ylab("Relative Interest")+ xlab("")+ ggtitle("Relative Interest mensual")
#Suavizacion
Media_Crafbeer<- rollmean(SA_CraftBeer, 3, fill = NA, align = "right")
autoplot(Media_Crafbeer)+ theme_minimal()+ ylab("Relative Interest")+ xlab("")+ ggtitle("Relative Interest mensual")
Con el fin de realizar un análisis descriptivo de la serie de tiempo se construyó el gráfico anterior, un gráfico de serie para tendencia de búsqueda de “CraftBeer”, en el gráfico de serie se puede evidenciar una fuerte tendencia positiva de las búsquedas de la palabra, también se puede ver un comportamiento estacional, ya que en cada año aproximadamente en el quinto mes del año que corresponde a mayo se evidencia un aumento significativo de búsquedas, de igual manera cada año se evidencia una caída de búsquedas en el onceavo mes que corresponde a noviembre.
El segndo gráfico anterior también es un gráfico de serie, sin embargo, este cuenta con una suavización, este se realizó con el motivo de clarificar la visualización de un comportamiento estacional, este definitivamente es más claro, ya que se muestra un aumento generalizado de búsquedas para los meses de mayo a julio.
ggseasonplot(SA_CraftBeer, year.labels = TRUE, year.labels.left = TRUE)+ ggtitle("Gráfico estacional de las relative interest")+ xlab("Mes")+ ylab("Relative Interest")+ theme_minimal()
Este gráfico de estacionalidad nos ayuda a identificar de manera más precisa el comportamiento de las búsquedas año a año, en este se puede ver de manera exacta que en todos los años hay un aumento en las búsquedas y sigue el mismo comportamiento, a aportar de febrero aumentan levemente las búsquedas, de abril a mayo se muestra un fuerte crecimiento en las búsquedas de CraftBeer, de mayo a junio las búsquedas caen hasta aproximadamente el mismo número de búsquedas de marzo, por último de octubre a noviembre se encuentra una nueva caída en las búsquedas.
Acf(SA_CraftBeer, main="Gráfico ACF")
Pacf(SA_CraftBeer, main="Gráfico PACF")
En el gráfico anterior se muestra el ACF, en este se puede evidenciar un decaimiento exponencial marcado, lo cual nos confirma que la serie de tiempo en cuestión es una serie con tendencia marcada.
Del gráfico anterior, el PACF se pueden identificar dos autocorrelaciones significativas, de la cual cabe resaltar que la primera en su gran magnitud indica nuevamente que la serie tiene tendencia, y también nos da a conocer que tiene autocorrelación positiva en la primera, si al mes pasado hubo un incremento en las búsquedas, en el siguiente mes se espera que también suceda, de igual manera con una disminución en las búsquedas, que para el siguiente mes también disminuirán.
yt_CraftBeer <- c(D_CraftBeer$`Relative Interest`)
yt_año_CraftBeer <- c()
for (i in 1:length(D_CraftBeer$`Relative Interest`)){
if (i<=12){
yt_año_CraftBeer[i]<- NA
}else{
yt_año_CraftBeer[i]<- yt_CraftBeer[i-12]
}
}
p_CraftBeer <- data.frame(yt_CraftBeer,yt_año_CraftBeer)
ggplot(p_CraftBeer, aes(yt_CraftBeer,yt_año_CraftBeer)) + geom_point()
## Warning: Removed 12 rows containing missing values (geom_point).
En este gráfico podemos evidenciar un comportamiento muy compacto entre el tiempo en años y las búsquedas de la palabra “CraftBeer”, de lo que podemos afirmar que la varianza de la serie no es constante, Homocedasticidad y con esto se incumple el requisito para ser considerada una serie estacionaria.
Serie_sin_estacionalidad_CraftBeer <- diff(SA_CraftBeer, lag=12, differences = 1)
Serie_sin_estacionalidad_tendencia_CraftBeer <- diff(Serie_sin_estacionalidad_CraftBeer, lag=1, differences = 1)
#Sin Tendencia
autoplot(Serie_sin_estacionalidad_tendencia_CraftBeer)+ theme_minimal()+ ylab("Relative Interest")+ xlab("")+ ggtitle("Relative Interest mensual")
#Sin Estacionalidad
ggseasonplot(Serie_sin_estacionalidad_tendencia_CraftBeer, year.labels = TRUE, year.labels.left = TRUE)+ ggtitle("Gráfico estacional de las relative interest")+ xlab("Mes")+ ylab("Relative Interest")+ theme_minimal()
#sin heterocedasticidad
yt_CraftBeer_Trans <- c(Serie_sin_estacionalidad_tendencia_CraftBeer)
yt_año_CraftBeer_Trans <- c()
for (i in 1:length(yt_CraftBeer_Trans)){
if (i<=12){
yt_año_CraftBeer_Trans[i]<- NA
}else{
yt_año_CraftBeer_Trans[i]<- yt_CraftBeer_Trans[i-12]
}
}
p_CraftBeer_Trans <- data.frame(yt_CraftBeer_Trans,yt_año_CraftBeer_Trans)
ggplot(p_CraftBeer_Trans, aes(yt_CraftBeer_Trans,yt_año_CraftBeer_Trans)) + geom_point()
## Warning: Removed 12 rows containing missing values (geom_point).
### Autocorrelacion
Acf(Serie_sin_estacionalidad_tendencia_CraftBeer, main="Gráfico ACF")
Pacf(Serie_sin_estacionalidad_tendencia_CraftBeer, main="Gráfico PACF")
Para obtener que las series de las búsquedas de “CraftBeer” y “Whiskey” sean estacionarias, se debe eliminar la tendencia, eliminar la estacionalidad y verificar que su varianza sea constante, para esto se procedió a realizar un proceso de diferenciación estacional, de s = 12, esto ya que se identifico en el punto anterior que en ambas ocasiones el periodo repite su comportamiento cada año, lo que equivale a un periodo de 12 meses, este proceso se realizó con una sola iteración, ya que el resultado ya cumplía con los requerimientos para ser considerada una serie estacionaria.
En el gráfico de líneas anterior perteneciente a la serie de búsquedas de “CraftBeer”, se puede evidenciar que el procedimiento realizado si cumplió con su cometido, la serie es estacionaria ya que no se evidencia tendencia ni comportamiento estacional, lo que permite la aplicación y desarrollo de modelos a esta serie diferenciada.
En el gráfico estacional se presenta un análisis estacional para los años pertenecientes a la muestra para el número de búsquedas de “CrafBeer”, en este podemos ver que para ninguno de los años se evidencia estacionalidad, ni comportamiento de tendencia positivo o negativo.
En el gráfico de dispersion podemos evidenciar un comportamiento disperso entre los años y las búsquedas de la palabra “CraftBeer”, heteroscedasticidad, de lo que podemos afirmar que la varianza de la serie es constante y con esto se cumple el requisito para ser considerada una serie estacionaria.
En este gráfico ACF podemos evidenciar un comportamiento disperso entre los años y las búsquedas de la palabra “CraftBeer”, heteroscedasticidad, de lo que podemos afirmar que la varianza de la serie es constante y con esto se cumple el requisito para ser considerada una serie estacionaria.
En el gráfico anterior se muestra el ACF, en este no se evidencia un decaimiento exponencial ni tendencia a cero, se observan dos correlaciones significativas, a partir de allí las correlaciones no tienen ningún patrón, lo que indica que ya la serie no tiene tendencia alguna.
Teniendo en cuenta los resultados en R, si se realiza la serie bajo la original, más no en la ajustada, se puede evidenciar que esta indica una serie no estacionaria presentando problemas de estacionalidad y tendencia, por lo tanto, podría modelarse como un proceso SARIMA más no ARIMA.
Así mismo, se pudo identificar que la serie contaba con indicios de estacionalidad, ya que esta tendía a repetirse cada 12 meses. Por esta razón, en el punto 2, se realizó una diferenciación (differences=d=1) con el fin de eliminar la estacionalidad contando con un Lag de 12, correspondiente a los 12 meses.
De tal manera, se puede referir a un posible modelo de ajuste para la parte estacional como:
Por otro lado, se logró identificar la tendencia presentada en la serie. Por esta razón, en el punto 2 se realizó una diferenciación (differences=d=1) con el fin de eliminar la tendencia presentada con un Lag 1.
De tal manera se puede referir a un posible modelo de ajuste para la tendencia como:
Como se pudo observar en el punto 2 hay autocorrelaciones significativas, en los correlogramas, tanto el ACF como el PACF tienen un comportamiento sinusoidal, por lo tanto, se podría hablar de un comportamiento semejante al ARMA, es decir, al modelo AR(p) se le incluye el modelo MA(q). Teniendo en cuenta lo anterior, se presentan 2 modelos de ajuste al comportamiento presente en la serie:
arima(SA_CraftBeer,order = c(1,1,0), seasonal = list(order= c(0,1,1), period(12)))
##
## Call:
## arima(x = SA_CraftBeer, order = c(1, 1, 0), seasonal = list(order = c(0, 1,
## 1), period(12)))
##
## Coefficients:
## ar1 sma1
## -0.3508 -0.5828
## s.e. 0.1347 0.2114
##
## sigma^2 estimated as 26.59: log likelihood = -146.3, aic = 298.61
Teniendo en cuenta esta alternativa, se obtuvo la siguiente ecuación:
Xt=yt-yt-12
Xt= 0 + (-0,5828)Yt-1+(-3,508)et-1+Et
arima(SA_CraftBeer,order = c(0,1,1), seasonal = list(order=c(1,1,0), period(12)))
##
## Call:
## arima(x = SA_CraftBeer, order = c(0, 1, 1), seasonal = list(order = c(1, 1,
## 0), period(12)))
##
## Coefficients:
## ma1 sar1
## -0.7329 -0.4767
## s.e. 0.0998 0.1517
##
## sigma^2 estimated as 22.69: log likelihood = -141.99, aic = 289.99
Teniendo en cuenta esta alternativa, se obtuvo la siguiente ecuación:
Xt=yt-yt-12
Xt= -0.7329et-1+(-0,4767)Xt-12
Modelo_Sarima <- auto.arima(SA_CraftBeer)
Modelo_Sarima
## Series: SA_CraftBeer
## ARIMA(0,1,1)(1,1,0)[12]
##
## Coefficients:
## ma1 sar1
## -0.7329 -0.4767
## s.e. 0.0998 0.1517
##
## sigma^2 estimated as 23.69: log likelihood=-141.99
## AIC=289.99 AICc=290.55 BIC=295.54
Teniendo en cuenta la herramienta estadística R y utilizando la función: auto.arima se obtiene cuál de los modelos planteados es el mejor para ajustarlo al modelo, en este caso, al programa estimó un modelo SARIMA con los siguientes parámetros
Teniendo en cuenta esta modelo, se obtuvo la siguiente ecuación:
Xt=yt-yt-12
Xt= -0.7329et-1+(-0,4767)Xt-12
Los supuestos del modelo SARIMA recaen sobre los residuos: deben formar un proceso de ruido blanco; es decir, deben ser estacionarios y no ser autocorrelacionados entre sí.
Residuales <- Modelo_Sarima$residuals
Acf(Residuales, main = "Gráfico ACF residuales SARIMA")
Pacf(Residuales, main = "Gráfico PACF residuales SARIMA")
Se demuestra un comportamiento de ruido blanco, ya que en ninguno de los graficos de los residuos (ACF y PACF) se tiene correlación significativa (esto basado en que ningún valor de correlación sobrepasa el umbral de decisión), lo que hace que la serie sea totalmente aleatoria y no se pueda modelar.
La prueba Dickey-Fuller H0: Los errores no son estacionarios H1: Los errores son estacionarios
adf.test(Residuales)
##
## Augmented Dickey-Fuller Test
##
## data: Residuales
## Dickey-Fuller = -3.9452, Lag order = 3, p-value = 0.01809
## alternative hypothesis: stationary
Prueba Ljung-Box H0: Los errores se distribuyen de forma independiente (las correlaciones son cero) H1: Los errores no se distribuyen de forma independiente (hay correlación)
Box.test(Residuales)
##
## Box-Pierce test
##
## data: Residuales
## X-squared = 0.17636, df = 1, p-value = 0.6745
La prueba Dickey-Fuller H0: Los errores no son estacionarios H1: Los errores son estacionarios
Teniendo en cuenta el p valor 0.018 y al ser menor que la significancia 0,05, es decir 0.018 < 0.05 se rechaza la hipótesis nula y se concluye que Los errores son estacionarios.
Teniendo en cuenta el p valor 0.6745 y al ser mayor que la significancia 0.05, es decir 0.6745 > 0.05 No se rechaza la hipótesis nula y se concluye que Los errores se distribuyen de forma independiente (las correlaciones son cero).
Es así que efectivamente los supuestos del modelo SARIMA, los cuales recaen sobre los residuos, si se cumplen al establecer el modelo obtenido automáticamente.
auto.arima(Residuales)
## Series: Residuales
## ARIMA(0,0,0) with non-zero mean
##
## Coefficients:
## mean
## 0.8268
## s.e. 0.5337
##
## sigma^2 estimated as 17.38: log likelihood=-170.29
## AIC=344.57 AICc=344.78 BIC=348.76
Para el Cafebar, en el punto 3 se propusieron dos modelos SARIMA (1,1,0)(0,1,1)[12]
y SARIMA (0,1,1)(1,1,0)[12] basados en que se pudo identificar una estacionalidad, ya que tiende a repetirse cada 12 meses, es por esto que en el segundo punto para quitar la estacionalidad se realizo una diferenciación (differences=d=1) para eliminar la estacionalidad con un Lag=12 correspondiente a 12 meses, por esto se puede referir a un posible modelo de ajuste para la parte estacional. Realizando el cálculo con el programa R, los resultados arrojan que el modelo para este caso en especifico se asemeja a un SARIMA(0,1,1)(1,1,0)[12] en este se encuentra que el modelo 2 propuesto con técnicas de visualización y el uso de un software automático coinciden. Demostrando asi que el modelo se asemeja efectivamente a un modelo que necesita de una diferenciación de 1 iteración, para eliminar la tendencia y que por otro lado se necesito de una diferenciación estacional para que la serie sea estacionaria.
Pron_Sarima <- forecast(SA_CraftBeer,12, model=Modelo_Sarima)
Pron_Sarima
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## Jan 2015 66.29507 60.05686 72.53328 56.75455 75.83559
## Feb 2015 74.15507 67.69812 80.61202 64.28002 84.03012
## Mar 2015 79.10840 72.43989 85.77692 68.90979 89.30702
## Apr 2015 87.34173 80.46816 94.21531 76.82951 97.85396
## May 2015 102.34173 95.26905 109.41442 91.52499 113.15848
## Jun 2015 81.63174 74.36539 88.89809 70.51881 92.74466
## Jul 2015 83.24840 75.79342 90.70338 71.84699 94.64981
## Aug 2015 83.81840 76.17944 91.45736 72.13562 95.50118
## Sep 2015 87.20174 79.38313 95.02034 75.24421 99.15926
## Oct 2015 84.24840 76.25419 92.24262 72.02230 96.47450
## Nov 2015 80.06174 71.89569 88.22779 67.57284 92.55064
## Dec 2015 85.58507 77.25072 93.91942 72.83879 98.33135
accuracy(Pron_Sarima,SP_CraftBeer)
## ME RMSE MAE MPE MAPE MASE
## Training set 0.826848 4.215553 2.845769 1.108172 7.400724 0.212107
## Test set -6.586458 9.689438 7.700590 -9.812508 11.162378 0.573957
## ACF1 Theil's U
## Training set 0.05421527 NA
## Test set 0.44414837 0.794519
autoplot(SA_CraftBeer)+ theme_minimal()+ ylab("haaaaaaaaaaaaaaaa")+ xlab("")+ ggtitle("Pronósticos modelo SARIMA")+ autolayer(Pron_Sarima, series = "SARIMA", PI = FALSE)
RMSE:
El error de pronóstico calculado por medio de un error cuadrático medio, arroja un valor de 4.2155 lo cual indica que existe un desface de 4.21 unidades de Cafebar con referencia al valor real de búsquedas para este producto. Este valor es de igual forma interpretable por exceso o por faltante indicando que el desfase con respecto al valor real puede ser en una sobreestimación como en subestimación pero dentro de este rango (+- 4.21). Dentro del contexto este valor es relativamente bajo, ya que una empresa productora de estos productos
MAE:
El error de pronóstico calculado por medio de un error cuadrático medio, arroja un valor de 6.228738 lo cual indica que existe un desface de 2.84 unidades de Cafebar con referencia al valor real de búsquedas para este producto. Este valor es de igual forma interpretable por exceso o por faltante indicando que el desfase con respecto al valor real puede ser en una sobreestimación como en subestimación pero dentro de este rango (+- 2.84).
MAPE:
El error de pronóstico tiene por medio del error absoluto porcentual promedio , tiene un valor de 7.4% lo que no es un valor que se pueda declarar mínimo. Estar cometiendo errores de tal magnitud puede generar perdidas más que considerables a nivel económico. En el caso de las búsquedas de Whisky se puede indicar una sobre o subestimación más que considerable para una posible perdida.
El método de cálculo del error de pronóstico que se recomienda como ideal es el MAE, ya que ofrece el error más pequeño, el cual evita caer en fallos en el pronóstico de búsquedas de Cafebar, tanto por encima como por debajo de las 2.84 unidades.
autoplot(SA_Whiskey)+ theme_minimal()+ ylab("Relative Interest")+ xlab("")+ ggtitle("Relative Interest mensual")
#Suavizacion
Media_Whiskey<- rollmean(SA_Whiskey, 3, fill = NA, align = "right")
autoplot(Media_Whiskey)+ theme_minimal()+ ylab("Relative Interest")+ xlab("")+ ggtitle("Relative Interest mensual")
Para realizar un análisis descriptivo de la serie de tiempo se muestra un gráfico de serie para tendencia de búsqueda de “Whiskey”, en el gráfico de serie se puede evidenciar una fuerte tendencia positiva de las búsquedas de la palabra, también se puede ver un comportamiento estacional, ya que en cada año aproximadamente en el sexto mes del año que corresponde a junio se evidencia una caída significativa de búsquedas, de igual manera cada año se evidencia una alza importante de búsquedas en el doceavo mes que corresponde a diciembre.
El gráfico anterior también es un gráfico de serie, sin embargo, este cuenta con una suavización, este se realizó con el motivo de evidenciar de un comportamiento estacional, en este se muestra un aumento generalizado de búsquedas para los meses de junio y julio.
ggseasonplot(SA_Whiskey, year.labels = TRUE, year.labels.left = TRUE)+ ggtitle("Gráfico estacional de las relative interest")+ xlab("Mes")+ ylab("Relative Interest")+ theme_minimal()
Este gráfico de estacionalidad nos permite identificar de manera más precisa el comportamiento de las búsquedas de los años evaluados, en este se puede ver de manera exacta que en todos los años hay un aumento en las búsquedas y sigue el mismo comportamiento, en este se evidencia que en febrero disminuyen las búsquedas, en marzo aumentan, de abril a septiembre se muestra un comportamiento no muy variable y a partir de septiembre se evidencia un fuerte crecimiento en las búsquedas de Whiskey.
Acf(SA_Whiskey, main="Gráfico ACF")
Pacf(SA_Whiskey, main="Gráfico PACF")
En el gráfico anterior se muestra el ACF, en este se puede evidenciar un decaimiento exponencial hasta la octava autocorrelación, lo cual nos confirma que la serie de tiempo en cuestión es una serie con tendencia, a partir de allí aumenta levemente la magnitud y nuevamente continua con el decaimiento.
Del gráfico anterior, el PACF se pueden identificar tres autocorrelaciones significativas, de la cual cabe resaltar que la primera en su gran magnitud indica nuevamente que la serie tiene tendencia, y también nos da a conocer que tiene autocorrelación positiva en la primera, esto quiere decir que si al mes pasado hubo un incremento en las búsquedas, en el siguiente mes se espera que también suceda, de igual manera con una disminución en las búsquedas, que para el siguiente mes también disminuirán.
yt_Whiskey <- c(D_Whiskey$`Relative Interest`)
yt_año_Whiskey <- c()
for (i in 1:length(D_Whiskey$`Relative Interest`)){
if (i<=12){
yt_año_Whiskey[i]<- NA
}else{
yt_año_Whiskey[i]<- yt_Whiskey[i-12]
}
}
p_Whiskey <- data.frame(yt_Whiskey,yt_año_Whiskey)
ggplot(p_CraftBeer, aes(yt_Whiskey,yt_año_Whiskey)) + geom_point()
## Warning: Removed 12 rows containing missing values (geom_point).
En este gráfico podemos evidenciar un comportamiento muy compacto entre el tiempo en años y las búsquedas de la palabra “Whiskey”, de lo que podemos afirmar que la varianza de la serie no es constante, Homocedasticidad y con esto se incumple el requisito para ser considerada una serie estacionaria.
Serie_sin_estacionalidad_Whiskey <- diff(SA_Whiskey, lag=12, differences = 1)
Serie_sin_estacionalidad_tendencia_Whiskey <- diff(Serie_sin_estacionalidad_Whiskey, lag=1, differences = 1)
#Sin Tendencia
autoplot(Serie_sin_estacionalidad_tendencia_Whiskey)+ theme_minimal()+ ylab("Relative Interest")+ xlab("")+ ggtitle("Relative Interest mensual")
#Sin Estacionalidad
ggseasonplot(Serie_sin_estacionalidad_tendencia_Whiskey, year.labels = TRUE, year.labels.left = TRUE)+ ggtitle("Gráfico estacional de las relative interest")+ xlab("Mes")+ ylab("Relative Interest")+ theme_minimal()
#sin heterocedasticidad
yt_Whiskey_Trans <- c(Serie_sin_estacionalidad_tendencia_Whiskey)
yt_año_Whiskey_Trans <- c()
for (i in 1:length(yt_Whiskey_Trans)){
if (i<=12){
yt_año_Whiskey_Trans[i]<- NA
}else{
yt_año_Whiskey_Trans[i]<- yt_Whiskey_Trans[i-12]
}
}
p_Whiskey_Trans <- data.frame(yt_Whiskey_Trans,yt_año_Whiskey_Trans)
ggplot(p_Whiskey_Trans, aes(yt_Whiskey_Trans,yt_año_Whiskey_Trans)) + geom_point()
## Warning: Removed 12 rows containing missing values (geom_point).
### Autocorrelacion
Acf(Serie_sin_estacionalidad_tendencia_Whiskey, main="Gráfico ACF")
Pacf(Serie_sin_estacionalidad_tendencia_Whiskey, main="Gráfico PACF")
En el gráfico de líneas anterior perteneciente a la serie de búsquedas de “Whiskey”, se puede evidenciar que el procedimiento realizado también cumplió con su cometido con esta serie, la serie es estacionaria ya que no se evidencia tendencia ni comportamiento estacional, lo que permite la aplicación y desarrollo de modelos a esta serie diferenciada.
En el gráfico estacional se presenta un análisis estacional para los años pertenecientes a la muestra para el número de búsquedas de “Whiskey”, en este podemos ver que para ninguno de los años se evidencia estacionalidad, ni comportamiento de tendencia positivo o negativo.
En el gráfico de dispersion podemos evidenciar un comportamiento disperso entre los años y las búsquedas de la palabra “Whiskey”, heteroscedasticidad, de lo que podemos afirmar que la varianza de la serie es constante y con esto se cumple el requisito para ser considerada una serie estacionaria.
En el gráfico anterior se muestra el ACF, en este no se evidencia un decaimiento exponencial ni tendencia a cero y tampoco comportamiento sinusoidal, se observan una correlación significativa la r_18, a partir de allí las correlaciones no tienen patrón alguno, lo que indica también que ya la serie no tiene tendencia alguna.
Del gráfico anterior, el PACF se pueden identificar dos auto correlaciones significativas, ambas ligeramente pasa de los límites de control, a partir de esta las correlaciones no presentan ningún comportamiento, una ligera tendencia a cero, comportamiento que puede ayudar identificar el tipo de modelo a utilizar.
Teniendo en cuenta los resultados en R, si se realiza la serie bajo la original, más no en la ajustada, se puede evidenciar que esta indica una serie no estacionaria presentando problemas de estacionalidad y tendencia, por lo tanto, podría modelarse como un proceso SARIMA más no ARIMA.
Así mismo, se pudo identificar que la serie contaba con indicios de estacionalidad, ya que esta tendía a repetirse cada 12 meses. Por esta razón, en el punto 2, se realizó una diferenciación (differences=d=1) con el fin de eliminar la estacionalidad contando con un Lag de 12, correspondiente a los 12 meses.
De tal manera, se puede referir a un posible modelo de ajuste para la parte estacional como:
Por otro lado, se logró identificar la tendencia presentada en la serie. Por esta razón, en el punto 2 se realizó una diferenciación (differences=d=1) con el fin de eliminar la tendencia presentada con un Lag 1.
De tal manera se puede referir a un posible modelo de ajuste para la tendencia como:
Como se pudo observar en el punto 2 hay autocorrelaciones significativas, en los correlogramas, tanto el ACF como el PACF tienen un comportamiento sinusoidal, por lo tanto, se podría hablar de un comportamiento semejante al ARMA, es decir, al modelo AR(p) se le incluye el modelo MA(q). Teniendo en cuenta lo anterior, se presentan 2 modelos de ajuste al comportamiento presente en la serie:
arima(SA_Whiskey,order = c(1,1,0), seasonal = list(order= c(0,1,1), period(12)))
##
## Call:
## arima(x = SA_Whiskey, order = c(1, 1, 0), seasonal = list(order = c(0, 1, 1),
## period(12)))
##
## Coefficients:
## ar1 sma1
## -0.2507 -0.2692
## s.e. 0.1429 0.1868
##
## sigma^2 estimated as 13.54: log likelihood = -128.4, aic = 262.81
Teniendo en cuenta esta alternativa, se obtuvo la siguiente ecuación:
Xt=yt-yt-12
Xt= 0 + (-2,507)Yt-1+(-’,2692)et-1+Et
arima(SA_Whiskey,order = c(0,1,1), seasonal = list(order=c(1,1,0), period(12)))
##
## Call:
## arima(x = SA_Whiskey, order = c(0, 1, 1), seasonal = list(order = c(1, 1, 0),
## period(12)))
##
## Coefficients:
## ma1 sar1
## -0.7857 -0.2390
## s.e. 0.1369 0.1715
##
## sigma^2 estimated as 12.02: log likelihood = -125.97, aic = 257.93
Teniendo en cuenta esta alternativa, se obtuvo la siguiente ecuación:
Xt=yt-yt-12
Xt= -0.7857et-1+(-0,20390)Xt-12
Modelo_Sarima <- auto.arima(SA_Whiskey)
Modelo_Sarima
## Series: SA_Whiskey
## ARIMA(1,1,1)(0,1,0)[12]
##
## Coefficients:
## ar1 ma1
## 0.3436 -0.9098
## s.e. 0.1589 0.0721
##
## sigma^2 estimated as 11.96: log likelihood=-124.56
## AIC=255.11 AICc=255.67 BIC=260.66
Teniendo en cuenta la herramienta estadística R y utilizando la función: auto.arima se obtiene cuál de los modelos planteados es el mejor para ajustarlo al modelo, en este caso, al programa estimó un modelo SARIMA con los siguientes parámetros
Teniendo en cuenta esta modelo, se obtuvo la siguiente ecuación:
Xt=yt-yt-12
Xt= et-1+(-0,4767)Xt-12
Los supuestos del modelo SARIMA recaen sobre los residuos: deben formar un proceso de ruido blanco; es decir, deben ser estacionarios y no ser autocorrelacionados entre sí.
Residuales <- Modelo_Sarima$residuals
Acf(Residuales, main = "Gráfico ACF residuales SARIMA")
Pacf(Residuales, main = "Gráfico PACF residuales SARIMA")
Se demuestra un comportamiento de ruido blanco, ya que en ninguno de los graficos de los residuos (ACF y PACF) se tiene correlación significativa (esto basado en que ningún valor de correlación sobrepasa el umbral de decisión), lo que hace que la serie sea totalmente aleatoria y no se pueda modelar.
La prueba Dickey-Fuller H0: Los errores no son estacionarios H1: Los errores son estacionarios
adf.test(Residuales)
## Warning in adf.test(Residuales): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: Residuales
## Dickey-Fuller = -4.838, Lag order = 3, p-value = 0.01
## alternative hypothesis: stationary
Prueba Ljung-Box H0: Los errores se distribuyen de forma independiente (las correlaciones son cero) H1: Los errores no se distribuyen de forma independiente (hay correlación)
Box.test(Residuales)
##
## Box-Pierce test
##
## data: Residuales
## X-squared = 0.00021274, df = 1, p-value = 0.9884
La prueba Dickey-Fuller H0: Los errores no son estacionarios H1: Los errores son estacionarios
Teniendo en cuenta el p valor 0.01 y al ser menor que la significancia 0,05, es decir 0.01 < 0.05 se rechaza la hipótesis nula y se concluye que Los errores son estacionarios.
Teniendo en cuenta el p valor 0.9884 y al ser mayor que la significancia 0.05, es decir 0.9884 > 0.05 No se rechaza la hipótesis nula y se concluye que Los errores se distribuyen de forma independiente (las correlaciones son cero).
Es así que efectivamente los supuestos del modelo SARIMA, los cuales recaen sobre los residuos, si se cumplen al establecer el modelo obtenido automáticamente.
auto.arima(Residuales)
## Series: Residuales
## ARIMA(0,0,0) with zero mean
##
## sigma^2 estimated as 8.972: log likelihood=-150.96
## AIC=303.92 AICc=303.99 BIC=306.01
Recapitulando, para el Whisky en el punto 3 se propusieron dos modelos SARIMA (p,d,q) (1,1,0)[12] y
SARIMA (p,d,q) (0,1,1)[12] basados en que se pudo identificar una estacionalidad, ya que tiende a repetirse cada 12 meses, es por esto que en el segundo punto para quitar la estacionalidad se realizo una diferenciacion (differences=d=1) para eliminar la estacionalidad con un Lag=12 correspondiente a 12 meses, por esto se puede referir a un posible modelo de ajuste para la parte estacional. Realizando el calculo con el programa R, los resultados arrojan que el modelo para este caso en especifico se asemeja a un SARIMA(1,1,1)(0,1,0)[12], en este se encuentra que no existe una grna variacion entre los modelos propuestos con tecnicas de visualizacion y el uso de un software automatico. Demostrando asi que el modelo se asemeja efectivamente a un modelo que necesita de una diferenciacion de 1 iteración, para eliminar la tendencia y que por otro lado se necesito de una diferenciacion estacional para que la serie sea estacionaria.
Pron_Sarima <- forecast(SA_Whiskey,12, model=Modelo_Sarima)
Pron_Sarima
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## Jan 2015 78.12744 73.69489 82.56000 71.34844 84.90645
## Feb 2015 72.82766 67.99608 77.65924 65.43840 80.21693
## Mar 2015 79.72467 74.77811 84.67123 72.15956 87.28978
## Apr 2015 72.68928 67.68407 77.69450 65.03446 80.34410
## May 2015 72.67713 67.62829 77.72596 64.95560 80.39866
## Jun 2015 71.67295 66.58530 76.76060 63.89206 79.45384
## Jul 2015 73.67151 68.54680 78.79623 65.83393 81.50909
## Aug 2015 72.67102 67.50999 77.83205 64.77791 80.56413
## Sep 2015 71.67085 66.47394 76.86777 63.72286 79.61885
## Oct 2015 78.67079 73.43829 83.90329 70.66838 86.67321
## Nov 2015 90.67077 85.40295 95.93860 82.61433 98.72722
## Dec 2015 104.67077 99.36786 109.97368 96.56067 112.78087
accuracy(Pron_Sarima,SP_Whiskey)
## ME RMSE MAE MPE MAPE MASE
## Training set -0.3881775 2.995332 2.133086 -0.8032252 3.524922 0.2836237
## Test set -6.2287381 6.817893 6.228738 -9.0184180 9.018418 0.8281979
## ACF1 Theil's U
## Training set -0.001882978 NA
## Test set -0.004426680 0.8508029
autoplot(SA_Whiskey)+ theme_minimal()+ ylab("haaaaaaaaaaaaaaaa")+ xlab("")+ ggtitle("Pronósticos modelo SARIMA")+ autolayer(Pron_Sarima, series = "SARIMA", PI = FALSE)
RMSE:
El error de pronóstico calculado por medio de un error cuadrático medio, arroja un valor de 6.817893 lo cual indica que existe un desface de 6.81 unidades de Whisky con referencia al valor real de búsquedas para este producto. Este valor es de igual forma interpretable por exceso o por faltante indicando que el desfase con respecto al valor real puede ser en una sobreestimación como en subestimación pero dentro de este rango (+- 6.81). Dentro del contexto este valor es relativamente bajo, ya que una empresa productora de estos productos
MAE:
El error de pronóstico calculado por medio de un error cuadrático medio, arroja un valor de 6.228738 lo cual indica que existe un desface de 6.22 unidades de Whisky con referencia al valor real de búsquedas para este producto. Este valor es de igual forma interpretable por exceso o por faltante indicando que el desfase con respecto al valor real puede ser en una sobreestimación como en subestimación pero dentro de este rango (+- 6.22).
MAPE:
El error de pronostico tiene por medio del error absoluto porcentual promedio , tiene un valor de 9.02% lo que no es un valor que se pueda declarar minimo. Estar cometiendo errores de tal magnitud puede generar perdidas mas que considerables a nivel economico. En el caso de las busquedas de Whisky se puede indicar una sobre o subestimacion mas que considerable para una posible perdida.
El metodo de calculo del error de pronostico que se recomienda como ideal es el MAE, ya que ofrece el error mas pequeño, el cual evita caer en fallos en el pronostico de busquedas de Whisky, tanto por encima como por debajo de las 6.22 unidades.
Los pronósticos les permiten a los administradores financieros que anticipe los hechos antes de que ocurra cualquier evento financiero en la empresa, la tendencia de los pronósticos ha hecho que las empresas tomen decisiones serias en el manejo de sus finanzas
Una demanda mal calculada o pronosticada puede llevar a que la empresa padezca de fierros problemas: a) Si los pronósticos se pronostica de más o de menos , la empresa incurre en sobre costos haciendo que su costo de oportunidad se incremente afectando en un futuro cercano las metas de rentabilidad y liquidez por obtener una baja rotación en cuanto a sus ventas.
Por un lado Skyrose tiene un fuerte deseo por mantener a Downtown Brew como cliente, el uso de herramientas de pronóstico puede ayudar a la empresa a dar información asertiva a la empresa de cervezas para una producción más precisa con base a la demanda que tienen actualmente, dándoles el crecimiento que ellos están buscando a lo largo de Canadá
Por otro lado una sobre estimación podría significar una sobre producción para la empresa de cervezas que puede resultar en que la empresa cervecera incurra en excedentes de inventario y finalmente terminar perdiendo el negocio por falta de precisión y confianza en los pronósticos de Skyrose
Por otro lado, On The rocks whisky tiene plena confianza en los pronósticos de Skyrose y esto le ha significado una gran ventaja para el contacto con los proveedores y así hacerse de los mejores insumos para su producción, lo que se puede traducir para ellos en un punto de valor agregado para sus productos
Una subestimación por parte de Skyrose puede ser fuertemente dañina para On the rocks, ya que los segundos pueden incurrir en pedidos que no cumplirán con el abastecimiento necesario para satisfacer la demanda de sus clientes, costándoles así el ascenso pronunciado que están teniendo en el mercado.