gpt4 book ai didi

r - 如何在 R 中使用 Monte Carlo 进行 ARIMA 模拟函数

转载 作者:行者123 更新时间:2023-12-03 09:26:32 33 4
gpt4 key购买 nike

这是我想用 R 做的算法:

  • 模拟来自 ARIMA 的 10 个时间序列数据集模型通arima.sim()功能
  • 将系列拆分为可能的子系列 2s , 3s , 4s , 5s , 6s , 7s , 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, 0)
    #n=10; phi <- 0.6; order <- c(1, 0, 0)
    bootstrap1 <- function(n, phi){
    ts <- arima.sim(n, model = list(ar=phi, order = c(1, 0, 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 1step 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 中得到类似下面的结果形式:
    ##           [,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
    但我收到以下错误消息:

    Error in MonteCarlo(func = bootstrap1, nrep = 3, param_list = param_list) :func has to return a list with named components. Each component has to be scalar.


    如何找到获得上述所需结果并使结果可重现的方法?

    最佳答案

    您收到此错误消息是因为 MonteCarlo 期望 bootstrap1()接受模拟的一个参数组合,并且每次复制只返回一个值 ( RMSE )。此处并非如此,因为块长度 ( lb ) 由 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, 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() .要并行执行模拟,您需要通过 ncpus 设置内核数. MonteCarlo 包使用 snowFall,因此它应该在 Windows 上运行。
    请注意,我还设置了 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)
    结果是一个数组。我认为最好通过 MakeFrame() 将其转换为 data.frame :
    Frame <- MakeFrame(MC_result)
    获取 reps x lb 很容易矩阵虽然:
    matrix(Frame$RMSE, ncol = length(lb), dimnames = list(1:reps, lb))

    关于r - 如何在 R 中使用 Monte Carlo 进行 ARIMA 模拟函数,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/64222355/

    33 4 0
    Copyright 2021 - 2024 cfsdn All Rights Reserved 蜀ICP备2022000587号
    广告合作:1813099741@qq.com 6ren.com