統計数理（2002）
第 50 巻 第 2 号 133–147
c 2002 統計数理研究所


特集「ファイナンス統計学」［研究詳解］

モンテカルロフィルタを用いた金利モデルの推定1

高橋 明彦 ・佐藤 整尚2（受付 2002 年 2 月 15 日；改訂 2002 年 7 月 5 日）

要旨

金利の期間構造の実証分析に関し，一般化状態空間モデルに基づくモンテカルロフィルタを用いた新しい推定方法を紹介する．さらに，日本国債市場の時系列データを用いた例を示した．
キーワード： 一般化状態空間モデル，モンテカルロフィルタ，金利モデル，アフィン型期間構造モデル，日本国債，スワップ金利．
1.

はじめに

本稿では，Takahashi and Sato（2001）により提案された金利の期間構造の実証分析に対する，一般化状態空間モデルに基づく新しい推定方法の枠組みを紹介する．特に，拡散過程によりその動態が記述される状態変数を基礎とした一般的な期間構造モデルの推定に対し，モンテカル
ロフィルタの具体的な適用法を新しい実証例を交えて解説する．
金利の期間構造モデルとは，複数の異なる満期の金利を共通の要素で統一的に記述すること等を通じ，金利変動を相互関連に重点を置いて説明しようとするモデルである．また，実務においては，債券及びスワップ，債券オプション，キャップ・フロア，スワップション等の金利の派生商品取引の値付け及びリスク管理の基盤として金利の期間構造モデルは重要である．代
表的な金利の期間構造モデルにおいては，まず状態変数のダイナミクスが（多次元を含む）
マル
コフ型確率微分方程式で表現され，瞬時的短期金利（the instantaneous short rate）
がそれらの状
態変数の関数で表される．さらに，無裁定条件及び一般均衡理論等に基づき，割引債価格従って利付き債券価格，スワップレートが確率微分方程式のパラメータ，リスクの市場価格に関するパラメータを含む状態変数の関数で表される．一方，現実においては，いくつかの満期の利付き債券価格，スワップレートが観測されるので，金利の期間構造モデルの推定に対する一つのアプローチは，利付き債券価格，スワップレートの時系列データを所与として仮定されたモデルの状態変数及びパラメータを推定することになる．
本稿では，金利の期間構造モデルが一般化状態空間モデルによって表現できることを示し，その例としてアフィン型期間構造モデル（ATSM）に属する 2 ファクターモデルを用いた国債金利への適用例を取り上げる．また，Takahashi and Sato（2001）でその詳細が示されている，Hulland White（1994, 1997）のクラスに属する割引債価格は陽に表現できないモデルに対する適用の概要も紹介する．

1 東京大学 大学院数理科学研究科：〒153–8914 東京都目黒区駒場 3–8–1; QZD12245@niftyserve.or.jp
2 統計数理研究所：〒106–8569 東京都港区南麻布 4–6–7; sato@ism.ac.jp

134

2.

統計数理

第 50 巻

第 2 号 2002

先行研究

市場価格の観測データから推定可能な金利の期間構造モデルの理論的枠組みは，Cox, Ingersoll（以下，CIR と略す）
により動学的一般均衡理論に基づき最初に与えらand Ross（1985a, 1985b）
れた．また，実証分析においては，単一の状態変数では期間構造の変動を説明するには不十分であり，通常複数の状態変数のモデルが用いられている（瞬時的短期金利が相互に独立な複数の状
態変数の和として表現される，所謂，マルチファクター CIR モデル Chen and Scott（1993），
Pearson and Sun（1994），Singh（1995），Duﬃe and Singleton（1997）
等）．特に，各満期の割引債

金利及び瞬時的短期金利が状態変数に関しアフィン関数（h : RN → R，h(y) = α + β y ，y ∈ RN
は状態変数，α ∈ R，β ∈ RN は定数）で表される場合アフィン型期間構造モデル（ATSM）と呼ばれ（Duﬃe（1996），Chapter 7），Dai and Singleton（2000）は，ATSM を統一的に考察し，このクラスに当てはまる既存のモデルを分類すると共に，米国のスワップレートの時系列データを用いて代表的なモデルを検証している．但し，観測される金利が名目金利であるということから，金利の期間構造モデルの必要要件として，金利が負にならない（瞬時的短期金利の非負条件）という条件を満たす必要がある．また，主に長期金利のレベルを表す状態変数と短期
金利と長期金利のスプレッドを表す状態変数の間において，負の相関が存在するという結果
が円，US スワップレート市場を対象とした実証分析等により報告されている（Takahashi and．この双方の性質を満たすことが保障できない点にSato（2001），Dai and Singleton（2000）等）
上記のクラスのモデルの限界がある（Dai and Singleton（2000），p. 1956）．この点は，ゼロ金利状態の続く円金利市場の期間構造モデルの構築において特に問題になろう．アフィン型期間
構造モデル（ATSM）は，Duﬃe and Kan（1996）により理論的に考察され，このクラスに含まれるモデルには，正規過程型モデル（Vasicek（1977），Langetieg（1980），Hull and White（1990,，CIR 型モデル（Cox et al.（1985b），Chen and Scott（1993），Pearson and Sun（1994），
1994））
Singh（1995），Duﬃe and Singleton（1997）），Longstaﬀ and Schwartz（1992）の 2 ファクターモ
デル，Chen（1996），Balduzzi et al.（1996）の 3 ファクターモデルなどがある．このクラスに属さない主なモデルには，Brennan and Schwartz（1982），Black and Karasinski（1991）の他，Hulland White（1994, 1997）で紹介されている瞬時的短期金利が，正規過程で表現される状態変数の非線形変換として表されるモデルがある．これらのモデルの中には，割引債価格は解析的に求まらないが，瞬時的短期金利の非負条件と状態変数間の負の相関の両立が可能となるモデルが存在する．
既存研究における実証分析手法は，一般化モーメント法（GMM（）Chan et al.（1992），Longstaﬀ），Simulated Moment Method（SMM（）Dai and Singleton（2000）），主成分分and Schwartz（1992）
析と三段階最小二乗法（three-stage least square method with the principal component analysis）（Singh（1995）），最尤推定法（Chen and Scott（1993），Pearson and Sun（1994）
および Duﬃe andなどであるが，割引債価格が簡単にもとめられない場合，これらの手法を適Singleton（1997））
用するのは極めて難しいことには注意を要する．
また，従来の実証研究では，瞬時的短期金利に代表される観測不可能な状態変数が観測可能な変数により代用されている場合が多く（Chan et al.（1992），Longstaﬀ and Schwartz（1992）
等），観測誤差（measurement error）の影響を受けている可能性が多分にある．さらに，観測誤差を明示的に考慮していても，その導入の仕方が不自然（Chen and Scott（1993））
或いは恣意的（Duﬃe and Singleton（1997））である場合も少なくない．本稿では，先行研究のこれらの問題点を克服するため，推定方法として一般化状態空間モデルに基づくモンテカルロフィルタを活用する．一般化状態空間モデルは，直接観測不可能な状態変数の変動を記述するシステムモデルと，観測データとして得られる債券価格，スワップ金利，期間金利と状態変数の関数関係を記

モンテカルロフィルタを用いた金利モデルの推定135

述する観測モデルにより構成され，状態変数を基礎とする金利の期間構造モデルの実証分析に
対し自然に適用することができる．観測誤差に関しても，観測モデルにおいて一般的にかつ，先入観を入れないかたちで導入することが可能である．また，モンテカルロフィルタを適用できる期間構造モデルの範囲は，先行研究で利用されている各手法に比べて極めて広く，割引債
価格が解析的に陽にもとまらない場合でも適用可能である．
3.

金利モデルの状態空間表現

本章では，状態空間表現に基づく金利の期間構造モデルを提示し，そのパラメータ推定と共通成分の抽出について説明する．
まず，マルコフ過程に基づく期間構造モデルについて簡単に説明しておく（より詳しくは Björk（1996），Duﬃe（1996）等を参照されたい）．ある T ∗ < ∞ に対する期間 [0, T ∗ ] においてフィル
タ付き確率空間 (Ω, F, {Ft }, P ) が与えられているとする．N 次元の状態ベクトル Yt が以下のような N 次元のマルコフ過程に従うとする．
(3.1)

dYt = µ(Yt , t)dt + S(Yt , t)dBt ,ただし，Bt は d 次元標準ブラウン運動であり，µ(Yt , t) と S(Yt , t) はそれぞれ RN ×[0, T ∗ ] → RN ，
RN × [0, T ∗ ] → RN×d となるような実数値関数である．さらに，t 時点における瞬間的な短期金利 r(Yt , t) と満期 T の t 時点における割引債の価格 P (Yt , t; T ) が，いずれも Yt の関数で表されるとする．ただし，t ∈ [0, T ∗ ]，T ∈ [t, T ∗ ] である．ここでいう割引債とはクーポンがなく，満期時点において表面価格 1 が支払われるものである．ここで仮定した割引債の価格の集合{P (Yt , t; T )}T ∈[t,T ∗ ] が t における金利の期間構造を示すという点を注意しておく．ここでの仮定は比較的少数の共通因子，Yt によって期間構造全体を説明できるというものであり，実際の共通因子が具体的にどのようなものであるかはモデルの立て方に依存する．
よって，ファイナンスにおける無裁定に関する議論をもとにすれば，P (Yt , t; T ) は以下の偏微分方程式
1
trace(SS  PY Y ) + [µ − φ] PY + Pt − rP = 0 ,22P∂Pを境界条件，P (Yt , T ; T ) = 1, のもとでの解となる．ただし，PY Y ≡ ∂Y∂ ∂Yであ , PY ≡ ∂Yる．ここで，RN ベクトル φ はいわゆる，リスクプレミアムとよばれるもので，Yt と t の関数
φ(Yt , t) であるとする．この関数は，例えば Cox et al.（1985a）により提案された資産価格の一
般均衡理論によって決定されるものである．さらに，よく知られているように，この確率微分方程式の解は適当な数学的正則条件の下で t 時点において与えられる情報で条件づけられた期待値によってあらわされる．
(3.2)

P (Yt , t; T ) = Q [e−ÊTtr(Yu ,u)du|Ft ] ,



ただし， Q [·|Ft ] はリスク中立確率測度 Q の下での t における情報で条件づけられた条件付期
待値をあらわす（詳細は Björk（1996）を見よ）．また，この確率測度 Q のもとで，状態ベクトル
Yt は
(3.3)

dYt = {µ(Yt , t) − φ(Yt , t)}dt + S(Yt , t)dBt∗という確率微分方程式に従うことが知られている．ここで，B ∗ は確率測度 Q のもとでの d 次
元標準ブラウン運動である（詳しくは Björk（1996），または Duﬃe（1996）の Chapter 7 およびAppendix E を参照のこと）
．

136

統計数理

第 50 巻

第 2 号 2002

さて，ここで，期間構造モデルに実際の観測データをあてはめるために，一般化状態空間モデルを導入することにする（一般化状態空間モデルについては Kitagawa and Gersh（1996）を参照）．ここでの状態空間モデルは以下のようなシステムモデルと観測モデルからなる．



(3.4)

Yt =
Zt =

 (Yt−∆t, vt ) システムモデル
 (Yt, ut ) 観測モデル

ただし，Yt ，Zt はそれぞれ，N 次元の状態ベクトル，M 次元の観測ベクトルをあらわす．また，∆t は観測の時間間隔である．さらに，vt ，ut はそれぞれ，N 次元のシステムノイズ，M 次元の観測ノイズであり，その密度はそれぞれ，q(v)，ψ(u) であらわされるとする．一般に (·, ·)と (·, ·) は非線形関数であり，状態の初期値ベクトル Y0 は密度関数 p0 (Y ) に従う確率変数である．
次にここで導入した状態空間モデルがどのようにして期間構造モデルの推定に適用されるかをみていく．∆t が十分に小さいとき，（3.1）式に対するオイラー近似によって，システムモデル
Yt = (Yt−∆t , vt ) が導出される．具体的には√(3.5)
Yt = Yt−∆t + µ(Yt−∆t , t − ∆t)∆t + S(Yt−∆t , t − ∆t)vt ∆t







となり，システムノイズ vt は N 次元の標準正規分布に従う．もちろん，（3.1）式の離散化については他の手法を適用することも可能である（たとえば，Shoji and Ozaki（1998）や Shoji（1998）
で提案されている方法もある）．さらに，線形確率微分方程式の場合のように Yt−∆t によってYt が解析的に解ける場合は，その表現を使うことも可能である．たとえば，（3.1）
式において Ytが
(3.6)

dYt = (AYt + β ∗ (t))dt + SdBtというような線形確率微分方程式に従う場合を考える．ただし，β ∗ (t) と A は，それぞれ，RN
値をとるような時間 t に依存する関数と N × N の定数行列であり，S は Σ = SS  が正定値符号となるような N × d の定数行列である．このとき，Yt−∆t があたえられたもとで Yt は
(3.7)

Yt = eA∆t Yt−∆t +

t

(∆t)

e(t−s)A β ∗ (s)ds + vt

t−∆t
(∆t)

= F Yt−∆t + β(t) + vt

teA(t−s) β ∗ (s)ds は
t−∆t
(∆t)
N × 1 のベクトル値をとる時間 t の関数である．さらに，vt
は平均 0 で以下のような分散共と表される．ただし，F ≡ e∆tA は N × N の定数行列であり，β(t) ≡分散行列 Σ∆t をもつような正規分布に従う．
(3.8)

Σ∆t ≡

 ∆t



esA ΣesA ds .

0

この場合，システムモデルは（3.7）式で与えられ，システムノイズ vt の密度関数 q(v) は（3.8）式で示されたような分散共分散行列をもつ正規分布になる．
一方，観測モデルにおいて，t 時点の観測ベクトル Zt は k(≥ 1) 種類の割引債の価格系列と観測ノイズベクトル ut の関数として表される．
(3.9)

Zt = h(P (Yt , t; t + T1 ) , . . . , P (Yt , t; t + Tk )) + ut

つまり，Zt の各要素を構成する観測された債券価格や金利は，異なる満期 (Ti , i = 1, . . . , k) の

モンテカルロフィルタを用いた金利モデルの推定137

割引債の価格と観測誤差によって表されると考える．さらに，おのおのの P (Yt , t; Ti ) が Yt の関数であることから，Zt も Yt の関数として表されることになる．
(3.10)

Zt = H(Yt ) + ut .

（3.3）のプロセスのもとで（3.2）式を計算することによって得h(·) において，P (Y, t; t + Ti ) はられる．また，ここで，誤差項 u の密度関数 ψ(u) が平均 0，分散共分散行列 Σu であるような多変量正規分布に従うと仮定する．ただし，Σu は M × M の対角行列でその対角要素は正であるとする．ここで，M は各時点の観測データの次元である．
t 時点の割引債金利および LIBOR（London Inter Bank Oﬀered Rates）
やスワップ金利は観測ベクトル Zt の典型的な例である．LIBOR は欧州通貨市場における銀行間預金のオファーレー
トであり，短期市場金利の代表的な指標である．また，ここで言うスワップとは金利スワップのことで，ある一定の期間にわたって，ある同額の想定元本に対する固定金利支払と変動金利支払とを，相対で交換するという契約をさす．通常，変動金利には LIBOR の金利が想定される．さて，この場合，h(·) は割引債金利はもちろん，LIBOR やスワップ金利と割引債の価格との理論的な関係により定式化される．つまり，満期が τn の割引債金利 (Rt (Yt , τn ))，LIBOR(Lt (Yt , τn )) やスワップ金利 (St (Yt , τn )) は以下のように割引債の価格によってあらわすことができる（詳しくは Björk（1996），Duﬃe（1996），Duﬃe and Singleton（1997）などを参照）．
(3.11)Rt (Yt , τn ) = −1
log P (Yt , t; t + τn )
τn(3.12)Lt (Yt , τn ) =1−1
P (Yt , t; t + τn )(3.13)St (Yt , τn ) =


δ

1
τn

1 − P (Yt , t; t + τn ),τ /δP (Y , t; t + iδ)
n

i=1

t

ただし，δ はキャッシュフローの発生間隔であり，たとえば，日本の円マーケットでは δ = 0.5（半年）とするのが通例である．
ここまでの話から，（3.1）（3.3）式で表される金利の期間構造モデルは（3.4）式の一般化状態空間モデルによって書き直されることが分かった．したがって，この状態空間モデルによって状態推定を行うことで期間構造モデルの Yt が推定できるのである．
4.

モンテカルロフィルタ

この節では具体的な推定手法について説明する．前節で説明された一般化状態空間モデル
は線形ガウスモデルであるとは限らないので，標準的なカルマンフィルタは適用できない．そこで，一般化状態空間モデルの推定方法として，モンテカルロフィルタを活用する．モンテカ
ルロフィルタに関してはいくつかの方法が提案されている（Doucet et al.（1995），Durbin and，Gordon et al.（1993），Tanizaki（1993）
等）
が本稿では Kitagawa（1996）の方法
Koopman（1997）に従う．その理論の詳細は Kitagawa（1996）に譲り，以下において，本稿で利用するモンテカル
ロフィルタの具体的なアルゴリズムの概要を説明する．まず Kitagawa（1996）の記述に従って，記号を整理する．∆t をデータ間隔とすると p(Yt |Zt−∆t ) は，Zt−∆t を与えたときの Yt の条件付
き密度関数で，「一期先予測」と呼ばれ，p(Yt|Zt ) は，Zt を与えたときの Yt の条件付き密度関数
(1)
(m)
(1)
(m)
を表し，
「フィルタ」と呼ばれる．{pt , . . . , pt }, {ft , . . . , ft } はそれぞれ p(Yt |Zt−∆t ) 及び(1)(m)
p(Yt |Zt ) に関する m 個のモンテカルロ試行の実現値ベクトルを表す．さらに，{f0 , . . . , f0 }を初期状態ベクトル Y0 の密度関数 p0 (Y ) に従う乱数により生成される実現値ベクトルとする

138

第 50 巻統計数理

第 2 号 2002と，アルゴリズムは以下のように与えられる．
[モンテカルロフィルタのアルゴリズム概要](1)(m)
1. 状態変数ベクトルの初期値ベクトル {f0 , . . . , f0 } を正規分布，N (µf0 , Σf0 ) に従うと仮定して発生させる．但し，Σf0 は対角要素が正であるような対角行列とする．
2. 以下のステップ（a）〜（d）を各時点 t = 0, ∆t, 2∆t, . . . , (T ∗ − ∆t), T ∗ に対して適用する．
∗
但し，T はデータの最終時点を表す．
・（a）システムノイズ vt(j) , j = 1, . . . , m を正規分布 N (0, Σ∆t ) (システムモデルが（3.7）式の場合），または N (0, I)（（3.5）式の場合）に従って発生させる．
・（b）
各 j = 1, . . . , m に対して
(j)

pt

=

(j)
(j)
 (ft−∆t, vt )

を計算する． の具体的な関数形は（3.7）式，または（3.5）式になる．
・（c）N (0, Σu ) の密度関数 n[x; 0, Σu ] を xt = Zt − H(p(j)
t ), j = 1, . . . , m において評価し，これを α(j)
,
j
=
1,
.
.
.
,
m
とする．ここでのH(·)
は，状態変数Y の関数として表現された
t（3.11）式，（3.12）式及び（3.13）式である．これらの式における割引債価格は（3.3）のプロセ
スのもとで（3.2）式を評価することによって得られるが，後に挙げる例のように解析的に解けない場合は，モンテカルロ・シミュレーションにより評価する必要がある．
(1)
(m)
(1)
(m)
・（d）
{ft , . . . , ft } を {pt , . . . , pt } から再抽出（リサンプリング）する．具体的には，各
(i)
(1)
(m)
ft , i = 1, . . . , m を {pt , . . . , pt } の中から，(i)Prob.(ft(j)= pt |Zt ) =(j)mαt α(j)j=1

i = 1, . . . , mtの確率に基づき再抽出する．
次に未知パラメータの推定は最尤推定法により行う．まず，θ を未知パラメータを総称したベクトルとすると，尤度 L(θ) は p(Z∆t |Z0 ) = p0 (Z∆t ) として，


T∗
∆t

L(θ) = p(Z∆t , . . . , ZT ∗ |θ) =p(Zk∆t |Z∆t , . . . , Z(k−1)∆t , θ)k=1と表される． この対数尤度 l(θ) は，モンテカルロフィルタの枠組みでは，


T∗
∆t

l(θ) =

k=1

log

m

(j)

αk∆t

j=1

−

T∗
log m
∆t

により近似的に計算され，これを θ に関して数値的最適化により最大化することで最尤推定量，θ̂ がもとめられる．最大化においてはグリッドサーチや自己組織化の手法が利用可能である（自己組織化については Kitagawa（1998）を参照）．最後に，複数のモデルの候補を比較する
場合は，モデルの あてはまり の評価基準として AIC（Akaike’s Information Criterion）を使うことができる（Akaike（1973））． これは
AIC = −2l(θ̂) + 2(the number of parameters)により容易に計算でき，この値が小さいほど相対的に良いモデルであると判断することができる．

モンテカルロフィルタを用いた金利モデルの推定

5.

139

例

（例 1）
は，Hull and White
本節では，モンテカルロフィルタの適用例を 2 つ紹介する．1 つ目（1994, 1997）のクラスに属し，状態変数が正規過程で記述され瞬時的短期金利がその非線形変換で表現されるモデルである．この場合，瞬時的短期金利の非負条件と状態変数間の負の相関が両立するが，割引債が解析的に解けないモデルとなる．円 LIBOR，スワップ市場を対象にした実証分析の詳細は Takahashi and Sato（2001）にあるので，本稿ではモンテカルロフィルタ
適用法の概要のみ説明する．2 つ目の例（例 2）は，アフィン型期間構造モデル（ATSM）に属す
る 2 ファクターモデルで，割引債価格が簡単な数値計算を伴うがほぼ解析的にもとめられるモデルである．このモデルに対しては，日本国債市場のデータを用いた適用例を紹介する．
例 1. 本例でとりあげる期間構造モデルは，状態変数が Hull and White（1994, 1997）に基づ
く 2 ファクターモデルで，その状態空間表現は以下のように書ける．
（システムモデル）
(∆t)

(5.1)

Yt = F Yt−∆t + β(t) + vt,

ただし，Yt = (Y1t , Y2t ) である．
（観測モデル）
Zt = (Z1,t , . . . , Z8,t )

(5.2)


(5.3)

Zn,t =(n = 1, 2 の場合）
(n = 3, . . . , 8 の場合) ,Lt (Yt , τn ) + un,t
St (Yt , τn ) + un,t

ただし，この例で使われた観測データは LIBOR の半年，1 年とスワップ金利の 2 年，3 年，4年，5 年，7 年，10 年の 8 次元データであるので，M = 8, τ1 = 0.5, τ2 = 1, τ3 = 2, τ4 = 3, τ5 =4, τ6 = 5, τ7 = 7, τ8 = 10 である．ここで，Lt (Yt , Tn ) と St (Yt , Tn ) の評価に必要な割引債価格
P (Yt , t, T ) は(5.4)

P (Yt , t; T ) = Q [e−ÊTtg(Yu )du|Yt ]と表せるが，g(·) が非線形関数の場合は解析的に求めることができないので，以下のようにモンテカルロシュミレーションによって求める（この例では g(·) は非線形な正の値をとる関数を
用いた）．
(i)

P (Yt , t; T ) 

(5.5)

J


1
exp
J j=1



T /∆t

−

(i ,j)

g(Yt+l∆t )∆t

,

l=0

(i)

(i,j)

ただし，Yt はモンテカルロフィルタにおける Yt の i 番目の粒子をしめす．また，Yt+l∆t は(i)Yt を初期点として
(i,j)

(i,j)

Yt+l∆t = F Yt+(l−1)∆t + β ∗ (t + l∆t) + v (i,j) (t + l∆t) ,
(i,j)

Yt

(i)

= Yt

により発生させた Yt+l∆t の j 番目のパスをしめす．なお，この式は一般論での（3.3）式に対応するものである．

140

統計数理

第 50 巻

第 2 号 2002

このような手法を使うことによって，明示的に割引債価格が求まらないようなモデルでも推定が可能となる．
例 2. 本例でとりあげる期間構造モデルは，Dai and Singleton（2000）のアフィン型期間構造モデルで，一般形は以下のように表される．
（理論モデル）√dYt = K(Θ − Yt )dt + Σ St dBt ,(5.6)

ただし，Σ は N × N 行列，S(t) はその対角要素が(S(t))ii = αi + γi Ytで表されるような対角行列である．さらに，金利モデルは以下のような式で表される．
(5.7)

rt = δ0 +N

δi Yit .

i=1

次に測度 Q のもとでは(5.8)√
dYt = K ∗ (Θ∗ − Yt )dt + Σ St dBt∗となる．ただし，K ∗ = K + ΣΦΘ∗ = K ∗−1 (KΘ − Σψ)Φ = (λ1 γ1 , λ2 γ2 , . . . , λN γN )
ψ = (λ1 α1 , λ2 α2 , . . . , λN αN )

である．ここでは，リスクの市場価格（market prices of risk）Λ(t) を以下のように仮定している．
Λ(t) =

S(t)λλ = (λ1 , λ2 , . . . , λN ) .

このときに，割引債価格は
(5.9)



P (Yt , t; t + τ ) = eb0(τ )−B(τ ) Ytと表され，その係数 b0(τ ), B(τ ) は以下の微分方程式を初期条件 b0(0) = 0, B(0) = 0 により解いた解として与えられる．
N





(5.10)1db0(τ )
2= −θ ∗ K ∗ B(τ ) +Σ  B(τ ) i αi − δ0 ,
dτ
2 i=1

(5.11)1dB(τ )
2= −K ∗ B(τ ) −Σ  B(τ ) i γi + δy .
dτ
2 i=1

N





ただし，δy = (δ1 , δ2 , . . . , δN ) ．上記の微分方程式は簡単な数値計算で各値を求めることがで
きる．

141

モンテカルロフィルタを用いた金利モデルの推定

今回，ここで扱う具体的な形は以下の 2 ファクターのケースである．
(5.12)

Yt = (Y1t , Y2t ) ,(5.13)

Σ=

1
σ21

0
1

k11
0

K=

,

√

0
k22
√

St =

Θ = (θ1 , 0) ,,β1 Y1t0√0α2 + β2 Y1t

rt = δ0 + Y1t + Y2t .

(5.14)

また，リスクの市場価格（the market price of risk）
に関しては，λ が 0 であると仮定した．このとき，割引債金利は以下のように与えられる．
1(5.15)Rt (Yt , τ ) = − log P (Yt , t; t + τ ) = −b0 (τ ) + b1 (τ )Y1t + b2 (τ )Y2t .
τ
ここで，b0(τ ), B(τ ) = (b1(τ ), b2(τ )) は（5.10），（5.11）式の微分方程式を解くことにより求められ，一部，数値計算が必要である．次に，（5.6）式を離散化し，状態空間モデルにおけるシステ
ムモデルを導出する．
（システムモデル）(5.16)

Y1,n+1 = Y1,n + k11 (θ1 − Y1,n )∆t +(5.17)

Y2,n+1 = Y2,n − k22 Y2,n ∆t + σ21

β1 Y1,n v1,nβ1 Y1,n v1,n +

α2 + β2 Y1,n v2,n .

ただし，v1,n , v2,n ∼ i.i.d. N (0, ∆t) であり，n は離散化された時点を示す．
実証分析に用いる日本国債のデータは 1996 年 5 月 29 日から 2000 年 6 月 20 日までの 1001
時点の日次データで，各時点の観測値 Z1,n , . . . , Z6,n は利付き国債価格より変換して得られた
1 年，2 年，5 年，7 年，10 年，20 年満期の割引債金利である．従って，観測モデルは以下の式で表現される．
（観測モデル）
Z(5.18)c(τ )Z6,n(5.19)

H (τ )

H (τ )1,n  1   1 1 2 1 Z2,n   c(τ2 )   H1 (τ2 ) H2 (τ2 ) 


 
 




Z3,n 
c(τ3 ) 
H1 (τ3 ) H2 (τ3 ) 






=
+

Z4,n 
c(τ4 ) 
H1 (τ4 ) H2 (τ4 ) 






 Z5,n   c(τ5)   H1(τ5 ) H2(τ5) c(τ6 )

τ1 = 1 ,τ2 = 2 ,H1 (τ6 )τ3 = 5 ,

u

Y1,n
Y2,n

1,n


u2,n 

u3,n 
+
u4,n  u5,n 

H2 (τ6 )τ4 = 7 ,τ5 = 10 ,u6,nτ6 = 20c(τ ) = −b0 (τ, k11 , k22 , θ1 , β1 , β2 , σ21 , α2 , δ0 )(5.20)

H1 (τ ) = b1 (τ, k11 , k22 , β1 , β2 , σ21 )H2 (τ ) = b2 (τ, k22 ) .

最後に，上記のモデルをモンテカルロフィルタによって推定した結果を示す．

142

統計数理

図 1.

第 50 巻

第 2 号 2002

日本国債の観測値（実線）とモンテカルロフィルタによる予測値（点線）．

モンテカルロフィルタを用いた金利モデルの推定

図 2.

予測誤差．

143

144

統計数理

第 50 巻

第 2 号 2002

図 3. 推定されたファクターの系列．

図 1 は，割引債金利とモデルによる一期先予測，図 2 はその予測誤差を示している．これ
らをみると，満期 20 年の割引債金利に対する予測誤差が他の満期のものに比べて，やや大きいという結果になっている．また，図 3 は推定されたファクターの時系列を示している．第 1
ファクター Y1 はほぼ 10，20 年物などの長期金利の動きと類似しており，第 1 ファクターと第2 ファクターの和（ Y1 + Y2 ）の変動は，ほぼ 1 年物の金利の変動に対応している．
6.

まとめ

本稿では，金利の期間構造モデルを一般化状態空間で表現し，状態変数やパラメータの推定法としてモンテカルロフィルタを用いる方法を紹介すると共に，日本国債市場を対象とした適用例を示した．
謝辞

本論文をまとめるにあたって電気通信大学の宮崎浩一講師に大変お世話になった．また，査読者の方からは貴重なコメントを頂いた．ここに謝辞を述べたい．

モンテカルロフィルタを用いた金利モデルの推定145

参 考 文 献Akaike, H.（1973）. Information theory and an extension of the maximum likelihood principle, Second, 267–281,
International Symposium on Information Theory（eds. B. N. Petrov and F. Csáki）
Akademiai Kiádo, Budapest.
Balduzzi, P., Das, S. R., Foresi, S. and Sundaram, R.（1996）. A simple approach to three factor aﬃneterm structure models, Journal of Fixed Income, 6, 43–53.
. Interest rate theory, Financial Mathematics Bressanone（ed. W. J. Runggaldier）,
Björk, T.（1996）
Springer, Berlin.
. Bond and option pricing when short rates are log-normal,
Black, F. and Karasinski, P.（1991）
Financial Analysts Journal, July-August, 52–59.
. An equilibrium model of bond pricing and a test of market
Brennan, M. J. and Schwartz, E. S.（1982）
eﬃciency, Journal of Financial and Quantitative Analysis, 17（3）, 301–329.
Chan, K. C., Karolyi, G. A., Longstaﬀ, F. A. and Sanders, A. B.（1992）. An empirical comparison ofalternative models of the short-term interest rate, Journal of Finance, 47, 1209–1227.
. Stochastic Mean and Stochastic Volatility—A Three-Factor Model of the Term
Chen, L.（1996）
Structure of Interest Rates and Its Application to the Pricing of Interest Rate Derivatives,
Blackwell, Oxford, U.K.
Chen, R. and Scott, L.（1993）. Maximum likelihood estimation for a multifactor equilibrium modelof the term structure of interest rates, Journal of Fixed Income, December, 14–31.
Cox, J. C., Ingersoll, J. E. and Ross, S. A.（1985a）. An intertemporal general equilibrium model asset
prices, Econometrica, 53, 363–384.
Cox, J. C., Ingersoll, J. E. and Ross, S. A.（1985b）. A theory of the term structure of interest rates,
Econometrica, 53, 385–408.
Dai, Q. and Singleton, K. J.（2000）. Speciﬁcation analysis of aﬃne term structure models, Journalof Finance, LV（5）, 1942–1978.
Doucet, A., Barat, E. and Duvaut, P.（1995）. A Monte Carlo approach to recursive Bayesian state
estimation, Proceedings IEEE Signal Processing/Athos Workshop on Higher Order Statistics,
Girona, Spain.
Duﬃe, D.（1996）. Dynamic Asset Pricing Theory, 2nd ed., Princeton University Press, Princeton,
New Jersey.
Duﬃe, D. and Kan, R.（1996）. A yield-factor model of interest rates, Math. Finance, 6, 379–406.
. An econometric model of the term structure of interest-rateDuﬃe, D. and Singleton, K. J.（1997）
swap yields, Journal of Finance, 52, 1287–1321.
Durbin, J. and Koopman, S. J.（1997）. Monte Carlo maximum likelihood estimation for non-Gaussianstate space models, Biometrika, 84, 669–684.
Gordon, N., Salmond, D. J. and Smith, A.F.M.（1993）. Novel approach to nonlinear/non-Gaussian
Bayesian state estimation, IEE Proceedings-F, 140（2）, 107–113.
. Pricing interest rate derivative securities, The Review of Financial
Hull, J. and White, A.（1990）
Studies, 3（4）, 573–592.
. Numerical procedures for implementing term structure models II:
Hull, J. and White, A.（1994）
Two-factor models, Journal of Derivatives, 2, 37–48.
Hull, J. and White, A.（1997）. Taking rates to the limits, Risk, December, 168–169.
Kitagawa, G.（1996）. Monte Carlo ﬁlter and smoother for non-Gaussian nonlinear state space models,J. Comput. Graph. Statist., 5（1）, 1–25.
Kitagawa, G.（1998）. A self-organizing state-space model, J. Amer. Statist. Assoc., 93, 1203–1215.

146

統計数理

第 50 巻

第 2 号 2002

. Smoothness Prior Analysis of Time Series, Lecture Notes in
Kitagawa, G. and Gersch, W.（1996）
Statist., No. 116, Springer, Berlin.
Langetieg, T. C.（1980）. A multivariate model of the term structure, Journal of Finance, 35, 71–91.
. Interest rate volatility and the term structure: A two-factor
Longstaﬀ, F. and Schwartz, E. S.（1992）
general equilibrium model, Journal of Finance, 47, 1259–1282.
Pearson, N. D. and Sun, T. S.（1994）. Exploiting the conditional density in estimating the term
structure: An application to the Cox, Ingersoll, and Ross model, Journal of Finance, 49,
1279–1304.
Shoji, I.（1998）. Approximation of continuous time stochastic processes by a local linearization
method, Math. Comp., 67, 287–298.
Shoji, I. and Ozaki, T.（1998）. A statistical method of estimation and simulation for systems ofstochastic diﬀerential equations, Biometrika, 85, 240–243.
Singh, M. K.（1995）. Estimation of multifactor Cox, Ingersoll, and Ross term structure model, Journalof Fixed Income, September, 8–28.
Takahashi, A. and Sato, S.（2001）. Monte Carlo ﬁltering approach for estimating the term structureof interest rates, Ann. Inst. Statist. Math., 53（1）, 50–62.
. Nonlinear Filters, Lecture Notes in Econom. and Math. Systems, No. 400,
Tanizaki, H.（1993）
Springer, Berlin.
Vasicek, O. A.（1977）. An equilibrium characterization of the term structure, Journal of Financial
Economics, 5, 177–188.

Proceedings of the Institute of Statistical Mathematics Vol. 50, No. 2, 133–147 (2002)

147

An Application of Monte Carlo Filter for Estimating
the Term Structure of Interest Rates
Akihiko Takahashi(Department of Mathematical Science, University of Tokyo)

Seisho Sato
(The Institute of Statistical Mathematics)

We have developed a new methodology for estimating the general class of term structure models based on a Monte Carlo ﬁltering approach. We utilize the generalized state
space model, which can be naturally applied to the estimation of term structure models
based on Markovian processes. It is also possible to introduce measurement errors in a
general way without any bias. Moreover, we illustrate our method using an aﬃne term
structure model and JGB data.

Key words: Generalized state space model, Monte Carlo ﬁlter, interest rate model, aﬃne term
structure model.