R中的蒙特卡罗模拟

Question

我试图模拟来自rho=0.7的AR（1）模型的数据（Y）。然后我将使用这些数据在截距上运行Y的回归（通过这样做，参数估计成为Y的平均值），然后使用鲁棒的标准错误。我想对这个假设运行一个蒙特卡罗模拟，使用2000次重复不同的滞后值。目的是显示当滞后变化时Newey West估计器的有限样本性能

A<-array(0, dim=c(2000,1))
for(i in 1:2000){
  y_new<-arima.sim(model=list(ar=0.7), n=50, mean=0,sd=1)
  reg<-lm(y_new~1)
  ad<-coeftest(reg, alternative="greater", vcov=NeweyWest(reg, lag=1, prewhite=FALSE))
  A[i]<-ad[,3]
}

我的问题是：上面的代码是进行这种模拟的正确方法吗？如果是，我如何得到一个代码来重复这个过程在HAC测试中的不同滞后值。我想运行测试每次增加1的滞后，因此我将这样做50次滞后1,2,3,4......,50，每次将2000个模拟测试统计信息存储在不同名称的向量中（sig. level=0,05，使用1.645的临界值），并根据各种滞后值绘制它们（拒绝概率）。请帮帮忙

Answer 1

因为你没有提到模拟的可能目的，所以很难说这是不是正确的方式。

通过为每个模拟样本计算50个测试统计信息，而不是为每个延迟重复模拟2000次（即模拟次数为2000*50），可以节省大量时间。

做模拟的更好的格式是

library(AER)
library(dplyr)
lags <- 1:50
nreps <- 2000

sim <- function (){
  ynew <- arima.sim(model = list(ar=0.7), n=50, mean=0, sd=1)
  reg <- lm(ynew ~ 1 )
  s <- rep(NA, 50)
  for(i in lags){    
    ad <- coeftest(reg, alternative="greater", vcov=NeweyWest(reg, lag = i, prewhite=FALSE))
    s[i] <- ad[ ,4]
  }
  s
}

以下代码将模拟结果存储在数据中。框架

 result <- lapply(1:nreps, function(i)data.frame(simulation = i, lag = lags, pvalues = sim())) %>%
 rbind_all

根据你模糊的描述，我推断出你想要的是什么样子

library(ggplot2)
result %>% 
  group_by(lag) %>% 
  summarize(rejectfreq = mean(pvalues > 0.05)) %>% 
  ggplot(., aes(lag, rejectfreq)) + geom_line()+
  coord_cartesian(ylim = c(0,1)) +
  scale_y_continuous(breaks=seq(0, 1, by=0.1))

虽然该图仅使用100次模拟创建，但很明显，当扰动项被确认时，Newey-West中滞后的选择并不重要。