M B K 48

もし貧乏人が経済学を学んだら

経済学モデルを解くための道具

Generalized Schur Method(シュール分解)―― The ABCs of RBCs 6章

(3具体例で行列Gが抜けていたので修正しました。)

 前回の記事で、
150321 (1)

   (1)

 という線形方程式を解く方法のひとつとして、Blandhard and Kahn method を取り上げましたが、Generalized Schur method (Generalized Schur decomposition シュール分解)も、同様な方法のひとつです。B が逆行列を持たない場合、Blandhard and Kahn method は使えませんが、Generalized Schur method ならそういう場合でも使えます。

 ここで、 xは状態変数(決定変数、t 期にすでに決定されてしまっている。例えば、資本ストック k など)で、yはジャンプ変数(未決定変数、t 期にはまだ決定されていない。例えば、消費 c )です。上の表記では、かく乱項( ε)がひとつだけですが、複数あっても(列ベクトルになる)この方法が使えます(たんに、εに[・] をつけるのを忘れただけ・・・)。この場合、A,B は正方行列でないといけません(行と列の数が同じ)。

 どちらかというと、Blandhard and Kahn method よりも、Generalized Schur method のほうが簡単で、かつ、Bが逆行列を持たなくても使えるので、こちらのほうが使いやすいと思います。

 まず、行列 B と行列 A を次のように分解します(AとBは正方行列でなければいけません)。
B=QTZ’
A=QSZ’
 ここで、行列 Q,行列 Z は、
QQ’=Q’Q=I=ZZ’=Z’Z
 になる性質をもっています(Q’,Z’のアポストロフィー(’)は、転置行列を表しています)。Matlab や Scilab などのプログラムに、この変換を行ってくれるコードがあるのでそれを利用します。そのとき、システムの固有値が小さい順に並ぶように変換される必要がありますが、上記のプログラムが自動でやってくれます。(ABCs of RBCs には、Scilab にはその機能がない、と書いてありますが、あります。)

img_JIRO

 「ただ」だからといって、なめてもらっちゃ困るぜ





 そして、S,T は上三角行列になっています(左下半分の要素が0)。最後に具体例を挙げますが、その例では、S,T がそれぞれ、
S=
    0    6.0682764     -1.7409307    -5.3804398    1.0028218 
    0.   5.9070823     -2.67548        -6.338396      1.6699841 
    0.           0.          0.8002334     -0.3662861    0.2318056 
    0.           0.              0.               1.1313905     0.3079086 
    0.           0.              0.                     0.           0.7914712

T=
1.689645     6.6220135      -2.7886191  -5.6935763    1.3700334 
    0.           6.1940275      -2.2035468   -5.2838763    1.8324852 
    0.                 0.             0.8423509     -0.3532317    0.9531660 
    0.                 0.                  0.             1.0681878    0.4692144 
    0.                 0.                  0.                 0.                   0.
 になります。

1 ショックがない場合
15032603
   
(2)

B=QTZ’
A=QSZ’
 に変換されているので、(2)式に代入すれば、
150323 (2)



 となります。左から Q’ をかけます。( Q’Q=I [単位行列]になるので)
150323 (3)

 

 次に、Z’ を以下のように分割します。
150323 (4)



 そのとき、システムの固有値のうち絶対値が1よりも小さい固有値の数がZ11(正方行列)の次数と同じようになるようにします(プログラムがその数を教えてくれます)。以下の具体例では3つです。したがって、その場合は、Z11が3×3行列になるように分割します。(ちなみに、一般的にこのような線形方程式の場合同じですが、その数が状態変数(state variable あるいは、決定変数 predetermined variable)の数と一致している場合のみ安定的な解をもちます。)
 同様に他の行列、A,B,S,T も同じように分割します。S,T の場合は、上三角行列なので、次のようになります。
150323 (5)



 そこで、(2)式を書き換えると、

150323 (6)



150323 (7)



 となります。これを展開し、下半分の方程式だけを取り出すと、
150323 (8)

 です。これはシステムの固有値のうち絶対値が1よりも大きいものに対応しています。したがって、発散します。システムが安定的な解を持つためには、左辺と右辺の[・・・]の中が0にならなければなりません。したがって(右辺の[・・・]に関して)
150323 (9)

 が成り立ちます。この式から、
150323 (10)

 が得られます。これで、ジャンプ変数 yを状態変数 xで表すことができました。
150323 (11)

 とおけば、
150323 (12)
   (3)
 です。この(3)式から、yt+1 =-N E[xt+1]  となります((3)式の t を t+1 に変えるだけ)。
 次に一番最初の(2)式に戻って、yにこの(3)式を、yt+1 に-N E[xt+1] を代入します。A,B を分割したかたちで表記すれば、

150323 (14)



 となります。これを展開し、今度は、上半分の方程式を取り出すと、
150323 (15)

 です。この式から( xt+1 の前の行列の逆行列を左からかければ)、
150323 (16)

 が得られます。これで、t+1期の状態変数 xt+1 も、t 期の状態変数 xで表されたことになります。

2 ショックがある場合
B=QTZ’
A=QSZ’
 と変換されているので、(1)式に代入すれば、
150323 (1)

   

 左から Q’をかけます。
150323 (17)



 それぞれの行列を分割したかたちで表せば、

150323 (18)



150323 (19)



 となります。G の分割は、上記の1より小さい固有値の数の行で分割します(ただし、列の数がその数以下だったら)。つまり、その固有値数が3で、Gの列が3以下だったら、上から3行までをG1にする、ということです。
 これを展開し、再び下半分の方程式を取り出し、それが発散しないための条件から、
150323 (20)

 が得られます(これは、下半分の方程式の右辺が 0 になる、という条件から得られています)。
 ここから、yを求めれば、
150323 (21)

150323 (22)

 (ここで、(S22Z’22-1=Z’22-122-1という関係を使っています。) これで、ジャンプ変数 yを状態変数 xで表すことができました。
150323 (11)

150323 (23)

 とおけば、
150323 (24)   (4)
 
 と表すことができます。

 次に、(4)式のtを将来に1期ずらし、つまり、t を t+1 に変え、yt+1 を表せば、
150323 (25)   (5)
 
 となります。t+1 期のかく乱項 εt+1 の t 期から見た期待値( E[εt+1] )は 0 になる、と想定できます(将来のショックは予測がつかない)。したがって、E[εt+1] の項は消えます。このように、将来のかく乱項の期待値が0になる、という性質を用いて方程式を解いているので、Blanchard and Kahn method の場合と同様に、ショックを表す方程式を連立方程式のシステムの中に入れておくことが必要です
  そこで、最初の(1)式に戻って、yと yt+1 とにそれぞれ、(4)式と(5)式を代入すれば、
150323 (26)



150323 (27)



 となります。これを展開し、今度は、上半分の方程式を取り出すと、
150323 (28)

 です。この式から( xt+1 の前の行列の逆行列を左からかければ)、
150323 (29)

150323 (30)

 が得られます。これで、t+1 期の状態変数 xt+1 も、t 期の状態変数 xで表されたことになります。
それぞれの変数の前の行列を別の文字で置き換えれば、
150321 (48)

150321 (38)   (6)

 と表すことができます。


3 具体例
 ABCs of RBCs では Hansen's Indivisible Labor model (労働時間が決まっていて、雇用されるとそれを調節できない)を例としています。そこで、ここでは、Hansen's (basic) model を使って計算してみます。
 対数化した後に導出される方程式は、次の6つです(101ページ)。

150323 (31)

150323 (32)

150323 (33)

150323 (34)



150323 (35)

150323 (36)

 最後の方程式はショックを表すものです。変数は5つ(K,Y,C,H,r)で、方程式は6つです(また、この場合、状態変数は、K,λ,Y,ジャンプ変数はC,rです)。A,Bは正方行列でなければいけないので、上記の方程式からHを消去し、方程式を5つにします。そして、x=[K,λ,Y]’,y=[C,r]’(どちらも、列ベクトルなので、転置を表す(’)を使って表記しています)です。

150323 (37)

150323 (36)

150323 (38)

150323 (39)

150323 (40)


 行列を使って表せば、
150323 (41)








150323 (42)









 となります。上の行列がBで、下の行列がAです。G=[1 0 0 0 0]’ です。

 そこで、この二つの行列 B,A を変換します。
B=QTZ’
A=QSZ’
 Scilab の場合、
[S,T,Z,dim] = schur(A,B,'d');
 というコードを使います(左の[・・・]の中は何でもよい。こちらで名前を与えればよい)。このシステムの固有値で、絶対値が1より小さい固有値の数は、dim で示されます。この例の場合、dim と入力すれば、3と答えを返してくるので、その数は3です。これは状態変数の数と一致しています(したがって、安定的な解をもちます)。

 しかし、Scilab では Q を返してくれません(私がコードを知らないだけかもしれませんが)。そこで、こちらで計算するのですが、この場合、B=QTZ’ というひとつの式だけでは求まりません(S,T)が逆行列を持たないため)。そこで、
B=QTZ’
A=QSZ’
 の2つの式を利用します。それぞれ、左から Z をかけます。
BZ=QT
AZ=QS
 両方の式を足します。
(B+A)Z=Q(T+S)
 そして、右からT+Sの逆行列をかけます。
Q=(A+B)*Z*inv(S+T)
(他の方法もあります。最初に両方の式を足して、(TZ’+SZ’)の逆行列を右からかける。)

 S,T はすでに上に示したとおりです。
 Z は、
Z=
  2.057D-17    0.7148765   - 0.2785001   - 0.6165709    0.1767189 
      0.            0.           - 0.6153090    0.0524966   - 0.7865360 
    - 1.        2.801D-17      2.552D-18   - 1.823D-17  - 3.214D-18  
- 8.974D-19    0.4068383    - 0.387616      0.7480079    0.3531580  
- 1.279D-17  - 0.5687127   - 0.6273646   - 0.2399345    0.4747746 

 となります。Z の転置行列 Z’は、
Z’=
    2.057D-17           0.               - 1.         - 8.974D-19   - 1.279D-17 
    0.7148765           0.            2.801D-17    0.4068383   - 0.5687127 
  - 0.2785001  - 0.6153090     2.552D-18    - 0.387616    - 0.6273646 
  - 0.6165709    0.0524966   - 1.823D-17      0.7480079  - 0.2399345  
    0.1767189  - 0.7865360   - 3.214D-18      0.3531580    0.4747746
 
 となります。固有値のうち1より小さい固有値の数が3(dim=3)なので、Z’11 が 3×3 になるように分割します。
15032602




 後は、同様に、S,T,A,Bを分割し、上記の方法で計算していくだけです。
 
 最終的に(6)式の C,D,N,L が
C  =
    0.9536739    0.1075239  - 3.109D-18 
    0.                  0.95               0.        
    0.2044602    1.3796686  - 1.594D-17
D=
    0.1131831 
    1.        
    1.4522827
N=
  - 0.5691029  - 0.3723670  - 2.143D-17 
    0.7955398  - 1.3796686    9.174D-18
L=
  - 0.3919653 
  - 1.4522827
 となります。


 Scilab のコード

AA=1.72;
beta=0.99;
theta=0.36;
delta=0.025;
gamma=0.95;

hbar=1/(1+AA*(1-beta*delta*theta/(1-beta*(1-delta)))/(1-theta));
kbar=hbar*(theta/((1/beta)-1+delta))^(1/(1-theta));
ybar=kbar^theta*hbar^(1-theta);
cbar=ybar-delta*kbar;
rbar=theta*ybar/kbar;

B=[kbar,0,-ybar,0,0;
0,1,0,0,0;
0,-1,1-(1-theta)*(1-hbar),0,0;
0,0,1,0,0;
0,0,0,1,-beta*rbar];

A=[(1-delta)*kbar,0,0,-cbar,0;
0,gamma,0,0,0;
theta,0,0,-(1-theta)*(1-hbar),0;
1,0,0,0,1;
0,0,0,1,0];

G=[0;1;0;0;0];

[S,T,Z,dim] = schur(A,B,'d');
Q=(A+B)*Z*inv(S+T);

Z2=Z';
Z21=Z2(4:5,1:3);
Z22=Z2(4:5,4:5);
N=inv(Z22)*Z21;
Q2=Q';
Q21=Q2(4:5,1:3);
Q22=Q2(4:5,4:5);
S22=S(4:5,4:5);
G1=G(1:3);
G2=G(4:5);L=inv(Z22)*inv(S22)*(Q21*G1+Q22*G2);
B11=B(1:3,1:3);
B12=B(1:3,4:5);
A11=A(1:3,1:3);
A12=A(1:3,4:5);
C=inv(B11-B12*N)*(A11-A12*N);
D=inv(B11-B12*N)*(G1-A12*L);


Blanchard and Kahn Method ―― The ABCs of RBCs 6章

 Blanchard and Kahn Method は、経済学モデルで次のような線形方程式を解く方法のひとつです(Bの逆行列を左からかけるだけでは解けない場合がよくある)。

$B\left[\begin{array}{l}
x_{t+1}\\
E_{t} y_{t+1}
\end{array}\right]= A\left[\begin{array}{l}
x_{t}\\
y_{t}
\end{array}\right]+G\varepsilon_{t}$
   (1)

 ここで、 $ x_{t}$ は状態変数(先決変数 predetermined variable、t 期にすでに決定されてしまっている。例えば、資本ストック $ k_{t}$ など)で、$ y_{t}$ はジャンプ変数(未決定変数 non-predetermined variable、t 期にはまだ決定されていない。例えば、消費 $ c_{t}$ )です。上の表記では、かく乱項( $\varepsilon_{t}$ )がひとつだけですが、複数あっても(列ベクトルになる)この方法が使えます(たんに、$\varepsilon_{t}$ に[・] をつけるのを忘れただけ・・・)。この場合、A,B は正方行列でないといけません(行と列の数が同じ)。

 まず、ショックがない( $G\varepsilon_{t}$ がない)場合を考えます。

1 ショックがない場合
 B の逆行列を左からかけます(したがって、Blanchard and Kahn Method は、B が逆行列を持たない場合には使えません。そのような場合には、シュール分解という方法があります。これについては別の記事で)。

$\left[\begin{array}{l}
x_{t+1}\\
E_{t} y_{t+1}
\end{array}\right]= B^{-1}A\left[\begin{array}{l}
x_{t}\\
y_{t}
\end{array}\right]+B^{-1}G\varepsilon_{t}$
   (2)

 B-1A を Z とおき、Z を以下のように対角化します。

$Z=B^{-1}A=M\Lambda M^{-1}$

 コンピューターの計算ソフトに行列を対角化するコードがあるのでそれを使います(後の例で見るように、Scilab の場合は、
[M,D]=spec(Z);
 です([ , ]の中はこちらで名前を与えればよい)。
 Λ は Z の固有値が並んだ対角行列です。以下に紹介する具体的な例では、その固有値は、
0                       
1.1319444 + 0.2196525i  
1.1319444 - 0.2196525i  
0.5 
 になります。これを絶対値が小さい順に並び変え、新たな対角行列をつくります。(ちなみに、絶対値が1より小さい固有値の数と状態変数の数が一致していないと上の方程式は安定的な均衡解をもちません)。

 上の例の場合、Λ を
0    0                      0                                      0                       
0    0.5                   0                                      0                       
0    0      1.1319444 - 0.2196525i                     0                       
0    0                     0                         1.1319444 + 0.2196525i 

 に変えます。この場合、元の Λ を、1列、4列、3列、2列の順番に並び変えています。そこで、行列 M も同じように並びかえます。元の M は、

0            0                           0                        0.8092927  
1.    0.0779335 - 0.4777554i    0.0779335 + 0.4777554i      0.3447199  
0       - 0.7786466                 - 0.7786466               - 0.4611485  
0     0.1936582 - 0.3491162i    0.1936582 + 0.3491162i    - 0.1164286

 です。これを、1列、4列、3列、2列の順番に並び変えます。したがって、新たな M は、

0     0.8092927           0                             0                       
1.    0.3447199    0.0779335 + 0.4777554i    0.0779335 - 0.4777554i  
0   - 0.4611485     - 0.7786466                  - 0.7786466               
0   - 0.1164286    0.1936582 + 0.3491162i    0.1936582 - 0.3491162i 

 になります。この新しい M を $\overline{M}$ と名付けておきます。この $\overline{M}$ の逆行列を求めます。そして、その $\overline{M}$ の逆行列を次のように4つに分割(ブロック化)します。

$\overline{M}^{-1} =\left[\begin{array}{l}
\hat{M}_{11} \hat{M}_{12}\\
\hat{M}_{21} \hat{M}_{22}
\end{array}\right]$

 そのとき、M11 行列の次数が、Zの固有値のうち絶対値が1より小さい固有値の数と等しくなるようにします。上の例では、その数が2なので、M11 を2×2行列に分割する、ということです。Mバーの逆行列が以下のようになるので、それを、

15032601



 のように分割します。他の行列、Λ(並び変えたもの)、A,B も同じように分割します。Λバー (元のΛを並び替えた行列)は対角行列なので、対角要素以外はゼロです。

$\overline{\Lambda} = \left[\begin{array}{l}
\overline{\Lambda}_{11}\hspace{7pt}0\\
\hspace{2pt}0\hspace{13pt}\overline{\Lambda}_{22}
\end{array}\right]$

  そこで、元の方程式、(1)式に戻って、B-1A が上記のように対角化できるので、(2)式は、

$\left[\begin{array}{l}
x_{t+1}\\
E_{t} y_{t+1}
\end{array}\right]= =\overline{M}^{ }\overline{\Lambda} \overline{M}^{-1}\left[\begin{array}{l}
x_{t}\\
y_{t}
\end{array}\right]$

 となります。Mバーの逆行列を左かからかければ、

$\overline{M}^{-1}\left[\begin{array}{l}
x_{t+1}\\
E_{t} y_{t+1}
\end{array}\right]= \overline{\Lambda} \overline{M}^{-1}\left[\begin{array}{l}
x_{t}\\
y_{t}
\end{array}\right]$

 となります。上記のように分割した形で表せば、

$\left[\begin{array}{l}
\hat{M}_{11} \hat{M}_{12}\\
\hat{M}_{21} \hat{M}_{22}
\end{array}\right]\left[\begin{array}{l}
x_{t+1}\\
E_{t} y_{t+1}
\end{array}\right]= \left[\begin{array}{l}
\overline{\Lambda}_{11}\hspace{7pt}0\\
 \hspace{2pt}0\hspace{13pt}\overline{\Lambda}_{22}
\end{array}\right]\left[\begin{array}{l}
\hat{M}_{11} \hat{M}_{12}\\
\hat{M}_{21} \hat{M}_{22}
\end{array}\right]\left[\begin{array}{l}
x_{t}\\
y_{t}
\end{array}\right]$

 です。これを計算し、上半分と下半分を別々に方程式で書けば、

$[\hat{M}_{11}x_{t+1}+\hat{M}_{12}E_{t}y_{t+1}]=\overline{\Lambda}_{11}[\hat{M}_{11}x_{t}+\hat{M}_{12}y_{t}]$   (3)

$[\hat{M}_{21}x_{t+1}+\hat{M}_{22}E_{t}y_{t+1}]=\overline{\Lambda}_{22}[\hat{M}_{21}x_{t}+\hat{M}_{22}y_{t}]$   (4)

 となります。下の(4)式の Λ22 は、1より大きい要素が対角に並んだものです。この式は発散することになります([・・・]の中を何か別の名前を与えて、再帰的に代入していけば、右辺の[・・・]の前の行列がΛ22の累乗になる)。したがって、このシステムが安定的な均衡解を持つためには、(4)式の左辺と右辺の[ ・・・ ]の中がそれぞれ 0 になる必要があります。

$\hat{M}_{21}x_{t}+\hat{M}_{22}y_{t}=0$   (5)

$\hat{M}_{21}x_{t+1}+\hat{M}_{22}E_{t}y_{t+1}=0$   (6)

 (5)式から、

$y_{t}= -\hat{M}_{22}^{ -1}\hat{M}_{21}x_{t}$     (7)

 となります。これで、ジャンプ変数 $y_{t}$ を状態変数 $x_{t}$ で表すことができました。

 次に、(7)式の t を t+1 に変えて(つまり、期間を t から t+1 にずらす)、(6)式に代入して、$y_{t+1}$ を $x_{t+1}$ で表せば、

$E_{t}y_{t+1}= - \hat{M}_{22}^{ -1}\hat{M}_{21}x_{t+1}$   (8)

 となります。今度は、これを(3)式(上半分の式)の左辺の $y_{t+1}$ に代入し、右辺の $y_{t}$ には(7)式を代入すれば、

$[\hat{M}_{11}-\hat{M}_{12}\hat{M}_{22}^{ -1}\hat{M}_{21}] x_{t+1}= \overline{\Lambda}_{11}[\hat{M}_{11}-\hat{M}_{12} \hat{M}_{22}^{ -1}\hat{M}_{21} ]x_{t}$

 となります。$x_{t+1}$ の前の行列の逆行列を左からかければ、

$x_{t+1}= [\hat{M}_{11}-\hat{M}_{12}\hat{M}_{22}^{ -1}\hat{M}_{21}]^{-1} \overline{\Lambda}_{11}[\hat{M}_{11}-\hat{M}_{12}\hat{M}_{22}^{ -1}\hat{M}_{21}]x_{t}$

 を得ます。これで、t+1 期の状態変数 $x_{t+1}$ も、t 期の状態変数 $x_{t}$ で表されたことになります。


2 ショックがある場合
 (1)式に戻って、左から B の逆行列をかければ、

$\left[\begin{array}{l}
x_{t+1}\\
E_{t} y_{t+1}
\end{array}\right]= B^{-1}A\left[\begin{array}{l}
x_{t}\\
y_{t}
\end{array}\right]+B^{-1}G\varepsilon_{t}$

 となります。BA は上記のように対角化できるので、

$\left[\begin{array}{l}
x_{t+1}\\
E_{t} y_{t+1}
\end{array}\right]= \overline{M}^{ }\overline{\Lambda} \overline{M}^{-1}\left[\begin{array}{l}
x_{t}\\
y_{t}
\end{array}\right]+B^{-1}G\varepsilon_{t}$

 となります。Mバーの逆行列を左からかければ、

$\overline{M}^{-1}\left[\begin{array}{l}
x_{t+1}\\
E_{t} y_{t+1}
\end{array}\right]= \overline{\Lambda} \overline{M}^{-1}\left[\begin{array}{l}
x_{t}\\
y_{t}
\end{array}\right]+\overline{M}^{-1}B^{-1}G$

 です。右辺の第2項を

$\left[\begin{array}{l}
\hat{G}_{1}\\
\hat{G}_{2}
\end{array}\right] = \overline{M}^{-1}B^{-1}G$

 と表しておきます。ここで、G1 と G2 の分割は、G1 の行が、Zの固有値のうち1より小さい固有値の数と等しくなるようします(上の例では、2行/2行で分割する。ただし、G1の列がその数より少ない場合)。ショックがない場合と同様に、それぞれの行列を分割します。

$\left[\begin{array}{l}
\hat{M}_{11} \hat{M}_{12}\\
\hat{M}_{21} \hat{M}_{22}
\end{array}\right]\left[\begin{array}{l}
x_{t+1}\\
E_{t} y_{t+1}
\end{array}\right]= \left[\begin{array}{l}
\overline{\Lambda}_{11}\hspace{7pt}0\\
 \hspace{2pt}0 \hspace{13pt}\overline{\Lambda}_{22}
\end{array}\right]\left[\begin{array}{l}
\hat{M}_{11} \hat{M}_{12}\\
\hat{M}_{21} \hat{M}_{22}
\end{array}\right]\left[\begin{array}{l}
x_{t}\\
y_{t}
\end{array}\right]+\left[\begin{array}{l}
\hat{G}_{1}\\
\hat{G}_{2}
\end{array}\right]\varepsilon_{t}$ 
(9)

 これを計算し、下の部分(発散する方程式)だけとりだすと、

$[\hat{M}_{21}x_{t+1}+\hat{M}_{22}E_{t}y_{t+1}]=\overline{\Lambda}_{22}[\hat{M}_{21}x_{t}+\hat{M}_{22}y_{t}]+\hat{G}_{2}\varepsilon_{t}$   (10)

 となります。計算をわかりやすくするてめに、左辺と右辺の[・・・]の中の $\hat{M}_{21}x_{t}+\hat{M}_{22}E_{t}y_{t}$ を $\lambda_{t}$ とおきます。

$\lambda_{t}=\hat{M}_{21}x_{t}+\hat{M}_{22}y_{t}$

 $\lambda_{t}$ をつかえば、(10)式は、

$E_{t}\lambda_{t+1}= \overline{\Lambda}_{22}\lambda_{t}+\hat{G}_{2}\varepsilon_{t}$

 と表されます。右辺の第2項 $\hat{G}_{2}\varepsilon_{t}$ を左辺に移行し、左から Λ22 の逆行列 Λ22-1 をかければ(つまり、$\lambda_{t}$ を表す式に変えれば)、

$\lambda_{t}=\overline{\Lambda}_{22}^{ -1}E_{t}\lambda_{t+1}- \overline{\Lambda}_{22}^{ -1}\hat{G}_{2}\varepsilon_{t}$  (11)

 となります。ここで、(11)式の期間 t をひとつ先送りし( tをt+1に変える)、右辺の $\lambda_{t+1}$ を $\lambda_{t+2}$ に変え、 $\varepsilon_{t}$ を  $\varepsilon_{t+1}$ に変え、 $\lambda_{t+1}$ を表す式をつくり、それを(11)式の右辺の  $\lambda_{t+1}$ に代入するというように、再帰的に代入を繰り返していくと、(11)式は、

$\displaystyle \lambda_{t}=-\sum_{i=0}^{\infty}\overline{\Lambda}_{22}^{ -i-1}\hat{G}_{2}E_{t}[\varepsilon_{t+i}]$   (12)

 になります。再帰的に代入していくときの右辺の最初の項、$(\overline{\Lambda}_{22}^{ -1})^{-i}\lambda_{t+1+k}$ は、i を無限大にすると、0 になります。なぜなら、$\overline{\Lambda}_{22}$ は、1より大きい要素が対角上に並んだ行列で、その逆行列(分数と考えるといい)を累乗していけば0になるからです。そのため、(11)式は(12)式のような形で表すことができるわけです。
 通常、経済学モデルでショックを考える場合、t+1 期のかく乱項 $\varepsilon_{t+1}$ の t 期から見た期待値 ( $E_{t}[\varepsilon_{t+1}]$ )は 0 になる、と想定できます(将来のショックは予測がつかない)。したがって、(11)式の総和記号の中で残るのは t 期( i=0 の場合)だけです(このように、Blanchard and Kahn method では、t+1期のかく乱項の期待値が0になる、という性質を使って連立方程式を解くので、ショックを表す方程式を連立方程式のシステムに入れておくことが必要です)。そこで、i=0 とすれば、(12)式は、

$\lambda_{t}=-\overline{\Lambda}_{22}^{ -1}\hat{G}_{2}\varepsilon_{t}$

 となります。$\lambda_{t}=\hat{M}_{21}x_{t}+\hat{M}_{22}y_{t}$ とおいているので、$\lambda_{t}$ に代入すれば、

$\hat{M}_{21}x_{t}+\hat{M}_{22}y_{t}=-\overline{\Lambda}_{22}^{ -1}\hat{G}_{2}\varepsilon_{t}$

 となります。ここから $y_{t}$ を求めれば、
$y_{t}=-\hat{M}_{22}^{ -1}\hat{M}_{21}x_{t}- \hat{M}_{22}^{ -1}\overline{\Lambda}_{22}^{ -1}\hat{G}_{2}\varepsilon_{t}$   (13)

 です。まず、$y_{t}$ を $x_{t}$ で表すことができました。

 次に、(13)式の期間をずらせば( t を t+1 に変えれば)、$y_{t+1}$ を表す式になります。このときかく乱項 $\varepsilon_{t+1}$ は消えます(t期から見た期待値は0なので)。

$E_{t}y_{t+1}=-\hat{M}_{22}^{ -1}\hat{M}_{21}x_{t+1}$   (14)

 ここで、(9)式に戻って、安定的な部分(上半分)をとりだせば、

$[\hat{M}_{11}x_{t+1}+\hat{M}_{12}E_{t}y_{t+1}]=\overline{\Lambda}_{11}[\hat{M}_{11}x_{t}+\hat{M}_{12}y_{t}]+\hat{G}_{1}\varepsilon_{t}$  (15)

 です。この(15)式の $y_{t}$,$y_{t+1}$ にそれぞれ(13)式、(14)式を代入します。

$[\hat{M}_{11}-\hat{M}_{12}\hat{M}_{22}^{ -1}\hat{M}_{21}] x_{t+1}=$

$\overline{\Lambda}_{11}[\hat{M}_{11}-\hat{M}_{12} \hat{M}_{22}^{ -1}\hat{M}_{21} ]x_{t}-[\overline{\Lambda}_{11}\hat{M}_{12}\hat{M}_{22}^{ -1}\overline{\Lambda}_{22}^{ -1}\hat{G}_{2}-\hat{ G}_{1}]\varepsilon_{t}$


 この式から $x_{t+1}$ を求めれば( $x_{t+1}$ の前の行列の逆行列を左からかける)、

$x_{t+1}=[\hat{M}_{11}-\hat{M}_{12}\hat{M}_{22}^{ -1}\hat{M}_{21}]^{-1}\overline{\Lambda}_{11}[\hat{M}_{11}-\hat{M}_{12} \hat{M}_{22}^{ -1}\hat{M}_{21} ]x_{t}$

$-[\hat{M}_{11}-\hat{M}_{12}\hat{M}_{22}^{ -1}\hat{M}_{21}]^{-1}[\overline{\Lambda}_{11}\hat{M}_{12}\hat{M}_{22}^{ -1}\overline{\Lambda}_{22}^{ -1}\hat{G}_{2}-\hat{ G}_{1}]\varepsilon_{t}$

 となります。これで、$x_{t+1}$ も $x_{t}$ で表すことができました。
 それぞれの変数の前の行列を別の文字で置き換えれば、

$y_{t}=-Nx_{t}- L\varepsilon_{t}$

$x_{t+1}=Cx_{t}+D\varepsilon_{t}$    (16)

 と表すことができます。


3 具体例(ニューケインジアンモデル)
 ABCs of RBCs 6章の Hansen model では、Bが逆行列を持たず、Blanchard and Kahn method が使えないので、Gali の Monetary Policy, Inflation, and the Business Cycle の3章のニューケインジアンモデルで計算してみます。

 金融政策当局は次のような金利ルールにしたがって、金融政策を行っています。

$i_{t}=\rho+\phi_{\pi}\pi_{t}+\phi_{y}\tilde{y}_{t}+v_{t}$  (17)

 $\rho=-\log\beta$(βは主観的時間割引因子)です。$v_{t}$ は金利ショックを表す項で(したがって、ショックがない場合は $v_{t}=0$ )、次のようなAR-1過程にしたがいます。

$v_{t}=\rho_{v}v_{t-1}+\varepsilon_{t}$   (18)

 (17)式の両辺から ρ を引き、名目金利を均衡値(ρ)からのかい離値で表せば、

$\hat{i}_{t}= \phi_{\pi}\pi_{t}+\phi_{y}\tilde{y}_{t}+v_{t}$   (19)

 となります。
 ニューケインジアンのIS曲線は、

$\displaystyle \tilde{y}_{t}= - \frac{1}{\sigma}(i_{t}-E_{t}[\pi_{t+1}]-r_{t}^{n}) + E_{t}[\tilde{y}_{t+1}]$

 です。$i_{t}$,$r_{t}^{n}$ (これは自然利子率)からそれぞれ ρ を引き、それぞれかい離値で表せば、

$\displaystyle \tilde{y}_{t}= - \frac{1}{\sigma}(\hat{i}_{t}-E_{t}[\pi_{t+1}]-\hat{r}_{t}^{n}) + E_{t}[\tilde{y}_{t+1}]$   (20)

 となります。ニューケインジアンのフィリップス曲線は、

$\pi_{t}=\beta E_{t}[\pi_{t+1}]+\kappa\tilde{y}_{t}$   (21)

 です。

 これらの方程式を解けば、未決定変数(nonpredetermined varialbe. このモデルでは、$y_{t}$ と $\pi_{t}$ )を先決変数(predetermined variable.  このモデルでは $v_{t}$ と $i_{t}$ )で表すことができます(結論から言うと、このモデルなら $v_{t}$ だけで表すことができる)。また、そのように変形すれば、ショックがあったときのインパルス応答を見ることができます。

 そこで、(18)、(19)、(20)、(21)式の連立方程式を Blanchard and Kahn method を使って解くのですが、このようなニューケインジアンモデルでは、通常、$y_{t}$ と $\pi_{t}$ が未決定変数です。行列を使って表したときに、上の一般的な説明と同じように、その二つの変数が下にくるようにしておきます。また、上の説明に書きましたが、ショックを表す方程式を入れておく必要があるので、

$x_{t+1}=\left[\begin{array}{l}
v_{t}\\
\hat{i}_{t}
\end{array}\right], \hspace{13pt}y_{t+1}=\left[\begin{array}{l}
E_{t}[\tilde{y}_{t+1}]\\
E_{t}[\pi_{t+1}]
\end{array}\right]$

 にします。
 (18)、(19)、(20)、(21)式を次のように変形します((18)式はそのまま。ここでは、ショックは金利ショックだけで、技術ショックはない、つまり、$\hat{r}_{t}^{n}=0$ と想定しているので、(20)式の $\hat{r}_{t}^{n}$ は 0 にします。)。

$v_{t}=\rho_{v}v_{t-1}+\varepsilon_{t}$

$-v_{t}+\hat{i}_{t}=\phi_{\pi}\pi_{t}+\phi_{y}\tilde{y}_{t}$

$-\displaystyle \frac{1}{\sigma}\hat{i}_{t}+E_{t}[\tilde{y}_{t+1}]+\frac{1}{\sigma}E_{t}[\pi_{t+1}]= \tilde{y}_{t}$

$\beta E_{t}[\pi_{t+1}]=-\kappa\tilde{y}_{t}+\pi_{t}$

 行列を使って表せば、

$\left[\begin{array}{l}
\hspace{4pt}1\hspace{17pt}0\hspace{17pt}0\hspace{17pt}0\\
-1\hspace{14pt}1\hspace{17pt}0\hspace{17pt}0\\
\hspace{4pt}0\hspace{4pt}-1/\sigma\hspace{8pt}1\hspace{8pt}1/\sigma\\
 \hspace{4pt}0\hspace{18pt}0\hspace{18pt}0\hspace{18pt}\beta
\end{array}\right]\left[\begin{array}{l}
v_{t}\\
\hat{i}_{t}\\
E_{t}[\tilde{y}_{t+1}]\\
E_{t}[\pi_{t+1}]
\end{array}\right]=$

$\left[\begin{array}{l}
\rho_{v}\hspace{7pt}0\hspace{10pt}0\hspace{10pt}0\\
\hspace{2pt}0\hspace{10pt}0\hspace{8pt}\phi_{y}\hspace{6pt}\phi_{\pi}\\
\hspace{2pt}0\hspace{10pt}0\hspace{10pt}1\hspace{10pt}0\\
\hspace{2pt}0\hspace{10pt}0\hspace{3pt} -\kappa\hspace{5pt}1
\end{array}\right]\left[\begin{array}{l}
v_{t-1}\\
\hat{i}_{t-1}\\
\tilde{y}_{t}\\
\pi_{t}
\end{array}\right] + \left[\begin{array}{l}
1\\
0\\
0\\
0
\end{array}\right]\varepsilon_{t}$

 となります。上の行列が B で下の行列が A です。

$B=\left[\begin{array}{l}
\hspace{4pt}1\hspace{17pt}0\hspace{17pt}0\hspace{17pt}0\\
-1\hspace{14pt}1\hspace{17pt}0\hspace{17pt}0\\
\hspace{4pt}0\hspace{4pt}-1/\sigma\hspace{8pt}1\hspace{8pt}1/\sigma\\
 \hspace{4pt}0\hspace{18pt}0\hspace{18pt}0\hspace{18pt}\beta
\end{array}\right]$

$A=\left[\begin{array}{l}
\rho_{v}\hspace{7pt}0\hspace{10pt}0\hspace{10pt}0\\
\hspace{2pt}0\hspace{10pt}0\hspace{8pt}\phi_{y}\hspace{6pt}\phi_{\pi}\\
\hspace{2pt}0\hspace{10pt}0\hspace{10pt}1\hspace{10pt}0\\
\hspace{2pt}0\hspace{10pt}0\hspace{3pt} -\kappa\hspace{5pt}1
\end{array}\right],\hspace{8pt}G=\left[\begin{array}{l}
1\\
0\\
0\\
0
\end{array}\right]$



 それぞれのパラメータ-の値は、$\sigma=1,\hspace{3pt}\beta=0.99,\hspace{3pt}\rho_{v}=0.5,\hspace{3pt}\phi_{y}=0.5/4,\hspace{3pt}\phi_{\pi}=1.5,\hspace{3pt}\kappa=0.1275$
 です。
 ここから、Blanchard and Kahn method を使って計算していきます。上の説明と同じ説明の繰り返しになりますが、まず、B-1A を対角化します。
 Scilab の場合は、
[M,D]=spec(Z);
 というコードで計算します([ , ]の中はこちらで名前を与えればよい)。
 Λ (上のコードでは D )は Z の固有値が並んだ対角行列です。その固有値は、
0                       
1.1319444 + 0.2196525i  
1.1319444 - 0.2196525i  
0.5 
 です。これを絶対値が小さい順に並び変え、新たな対角行列をつくります。
 したがって、新たな対角行列 Λ は

0    0                      0                                      0                       
0    0.5                   0                                      0                       
0    0      1.1319444 - 0.2196525i                     0                       
0    0                      0                         1.1319444 + 0.2196525i 

 です。この場合、元の Λ を、1列、4列、3列、2列の順番に並び変えています。
 そこで、行列M も同じように並びかえます。Mを1列、4列、3列、2列の順番に並び変える、ということです。元の M は、

0            0                           0                        0.8092927  
1.    0.0779335 - 0.4777554i    0.0779335 + 0.4777554i      0.3447199  
0       - 0.7786466                 - 0.7786466               - 0.4611485  
0     0.1936582 - 0.3491162i    0.1936582 + 0.3491162i    - 0.1164286

 でした。これを、1列、4列、3列、2列の順番に並び変えればいいので、新たな M は、

0     0.8092927           0                             0                       
1.    0.3447199    0.0779335 + 0.4777554i    0.0779335 - 0.4777554i  
0   - 0.4611485     - 0.7786466                  - 0.7786466               
0   - 0.1164286    0.1936582 + 0.3491162i    0.1936582 - 0.3491162i 

 となります。この新しい M の逆行列を求めます。そして、その逆行列を4つに分割します。

$\overline{M}^{-1} =\left[\begin{array}{l}
\hat{M}_{11} \hat{M}_{12}\\
\hat{M}_{21} \hat{M}_{22}
\end{array}\right]$


 そのとき、M11 の次数が、Zの固有値のうち絶対値が1より小さい固有値の数と等しくなるようにします。この例では、その数が2なので、M11 を2×2行列に分割する、ということです。Mバーの逆行列が以下のようになるので、それを、

15032601



 のように分割します。他の行列、Λ、A,B も同じように(それぞれ左上の分割行列が2×2になるようにする)分割します。

 続けて、上記の説明と同じように計算していくと、(16)式の表記の C,D,N,L が

C= 0.5              0.        
     0.2129760    0
D= 1.         
       0.4259520
N= 0.5698166    0  
       0.1438646    0
L= 1.1396333 
       0.2877292

 となります。

 したがって、

$v_{t}=0.5v_{t-1}+ \varepsilon_{t}$

$\hat{i}_{t}=0.2129769 v_{t-1}+ 0.4259520 \varepsilon_{t}$

$\tilde{y}_{t}=-0.5698166 v_{t-1}- 1.139633 \varepsilon_{t}$

$\pi_{t}=-0.1438646 v_{t-1}- 0.287792 \varepsilon_{t}$

 と表すことができます。最初の式は、もともとのショックを表す式と同じです。

Scilab のコード(Matlab では、固有値を求める[対角化する]コードを変える必要があります。spec[ ]→eig[  ]に変える)。

theta=2/3;
beta=0.99;
sigma=1;
rhov=0.5;
phiy=0.5/4;
phipi=1.5;
alpha=1/3;
eps=6;
phi=1;
theta2=(1-alpha)/(1-alpha+alpha*eps);
lambda=(1-theta)*(1-beta*theta)*theta2/theta;
kap=lambda*(sigma+(phi+alpha)/(1-alpha));

B=[1,0,0,0;
-1,1,0,0;
0,-1/sigma,1,1/sigma;
0,0,0,beta];

A=[rhov,0,0,0;
0,0,phiy,phipi;
0,0,1,0;
0,0,-kap,1];

G=[1;0;0;0];

Z=inv(B)*A;
[R,D1]=spec(Z);
v1=spec(Z);
v2=[v1(1),v1(4),v1(3),v1(2)];
D2=diag(v2);
M=[R(1:4,1),R(1:4,4),R(1:4,3),R(1:4,2)];
M1=inv(M);
M11=M1(1:2,1:2);
M12=M1(1:2,3:4);
M21=M1(3:4,1:2);
M22=M1(3:4,3:4);
D11=D2(1:2,1:2);
D22=D2(3:4,3:4);

N=inv(M22)*M21;
C=inv(M11-M12*inv(M22)*M21)*D11*(M11-M12*inv(M22)*M21);

GG=inv(M)*inv(B)*G;
G1=GG(1:2);
G2=GG(3:4);

L=inv(M22)*inv(D22)*G2;

D=-inv(M11-M12*inv(M22)*M21)*(D11*M12*inv(M22)*inv(D22)*G2-G1);


プロフィール

沢ひかる

貧乏人。

経済学「部」とは無縁です。

QRコード
QRコード
  • ライブドアブログ