如何解决如何在Windows上的R中将蒙特卡洛用于ARIMA仿真功能
这是我想对R执行的算法：
- 通过
ARIMA函数从arima.sim()模型中模拟10个时间序列数据集 - 将系列划分为可能的
2s，3s，4s，5s，6s，7s，{{ 1}}和8s。 - 对于每种大小，对块进行重新采样并进行替换，以得到新的序列，并通过
9s函数从每个块大小的子序列中获得最佳的ARIMA模型。 - 获取每个块大小
auto.arima()的每个子系列。
下面的RMSE函数可以完成此任务。
R调用函数
## Load packages and prepare multicore process
library(forecast)
library(future.apply)
plan(multisession)
library(parallel)
library(foreach)
library(doParallel)
n_cores <- detectCores()
cl <- makeCluster(n_cores)
registerDoParallel(cores = detectCores())
## simulate ARIMA(1,0)
#n=10; phi <- 0.6; order <- c(1,0)
bootstrap1 <- function(n,phi){
ts <- arima.sim(n,model = list(ar=phi,order = c(1,0)),sd = 1)
########################################################
## create a vector of block sizes
t <- length(ts) # the length of the time series
lb <- seq(n-2)+1 # vector of block sizes to be 1 < l < n (i.e to be between 1 and n exclusively)
########################################################
## This section create matrix to store block means
BOOTSTRAP <- matrix(nrow = 1,ncol = length(lb))
colnames(BOOTSTRAP) <-lb
########################################################
## This section use foreach function to do detail in the brace
BOOTSTRAP <- foreach(b = 1:length(lb),.combine = 'cbind') %do%{
l <- lb[b]# block size at each instance
m <- ceiling(t / l) # number of blocks
blk <- split(ts,rep(1:m,each=l,length.out = t)) # divides the series into blocks
######################################################
res<-sample(blk,replace=T,10) # resamples the blocks
res.unlist <- unlist(res,use.names = FALSE) # unlist the bootstrap series
train <- head(res.unlist,round(length(res.unlist) - 10)) # Train set
test <- tail(res.unlist,length(res.unlist) - length(train)) # Test set
nfuture <- forecast::forecast(train,model = forecast::auto.arima(train),lambda=0,biasadj=TRUE,h = length(test))$mean # makes the `forecast of test set
RMSE <- Metrics::rmse(test,nfuture) # RETURN RMSE
BOOTSTRAP[b] <- RMSE
}
BOOTSTRAPS <- matrix(BOOTSTRAP,nrow = 1,ncol = length(lb))
colnames(BOOTSTRAPS) <- lb
BOOTSTRAPS
return(list(BOOTSTRAPS))
}
我得到以下结果：
bootstrap1(10,0.6)
我想按时间顺序重复上述## 2 3 4 5 6 7 8 9
## [1,] 0.8920703 0.703974 0.6990448 0.714255 1.308236 0.809914 0.5315476 0.8175382
至step 1，然后想到step 4中的Monte Carlo技术。因此，我加载了它的包并运行以下功能：
R期望以param_list=list("n"=10,"phi"=0.6)
library(MonteCarlo)
MC_result<-MonteCarlo(func = bootstrap1,nrep=3,param_list = param_list)
形式获得以下类似结果：
matrix但是我收到以下错误消息：
MonteCarlo中的错误（func = bootstrap1，nrep = 3，param_list = param_list）： func必须返回包含命名组件的列表。每个分量都必须是标量的。
我如何找到获得上述期望结果并使结果可再现的方式？
编辑
我希望预期的## [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
## [1,] 0.8920703 0.703974 0.6990448 0.714255 1.308236 0.809914 0.5315476 0.8175382
## [2,] 0.8909836 0.8457537 1.095148 0.8918468 0.8913282 0.7894167 0.8911484 0.8694729
## [3,] 1.586785 1.224003 1.375026 1.292847 1.437359 1.418744 1.550254 1.30784
将在 Windows
解决方法
您收到此错误消息是因为MonteCarlo希望bootstrap1()接受模拟的一个参数组合，并且它仅返回一个值（{{1} }）。这里不是这种情况，因为块长度（RMSE）是由模拟时间序列（lb）的长度确定的（ in n，所以您将为每个调用获取bootstrap1个块长度的结果。
一种解决方案是将块长度作为参数传递，并适当地重写n - 2：
bootstrap1()要运行模拟，请将参数以及library(MonteCarlo)
library(forecast)
library(Metrics)
# parameter grids
n <- 10 # length of time series
lb <- seq(n-2) + 1 # vector of block sizes
phi <- 0.6 # autoregressive parameter
reps <- 3 # monte carlo replications
# simulation function
bootstrap1 <- function(n,lb,phi) {
#### simulate ####
ts <- arima.sim(n,model = list(ar = phi,order = c(1,0)),sd = 1)
#### devide ####
m <- ceiling(n / lb) # number of blocks
blk <- split(ts,rep(1:m,each = lb,length.out = n)) # divide into blocks
#### resample ####
res <- sample(blk,replace = TRUE,10) # resamples the blocks
res.unlist <- unlist(res,use.names = FALSE) # unlist the bootstrap series
#### train,forecast ####
train <- head(res.unlist,round(length(res.unlist) - 10)) # train set
test <- tail(res.unlist,length(res.unlist) - length(train)) # test set
nfuture <- forecast(train,# forecast
model = auto.arima(train),lambda = 0,biasadj = TRUE,h = length(test))$mean
### metric ####
RMSE <- rmse(test,nfuture) # return RMSE
return(
list("RMSE" = RMSE)
)
}
param_list = list("n" = n,"lb" = lb,"phi" = phi)
传递到bootstrap1()。为了并行进行仿真，您需要通过MonteCarlo()设置内核数。 MonteCarlo软件包使用snowFall，因此应在Windows上运行。
请注意，我还设置了ncpus（否则结果将是所有复制的平均值）。事先设置种子将使结果可重复。
raw = T结果是一个数组。我认为最好通过set.seed(123)
MC_result <- MonteCarlo(func = bootstrap1,nrep = reps,ncpus = parallel::detectCores() - 1,param_list = param_list,export_also = list(
"packages" = c("forecast","Metrics")
),raw = T)
将其转换为data.frame：
MakeFrame()尽管很容易获得Frame <- MakeFrame(MC_result)
矩阵：
reps x lb版权声明：本文内容由互联网用户自发贡献，该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容， 请发送邮件至 dio@foxmail.com 举报，一经查实，本站将立刻删除。