Using an arima.sim()
function to simulate time series data that follows a particular ARIMA
model requires a lot of trials (especially for a ridiculously small sample size) of changing set.seed
value like the example bellow:
library(forecast)
set.seed(1)
ar1 <- arima.sim(n = 15, model=list(ar=0.2, order = c(1, 0, 0)), sd = 1)
(ar2 <- auto.arima(ar1, ic ="aicc"))
set.seed(3)
ar1 <- arima.sim(n = 15, model=list(ar=0.2, order = c(1, 0, 0)), sd = 1)
(ar2 <- auto.arima(ar1, ic ="aicc"))
set.seed(3)
ar1 <- arima.sim(n = 15, model=list(ar=0.2, order = c(1, 0, 0)), sd = 1)
(ar2 <- auto.arima(ar1, ic ="aicc"))
until I get my desired result
I have found an R
code that prints out the right seed to set for AR(1)
:
library(future.apply)
FUN <- function(i) {
set.seed(i)
ar1 <- arima.sim(n=15, model=list(ar=0.6, order=c(1, 0, 0)), sd=1)
ar2 <- auto.arima(ar1, ic="aicc")
(cf <- ar2$coef)
if (length(cf) == 0) {
rep(NA, 2)
}
else if (all(grepl(c("ar1|intercept"), names(cf))) &
substr(cf["ar1"], 1, 6) %in% "0.6000") {
c(cf, seed=i)
}
else {
rep(NA, 2)
}
}
R <- 1e4
seedv <- 1:R
library(parallel)
cl <- makeCluster(detectCores() - 1 + 1)
clusterExport(cl, c("FUN"), envir=environment())
clusterEvalQ(cl, suppressPackageStartupMessages(library(forecast)))
res <- parLapply(cl, seedv, "FUN")
(res1 <- res[!sapply(res, anyNA)])
stopCluster(cl)
res2 <- Reduce(function(...) merge(..., all=T), lapply(res1,
function(x) as.data.frame(t(x))))
res2[order(res2$seed), ]
which produces this result:
#ar1 seed
#1 0.8000417 4195
which when I check as bellow:
ar1 <- arima.sim(n=15, model=list(ar=0.6, order=c(1, 0, 0)), sd=1)
(ar2 <- forecast::auto.arima(ar1, ic="aicc"))
it will be true
Series: ar1
ARIMA(1,0,0) with zero mean
Coefficients:
ar1
0.6000
s.e. 0.1923
sigma^2 estimated as 2.051: log likelihood=-26.38
AIC=56.75 AICc=57.75 BIC=58.17
What I want
I want a case of AR(2)
Bellow is my trial:
FUN <- function(i) {
set.seed(i)
ar1 <- arima.sim(n=15, model=list(ar=c(-0.9, -0.9), order=c(2, 0, 0)), sd=1)
ar2 <- auto.arima(ar1, ic="aicc")
(cf <- ar2$coef)
if (length(cf) == 0) {
rep(NA, 2)
}
else if (all(grepl(c("ar1|ar2|intercept"), names(cf))) &
substr(cf["ar1"], 1, 5) %in% "-0.90" & substr(cf["ar2"], 1, 5) %in% "-0.90") {
c(cf, seed=i)
}
else {
rep(NA, 2)
}
}
R <- 2e5
seedv <- 1e5:R
library(parallel)
cl <- makeCluster(detectCores() - 1 + 1)
clusterExport(cl, c("FUN"), envir=environment())
clusterEvalQ(cl, suppressPackageStartupMessages(library(forecast)))
res <- parLapply(cl, seedv, "FUN")
(res1 <- res[!sapply(res, anyNA)])
stopCluster(cl)
res2 <- Reduce(function(...) merge(..., all=T), lapply(res1,
function(x) as.data.frame(t(x))))
res2[order(res2$seed), ]
The search stops here
######################################################################
I confirm my result as bellow:
set.seed(1509)
ar1 <- arima.sim(n=20, model=list(ar=c(-0.9, -0.9), order=c(2, 0, 0)), sd=1)
library(forecast)
(ar2 <-auto.arima(ar1, ic="aicc"))
My output on coefficients of AR(2) are not exactly -0.90, -0.90
:
Series: ar1
ARIMA(2,0,0) with zero mean
Coefficients:
ar1 ar2
-0.9369 -0.8997
s.e. 0.0940 0.0752
sigma^2 estimated as 0.9048: log likelihood=-28.12
AIC=62.24 AICc=63.74 BIC=65.23
question from:
https://stackoverflow.com/questions/65894592/searching-for-the-right-ar2-seed-like-i-did-for-ar1-in-r