Welcome to WuJiGu Developer Q&A Community for programmer and developer-Open, Learning and Share
Welcome To Ask or Share your Answers For Others

Categories

0 votes
240 views
in Technique[技术] by (71.8m points)

R: search for the Right Seed(s) that Simulate AR(2) with its Coefficients(phi_1 and phi_2)

I want R to search for the seed(s) that will simulate `AR(2):

ar1.2 <- arima.sim(n = 100, model=list(ar=c(0.5, 0.4), order = c(2, 0, 0)), sd = 1)

Such that when I run the following:

set.seed()
ar1.2 <- arima.sim(n = 100, model=list(ar=c(0.5, 0.4), order = c(2, 0, 0)), sd = 1)
auto.arima(ar1.2)

It will confirm the AR- coefficients as phi_1 = 0.5... and phi_2 = 0.4...

I have two conditions to be met chronologically as follows:

  1. Make sure that the ARIMA order is (2, 0, 0)
  2. Out of the seeds that satisfied condition 1 above, search for the one that its coefficients as 0.5... and 0.4...

Edit

It is done for `AR(1) as follows:

FUN <- function(i) {
  set.seed(i)
  ar1 <- arima.sim(n=100, model=list(ar=0.4, order=c(1, 0, 0)), sd=1)
  ar2 <- auto.arima(ar1, ic="aicc")
  c(arimaorder(ar2), seed=i)
}

R <- 24000  ## this would be your 1e5
seedv <- 23000:R

library(parallel)
cl <- makeCluster(detectCores() - 1 + 1)
clusterExport(cl, c("FUN"), envir=environment())
clusterEvalQ(cl, suppressPackageStartupMessages(library(forecast)))

res <- parSapply(cl, seedv, "FUN")

stopCluster(cl)

seed_out <- res["seed", which(apply(res, 2, function(x) all(x[1:3] == c(1, 0, 0))))]


#########################################################################
sink("ARIMA.SIM_n20_ar0.4.txt")
##########################################################################
arima_order_results = data.frame()
for (my_seed in seed_out){
   set.seed(my_seed)
  ar1 <- arima.sim(n = 100, model=list(ar=0.4, order = c(1, 0, 0)), sd = 1)
  ar2 <- auto.arima(ar1, ic = "aicc")
  arr <- as.data.frame(t(ar2$coef))
  if(substr(as.character(arr[1]), 1, 5) == "0.400") {

    arr <- cbind(data.frame(seed=my_seed),arr)
    print(arr)

    arima_order_results = bind_rows(arima_order_results,arr)
    # write.csv(my_seed, paste0(arr, ".csv"), row.names = FALSE)

  }
}
##########################################################################
sink()
#######################################################################

Get the output in your Working Directory

   #seed       ar1
#1 23027 0.4009039
   #seed       ar1 intercept
#1 23305 0.4005298 0.4055362
    #seed       ar1
#1 23443 0.4004223
    #seed      ar1
#1 23845 0.400621

I confirm the output with a seed of 23027

set.seed(23027)
ar1 <- arima.sim(n=100, model=list(ar=0.4, order=c(1, 0, 0)), sd=1)
library(forecast)
auto.arima(ar1, ic="aicc")

A parallel and better Answer for the case of AR(1) is here

I want such a solution for AR(2)


与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome To Ask or Share your Answers For Others

1 Answer

0 votes
by (71.8m points)

I used the approach from the solution of the link in the question above. I also use ar1 = 0.3 while ar2 = 0.4

FUN2 <- function(i) { 
  set.seed(i) 
  ar1 <- arima.sim(n=20, model=list(ar=c(0.3, 0.4), 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, 3) %in% "0.3" & substr(cf["ar2"], 1, 3) %in% "0.4") { 
    c(cf, seed=i) 
  } 
  else { 
    rep(NA, 2) 
  } 
}

R <- 1e4 
seedv <- 1:R 

library(parallel) 
cl <- makeCluster(detectCores() - 1 + 1) 
clusterExport(cl, c("FUN2"), envir=environment()) 
clusterEvalQ(cl, suppressPackageStartupMessages(library(forecast))) 

res <- parLapply(cl, seedv, "FUN2") 

(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), ]

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome to WuJiGu Developer Q&A Community for programmer and developer-Open, Learning and Share
...