# We'll use generic data stored at R data basis:
library(KFAS)
## Please cite KFAS in publications by using:
##
## Jouni Helske (2017). KFAS: Exponential Family State Space Models in R. Journal of Statistical Software, 78(10), 1-39. doi:10.18637/jss.v078.i10.
data <- AirPassengers
plot(data)
# As we can see, the series presents 1) Seasonality, and 2) Trending.
# We proceed by subtracting 1) and then 2):
# Decomposing data:
decomp <- decompose(data, type = c("multiplicative"))
# Unseasonalization:
data_detr<- data/decomp$seasonal
#data_detr <- data_unse/decomp$trend
plot(data_detr, type = "l", main = "Data Unsesionalized",
xlab = "Time", ylab = " ", col = "blue")
inits_manual <- c(1, 1, 1)
modelo <- SSModel(data_detr ~ SSMtrend(1, Q = 0.01)+
SSMcycle(period = 10, Q = 0.01), H = 0.01)
fit <- fitSSM(modelo, inits = inits_manual, method = "BFGS")
results <- KFS(fit$model, smoothing = c("state", "signal","disturbance"))
plot(data_detr, type = "l", col = "blue", lwd = 1,
ylab = " ", xlab = "Time", main = "Unsesionalized vs Filtered-Smoothed Series")
lines(results$muhat, col = "red", lty = 2,
lwd = 2)
lines(data, col = "purple", lty = 3)
legend("topleft", legend = c("Unsesionalized", "Filtered-Smoothed", "Original"),
col = c("blue", "red","purple"), lty = c(1, 2, 3), lwd = 1,cex = 0.5)
results_manual <- KFS(fit$model, smoothing = "disturbance")
# Trend residuals
plot(results_manual$etahat[, 1], type = "l", col = "blue",
main = "Trend Residuals",
xlab = "Tiempo", ylab = "Residuals")
# Cycle residuals
lines(results_manual$etahat[, 2], type = "l", col = "red",
main = "Cycle Residuals",
xlab = "Tiempo", ylab = "Residuals")
abline(a=0,b=0, h=1965, v = 0, lty = 2)
legend("topleft", legend = c("Trend Residuals", "Cycle Residuals"),
col = c("blue", "red"), lty = c(1, 1),cex = 0.8)
library(seasonal)
m <- seas(data, x11 = "")
plot(m)
series <- final(m)
inits_manual <- c(1, 1, 1)
modelo_1 <- SSModel(series ~ SSMtrend(1, Q = 0.01)+
SSMcycle(period = 10, Q = 0.01), H = 0.01)
fit_1 <- fitSSM(modelo_1, inits = inits_manual, method = "BFGS")
results_1 <- KFS(fit_1$model, smoothing = c("state", "signal","disturbance"))
par(mfrow = c(1,2))
plot(data_detr, type = "l", col = "blue", lwd = 1,
ylab = " ", xlab = "Time", main = "Multiplicative desagriggation")
lines(results$muhat, col = "red", lty = 2,
lwd = 2)
lines(data, col = "purple", lty = 3)
legend("topleft", legend = c("Unsesionalized", "Filtered-Smoothed", "Original"),
col = c("blue", "red","purple"), lty = c(1, 2, 3), lwd = 1,cex = 0.5)
mtext("Unsesionalized and Filtered-Smoothed Series Comparison", outer = TRUE, cex = 1.5)
plot(series, type = "l", col = "blue", lwd = 1,
ylab = " ", xlab = "Time", main = "ARIMA X-13")
lines(results_1$muhat, col = "red", lty = 2,
lwd = 2)
lines(data, col = "purple", lty = 3)
legend("topleft", legend = c("Unsesionalized", "Filtered-Smoothed", "Original"),
col = c("blue", "red","purple"), lty = c(1, 2, 3), lwd = 1,cex = 0.5)
results_manual <- KFS(fit$model, smoothing = "disturbance")
results_manual_1 <- KFS(fit_1$model, smoothing = "disturbance")
par(mfrow = c(2,2))
# Trend Residuals
plot(results_manual$etahat[,1], type = "l", col = "blue",
main = "Multiplicative Desagriggation \n Trend Residuals",
xlab = "Time", ylab = "Residuals")
abline(a=0,b=0, h=1965, v = 0, lty = 2)
plot(results_manual_1$etahat[,1], type = "l", col = "purple",
main = "ARIMA X-13 \n Trend Residuals",
xlab = "Time", ylab = "")
abline(a=0,b=0, h=1965, v = 0, lty = 2)
# Cycle Residuals
plot(results_manual$etahat[,2], type = "l", col = "red",
main = "Cycle Residuals",
xlab = "Time", ylab = "Residuals")
abline(a=0,b=0, h=1965, v = 0, lty = 2)
plot(results_manual_1$etahat[,2], type = "l", col = "green",
main = "Cycle Residuals",
xlab = "Time", ylab = "")
abline(a=0,b=0, h=1965, v = 0, lty = 2)
