使用R做蒙特卡罗模拟，解决以下问题

Question

我的问题是：如何生成两个时间序列从下列数据生成过程中产生的100长度的：

xt = xt−1 + et, 
yt = yt−1 + ut, 

其中et和ut是相互独立的标准正常的iid系列。在包含截距的上的yt的回归中，模拟拒绝零假设的概率，即斜率系数α为零，公称大小为5％。

1. 对α= 0,0.5,0.9,0.95,0.99的不同值进行模拟。
2. 讨论你的结果。随着阿尔法值的变化，拒绝概率如何变化？

我的想法是，我通过编码

y <- arima.sim(100,model=list(ar=0)) 
x <- arima.sim(100,model=list(ar=0)) 

然后做回归Ÿx和存储的α的值，并重复上述步骤，为1000倍生成两个华宇系列，并找到阿尔法的分布。 但我的问题是：

1. 我不知道如何存储阿尔法值
2. 如何重复回归1000次

我R的新学员，希望有人可以通过编写R代码告诉我如何解决这个问题。谢谢！

Answer 1

一个tidyverse解决方案：

library(dplyr) 
library(broom) 
library(purrr) 
library(ggplot2) 

# define a function the returns the alpha -- its point estimate, standard error, etc. -- from the OLS 
iteration <- function() { 
    y <- arima.sim(100,model=list(ar=0)) 
    x <- arima.sim(100,model=list(ar=0)) 
    lm(y~x) %>% 
    broom::tidy() %>% 
    filter(term == 'x') 
} 

# 1000 iterations of the above simulation 
alphas <- map_df(1:1000, ~iteration()) 

# plot the results 
alphas %>% 
    ggplot(aes(estimate)) + geom_density()