作者：許夢潔 (編譯) (知乎 | 簡書 | 碼雲)
 

Stata連享會 「Stata 現場班報名中……」

Source: Ashish Rajbhandari → Vector autoregression—simulation, estimation, and inference in Stata

文章目錄

• Stata連享會 [「Stata 現場班報名中……」](https://gitee.com/arlionn/stata_training/blob/master/README.md)
• 1. 模擬
• 2. 設定係數值
• 3. 由多維正態分佈中生成擾動項
• 3. 生成觀測序列
• 4. 估計
• 5. 推斷：脈衝響應函式
• 6. 正交化的脈衝響應函式
• 結論
• 參考文獻
• 關於我們
• 聯絡我們
• 往期精彩推文

1. 模擬

首先使用如下設定模擬雙變數 VAR(2) :

$\begin{array}{}\text{(1)}\end{array}$

其中 $y_{ 1,t }$ 是在時間 $t$ 的觀測變數， $\mu$ 是一個 $2\times 1$ 的截距項向量， $A_1$ 和 $A_2$ 均為 $2\times 2$ 的引數矩陣， $\epsilon _{ t }$ 是不與時間相關的擾動項。我假設 $\epsilon _{ t }$ 服從分佈 $N(0,\Sigma)$ ，其中 $\Sigma$ 是一個 $2\times 2$ 的協方差矩陣。

設定樣本量為 1000 ，並在 Stata 中生成需要的變數：

 clear all

. set seed 2016

. local T = 1100

. set obs `T'
number of observations (_N) was 0, now 1,100

. gen time = _n

. tsset time
        time variable:  time, 1 to 1100
                delta:  1 unit

. generate y1 = .
(1,100 missing values generated)

. generate y2 = .
(1,100 missing values generated)

. generate eps1 = .
(1,100 missing values generated)

. generate eps2 = .
(1,100 missing values generated)

在第 1-6 行，我設定了隨機種子，將樣本量設為 1000，並且生成了一個時間變數 time 。接下來，我生成了變數 y1 ，y2 ， eps1 和 eps2 來存放觀測序列和擾動項。

2. 設定係數值

我為本文的 VAR(2) 模型選擇瞭如下的引數值：
$A_{ 1 }= \begin{bmatrix} 0.6 &-0.3 \\ 0.4& 0.2 \end{bmatrix},\quad A_{ 2 }= \begin{bmatrix} 0.2 &0.3 \\ -0.1& 0.1 \end{bmatrix} \tag{2}$

. mata:
------------------------------------------------- mata (type end to exit) -----
: mu = (0.1\0.4)

: A1 = (0.6,-0.3\0.4,0.2)

: A2 = (0.2,0.3\-0.1,0.1)

: Sigma = (1,0.5\0.5,1)

: end
-------------------------------------------------------------------------------

在 Mata 中，我分別建立了矩陣 $\mu$ ， $A_1$ ， $A_2$ 和 $\Sigma$ 來放置引數值。在模擬得到樣本之前，我檢查了由這些引數值是否會得到一個平穩的 VAR(2) 模型。令：

$F= \begin{bmatrix} A_{ 1 }&A_{ 2 } \\ I_{ 2 } & 0 \end{bmatrix}$

其中 $I_2$ 是 $2\times 2$ 的識別矩陣， $0$ 是一個每個元素均為 0 的 $2\times 2$ 矩陣。如果矩陣 F 的所有特徵根均小於 1 ，則此 VAR(2) 過程即為一個平穩過程。下面放上計算特徵根的程式碼：

. mata:
------------------------------------------------- mata (type end to exit) -----
: K = p = 2               // K = number of variables; p = number of lags

: F = J(K*p,K*p,0)

: F[1..2,1..2] = A1

: F[1..2,3..4] = A2

: F[3..4,1..2] = I(K)

: X = L = .

: eigensystem(F,X,L)

: L'
                              1
    +----------------------------+
  1 |                .858715598  |
  2 |  -.217760515 + .32727213i  |
  3 |  -.217760515 - .32727213i  |
  4 |                .376805431  |
    +----------------------------+

: end
------------------------------------------------------------------------------

我根據之前的設定建立了矩陣 F 並使用函式 eigensystem( ) 計算其特徵根。 矩陣 X 中儲存了特徵向量，L 中儲存特徵根。所有 L 中的特徵根均小於 1 ，因此該 VAR(2) 過程是平穩的。在檢驗過是否平穩之後，接下來生成 VAR(2) 模型的擾動項。

3. 由多維正態分佈中生成擾動項

我從分佈 $N(0,\Sigma)$ 中生成兩個隨機正態變數，並將它們分別賦值給變數 eps1 和 eps2。

. mata:
------------------------------------------------- mata (type end to exit) -----
: T = strtoreal(st_local("T"))

: u = rnormal(T,2,0,1)*cholesky(Sigma)

: epsmat = .

: st_view(epsmat,.,"eps1 eps2")

: epsmat[1..T,.] = u

: end
-------------------------------------------------------------------------------

我將樣本規模(在 Stata 中定義為一個區域性巨集變數 T )賦值給一個 Mata 數值變數，這步將可以簡化之後的工作。在 Mata 中，我使用了兩個函式：st_local( ) 和 strtoreal( ) 來儲存樣本大小。第一個函式可以從 Stata 巨集中獲取文字值，第二個函式則可以將文字值轉變為實數值。

第二行生成了一個 $1100\times 2$ 的正態擾動項矩陣，其每個元素都服從分佈 $N(0,\Sigma)$ 。我使用函式 st_view( ) 將生成的正態擾動項存放