HITOTSUBASHI

UNIVERSITY

Center for Financial Engineering Education (CFEE)CFEE Discussion Paper Series No. 2009-1

信用リスク計測に関する一考察

Implied Data にもとづく統計的推論と格付けごとの倒産確率の推定と推定誤差の計測方法

高橋 一

Center for Financial Engineering Education (CFEE)
Department of Economics, Hitotsubashi University
2-1 Naka, Kunitachi Tokyo, 186-8601 Japan
http://www.econ.hit-u.ac.jp/~finmodel信用リスク計測に関する一考察
Implied Data にもとづく統計的推論と格付けごとの倒産確率の推定と推定誤差の計測方法
2009 年 6 月 20 日
2009 年 9 月 15 日改訂
2009 年 9 月 30 日改訂
一橋大学大学院経済学研究科
CFEE
高橋 一

本稿は日本商品先物振興協会 2004 年度研究助成及び CFEE からの援助により遂行された。同協会と CFEE のご厚意に深く感謝致します。

-1-

始めに
ある企業が倒産する確率はその企業の信用度により異なる。逆に信用度はその確率により決まるといっても良い。例えばムーディースは企業の信用度を Aaa, Aa, A, Baa, Ba, B,
Caa の７グループに分類、倒産状態を入れると全体では８のクラスに分けている。1970 年
から 2003 年の間でそれらの企業が６年以内に倒産する確率は 0.12%(Aaa), 0.24%(Aa),
0.54%(A), 2.16%(Baa), 11.17%(Ba), 31.9%(B), 60.36%(CCC)となっている。ムーディース
がどのような方法でヒストリカルなデータより倒産確率を算出しているか、少なくとも筆者には不明であるが本論文では各クラス毎の倒産確率の推定方法について考察する。又、本稿では現実の確率の元での倒産確率と所謂マルチンゲール測度に於ける倒産確率の双方での推定法を提示・提案し、更に推定精度を求める。マルチンゲール測度のもとでの倒産確率は派生証券価格決定の際に必要とされる、一方現実の確率のもとでの倒産確率の必要性は言うまでもない。
両者はリスクの市場価格を推定することにより統合可能ではあるが、本論文ではその部分には立ち入らないこととする。何れのケースでも推定精度を与える事によりリスクを加味した価格設定、ポートフォリオ評価等が可能となる。
本論文は、第１部”既存理論のレビュー”。新たな推定法の提案は”第２部”Implied 倒
産確率の推定と統計学的考察”、そして第３部”RealProbability のもとでの倒産確率の推定”より構成されている。第１部では、現代ファイナンス理論の標準理論というべき
Harrison-Pliska による派生証券価格決定の基本定理とそれに基ずくゼロクーポン債の価格決定論。又、過去３０年の間、主に生物・医学統計の分野で著しく発展してきた生存時間
解析についての簡単な解説を行う。第 2 部では Implied 倒産確率をデータと見立てた統計解析を、第３部では Bootstrap 方を中心に議論を進めて行く。

-2-

第１部 既存理論のレビュー
以下、(Ω, F, Pr) を基本となる確率空間とする。測度空間 (Ω, F)はブラウン運動が定義できるだけの豊かさを持ち、F は Pr に関し完備であると仮定する。{W(t), t ≧0}を(Ω,F, Pr)で定義される標準ブラウン運動、{Ft＝σ(W(s), s ≦t),t > 0} （但し、Ft⊂F ）をブラウン運動に関する情報集合列(Brownian Filtration)とする、但し、σ(W(s), s ≦t) はブ
ラウン運動 (W(s), s ≦t) より生成されたσ加法属である。一方、離散時間のフレームワ
ークでの議論の祭には ε1, ...,εnσ(ε1, ...,εn ),第１-1 節を独立で標準正規分布に従う確率変数列とし、 {F#n＝n≧1}, (F#n⊂F ) を正規ランダムウォークに関する情報集合列とする。

金利の期間構造とスプレッド

市場には各種の利子が存在する。銀行預金の利子が一番判りやすいが、それ以外にも
LIBOR、国債や社債に代表されるクーポン債やゼロクーポン債、さらにはスワップ・レー
が市場利子率を決定する。さて、借り入れや貸し出し期間の長さ、借り手・貸し手の信用力や流動性によって利子率は異なる。一般に国債には倒産の危険がなく(Default Free)流動性も十分あるとされている為、国債の利子率は一番低い。これは、市場では日常的に観察されるし、直感的にもわかりやすい事実であろう。一方、社債等は発行主体の倒産可能性（Defaultable)がある為、利子率は高くなる。そこで倒産リスクを持った債券と無リスクな国債の利回りとの間に金利スプレッド(Spread)と呼ばれる差が生る。金利スプレッドの存在を説明し、更にそれを定量的に算出する理論の構築とは過去 30 年以上に渡り金融工学・
数理ファイナンス世界で大きな問題の一つとなっている。本稿では取りあえず、倒産等の信用リスクがスプレッドの唯一の要因と仮定し、以上の問題に答えていこう* 1。
利子率を現時点から満期日までの時間の関数と見たものを金利の期間構造(Term Structure of Interest Rate)とよび、そのグラフがイールドカーブ(Yield Curve)である。期間構造は日々変化する、その変化を分析する代表的な方法として（１）局所無裁定法(Local Method)と（２）マルチンゲール法(Martingale Method)がある。局所無裁定法に基づく派生証券価格決定理論の代表的な論文は Vasicek[1977]、マルチンゲール法では Harrison andKreps[1979], Harrison and Pliska[1981], Heath, Jarrow and Morton [1992]等であるが、本稿では、次節以下、マルチンゲール法に基づいた議論を紹介する。

第１－２節

マルチンゲール法に基づく無リスクな利子率の決定

本節ではまず、利子率決定のメカニズムを国債等に代表されるクーポン支払いがなく（Zero Coupone）Default Freeな割引債券の価格（安全債券価格）の決定を解説する。又、簡単化のため我々の考える経済では連続時間での取引が有限な時間内[0,Τ]で行われるも*1流動性もスプレッド存在の要素となっている。

-3-

のとし、P*(t,T)をTを満期日とするZero Coupon Default Free債のt時点での価格とする（0≦t≦T≦Τ）。以下、各債権を満期日により区別するために記号の乱用ではあるがP*(t,T)で債権そのものをも表すことにし、それをT債権と呼ぶことにする *2。ξ=T-t を当該債権の残存期間、R*(t,ξ)で利回り(Yield to Maturity)を表すことにしよう。上でも述べたようにR*(t,ξ)をξの関数と見たものが金利の期間構造である。金利の期間構造を分析する際の、もう一方の主役がスポット金利（Spot Rate）である。スポット金利r(t)は時点tからt+dtの間の利子率である（離散時間モデルでは次時点までの間の利子率）。この金利は時点tで決定され、当該区間では既知であるという意味で安全利子率(Risk Free Rate)と呼ばれる。利子率r(u)で１円を[0,t]の間連続複利で運用したときの元利合計を基準資産とよびB(t)と書く。ここで、以下の関係式が成立する、(1-1)

B(t) =

exp{∫0t r(u)du }次に P*(t,T)に対する基本的な仮定を述べ、その後に無裁定条件下での P*(t,T), R*(t,ξ),r(t)の関係を考える。特に P*(t,T)と r(t)との関係は基本的である。
[仮定１－１] 時点ｔにおける割引債価格 P*(t,T)は以下の条件を満たすものとする。
(ⅰ）任意の 0≦t≦Τ に対し(1-2)

P*(t,t) = 1(ⅱ）任意の 0≦t≦T≦Τ に対し、
(1-3)

P*(t,T) > 0(iii）一般に P*(s,T), t≦s≦T は満期日以外ではランダムに変動する確率過程である。即ち、時点ｔで全ての T ≧ 0 に対し、P*(t,T)は既知であるが、P*(s,T)は不確定要素を含む確率変数である。//
T 債権の時点 t 価値は P*(t,T)円、時点 T では１円の価値がある。利回り R*(t,ξ)との間に裁定関係が存在しないためには以下の式が成立していなければならない。
(1-4)

P*(t,T)・exp{ξR*(t,ξ)} = P*(T,T) = 1、

但し、ξ=T-t

又は、
(1-5)

R*(t,ξ) =−1
ξlog P*(t,T)*2 本稿では満期日により債券を区別している。即ち、満期日が同じであれば同じ債券とし
ている。

-4-

次に[仮定１－１（iii)]に関連し確率的な期間構造モデルを導入する。市場には一般に複数の確率変動の要素が存在するが、本節では簡単のため確率変動のソースは一つ、１－ファクターモデルを仮定する[仮定１－２：確率モデル]

W*(t)を(Ω, F, Pr)上で定義された標準ブラウン運動とする。

任意の 0 ≦t≦T≦Τ に対し、P*(t,T), r(t) の確率変動は以下の（伊藤の）確率微分方程式で与えられる。

(1-6)

(1-7)

dP ∗(t,T)= μ*(t,T)dt + σ(t,T)dW*(t),P∗(t,T)dr(t)= α(t)dt + β(t)dW*(t)//

さて、もしもスポットレート r(u)と債券価格 P*(u,T) が全ての 0≦u≦T に対し既知であったとすると（即ち、非確率的モデル）、市場に裁定機会が存在しない時にはP*(t,T)exp{∫tT r(u)du } = 1(1-8)

が成立する。P*(T,T) ＝ 1 より、(1-8)を書きかえると、(1-9)

P*(t,T) = exp{∫0t r(u) du }・exp{-∫0T r(u) du }・P*(T,T)

B(t) の定義より、(1-10)

P∗(t,T)
B(t)＝P∗(T,T)

for all 0≦t ≦T

B(T)が成立する。さらに(1-10)は無裁定条件と同値な条件となる事も証明される(cf. 高橋[2009]）。ここで重要なことは一般の確率モデルに於いても、或る確率測度のもとの期待値の世界で同様の結果が証明される。それが以下で紹介する Harrison-Pliska の定理”無裁定条件下における派生証券価格決定の基本定理”でマルチンゲール法の基礎である。
[定理１－４：派生証券価格決定の基本定理 Harrison and Pliska [1982]]市場に裁定機会が存在しない為の必要十分条件は Pr と同値な確率測度 Q が存在し、そのもとでP*(t,T)/B(t) が Ft マルチンゲールとなることである。即ち、任意の 0≦t≦T≦Τ に対し、

(1-11)-5-P∗(t,T)
B(t)P∗(T,T)= EQ{ B(T) |Ft}=

1

EQ{B(T)|Ft}

が成立する。
但し、同値な確率測度とは両確率の元での確率ゼロで起こる事象(不可能事象）を共有するという意味での同値関係である。また、EQ は測度 Q のもとでの期待値である。
//
これより、T 債権の時点ｔ価格は(1-12)

EQ{exp{-∫tT r(u)du |Ft}

P*(t,T) =となる。勿論この公式を実際に使うためには Q 測度の元でのスポット金利プロセス {r(t),t ≧0}とその同時分布が与えられなければならないが、幾つかの重要な特殊ケースにおいて実際に計算可能である。
これまでは債券の発行主体は倒産しないという仮定の下での議論であったが、倒産の可能性をもった一般企業が発行するゼロクーポン社債の時点 t 価格 P(t,T) はどのように決定されかも同様に重要な問題である。詳しくは、本論文の主要テーマとも関連するため第２部で議論するが、先ずは、ここでも P(T,T)=1,

P(t,T) ≧0（0≦t≦T)を仮定しよう。その上で安全債権と危険債権の時点ｔから満期日 T までの利回りをそれぞれ R*(t,ξ), R(t,ξ)とした時の両者の差が所謂スプレッド(Spread)である。
[定義１－５：金利スプレッド] 安全債権の利回りと危険債権の利回りの間のスプレッド
Υ(t) =Υ(t:P*, P)は(1-13)

Υ(t) =

Υ(t: P*, P) = R(t,ξ)- R*(t,ξ)=−1
ξP (t,T)

log[P ∗

(t,T)]

で定義される。//
倒産リスクを加味すれば P(t,T)≦P*(t,T)、即ち、R(t,ξ)≧R*(t,ξ)と予想され、一般にスプレッドは正値を取ると思われる。第２部では P(t,T)の価格決定モデルとスプレッドに基づいた倒産確率の推定について考察する。その前に次節に於いてハザードモデルを紹介
しよう。

第１-３節

ハザード(Hazard)モデルとセンサーされた(Censored)データ

生存時間分析は品質管理や医学統計の世界で重要な位置を占め続けてきている。そこで研究・開発されてきた手法は信用リスク問題へ応用可能である。企業倒産時刻と死亡時刻を対応させるのみならず、信用リスク問題における他の信用クラス（各付け）への移動は-6-

医学統計における観測個体の退出（センサーリング）に相当する。また、細かな格付けを行うことによる格付け毎のデータ数の相対的な減少は Cox の比例ハザードモデルにより改善される。以下本節ではハザード関数、比例ハザードモデルについて簡単な解説を行う。
[ハザード関数：記号と基本用語] すべての企業がある格付けクラスに入り、倒産するまで同じクラスに居つづけるならば問題は比較的単純である。しかし、実際にはそのまま倒産する場合もあれば､他の格付けへ移行するケースもある。実際、Aaa に格付けされた企業がそのまま倒産することは殆ど無いだろう。倒産する場合も、先ず下位のクラスに格付けされ、それから倒産となる場合が多い。そこで問題となるのは倒産前に他のクラスへ移行した企業の持つ情報、所謂センサーされたデータが持つ情報をどのように活用するかである。さて、センサーされたデータとは当該事象（倒産）の生起・不生起以前に対象が観測不能となる状態のことである。少なくともセンサーされた時点までは倒産していないという情報をデータは持っているが、それをどのように活用するかが生存解析の主題の一つである。 正確には倒産時点が不明という意味に於いて､これを右センサーされたデータと呼ぶことがある。一方、当該格付けクラスに参入した時点が不明であるようなデータは左センサーされたデータと呼ばれる。
まず、一つの格付けクラスを選び固定する。ある企業が当該格付けに入ってから倒産する迄の時間（生存時間）を確率変数τで表す。一方センサーされる迄の時間は C とする。
実際には、大部分の企業は観測が始まった時点では既に当該クラスに格付けされているという意味に於いて左センサーされたものと考えられるが、本稿では単純化のためこの立場は取らない。その為、例えば現時点を基準時 0 とし、C をセンサー時刻（時点）と言うこともある。τの分布関数と密度関数をそれぞれ F(t), f(t) とする。また、センサー時点 C
は生存時間τとは独立で分布関数 G(t)、密度関数 g(t) に従うと仮定する。但し、観測終了時点 T については自動的に定数値のセンサー時点とする（C は T で切断されている）。
さて、生存関数 S(t)とハザード関数 λ(t) はそれぞれ、
S(t) = 1 - F(t) = P{τ> t}(1-16)
λ(t) =f(t)S(t)

で定義する。生存関数は、時点 t を超えて企業が生き残る確率、また、ハザード関数はλ(t)dt ≒ P{ t <τ≦t+dt | τ>t } と書き直してみると分かりやすい。即ち、ハザード関数とは時点 t まで倒産しなかった企業が次の瞬間（微少時間内）に倒産するという条件付き確率である。簡単な積分計算より(1-17)-7-

S(t) = exp{-∫0tλ(u)du},

F(∞)=1 ⇔∫0∞λ(u)du = ∞が得られる。次に観測時点 Y とセンサー指示関数δを以下のように定義する。
(1-18)

Y = min{τ, C} = τ ∧C ～H(t)(1-19)δ=I[τ≦C] ＝�

1 if uncensored0 if censored

τの分布、又は生存関数やハザード関数に対しパラメトリック、ノンパラメトリック等
色々なモデルが考えられる。以下では代表的なパラメトリックモデルを紹介するが、これらは経験的合理性を伴ったハザード関数の形から考えられたものである。
[例１－６:指数(Exponential)モデル] 一番簡単なパラメトリックモデルはハザード関数
が時間に対して、一定となる指数分布型のモデルであろう。すべての t > 0 に対して、ハザード関数λ(t)と生存関数 S(t)は、それぞれ、
λ(t) = θ > 0(定数)S(t) =t≧0(1-20)
exp{-θt}である。分布関数 F(t)= 1-exp{-θt} t≧0, 密度関数 f(t) =θexp{-θt} t≧0, である。勿論、一般にθは未知で推定の対象となる。指数分布は簡単ではあるが、その単純性から多くのケースで活用されている。指数分布の拡張として以下のガンマモデルとワイバルモデ
ルがよく知られている[例１－７:ガンマ(Ganma)モデル]

指数分布の自然な拡張がガンマ分布である。θを尺度(Scale)パラメータとし、それに加え形状(Shape)パラメータαが導入される。密度関数は(1-21)

f(t) =

θαtα-1exp{-θt}

Γ(α)

t ≧0

ただし、Γはガンマ関数である。ガンマ分布関数は不完全ガンマ関数を含むため生存関数、
ハザード関数、分布関数の明示的な表現は得られない。
[例１－８：ワイバル (Weibull)モデル](1-22)

指数モデルにおいて、生存関数を

S(t) = exp{-(θt)α}と拡張したものがワイバルモデルである。簡単な計算より密度関数並びにハザード関数は-8-

f(t ) = αθ(θt)α -1 exp{-(θt)α}(1-23)
λ(t) = αθ(θt)α-1となることも直ちに示せる。//
以上、代表的な３種類のパラメトリックモデルを紹介したが、これらの他にも、ラレイ（Rayleigh)モデル、対数正規（Lognormal)モデル、パレート（Pareto)モデル等がよく知られている。比例ハザードモデル(cf. (1-24))に於けるベースラインハザードとしてどの様なモデルを選ぶべきかについて、信用リスク問題でも、信頼性工学や医療統計同様、外的に与えられる幾つかの経済変数の影響を無視できない。特にロングランの信用リスク問題では経済の構造変化等も考えられる。一方、実装の際には指数分布型モデルの単純さも魅力的である。

勿論データによる丁寧な検証は必要であるが、比較的短期期間に於いて指数分布を仮定することは一定の合理性を持つであろう。
[ランダム・センサーリングモデル(Random Censoring Model)] センサー時点について考えると、例えば観測終了時 T のように予め決まっている時点もあれば、何らかの外的理由による退場もある。ここでは、後者によるセンサー時刻 C を確率モデル化してみる (cf.
Efron[1980])。ランダムセンサーリングモデルは第２節で展開するセンサーされたデータに対するパラメトリックブーツトラップを実装するために用いられる。
(1-18)で倒産時刻τの代わりにセンサー時刻 C に注目し、その生存関数 1 - G(t) の推定問題を考える（但し、ここでは T で切断されてはいないと仮定しよう）。観測終了時刻に対応するセンサー時刻を除き C とτの役割を入れ替えてみると、C に取ってτはセンサ
ー時点に他ならない。τを T で切断したものにすると、ここでも観測値は上と同様 Y =min{τ, C} であるが、センサー指示関数は δ' = 1 - δ となる。[例１－６]～[例１－８]のパラメトリックモデルは C に関するモデルとしても有効である。//
勿論、格付けクラス毎に生存時間の分布も異なる、そこで本稿では以下のモデル、比例ハ
ザードモデルを導入する。ポイントはハザード関数を通じて異なるモデルを特徴付けるこ
とである。

第１－４節

比例ハザードモデル(Proportional Hazard Model)

本節では Cox [1974]が提唱した比例ハザードモデルを概観する。x = (x1,x2, ..., xp) を独
立変数、β=(β1,β2, ...,βp)を(未知)パラメータとし、比例ハザード関数を以下のように定義する：(1-24)-9-λ(t:x,β) = λ0(t)exp{β'x}以下λ0(t)をベースライン(Base Line)ハザード(関数)と呼ぶ。生存関数はS(t:x,β) = exp{-∫0t λ(u:x,β)du} = exp{-eβ 'x∫0tλ0(u)du}(1-25)= exp{-∫0tλ0(u)du}exp{eβ 'x} =S0(t)exp(β'x )また、分布関数は F(t:x,β) = 1- S0(t)exp(β'x )となる。例えば p = 1, x = 0 or 1 であれば、S(t:0,β)と S(t:1,β)により 2 標本問題がカバーされる。勿論これ以外にも幾つかの方法が考えられる。p=1 で同一のβに対し異なる x の値でクラスを区別する方法や p ="格付けクラスの数"とし x をダミー変数とする等、色々考えられる。何れにせよデータに基づいたモデル
選択とλ0(t)とβの推定が必要となる。この問題については第３部で詳しく論じる。第３－２節で展開する Parametric Bootstrap の場合ベースラインハザードについてはパラメトリックモデルを仮定する。

第１－５節

ハザード関数のパラメトリックモデルと推定論

本節ではパラメトリックな生存時間モデル{H(t,θ):θ∈Θ}を仮定し未知パラメータθ
の最尤推定法について議論する。(1-18),(1-19）にあるように、我々が用いることの出来る
データは(Y1,δ 1),(Y2,δ2), ... , (Yn,δn)である。ここで Y1, Y2, ..., Yn は独立で同一な分布
H(t:θ),θ∈Θ からの大きさ n のランダムサンプルである。先ず尤度関数の定義から始め
よう、(Yi,δi) の実現値 (yi,δi)の尤度関数 Lτi(θ) = Lτθ(yi,δi) は(1-26)

Lτi(θ) = Lτθ(yi,δi) = �=f(yi ) = f(yi ∶ θ)S(yi ) = S(yi ∶ θ)f(yi)δiS(yi)1-δiif�δi = 1 （uncensored）
δi = 0 （censored)(i = 1, ..., n)

従って、標本全体の尤度関数 Lτ(θ)はセンサー時刻と倒産時刻の独立性の仮定より以下のように表しても一般性は失われない（cf. Miller [1981])
Lτ(θ)= Lτ(y1, yi2, ... , yn ;δ1,δ2, ... ,δn :θ)= Πi=1n Lτθ(yi,δi)(1-27)=[Πuf(yi)][ΠcS(yi)]但し、Πu はセンサーされていないデータに対応する i に関する積、Πc はセンサーされたデータに関する積である。この様に、尤度関数が求められたなら、一定の正則条件下で未知パラメータは対数尤度関数を微分し、それをゼロとおくことにより求められる：- 10 -

0=(1-28)∂∂θjlog Lτ(θ) = Σi=1n

= Σu

∂

∂θj

∂

∂θj

log fθ (yi) + Σclog Lτθ (yi,δi)∂∂θjlog Sθ (yi),j = 1, ..., p�=(θ�1, ... , θ�p)は連立方程式
そして、最尤推定量（ベクトル） 𝛉𝛉(1-29)

∂

0=

∂θj

�)
log Lτ(𝛉𝛉j = 1, ..., pの解として与えられる。勿論、多くの場合上記方程式は解析的には解けず、数値的に解くしかない。数値解析の代表的方法は以下で説明する Newton-Raphson(NR)法である。
[Newton Raphson 法]

N-R 法は以下の手順で方程式(1-29)を解いていく：
�(1)とする。
1 Step: まず最尤推定値（(1-29)の解）を予想しそれを𝛉𝛉st∂

2nd Step:

∂θj

0=

∂

�の周りでテイラー展開、(1-29)式よりlog Lτ(θ)を𝛉𝛉

∂θj

�) =
log L(𝛉𝛉

∂

∂θj

2

�(1)) + (𝛉𝛉� - 𝛉𝛉
�(1) ) ∂ 2 log Lτ(𝛉𝛉�(1)) + ・・・・
log L(𝛉𝛉∂θ j( j = 1, …,p)�(2) を得る。以下の形
上式で二次以上の高次項を無視することより、更新された推定値 𝛉𝛉

でのベクトル表示にすると、(1-30)

�(2) = 𝛉𝛉�(1) + [i(𝛉𝛉�(1))]-1 ∂ log Lτ(𝛉𝛉
�(1))
𝛉𝛉
∂ 𝛉𝛉

但し、
(1-31)

2

�(1)) = − ∂2 log Lτ(𝛉𝛉
�(1))
i(𝛉𝛉
∂θ

は標本情報行列と呼ばれる。
3rd Step:(1-32)

2nd step を繰り返し、一般には以下の推定値を得る：
�(k+1) = 𝛉𝛉�(k) + [i(𝛉𝛉�(k))]-1 ∂ log L(𝛉𝛉
�(k))
𝛉𝛉
∂𝛉𝛉
P

k≧1

�)に対する適当な正則条件下で 𝛉𝛉�(k) → 𝛉𝛉� となることが示されるため、あらか
一方、Lτ(𝛉𝛉
じめ定められた限界値 a > 0 に対し(1-33)- 11 -

�(k) || < a�(k-1) － 𝛉𝛉|| 𝛉𝛉�(k) を最尤推定値とする、但し、||・|| はユークリとなるところで繰り返し計算を中止し、𝛉𝛉ッドノルムである。

もしも、Fisher 情報行列 I(θ)= E{i(θ)}が利用できる場合は、標本情報行列の代わりに
Fisher 情報行列を用い(1-32)の一般項の代わりに

�(k+1) = 𝛉𝛉�(k) + I(𝛉𝛉�(k))-1𝛉𝛉(1-34)∂

�(k))
log Lτ(𝛉𝛉

k≧1

∂𝛉𝛉

と求めることも出来る。//
一定の正則条件下で、センサーされたデータに基づく最尤推定量も一致性と漸近正規性を持つことが知られている：例えば、d�－θ)I(θ)1/2 →(𝛉𝛉(1-35)

N(0,Ip ) as n →∞ただし、 Ip は p× pの単位行列である。これより、θ0 を真のパラメータとするとき以下の結果が証明される：

(1)

Wald
�－θ0 )’I(θ0)(𝛉𝛉
�－θ0)
(𝛉𝛉

(2) Neyman-Pearson

d

→

χ2p

d

(3) Rao�)] → χ2p
-2 log[L(θ0 )/ Lτ(𝛉𝛉∂∂log Lτ(θ0 )T I(𝛉𝛉𝟎𝟎 )−1 ∂𝛉𝛉 log Lτ(θ0)∂𝛉𝛉as n →∞

as n → ∞d→ χ2p

as n →∞//

この様に、生存関数の推定に関してはセンサーされたデータに基づく場合でも、通常のデータ抽出の場合とほぼ同様な統計的学的な結果が示される。また、比例ハザードモデルの推定等も現在では標準的な統計学手法であり、その統計学的な性質について標準的な結果が数多く示されてきている。
しかしながら本論文では若干特異な立場から論じていく（パラメトリックモデルとノンパラメトリックモデルの融合）ため、この件については第３部
で今一度論じることとする。
以下は Parametric Bootstrap 法適用の際に必要となる。死亡時刻とセンサー時刻を入れ替えた時の尤度関数を与えておく。

- 12 -

ランダムセンサーリングモデルに於いても尤度関数は(1-26),(1-27)と同様に定義できる。(C,δc)に関する尤度関数を LC と置くと尤度はLC(θ)= LC(y1, yi2, ... , yn ;δc1,δc2, ... ,δcn :θ)= Πi=1n LCθ(yi,δci)(1-36)=[Πug(yi)][ΠcR(yi)]δc=(1-37)�I[τ≧C] ＝

1
0

if

uncensored�
censored

但し、Πu はセンサーされていないデータ（変数 C)に対応する i に関する積、Πc はセンサーされたデータ（変数 τ）に関する積である。ここでも上と同様、尤度関数に基づいた
統計学的推論が可能となる。
第１-６節

生存関数の Non-Parametric 推定

生存関数のノンパラメトリック推定法は古くから研究されてきている。代表的なものとして、例えば、Reduced Sample Method, Actuarial Method, そして Product-Limit Estimator、
即ち Kaplan-Meier（K-M)推定量等がある。以下 K-M 推定量と、その導出のアイディアを簡単に紹介する。先ず、時間軸を区間 Ii = (Y(i-1), Y( i)], i=1,…,n と分割する。ただし、（Y(1),Y(2), ..., Y( n))は観測時点 (Y1, Y2,..., Yn)を順番に並べたものでδ(i) は Y(i) に対応するセンサー指示関数とする。

( I ]( I ]( I ]

( I

0

1

2

3

n

]

Y(1)

Y( 2)

Y(3)

Y( n-1)

Y( n)

x

0

0

x

0

[図１－１]
但し、変数 Y の下の" 0" と"x" は "0 = センサー、x = センサーされていない"データということを意味する。ここで、R(t)を時点 t におけるリスク集合(Risk Set)（時点 t-まで生き残っている企業の集合）とし、さらに以下の記号を導入する。
ni =R(Y(i))に含まれる個体数 = 時点 Y(i)- で倒産していない企業数

di = 時点 Y(i)で倒産した企業数
pi = Pr{τ > Y(i) | τ > Y(i-1)}
qi = 1 - pi
時点 t における生存関数の値を求めるために、区間(0,t]を上記のような小区間に細分する。

- 13 -

すると、確率の乗法定理より、
S(Y(k)) = Pr{τ > Y(k) }(1-38)= Pr{τ > Y(1) }Pr{τ > Y(2) |τ > Y(1)}・・・Pr{τ > Y(k) |τ > Y(k-1)}= p1p2 ・・・・pkとなる。従って S(Y(k)) を推定するには(1-38)式中の各 pi = Pr{τ>Y(i+1)|τ>Y(i)}を推定すればよい。そこで、 pi と qi をdq� i = n i(1-39)
そして、(1-40)p� i = 1 - q� i =

i

�

1

1 − n
1

�

if

i

δ(i) = 1
δ(i) = 0で推定する。ここで、タイのない（No Tied Data, 同時刻で観測されたという意味での）データの場合、生存関数の K-M 推定量 S� (t)を(1-41)

で定義する。

- 14 -S� (t) =∏y(i)≦t p� i=∏y(i)≦t [

n−i

n−i+1

]δ (i)

第２部 Implied 倒産確率の推定と統計学的考察
第 1 部は準備として、本論文で用いられるファイナンスと統計学の標準的な手法を概観してきた。
第 2 部では金利スプレッドより得られる情報を用い倒産確率を推定する方法と、推定誤差について考えていく。第２部の主目的は推定値が必ず持っている誤差限界をキチンと求めていくことにある。

第２－１節：Duffie-Singleton の信用リスクモデル
本節では信用リスクを加味した金利の期間構造モデルの解説・構築、およびモデルの推定法について考えていく。1970 年代に Merton[1974]により提唱された構造型モデル(Structual Model)とその後、90 年代に主に Duffie-Singleton[1997]等により考察された誘導型(Reduced Model)が代表的な方法である。本稿では先ず Duffie-Singleton（D-S)モデルを紹介し、それに基ずく倒産確率の推定法”Implied 倒産確率" について議論していく：[モデルと価格式] 倒産リスクを持つ(Defaultable)企業の社債価格を求めるために、本節では先ず離散時間モデルで基本的なアイディアを解説し、その後連続時間モデルにおける価格式を提示する。以下、離散モデルではブラウン運動の代わりに正規ランダムウォークを仮定し、フィルトレーションも正規ランダムウォークに対応したものとする。又、離散時間モデルでは記号に# を付けることとする。
本節ではリスク中立的世界（同値マルチンゲール測度下）で議論を進めていく。以下 Q
を原確率 Pr と同値なマルチンゲール確率測度とする。本章でもこれまで同様 Defaultable
な企業 の発 行す る債券 価格 等を 表す 時は *無 し の記号 を用 いる ことと する 。例 えば
Defaultable 債券の時点 t 価格は離散時間の元では P# (t,T)、連続時間では単に P(t,T)と記す。
現時点は t、満期時点 T 等の記号はこれまで通りである。本節の前半部は一部を除き基本
的には Duffie-Singleton [1999]からの抜粋である。
ある特定の企業を企業倒産に関連しς以下の記号を導入する。まず、当該企業が倒産する
時刻を離散型確率変数τ#としその分布関数を F# とする：
(2-1)

F#(t) = Q{τ#≦t}

また、離散型ハザード関数を(2-2)

λ#(t) = Q{τ# = t | τ#≧t}で定義する。次にψ#(t)を倒産時における資金の回収率、r#(t) = r#t を区間 [t, t+1) における安全利子率とする。

- 15 -

以上の準備の元、満期時点 T で XT 円償還する Defaultable な債権の時点 t 価格 V#(t)を求めてみよう。時点 t 迄倒産が起こらないという条件下で債権の時点 t 価格 V#(t)を求める
ために、次の時点 t+1 で倒産する場合と生き残る場合に分けて考える。EQ を確率測度 Q のもとでの期待値とするとき、(2-3)

V#(t) = λ#(t)exp{-r#t}EQ{ψ#(t+1)|F#t}+ [1-λ#(t)]exp{-r#t}EQ{V#(t+1)|F#t},

但し、債務回収は t+1 時点で行われるものとする。上式を１期ずつ繰り返し t を T 迄増加させると、最終的にはV#(t) = EQ[Σj=0T-t-1 exp{-Σk=0j r#(t+k)}ψ#(t+j+1)Πk=0j{1-λ#(t+k-1)}|F#t](2-4)+

EQ[exp{-Σk=0jT-t-1 r#(t+k)}XTΠm=0T-t[1-λ#(t+m-1)]|F#t]を得るが、実際に使うには少々複雑すぎる。そこで簡単化のため時点 t に於ける回収率
ψ#(t+1)は V#(t+1)の一定割合で F#t－可予測関数 δ#t = (1-L#t)と仮定し、以下の関係式を得る：
EQ{ψ#(t+1)|F#t}=δ#tEQ{V#(t+1)|F#t}(2-5)

上式で導入された δ#t = (１-L#t)は市場価値回収率（RMV: Recovery of Market Value）と呼ばれる時間の関数で D-S モデルを特徴づけるポイントである。(2-5)より(2-3)式は
δ#tλ#(t)exp{-r#t}EQ{V#(t+1)|F#t}+ [1-λ#(t)]exp{-r#t}EQ{V#(t+1)|F#t}(2-6)となる。上と同様に t を T まで増やしていくことにより、(2-7)

V#(t) = EQ{exp{-Σj=0 T-t-1R#t+j}XT | F#t}を得る、但し(2-8)

exp{-R#t }=(1-λ#(t))exp{- r#(t)}+λ#(t)exp{-r#(t)}δ# tである。Taylor の定理より、(2-9)

exp{-R#t }= exp{-(r#(t)+λ#(t)(1-δ#t)}[1+o(r#t)]

as r#t→0

が 成 立 す る 。 即 ち 、 ハ ザ ー ド 関 数 λ #(t) と 回 収 率 δ #t=(1-L#(t)) で 特 徴 づ け ら れ る
Defaultable な企業の社債価格は無リスク債券価格と同様な式より求められる。RiskFree の場合との唯一の違いは割引率が Risk Free レート r#(t) でなく R#(t) = r#(t) +λ#(t)(1-δ#t)- 16 -

になる点である。
以下、この結果を連続時間のフレームワークで書き直す。記号は主な関数から#を外すことにする。離散時間モデルの結果は、幾つかの数学的な問題は生じるが基本的には観測時点間をゼロに”収束”させると、結果的にはその間における利子率もゼロの収束すること
より連続時間モデルに変換される。これより証明は省くが、以下の連続時間における価格決定式が求められるのである。
V(t) = EQ{exp{-∫t T R(u)du}XT | Ft}(2-10)

勿論、R#(t) = r(t) +λ(t)(1-δt)となることは言うまでもない。

第２－２節

Implied 倒産確率の推定： ノンパラメトリック法

本節では主として連続時間モデルを仮定する。仮に r(t)と L(t)が既知ならば、市場で観測された実際の価格 V(t)を(2-10)式の左辺に代入することにより”理論的には”ハザード率を逆算することが可能である（離散時間の場合は V#(t)を(2-7)式の左辺に代入）。この様にして求められたハザード関数をインプライドハザード関数（Implied hazard function)と呼ぶ。勿論、解析的に求めることはもとより、このままの形では数値的に計算することも非常に難しい。そこで、楠岡・青沼・中川 [2001]と同様な単純化をはかろう。
先ず、金融資産としてゼロクーポン割引債を仮定；XT=１（円）（即ち V(t)= P(t,T)

である）、RMV δ(t)は定数δとする。更なる簡単化のため t = 0 とし、Q{τ>s}＝exp{-∫0sλ(u)du}に注意し、改めて T = t と書くと、EQ{exp{-∫0 tR(u)du}}P(0,t) =(2-11)=

EQ{exp{-∫0 tr(u)du - (1-δ)∫0 tλ(u)du}}=EQ{exp{-∫0 tr(u)du}・exp{(1-δ)∫0 tλ(u)du}} ＝ P*(0,t)S(t)( 1-δ)

これと(1-13)式より、(2-12)

S(t) = exp{-Υ(t : P*, P)t}(1/1-δ)

ここでも、勿論回収率δの推定のみならず、最終的な結果は回収の方法と時期にも依存する。実際、D-S モデルでは Default の次の時点で債務回収が行われると仮定されているが、楠岡 et. al.[2001]では債務回収は満期日に支払われると仮定している為、彼等の算式は若干異なってきている；- 17 -

(2-13)S(t) = [exp{-Υ(t; P*, P)t}-δ]/(1-δ)(2-12)式に基づくにせよ(2-13)式から求めるにせよ、生存関数は満期日に対応するｔの値でのみ計算される。従って、求められた関数の平滑化と云う問題が次に出てくるだろう。
また、この様に計算された生存関数の精度等の統計学的な性質は一体どうなっているのか。
次節以降ではこれらの問題について考えていく。

第２－３節 Implied 倒産確率の推定：パラメトリックモデルとその推定生存関数 S(t) の推定問題を考えたとき、上記のようなノンパラメトリック法ではすべての t に対し S(t)を求めることは難しい。一方、生存関数の形態に関して何らかの事前情報が存在する場合もあるだろう。そのような時には生存関数（又はハザード関数）に以下のようなパラメトリックモデルを仮定することが考えられる。
[仮定２－１

確率モデル]

Θを Rp の部分集合とする、但し p はパラメータの次元。θ

∈Θを未知パラメータとし、G(t:θ) = Q{τ> t：θ} をパラメトリックな生存関数とする。
又、我々がスプレッド Υ(t) より求めることが出来る観測値 S(t:θ)は以下のモデルに従う

確率的な誤差を含むと仮定する(即ち、確率変数である）；
S(t i )1−δ = S(t i : θ)1−δ = G(t i : 𝛉𝛉)1−δ + σρt i ε ,(2-14)ti

但しεt は i.i.d.で平均と分散は(0,1) (i=1,2,…,n,)、σ> 0 は未知の標準偏差，ρt は ρ0 =i

0 を満足する t の“既知の”関数とする。市場で観測されるスプレッドΥ(t)より計算される S(t)を s(t)と書き実現値（”Implied 生存関数”）と呼ぶ。また、G(t:θ)は任意の t に対して、θに関し２回微分可能であると仮定する。//
1−P ∗(0,t)

以下、ρt i は既知とし議論を進める。応用の場面ではρt i = [1 − P∗ (0, t)]又はρt i = � P ∗(0,t) �等を仮定することが合理的であろう 2。

所で,モデル(2-14)は可変分散モデルであるが本論文ではρt i を既知と仮定しているため、変換

X(t i ) =き、(2-14)は(2-15)

S(t i )1−δρti

並びに Ξ(t i : 𝛉𝛉) =G(t i :θ)1−δρtを施すことにより等分散モデルへ変換でi

X(t i ) = Ξ(t i : 𝛉𝛉) + σεt εt are i.i.d. with (0,1) for i=1,2,…,n,iiとなる。上と同様 X(t)の実現値を x(t)=xt と書くことにし、それを擬似データと呼ぶ。モデル(2-15)のもとでのθの推定は、インプライドな生存確率より求められた擬似データ列2

このモデルの導出については (Takahashi[2009]) を参照。

- 18 -

� を以下の判定基準{x(ti) = xt i i =1,…n}を観測データ値としθの"擬似最小二乗"推定値𝛉𝛉
(2-16)のもとで（非線形最小二乗法を用い）求めることから始まる：
� = argθ∈Θ[minΣi=1n{Ξ(ti:θ) 𝛉𝛉

(2-16)

x(ti)}2]

�の統計学的性質については以下で展開する“直感的”な議論より、一定の正則条
推定量𝛉𝛉

件下で漸近正規性が成立する。

簡単化のため p=1, 即ち１-パラメータの場合を考える。まず、”擬似”スコアー関数(“Pseudo” score function)を以下のように定義する、
1 ∂

ψ(t,x(t),θ) = 2 ∂θ [Ξ(t:θ) - x(t)]2(2-17)∂= [ ∂θ Ξ(t:θ)][Ξ(t:θ) - x(t)]確率変数としての疑似スコアー関数ψ(t,X(t),θ)の平均と分散はそれぞれ、
E{ψ(t,X(t),θ)} = 0∂

Var{ψ(t,X(t),θ)} = σ2[∂θ Ξ(t:θ)]2となることは(2-15) から直ちに判る。
さて、”データ” {x(ti) =xt i i =1,…n} が与えられたとき、θ� がΘの内点であればθ�は以下の方程式の解となる；(2-18)

Σi=1nψ(ti,x(ti)θ�) = 0

また、推定量 �ϑn をΣi=1n ψ(ti,X(ti),ϑ�n) = 0 で定義し、以下、 �ϑn の漸近正規性を示そう。
その為に先ず、ψ(ti,X(ti), ϑ�n)を真のパラメータ値θの周りで Taylor 展開すると、(2-19)

0 = Σi=1nψ(ti,X(ti) ϑ�n)∂= Σi=1n [ψ(ti,X(ti),θ)+(ϑ�n - θ) ∂θ ψ(ti,X(ti),θ) + (1/2)(ϑ�n - θ)2

但し、|θ* -θ| < |ϑ�n - θ| 。これより、(2-20)- 19 -√n(ϑ�n - θ) =1√n∑ni=1 ψ(t i ,X(t i ),θ)1 n∂∑
ψ(t i ,X(t i ),θ)n i=1 ∂θ

+ Op(√n(ϑ�n - θ)2) as

∂2

∂θ 2

ψ(ti,X(ti),θ*)p| ϑ�n - θ| → 0

次なるステップのために以下の正則条件を考える。
[仮定 ２－２：漸近正規性のための正則条件](a)． {ψ(ti,X(ti),θ), i = 1, ... ,n}は Lindberg-Feller の中心極限定理に従う。
(b)． �ϑn → θ in probability (as n → ∞)12∂(c)．n ∑ni=1 �∂θ Ξ(t i , θ)� ⟶ A > 0
(d).

(e)．

1

∑ni=1 �n∂∂θ

(as n → ∞)∂2

Ξ(t i , θ)� �∂θ 2 Ξ(t i , θ)� ⟶ B ≠ 0 (as n → ∞)2√n( �ϑn - θ) → 0

in probability (as n → ∞)。//

ここで上記の正則条件下で次の定理が成立する。
[定理２－３][仮定 2-1,2-2]のもとで2√n( �ϑn - θ) は平均 0、分散σ A−1 の正規分布に分布収束する。 //[証明]以下、証明をスケッチしていく。先ず、(2-14)より

∂

2

∂

E{∂θ ψ�ti , X(t i ), θ�} = �∂θ Ξ�t i , θ��

(2-21)(2-22)

E{

∂2

∂θ 2

∂

ψ(ti,x(ti),θ*)} = �∂θ Ξ(t i , θ)� �

∂2

∂θ 2

Ξ(t i , θ)�
1∂

大数の法則と[仮定２－２] (c),(d),(e)より n→ ∞ の時、n ∑ni=1 ∂θ ψ�ti , X(t i ), θ� は漸近的に
1

n

1

n

∑ni=1 �∑ni=1(2-23)(2-24)

∂

∂θ
∂2

∂θ 2

2

Ξ(t i , θ)� と等しくなり A に確率収束する。議論は若干ややこしくなるが，
ψ(ti,x(ti),θ*)も同様にして B へ確率収束する。即ち、
1

n

∑ni=1
1

n

∂

1

ψ�ti , X(t i ), θ� ~
∂θ

∑ni=1

∂2

∂θ 2

n

∑ni=1 �Ξ(t i , θ)�
∂θ1∑ni=1 ψ(t i ,X(t i ),θ)2∂
Ξ(t i ,θ)�∂θσ �∑ni=1 �- 20 -∂→A∂2

ψ(ti,x(ti),θ*)~ n ∑ni=1 �∂θ Ξ(t i , θ)� �∂θ 2 Ξ(t i , θ)� → B

一方、Lindberg-Feller の中心極限定理より(2-25)2∂d→ N(01)

以上を組み合わせ、(2-20)より(2-26)√n(ϑ�n - θ) =1∑n ψ(t i ,X(t i ),θ)√n i=11 n∂∑
ψ(t i ,X (t i ),θ)n i=1 ∂θ=

+ Op(√n(ϑ�n - θ)2)

∑ni=1 ψ(t i ,X(t i ),θ)�2∂
σ �∑ni=1 � Ξ(t i ,θ)�1 n∂∑
ψ(t i ,X(t i ),θ)n i =1 ∂θ

∂θ

�
2
1
∂
σ � ∑ni=1 � Ξ(t i ,θ)�1

+2 √n�ϑ�n – θ�

2

n

∂θ

1 n
∂2
∑
ψ�t i ,X(t i ),θ∗�
n i=1 ∂θ 2
1 n∂∑
ψ�t i ,X(t i ),θ�n i =1 ∂θ

�1
n

�1
n

∂
∂2
Ξ�t i ,θ��� 2 Ξ�t i ,θ��∂θ∂θ

∑ni=1 �∂∂2
Ξ�t i ,θ��� 2 Ξ�t i ,θ��∂θ∂θ

∑ni=1 �上式最右辺第一項の分子は Lindberg-Feller の中心極限定理により標準正規分布に分布収束2する。分母は(2-23)より√A に収束。一方、第二項は Op�O√n�ϑ�n – θ� � となり、仮定よりゼロに確率収束する。以上より、(2-27)
が導かれる。//√n(ϑ�n - θ)d→N(0. σ2 A−1 )この様に定理自体は標準的な手法で証明される。実際に応用する場面では[仮定２－２]のチェックが問題となる。[定理２－３]よりθの(1-ε)x100%(2-28)

1

信頼区間は近似的に1( �ϑn -zε/2σ�nA , �ϑn + zε/2 zε/2σ�nA ）εで与えられる、ここで zε/2 は標準正規分布の上側 点である。σの推定値は
２

1

n−d

Σi=1n{Ξ(ti: �ϑn) -x(ti)}2 で与える事が合理的であろう。

これより、観測された Implied 倒産確率を用いたパラメトリックモデルの推定と推定量の漸近正規性が示され、その意味において推定値(2-15)は正当性を主張できる。勿論、仮定２－２の十分条件を求めることが現実には必要である。

- 21 -

第２－４ Implied 倒産確率の推定 :Bootstrap 法による推測

前節では中心極限定理を用い未知パラメータの信頼区間を求めたが、本節では Bootstrap（標本からの再標本法）を用い同様のことを行う（cf. Efron[1979],[1981]）。手順は非常に単純である。ここでも前節と同様パラメトリックな確率モデルを仮定する。パラメトリックモデルへの限定は単に以下の形で Bootstrap 法を用い、未知パラメータの推定と精度の推定の議論をするためである。本節における Bootstrap 法は標準的な回帰モデルに於けるBootstrap 法の準用である。その手順は次のステップから成る：[Bootstrap 法」
1st Step:

先ず、与えられた疑似データ {s(ti) = st i i =1,…n} よりθの推定値 θ を判定基準(2-15)のもとで（非線形最小二乗法を用い）推定する：参考のため(2-15)を再度示す。
θ� =(2-29)argθ∈Θ[min∑ni=1{Ξ (t i :θ)- x(ti )}2]

2nd Step: 次に、[仮定 2-1] (2-14)式中にある εt i を ei ( = e t i ) = x(ti ) - Ξ(ti:θ�)で推定、残差ベクトル E = {e1, ..., en}を構築する。これからが所謂 Bootstrap 法のメインパートである。
3rd Step:

E の要素から繰り返しを許し無作為にｎ個の ei を選び､それをE*( 1) = {e*(1) 1, ..., e*(1) n}(2-30)と書く。｛e*(1) 1, ... ,e*(1) n}はトータルでは nn 個の組み合わせの中から選ばれる。その中での極端なケースとしては全て異なる場合（順列を除けば元のデータと同じ中身）もあれば、他方では全ての要素が等しくなる場合があり得る。
4th Step: Bootstrap 標本を以下のように構築する：(2-31)

Q{τ> t i | θ�} + e*(1)ii = 1, … , n
Ξ(t i , θ�) + e*( 1)ix*(1)(ti) = x*(1)i ==5th Step: 疑似データ(2-31)に基づくθの推定値を(2-32)
で求める。

θ�*(1)=argθ∈Θ[minΣi=1n{Ξ(ti:θ) -x*( 1)(ti)}2]

6th Step: Step 3 ～Step５(2-30)～ (2-32)を M 回繰り返しθ の推定値の集合(2-33)- 22 -{θ�*(1)1 , ... , θ�*(M)n}を得る（M=10000 程度とれば、可成り精密なヒストグラムを描くことが出来る。計算時間については非線形最小自乗法の数値計算に大きく依存する）。
ここで、θの Bootstrap 推定値を以下の式で定義する：(2-34)1� ∗(k)
θ�(Boots) = M ∑M
k=1 θ

Bootstrap バイアス推定値を(2-35)

�
Bias(θ�)(Boots) =θ�(Boots) – θ�

一方、Bootstrap 分散推定値は(2-36)1� ∗(k) − θ�(Boots ) )2
Var(θ�)(Boots) = M−1 ∑M
k=1(θそして、推定量の Bootstrap 推定された分布関数は(2-37)

1
Dist θ� (t)(Boots) = M {# of k such that θ�∗(k) ≤ t}で求められる。

しかし、θ の信頼区間を構築するには {θ�*(1), …, θ�*(M)}のヒストグラムから直接求める方が簡単であろう。

- 23 -

第３部 Real Probability のもとでの倒産確率
第３部では現実の確率のもとでの倒産確率の推定について考える。倒産確率は当該企業の信用レベル、即ち”格付け”により異なる為、各格付けクラス毎の倒産確率を求める必要がある。又、例えば AAA に格付けされている企業が AAA のまま倒産することは通常無い。AAA に限らず、多くの高い格付けクラスの企業は同じクラスに居る時に倒産するより、先ず下のクラスへ降格し､その後倒産することが多い（勿論、逆に格上げされ倒産しないケースもある）。従って、我々の目の前にある倒産データは、所謂”センサーされたデータ（Censored Data）”となる。即ち、得られるデータは、観測期間中に倒産する場合、倒産しない場合だけでなく、他の格付けに移るという場合も考慮しなければならない。この様な、データを扱う一般的な方法として近年盛んに利用されているのが生存分析からの各種方法である。本論文では、その中でも D.R. Cox[1972]により提唱された比例ハーザードモ
デルとその変形を中心に論じていく。第１－３節で論じた様に、比例ハザードモデルでは
ハザード関数をベースラインと各クラスに特徴的な部分との積に分解する。以下ではベースラインについては主にパラメトリックなモデルを仮定する。そこで、次の小節では先ず始めに、パラメトリックモデルと推定量の統計学的性質について解説し、その後所謂 Cox
モデルに於ける統計学的諸性質を Bootstrap 法に基づき求める。
第 3 部でもデータ等に関する記号はこれまで通りとする。以下(1-24,25)式で定義された
比例ハザードモデルの推定について論じていく。参考のためにモデル(1-24),(1-25)を再度示しておく：[比例ハザードモデル]

ある共変量(Covariate)ベクトル x と回帰係数ベクトルβ に対し、

ハザード関数が(1-24)

λ(t:x,β) = λ0(t)exp{β'x},

β= (β1, ...., βp)',x = (x1, ...., xp)'である生存時間モデルを比例ハザードモデルと呼ぶ。ここで、λ0(t)はベースラインハザードと呼ばれる共通項である。簡単な積分計算より生存関数は
S(t,x) = exp{-∫0t λ(u;x)du} = exp{-eβ 'x∫0t λ0(u)du}(1-25)= exp{-∫0t λ0(u)du}exp{eβ 'x } =で与えられる。//- 24 -

S0(t)exp(β'x )第３－１節：回帰係数βの推定：ノンパラメトリック法
Cox の興味はβに関する検定にあるようだが本稿では推定問題に焦点を合わせていく。
この問題と解法に対する Cox の基本的な姿勢は「βに関する検定問題では（尤度比検定等を用いる場合）ベースラインハザード関数λ0 は任意である。極端に言えば倒産時点（failuretime)以外ではそれが恒等的にゼロであっても構わない。倒産時点やセンサーされた時刻のみがβに関する情報を持っている。従ってβに関する推論は倒産時点やセンサー時刻が与えられたという条件のもとで行えばよい」である（cf. Cox [1972])。そこで、以下の”条
件付き尤度関数(Conditional Likelihood Function)”が定義され、それに基づいた推測を行うことが正当化される。推定問題でもベースラインハザード関数とβを別々に推定するならば同様の議論が成立すると考えられる為、本稿でも条件付き尤度関数に基づいた議論を展開する。
[条件付き尤度関数]以下簡単化のため、データ Y = {Y1,Y2, ...., Yn }は全て異なる値を取るものとし（No tied data）、その順序統計量を {Y(1), Y(2), ...., Y(n)}と書く。又、これまでと同様に小文字表現で観測値を表すこととする。δ( i)を時点 Y(i) に対応するセンサー指示関Y(i) に対応するリスク集合、
R(Y(i))={時点 Y(i) 数、x(i) は対応する共変量ベクトルとする。
まで倒産していない（生き残っている）企業の集合} を簡単のため R(i) と書く。この時、任意のセンサーされていない Y(i)に対して
Pr{区間 [Y(i),Y(i) + Δy]で倒産する |R(i)}(3-1)≒ Σj∈ R(i) exp{β'xj}λ0(Y(i) )Δyが成立している事が示される。これより、Pr{会社(i)が時点 y(i)で倒産する | 時点 y( i)で R(i)内の１社が倒産する}(3-2)＝ exp{β'x(i)}/Σj∈R(i) exp{β'xj }そして、センサーされていない全てのデータに対する上式の積 Lc(β)を条件付き尤度関数
と呼ぶ。
(3-3)

Lc(β) ＝Πu exp{β'x(i)}/ Σj∈R(i) exp{β'xj }ただし、Πu はセンサーされていないデータに関する積である。//
一方、センサーされたデータを含む場合はどうなるか？
の側面から考えてみる。

- 25 -

そのために(3-3) の意味を別Y1, Y2, ...., Yn は独立で Yi～Fi 、Ri を Yi のランク（i=1,...,n)とする。ここでもY =(Y1, Y2, ...., Yn), R=(R1, R2, ...., Rn)とし、小文字で観測値（ベクトル）を表す。fi を Fiの密度関数とすれば(3-4)∫･･･∫ Πn
i=1 f(i)(ui)dui

p(r) = Pr{R = r }＝u1 <u2 <... <un
ただし、f(i) は順序統計量 Y(i)の密度関数である。ここで以下の結果が Kalbfleisch andPrentice[1973]により与えられている{} であれば、[定理３－１] もしも、Fi(t) = 1 - exp -exp{β'xi ∫0t λ0(u)du}p(r) ＝Πi=1n exp{β'x(i)}/ Σj∈R(i)exp{βxj }

(3-5)

//

即ち、Cox の条件付き尤度関数はランク統計量の密度（確率）関数に他ならない。そこで、この考えをセンサーされたデータに拡張するために、先ず、センサーされたデータセ
ットにおける順序統計量を以下のように定義する：[定義３－２：センサーされた観測値のランク](Y1,δ1), (Y2,δ2), .... , (Yn,δn)をセンサーされたデータセット、Yi は分布関数 Fi（密度関数 fi）に従うとする。ここで、センサーされたデータに対するランクベクトルを以下のように定義する：
Ru/c = (Ru/c 1, ..., Ru/c n)

(3-6)
ただし、
(3-7)

Ru/c i=�

センサーされていないデータ内の Y のランク δi = 1
直前のセンサーされていないデータのランク δi = 0

さて、δ=(δ1,δ2, .... ,δn)と ru/c = (ru/c 1, ..., ru/c n)を実現値のベクトルとするとき
p(ru/c,δ) = Pr{Ru/c = ru/c, δ=δ}(3-8))=∫･･･∫ [{Πn( u)fu( i)(ui)dui
i=1

u1 <u2 <... <un}{ Π [1 - F (u )]}]du ･･･du
j

j

1

n(u)

j∈Ci,i+1

ただし、u(i) は i 番目のセンサーされていないデータ、n(u)はセンサーされていないデー
タの総数。Ci,i+1 はセンサーされていない i 番目と i+1 番目のデータの間に含まれるセンサ
ーされたデータの集合である。ここでも上と同様の結果が Kalbfleisch and Prentice[1973]- 26 -

により証明されている：{}[定理３－３] もしも、Fi(t) = 1 - exp -exp{β'xi ∫0t λ0(u)du} であれば、(3-9)

p(ru/c,δ) ＝Πu exp{β'x(i)}/ Σj∈R(i) exp{βxj }=Lc (β) //即ち(3-3) で定義された条件付き尤度関数はセンサーされたデータに関してもそのまま用いることが出来るのである！Cox[1972] は上記 Lc(β)を、あたかも通常の尤度関数と同じに扱って良いと予想し、βの MLE も標準的な手法で求められると予想した。実際、この予想は後日、Prentice and Kalibfleisch[1979]等により厳密な形で証明されている。従ってβの最尤推定量を求めるには、スコアーベクトルと標本情報行列を用いた数値計算を行え∂ばよい。具体的には先ず、スコアーベクトル log Lc(β)と標本情報行列 i(β)を定義する。
∂𝛃𝛃(3-10)(3-11)

∂

∂𝛃𝛃

∂

∂

∂

log Lc(β) = [∂β log Lc(β), ∂β log Lc(β), ... , ∂β log Lc(β)]i(β)＝ -

∂21

⎡ ∂β 1 ∂β 1 log Lc(𝛃𝛃) ⋯
⎢
⋮
⋱
⎢ 2
⎢ ∂
log Lc(𝛃𝛃), ⋯⎣∂β 1 ∂p2

p∂2

log Lc(𝛃𝛃) ⎤
⎥
⋮
⎥
2
∂
log Lc(𝛃𝛃)⎥∂β p ∂β p⎦∂β 1 ∂β p∂

第 1 部 5 節のパラメトリックモデルの時と同様 ∂𝛃𝛃log Lc(β) = 0 を解くことでβの最尤推定量を求めるが、ここでも一般的には Newton-Raphson 法が用いられる。
� (1) を初期値（初期予測ベクトル）とする。そして、順次、
先ず 𝛃𝛃(3-12)

� (k) = 𝛃𝛃� (k-1) + i (𝛃𝛃� (k-1)) ∂ log Lc(𝛃𝛃� (k-1))
𝛃𝛃
∂𝛃𝛃

k ≧2

で次の推定値を順次求め、安定的な推定値が得られるまで続ける。
βの推定のあと、ベースラインハザード関数の推定を行うが、本稿では利便性の観点よりパラメトリックモデルを仮定する。従って問題は未知パラメータの推定となり、基本的には第１－５節の方法で推定される。ただし、ここで利用するデータは共変量を含む(Y1,x1,δ1),(Y2,x2,δ2), ..., (Yn,xn,δn)である。(Yi,xi,δi)の実現値 (yi,xi,δi)の尤度関数
Li(θ) = Lθ(yi,xi,δi) は(3-13)- 27 -f(yi , 𝐱𝐱𝐢𝐢 : 𝛉𝛉)Li(θ) = Lθ(yi,xi,δi) = �S(yi , 𝐱𝐱i : 𝛉𝛉)𝐢𝐢𝐢𝐢 𝛅𝛅𝐢𝐢 = 1（uncensored）

if 𝛅𝛅𝐢𝐢 = 0（censored）

=f(yi,x1:θ)δiS(yi,x1:θ)1-δi

但し、(3-14)� 'x)}
S(y,x:θ)=exp{-∫0yλ0(u:θ)du・exp(𝛃𝛃� 'x)・exp{-∫0yλ0(u:θ)du・exp(𝛃𝛃� 'x)}
f(y,x:θ)= λ0(y:θ)・exp(𝛃𝛃� を用いる。仮に、ベースラインとしてワイブルモ
又ここではβの代わりに(最尤)推定値𝛃𝛃デルを仮定すると、S0(t,x1:θ)=exp{-(ϑt)α}、f(t ,x1:θ)=αϑ (ϑt)α-1 exp{-(ϑt)α}、λ(t) = αϑ (ϑt)α-1 より、Li(θ) =

f(yi,x1:θ)δiS(yi,x1:θ)1-δi(3-15)� 'x)・exp{-(ϑyi)α}・exp(𝛃𝛃� 'x)}]δi= [αϑ(ϑt)α-1exp{-(ϑyi)α}・exp(𝛃𝛃�[{exp{-(ϑt)α}}exp ⁡(β ′x )]1-δi

ただし、𝛉𝛉 = (ϑ, α)である。(3-13)よりθに関する尤度関数はL (θ)= L(y1, yi2, ... , yn ;x1,x2, ... ,xn,δ1,δ2, ... ,δn:θ)(3-16)= Πi=1n Lθ(yi,xi,δi) = [Πuf(yi)δi][ΠcS(yi)1-δi]� )を求めることが可能である。
�=θ(𝛃𝛃
これより、１-５節の方法でθの最尤推定量𝛉𝛉- 28 -

第３－２節Bootstrap 法

モデル内のパラメータを推定した後必要なことは推定値の精度をどのように評価するかである。最尤推定量は､一定の正則性の条件下で漸近正規性を持っているが､前節の様な方法でβとベースラインハザード関数のパラメータθを求めた場合、標準的な方法による誤差評価の正当性は疑問である。そこで、ここでも第２－４節同様 Bootstrap 法を適用する。比例ハザードモデルに於ける Bootstrap 法はこれまでにも幾つかの方法が提案されてきているが本節では、①パラメトリック Bootstrap (Cf. Casellas, Tarres and Varona[2006])の拡
張と②Zelterman, Le and Louis [1996] による Beysian Bootstrap 法の変形、の二方法を紹介する。しかしながら、Casellas et. al の方法ではセンサーされたデータに対する Bootstrap は扱えない。そこで本節では Efron[1981]によるランダムセンサーリングの考えを利用する
ことにより、これを可能とした。

第３－２－１節

パラメトリック Bootstrap 法

パラメトリック Bootstrap 法とは、先ず原データより未知パラメータθを推定、推定値を�と記す。次に推定されたパラメータの元で推定された分布関数 F(t)=F(t|𝛉𝛉�)から再標本と
𝛉𝛉なる確率変数を発生させる方法の通称である。通常の Bootstrap が原データから直接再標本するのに対し、ここでは推定された分布関数より新たなデータを生み出すところ�)からの大きさｎのランダムサンプルで t1, t2, ...,が異なる。形式的には、T1, T2, ..., Tｎを F(t|𝛉𝛉

�=θ(t1, t2, ..., tｎ)をθの推定値とする。ここで、tｎをその実現値とする。又、𝛉𝛉�)からの大きさｎのサンプルである。具体的にはｎ再標本列 (T*1, T*2, ..., T*ｎ)は F(t|𝛉𝛉
個の独立な(0,1)の一様乱数(U1, U2, ..., Uｎ)を発生させ

T*i = F-1(Ui｜θ)とすればよい。後の手続きは２-４節と同様である。第１番目の再標本列 (T*( 1)1, T*( 1)2, ...,T*( 1) ｎ )から推
�*( 1)を求める。
定値𝛉𝛉
�*( k)を順次求め(k =1, ... ,以下同様に第ｋ番再標本列(T*( k)1, T*( k) 2, ..., T*(k) ｎ)から推定値𝛉𝛉
�*(1), ... ,θM)、bootstrap 推定値の列{𝛉𝛉}を得る。ここでも上と同様θの Bootstrap 推定

値を以下の式で定義する：(3-17)1� ∗(k)
θ�(Boots) = M ∑M
k=1 θ

Bootstrap バイアス推定値(3-18)

�
Bias(θ�)(Boots) =Bootstrap 分散推定値(3-19)

θ�(Boots) – θ�1
� ∗(k) − θ�(Boots ) )2
Var(θ�)(Boots) = M−1 ∑M
k=1(θ推定量の Bootstrap 推定された分布関数- 29 -

(3-20)1Dist θ� (t)(Boots) = M {# of k such that θ�∗(k) ≤ t}で求められる。

比例ハザードモデルのパラメトリック Bootstrap を行う際に、例えばワイブルハザード関数を仮定すると、倒産時刻 T の分布関数はF(t:x,β|ϑ,α) = [1- exp{-( ϑt)α}exp{β 'x}]

(3-21)となる。U を一様確率変数とするとき倒産確率変数 T* は(3-22)

T* = F-1(U:x,β| �ϑ,α
�)

(

̂

̂

)

1

= −log[1 − U]/[exp{αlog( ϑ) + 𝛃𝛃′𝐱𝐱}] ϑ̂で与えられる。従って、独立な(0,1)上の一様乱数(U1, U2, ..., Uｎ)より、倒産時刻の再標本列(T*( k)1, T*(k) 2, ..., T*(k)ｎ), k = 1,... ,M を発生させればよい。
一方、センサーされた再標本データの作り方は Efron のランダムセンサーリングの考え方を用いる。先ず、センサー時点 C（これは倒産時刻と観測最終時点でセンサーされている！）の分布を推定。ここでもベースラインハザードにはパラメトリックモデルを用いた
比例ハザードモデルを用いると良いだろう。即ち、共変量ベクトル xi を伴った、データ列(Ci, xi, δ'i) i = 1,...., n、に基づき C の分布関数 G(c :x,ψ)の未知パラメータψを推定する、但しδ' は観測最終時点を除きδ' = 1 - δである。次に、上と同様の方法で xi を与えたうえで C*i を発生する。そして、上記{(T*( k)1, T*(k)2, ..., T*(k)ｎ),k=1, ...,M}と組み合わせ
ることによりセンサーされたデータ列{(Y*(k)i ,δ*(k)i ), i = 1, ..., n},
させる、但し、(3-23)- 30 -

Y*( k)i = min{T*( k)i, C*(k)i },δ*( k)i= �∗(k)

1 if Ti
0

∗(k)

if Ti

∗(k)

≦ Ci

∗(k)

＞ Ci

k =1, ..., M を発生

第３－２－２節ノンパラメトリック Bootstrap
ノンパラメトリック比例ハザードモデルでの Bootstrap 標本発生方法は Zelterman,

Le

and Louis [1996]により与えられている。彼等の考え方は以下の通りである。
先ず、ここでは共変量ベクトル x もランダムベクトルと考え、(Y,δ,X)の同時分布を想定する。Bootstrap 標本についてはこれまで通り(Y*,δ*,X*)と記すことにする。再標本変量を発生する際、(Y*,δ*)を先ず求め、次に X*を(Y*,δ*)が与えられた時の条件付き分布より発生させる。
我々の手元にはｎ個のデータ｛(y( i),δi, xi), i=1, ..., n}がある。オリジナルのデータに関してタイ(Tei)を排除することは可能であったが、Bootstrap でタイは必然であるため、以下の議論ではタイを許容する。そこで、先ず、互いに異なるセンサーされていない y の値を
大きさの順に並べたものをy(1) < y( 2) < ・・・ < y(N)(3-24）とする、ここで N は互いに異なる y の数である。さらに、これまでと同様、
di = 時点 y(i) で観測されたセンサーされていない観測値の数、Ri = R(y(i))= 時点 y(i)におけるリスク集合を定義しておく。また、ここでは共変量ベクトル xi も確率ベクトルと仮定、その周辺分布を
ρi = Pr{x* = xi } （ここでは(3-25)

ρi =1/n である）とする。Bootstrap の第一ステップは(y*,δ*)を n 個のペアー {(y1,δ1), ... ,(yn,δn)}から抽出、ただし、当初の仮定よりここではδ*=1 である。次いで共変量ベクトル x を条件付き確率(3-26)

πij =πij(β) = Pr{x* = xj |(y*,δ*)=(y(i),1)}

に従って選択する。詳しい導出プロセスは Zelterman,

Le and Louis[1996]にあるが、ベイズの定理により
πij =πij(β)(3-27)= [exp{β'xj}(δ+ρj - ∑l=1i-1 πljdl)]/[∑k=1nexp{β'xj}(δ+ρk - ∑l=1i-1 πlkdl)]となる、但し、 ρ＋=∑δi = ∑di である。又ここでは ∑i=10 = 0 とする。実装の際には� ）を代入し、i=1
πij(β)は、βにその MLE（乃至は Maximal Partiaal Likelihood Estimate）𝛃𝛃から順次求められる；ρi =1/n に注意すると、例えば、- 31 -

(3-28)�)
π1j =π1j(𝛃𝛃� 'xj}/∑k=1n 𝛃𝛃
� 'xk}= exp{𝛃𝛃�)
π2j =π2j(𝛃𝛃� 'xk}]/[∑k=1nexp{𝛃𝛃� 'xj}(δ+ρk - d1π1j)]
� 'xj}(δ+ρj - d1π1j)]/∑k=1n 𝛃𝛃= [exp{𝛃𝛃となる。

一方、センサーされた標本 (y*(i),0) については
πCij = Pr{x* = xj | (y*,δ*)= (y(i),0)}(3-29)を求めなければならない。さて、(y*,δ*)= (y(i),0)が意味することは、時点 y(i)では倒産が起きたか否かは定かでない。唯一確定的なことは、倒産時刻は時点 y(i)より先と云うことである。従って、倒産時刻を T とすると、
πCij =

Pr{x* = xj | T≧y(i)}(3-30)= ∑l=uNπljdl/∑l=uN dl
但し、u = u(i) はmin yu ≧y(i) となる最小の u である。

最後に２点ほどコメントを付け加えておく。まず、y(n)がセンサーされている場合(3-30)
式は定義出来ない。そこで、πCij をそれまでで一番生存時間の長いデータの中で、センサ
ーされていない観測値と同じ πij とする。又、y(N)がセンサーされていない場合は対応するハザード関数値が１となる。結果として(3-27)のπij が負又は１以上の値を取る可能性が出てくる。この場合も、上と同様の工夫でπij を定義する。

- 32 -

[参考文献]Casellas, J., J. Tarres and L. Varona [2006] "Parmetric Bootstrap for Testing Modrl Fitting in
the Proportional Hazards Framwork: An Application to the Survival Analysis of Bruna Dels
Pirineus Beef Calves" American Society of Anmal Sciences 84, 2609-16
Cox, D.R. [1972] "Regression Models and Life-Tables" Journal ofthe Royal Statistical Society,
Ser. B 34, 187-202
Duffie, D. and K. Singleton [1997] "An Economic Model of the Term Structure of Intereest Rate
Swap Yields," Journal of Finance 52, 1287-1321
Duffie, D. and K. Singleton [1999] "Modseling Term Structures of Defaultable Bonds" Revie ofFinancial Studies 12, 678-720
Efron, B. [1979] "Bootstrap Methods: Another Look at the Jackknife" Annals of Statistics
7, 1-26
Efron, B. [1980] "Censored Data and the Bootstrap" Journal ofthe American Statistical
Association, 76, 312-319
Efron, B. [1982] The Jackknife, the Bootstrap and Other Resampling Planes, SIAM, Philadelphia
Harrison, M. and D. Kreps [1979] "Martingales and and Arbitrage in Multiperiod Security
Markets" Journal of Economic Theory, 20, 381-408.
Harrison, M. and S. PLiska [1981] "Martingales and Stochastic Integrals in the Theory ofContinuous Trading" Stochastic Processes and Their Applications, 11, 215-260.
Heath, D., R. Jarrow, and A. Morton [1992] "Bond Pricing and the Term Structure of Interest
Ratres: A New Methodology for Contingent Claims Valuation" Econometrica, 60, 77-106
Kalbfleisch, J.D. and R.L. Prentice [1973] "Marginal Likelihoods Based on Cox's Regression andLife Model" Biometrika 60, 267-278
Kaplan, E.L. and Meier, P. [1958] "Nonparametric Estimation from Incomplete Observations"
Journal of the American Statistical Association, 53, 457-481
Merton, R. [1974] "On the Pricing of Corporate Dept: The Risk Structure of Interest Rates" The- 33 -

Journal of Finance 29, 449-470
Miller Jr, R. G. [1981] Survival Analysis, John Wiley & Sons New York
Prentice, R.L. and J.D. Kalbfleisch [1979] "Hazard Rate and Models with Covariates,"
Biometrics 35, 25-39
Takahashi H. [2009] “Statistical Analysis of Implied Data” to appear
D. Zelterman, Le, C.T. and T.A. Louis [1996] "Bootstrap Techniques for Proportional Hazards
Models with Censored Observations" Statistics and Computing 6, 191-199
楠岡成雄、青沼君明、中川秀敏 [2001] クレジットリスクモデル, 金融財政事情研究会- 34 -