library(ggplot2)
library(forecast)
library(tseries)
library(readxl)
SNG_GDP <- read_excel("C:/Users/samorapa/OneDrive/Estudios/universidad/Economía Internacional/Proyecto de Aula 2/Anexos y Consolidados/Consolidados Fianles/data segundo proyecto de aula.xlsx")
names(SNG_GDP)
## [1] "year" "SGP" "CHN" "JPN" "KOR" "MYS" "HKG" "IDN"
summary(SNG_GDP$SGP)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -3.680 4.903 7.867 7.428 10.308 15.240
library(ggplot2)
plot_GDP <- ggplot(SNG_GDP, aes(y=SGP, x=year, colour="PIB de Singapur"))+geom_line()+ggtitle("Evolución del PIB de Singapure")+ylab("PIB")+xlab("Años")
plot_GDP

#Limpiando datos
library(forecast)
SNG_GDP$SGP_ts <- ts(SNG_GDP[, c("SGP")])
library(tseries)
library(forecast)
library(ggplot2)
SNG_GDP$ma5 <- ma(SNG_GDP$SGP_ts, order=5)
SNG_GDP$ma15 <- ma(SNG_GDP$SGP_ts, order = 15)
ggplot()+geom_line(data = SNG_GDP, aes(x=year, y=SGP, colour="PIB de Singapur"))+geom_line(data = SNG_GDP, aes(x=year, y=ma5, colour="Promedio Movil de 5 periodos"))+geom_line(data = SNG_GDP, aes(x=year, y=ma15, colour="Promedio Movil de 15 periodos"))+ggtitle("Tendencia del PIB de Singapure")+xlab("Año")+ylab("PIB de Singapur")
## Warning: Removed 4 rows containing missing values (geom_path).
## Warning: Removed 14 rows containing missing values (geom_path).

library(tseries)
library(forecast)
GDP_ma5 <- ts(na.omit(SNG_GDP$ma5), frequency=12)
decomp <- stl(GDP_ma5, s.window = "periodic")
deseasonal_GDP <- seasadj(decomp)
plot(decomp)

adf.test(GDP_ma5, alternative = "stationary")
## Warning in adf.test(GDP_ma5, alternative = "stationary"): p-value smaller
## than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: GDP_ma5
## Dickey-Fuller = -6.4567, Lag order = 3, p-value = 0.01
## alternative hypothesis: stationary
Acf(GDP_ma5, main="")

Pacf(GDP_ma5, main="")

GDP_D1 <- diff(deseasonal_GDP, differences = 1)
plot(GDP_D1)

#Test de Dickey-Fuller
adf.test(GDP_D1, alternative = "stationary")
##
## Augmented Dickey-Fuller Test
##
## data: GDP_D1
## Dickey-Fuller = -3.2453, Lag order = 3, p-value = 0.0898
## alternative hypothesis: stationary
Acf(GDP_D1, main="ACF para Primeras Diferencias")

Pacf(GDP_D1, main="Pafc para Primeras Diferencias")

auto.arima(deseasonal_GDP, seasonal = FALSE)
## Series: deseasonal_GDP
## ARIMA(2,1,2)
##
## Coefficients:
## ar1 ar2 ma1 ma2
## -0.6199 -0.9230 0.9031 0.8553
## s.e. 0.0722 0.0636 0.1093 0.0983
##
## sigma^2 estimated as 0.7268: log likelihood=-63.47
## AIC=136.95 AICc=138.28 BIC=146.61
fit <- auto.arima(deseasonal_GDP, seasonal = FALSE)
tsdisplay(residuals(fit), lag.max = 18, main = "(2,1,2) Residuos del Modelo")

fcast <- forecast(fit, h=5)
plot(fcast)

#Probando el modelo sin estacionalidad
hold <- window(ts(deseasonal_GDP), start=50)
fit_no_hold <- arima(ts(deseasonal_GDP[-c(50:55)]), order=c(2,1,2))
fcast_no_hold <- forecast(fit_no_hold, h=5)
plot(fcast_no_hold)
lines(ts(deseasonal_GDP))

fit_w_seasonality <- auto.arima(deseasonal_GDP, seasonal = TRUE)
fit_w_seasonality
## Series: deseasonal_GDP
## ARIMA(0,1,2)(0,0,1)[12]
##
## Coefficients:
## ma1 ma2 sma1
## 0.2183 -0.3669 -0.8234
## s.e. 0.1395 0.1627 0.5500
##
## sigma^2 estimated as 0.6806: log likelihood=-67.18
## AIC=142.36 AICc=143.23 BIC=150.09
seas_fcast <- forecast(fit_w_seasonality, h=5)
plot(seas_fcast)

tsdisplay(residuals(fit_w_seasonality), lag.max = 18, main = "(0,1,2)(0,0,1)[12] Residuos del Modelo Estacional")

#Probando el modelo estacional
hold <- window(ts(deseasonal_GDP), start=50)
fit2_no_hold <- arima(ts(deseasonal_GDP[-c(50:55)]), order=c(0,1,2),seasonal = list(order=c(0,0,1), period=12))
fcast2_no_hold <- forecast(fit2_no_hold, h=5)
plot(fcast2_no_hold)
lines(ts(deseasonal_GDP))

fit_w_seasonality2 <- arima(deseasonal_GDP,order = c(0,1,5), seasonal = list(order=c(0,0,1), period=12))
fit_w_seasonality2
##
## Call:
## arima(x = deseasonal_GDP, order = c(0, 1, 5), seasonal = list(order = c(0, 0,
## 1), period = 12))
##
## Coefficients:
## ma1 ma2 ma3 ma4 ma5 sma1
## 0.1085 0.0967 0.3114 -0.2953 -0.5991 -0.2681
## s.e. 0.1088 0.1471 0.1648 0.1550 0.1343 0.2191
##
## sigma^2 estimated as 0.6525: log likelihood = -64.45, aic = 142.9
seas_fcast2 <- forecast(fit_w_seasonality2, h=5)
plot(seas_fcast2)

tsdisplay(residuals(fit_w_seasonality2), lag.max = 18, main = "(0,1,5)(0,0,1)[12] Residuos del Modelo Estacional Ajustando Promedio Móvil")

fit_w_seasonality3 <- arima(deseasonal_GDP,order = c(8,1,5), seasonal = list(order=c(0,0,1), period=12))
## Warning in arima(deseasonal_GDP, order = c(8, 1, 5), seasonal = list(order
## = c(0, : possible convergence problem: optim gave code = 1
fit_w_seasonality3
##
## Call:
## arima(x = deseasonal_GDP, order = c(8, 1, 5), seasonal = list(order = c(0, 0,
## 1), period = 12))
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 ar7 ar8
## 0.7845 0.5484 -0.6352 0.0461 -0.6726 0.5098 0.7706 -0.9503
## s.e. 0.0843 0.1010 0.2014 0.1101 0.1721 0.1447 0.1342 0.2729
## ma1 ma2 ma3 ma4 ma5 sma1
## -0.7261 -0.9143 0.7954 0.0696 0.0850 0.9391
## s.e. 0.3262 0.2716 0.3524 0.2514 0.1956 1.9185
##
## sigma^2 estimated as 0.2867: log likelihood = -53.71, aic = 137.42
seas_fcast3 <- forecast(fit_w_seasonality3, h=5)
plot(seas_fcast3)

tsdisplay(residuals(fit_w_seasonality3), lag.max = 18, main = "(8,1,5)(0,0,1)[12] Residuos del Modelo Estacional Ajustando el Autoregresor")

#Probando el modelo estacional con Promedio Móvil y Autoregresor ajustado
hold2 <- window(ts(deseasonal_GDP), start=50)
fit2_no_hold2 <- arima(ts(deseasonal_GDP[-c(50:55)]), order=c(8,1,5),seasonal = list(order=c(0,0,1), period=12))
## Warning in arima(ts(deseasonal_GDP[-c(50:55)]), order = c(8, 1, 5),
## seasonal = list(order = c(0, : possible convergence problem: optim gave
## code = 1
fcast2_no_hold2 <- forecast(fit2_no_hold2, h=5)
plot(fcast2_no_hold2)
lines(ts(deseasonal_GDP))
