## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.width = 6,
  fig.height = 4
)

## ----eval = FALSE-------------------------------------------------------------
# remotes::install_github("lhvanegasp/mtarm")

## ----eval = FALSE-------------------------------------------------------------
# install.packages("mtarm")

## ----fig.width=9, fig.height=7------------------------------------------------
library(mtarm)
data(iceland.rf)       
str(iceland.rf)     

## ----fig.width=9, fig.height=5------------------------------------------------
summary(iceland.rf[,-5])

## ----fig.width=9, fig.height=5------------------------------------------------
plot(ts(as.matrix(iceland.rf[,-5])), main="Iceland")

## -----------------------------------------------------------------------------
set.seed(09102)
fits <- mtar_grid(~ Jokulsa + Vatnsdalsa | Temperature | Precipitation,
                  data=iceland.rf, subset={Date<="1974-11-06"}, 					      
                  row.names=Date, nregim.min=2, nregim.max=2, p.min=15,  				
                  p.max=15, q.min=4, q.max=4, d.min=2, d.max=2,  					      
                  n.burnin=500, n.sim=400, n.thin=2, ssvs=TRUE,
                  dist=c("Gaussian","Student-t","Laplace"),
                  plan_strategy="multisession")

fits

## -----------------------------------------------------------------------------
DIC(fits)

WAIC(fits)

## -----------------------------------------------------------------------------
newdata <- subset(iceland.rf, Date>"1974-11-06") 
set.seed(09102)

oos <- out_of_sample(fits, newdata=newdata, n.ahead=nrow(newdata), FUN=median) 
oos[,c(1,2,5,6)]

## -----------------------------------------------------------------------------
set.seed(09102)
oos2 <- out_of_sample(fits, newdata=newdata, n.ahead=nrow(newdata), 
                      rolling=5, FUN=median) 

for(i in 1:length(oos2)){                   
   cat("\n",i,"-step-ahead\n")   
   print(oos2[[i]][,c(1,2,5,6)]) 
}                            

## -----------------------------------------------------------------------------
summary(fits[["Laplace.2.15.4.2"]])

## ----fig.width=9, fig.height=5------------------------------------------------
res <- residuals(fits[["Laplace.2.15.4.2"]])  

## ----fig.width=9, fig.height=5------------------------------------------------
par(mfrow=c(1,2)) 
qqnorm(res[["full"]], pch=20, col="blue", main="") 
abline(0, 1, lty=3) 
hist(res[["full"]], freq=FALSE, xlab="Quantile-type residual",
     ylab="Density", main="") 
curve(dnorm(x), col="blue", add=TRUE)

## ----fig.width=9, fig.height=5------------------------------------------------
par(mfrow=c(1,2))  
acf(res[["full"]], col="blue", main="")
pacf(res[["full"]], col="blue", main="")

## -----------------------------------------------------------------------------
pred <- predict(fits[["Laplace.2.15.4.2"]], newdata=newdata,
                n.ahead=nrow(newdata), row.names=Date, credible=0.8)

head(pred$summary)
tail(pred$summary)

## -----------------------------------------------------------------------------
fitmcmc <- coda::as.mcmc(fits[["Laplace.2.15.4.2"]])
summary(fitmcmc)

## -----------------------------------------------------------------------------
geweke_diagTAR(fits[["Laplace.2.15.4.2"]])


## -----------------------------------------------------------------------------
effectiveSize_TAR(fits[["Laplace.2.15.4.2"]])