## Not run:
# #########################################################################
# #### ####
# #### Monte-Carlo Portmanteau Tests ####
# #### ####
# #########################################################################
# ## Monte-Carlo test for randomness series ##
# #########################################################################
# data("DEXCAUS")
# returns <- log(DEXCAUS[-1]/DEXCAUS[-length(DEXCAUS)])
# portest(returns) ## MC using one CPU takes about 25.16 seconds
# portest(returns, nslaves=4) ## MC using 4 CPUs takes about 9.51 seconds
# portest(returns, MonteCarlo=FALSE) ## asymptotic gvtest
# portest(returns,test="LjungBox", MonteCarlo=FALSE) ## asymptotic LjungBox
# #########################################################################
# ## Monte-Carlo goodness-of-fit arima test using 4 CPUs ##
# #########################################################################
# ## arima() or Arima() function takes about 14.32 seconds
# ans1 <- arima(WWWusage,order=c(3,1,0))
# portest(ans1, nslaves = 4)
# #
# ## arima0() function takes about 15.26 seconds
# ans2 <- arima0(WWWusage,order=c(3,1,0))
# portest(ans2, nslaves = 4)
# #
# ## auto.arima() function from package forecast takes about 13.59 seconds
# library("forecast")
# ans3 <- auto.arima(WWWusage)
# portest(ans3, nslaves = 4)
# #
# ## ar() function takes about 9.39 seconds
# ans4 <- ar(Nile,order.max=2)
# portest(ans4, nslaves = 4)
# #
# ## FitAR() function takes about 10.78 seconds
# library("FitAR")
# ans5 <- FitAR(Nile, p=2)
# portest(ans5, nslaves = 4)
# #########################################################################
# ## Monte-Carlo goodness-of-fit VAR test - Multivariate series ##
# #########################################################################
# data("IbmSp500")
# ibm <- log(IbmSp500[,2]+1)*100
# sp500 <- log(IbmSp500[,3]+1)*100
# IBMSP500 <- data.frame(cbind(ibm,sp500))
# ## ar.ols() function takes about 9.11 seconds
# ans6 <- ar.ols(IBMSP500, aic=FALSE, intercept=TRUE, order.max=5)
# portest(ans6, NREP=100, test="gvtest", nslaves=4)
# ## VAR() function takes about 11.55 seconds
# library("vars")
# ans7 <- VAR(IBMSP500, p=5)
# portest(ans7, NREP=100, test="gvtest", nslaves=4)
# portest(ans7,test="Hosking", MonteCarlo=FALSE) ## asymptotic Hosking test
# #########################################################################
# ## Monte-Carlo test for GARCH effects using 4 CPUs ##
# #########################################################################
# ## Example 1
# ## Test for ARCH effects on returns series takes about 14.65 seconds
# data("monthintel")
# returns <- as.ts(monthintel)
# lags <- c(5, 10, 20, 40)
# portest(returns, lags = lags, nslaves = 4, SquaredQ = TRUE)
# #
# ## Example 2
# library("fGarch")
# library("tseries")
# data("GNPDEF")
# z<-ts(GNPDEF[,2], start=1947, freq=4)
# r <- 100*diff(log(z))
# ## use garch() function takes about 6.75 seconds
# FitGarch1 <- garch(r, order = c(1,1))
# portest(FitGarch1,NREP=100,nslaves = 4)
# portest(FitGarch1,NREP=100,nslaves = 4,SquaredQ=TRUE)
# #
# ## use garchFit() function takes about 13.56 seconds
# GarchFit2 <- garchFit(formula = ~arma(4,0)+garch(1,1), data=r, trace=FALSE)
# portest(GarchFit2, NREP=100, nslaves = 4, SquaredQ = FALSE)
# portest(GarchFit2, NREP=100, nslaves = 4, SquaredQ = TRUE)
# #########################################################################
# ## Monte-Carlo test on residuals with infinite variances ##
# #########################################################################
# ## It takes about 32.7 seconds on personal PC with 4 CPUs
# data("CRSP")
# CRSP.AR5<- arima(CRSP, c(5, 0, 0))
# lags <- c(10, 20, 30)
# portest(CRSP.AR5,lags=lags,nslaves=4,NREP=1000,InfiniteVarianceQ=TRUE)
# #########################################################################
# ## Monte-Carlo test for Fractional Gaussian Noise, FGN. ##
# #########################################################################
# ## It takes about 55.06 seconds on personal PC with 4 CPUs
# library("FGN")
# data("NileMin")
# NILE.FGN <- FitFGN(NileMin)
# lags <- c(5, 10, 20)
# portest(NILE.FGN, lags = lags, nslaves = 4)
# portest(NILE.FGN, MonteCarlo=FALSE) ## asymptotic distribution method
# ##############################################################
# ## Write two functions to fit a model and simulate results
# ## Apply Monte-Carlo test on fitted obj with class "list"
# ##############################################################
# ## Example 1
# ## Threshold Autoregressive (TAR) Model example from TSA package
# ## It takes about 64.27 seconds on personal PC with 4 CPUs
# library("TSA")
# FitModel <- function(data){
# fit <- TSA::tar(y=log(data),p1=4,p2=4,d=3,a=0.1,b=0.9,print=FALSE)
# res <- ts(fit$std.res)
# parSpec <- list(res=res,fit=fit)
# parSpec
# }
# SimModel <- function(parSpec){
# fit <- parSpec$fit
# exp(tar.sim(fit)$y)
# }
# data(prey.eq)
# portest(FitModel(prey.eq),nslaves=4,func=list(SimModel,FitModel),pkg="TSA")
# #
# ## Example 2
# ## It takes about 10.75 seconds on personal PC with 4 CPUs
# FitModel <- function(data){
# fit <- ar(data,aic = FALSE, order.max=2)
# order <- 2
# res <- ts(fit$resid[-(1:order)])
# phi <- fit$ar
# theta <- NULL
# sigma <- fit$var.pred
# demean <- fit$x.mean
# list(res=res,phi=phi,theta=theta,sigma=sigma,demean=demean)
# }
# SimModel <- function(parSpec){
# res <- parSpec$res
# n <- length(res)
# innov <- sample(x=res,size=n,replace = TRUE)
# phi <- parSpec$phi
# theta <- parSpec$theta
# sigma <- parSpec$sigma
# demean <- parSpec$demean
# arima.sim(n = n, list(ar = phi, ma = theta), innov = innov,
# sd = sqrt(sigma), mean = demean)
# }
# Fit <- FitModel(Nile)
# portest(Fit,nslaves=4,func=list(SimModel=SimModel,FitModel=FitModel),pkg="stats")
# #########################################################################
# #### ####
# #### Simulation using varima.sim Function ####
# #### ####
# #########################################################################
# # Simulate MA(1) where innovation series is provided via argument innov
# ########################################################################
# set.seed(1234)
# n <- 200
# phi <- NULL
# theta <- 0.6
# d <- NA
# sigma <- 1.9
# Z <- varima.sim(phi, theta, d, sigma, n,innov=rnorm(n))
# plot(Z)
# ########################################################################
# # Simulate ARIMA(2,1,0) process with phi=c(1.3,-0.35), Gaussian innovations
# # The series is truncated at lag 50
# ########################################################################
# set.seed(1234)
# Trunc.Series <- 40
# n <- 1000
# phi <- c(1.3, -0.35)
# theta <- NULL
# d <- 1
# sigma <- 1
# Z <- varima.sim(phi,theta,d,sigma,n,Trunc.Series=Trunc.Series)
# coef(arima(Z,order=c(2,1,0)))
# ########################################################################
# # Simulate MA(1) process with theta = 0.5, t5-distribution innovations
# ########################################################################
# set.seed(1234)
# n <- 200
# phi <- NULL
# theta <- 0.5
# Z <- varima.sim(phi, theta, sigma=1, n=n, innov.dist="t", df=5)
# plot(Z)
# ########################################################################
# # Simulate univariate ARMA(2,1) process with length 500,
# # phi = c(1.3, -0.35), theta = 0.1. Drift equation is 8 + 0.05*t
# # Stable innovations with: ALPHA = 1.75, BETA = 0, GAMMA = 1, DELTA = 0
# ########################################################################
# set.seed(1234)
# n <- 500
# phi <- c(1.3, -0.35)
# theta <- 0.1
# constant <- 8
# trend <- 0.05
# demean <- 0
# d <- 0
# sigma <- 0.7
# ALPHA <- 1.75
# BETA <- 0
# GAMMA <- 1
# DELTA <- 0
# Stable <- c(ALPHA,BETA,GAMMA,DELTA)
# Z <- varima.sim(phi,theta,d,sigma,n,constant,trend,demean,
# innov.dist="stable",StableParameters=Stable)
# plot(Z)
# ########################################################################
# # Simulate a bivariate white noise series from a multivariate t4-distribution
# ########################################################################
# set.seed(1234)
# Z <- varima.sim(sigma=diag(2),n=200,innov.dist="t",df=4)
# plot(Z)
# ########################################################################
# # Simulate a trivariate VARMA(1,1) process with length 300.
# # phi = array(c(0.5,0.4,0.1,0.5,0,0.3,0,0,0.1), dim=c(k,k,1)), where k =3
# # theta = array(c(0,0.25,0,0.5,0.1,0.4,0,0.25,0.6), dim=c(k,k,1)).
# # innovations are generated from multivariate normal distribution
# # The process have mean c(10, 0, 12),
# # Drift equation a + b * t, where a = c(2,1,5), and b = c(0.01,0.06,0)
# # The series is truncated at default value: Trunc.Series=ceiling(100/3)=34
# ########################################################################
# set.seed(1234)
# k <- 3
# n <- 300
# Trunc.Series <- 50
# phi <- array(c(0.5,0.4,0.1,0.5,0,0.3,0,0,0.1),dim=c(k,k,1))
# theta <- array(c(0,0.25,0,0.5,0.1,0.4,0,0.25,0.6),dim=c(k,k,1))
# sigma <- diag(k)
# constant <- c(2,1,5)
# trend <- c(0.01,0.06,0)
# demean <- c(10,0,12)
# Z <- varima.sim(phi, theta, d = 0,sigma, n, constant,trend,demean)
# plot(Z)
# ########################################################################
# # Simulate a bivariate VARIMA(2,d,1) process with n=300, where d=(1,2).
# # phi = array(c(0.5,0.4,0.1,0.5,0,0.3,0,0),dim=c(k,k,2)),
# # theta = array(c(0,0.25,0,0), dim=c(k,k,1)).
# # innovations are generated from multivariate normal
# # The process have mean zero and no deterministic terms.
# # The variance covariance is sigma = matrix(c(1,0.71,0.71,2),2,2).
# # The series is truncated at default value: Trunc.Series=ceiling(100/3)=34
# ########################################################################
# set.seed(1234)
# k <- 2
# n <- 300
# Trunc.Series <- 50
# phi <- array(c(0.5,0.4,0.1,0.5,0,0.3,0,0),dim=c(k,k,2))
# theta <- array(c(0,0.25,0,0),dim=c(k,k,1))
# d <- c(1,2)
# sigma <- matrix(c(1,0.71,0.71,2),k,k)
# Z <- varima.sim(phi, theta, d, sigma, n)
# plot(Z)
# ########################################################################
# # Simulate a bivariate VAR(1) process with length 600.
# # Stable distribution: ALPHA=(1.3,1.6), BETA=(0,0.2), GAMMA=(1,1), DELTA=(0,0.2)
# # The series is truncated at default value: Trunc.Series=min(100,200)=100
# ########################################################################
# set.seed(1234)
# k <- 2
# n <- 600
# phi <- array(c(-0.2,-0.6,0.3,1.1),dim=c(k,k,1))
# theta <- NULL
# d <- NA
# sigma <- matrix(c(1,0.71,0.71,2),k,k)
# ALPHA <- c(1.3,1.6)
# BETA <- c(0,0.2)
# GAMMA <-c(1,1)
# DELTA <-c(0,0.2)
# Stable <- c(ALPHA,BETA,GAMMA,DELTA)
# Z <- varima.sim(phi,theta,d,sigma,n,innov.dist="stable",StableParameters=Stable)
# plot(Z)
# ## End(Not run)
Run the code above in your browser using DataLab