AR-Model-模擬教學

假設

1. PCA已完成。
2. 模型估計已完成—選到AR(2)，參數c=0.1; alpha1 = -0.72; alpha2 = -0.68。

• Model: y_t+1 = c + alpha1＊y_t + alpha2＊y_t-1 + e_t+1。
• initial_value: y_0 = 0.55; y_1 = 0.31。

3. 想模擬1000組60年情境。

> ## 亂數設定
> N = 1000
> t = 60
> set.seed(9527)
> e_ <- matrix(rnorm(N*t), N,t)
> ## 參數估計結果
> c=0.1; alpha_1 = -0.72; alpha_2 = -0.68
> param = c(c, alpha_1, alpha_2)
> ## 用AR(2)舉例，一開始會需要歷史前兩期的起始值... y_0 = 0.55; y_1 = 0.31，前兩期的模擬會受的起始值影響，所以另外處理。
> y_t = list()
> y_t[[1]] = sum(param*c(1, 0.55, 0.31)) +e_[,1]
> y_t[[2]] = param[1] + param[2]*y_t[[1]] + param[3]*0.31 + e_[,2]
> for ( i in 3:t){
+   y_t[[i]] = param[1] + param[2]*y_t[[i-1]]+ param[3]*y_t[[i-2]] +e_[,i]
+ }
> ### 把結果合併成一個data.frame()，列=情境，行=時間。
> result = do.call('cbind',  y_t) %>% data.frame()
> # write.csv(result, ...)

• 上面是1個PC的1000組模擬，如果想要模擬多個PC可以把上面程式碼包成function，搭配apply或是for迴圈執行。
• 下面挑幾舉模擬情境來看看。

第1組情境

> plot_scen <- function(x){
+ df = data.frame('Date'= 2018:(2018+t-1), result[x,] %>% t())
+ colnames(df)[2] = 'Scenario'
+ ggplot(data=df, aes(x=Date, y=Scenario)) +
+   geom_line(color = '#4d79ff')+
+   geom_point()+
+   scale_x_continuous(breaks = seq(2018, (2018+t-1), 10))
+ }
> 
> plot_scen(1)

第807組情境

> plot_scen(807)

第123組情境

> plot_scen(123)