現代数理統計学の基礎 ktatsuya77 MathStat typo

このファイルでは，「現代数理統計学の基礎」（共立出版）の誤植訂正及び補足説明を与 えています。

「現代数理統計学の基礎」訂正箇所

下の記述において，「-5 行」は，「下から 5 行目」を意味する。

○ 第１章 p.7, -2 行：P (A_{k−1} を P (Ak−1})_{} に訂正}

○ 第１章 p.7, -1 行：(5) を (1) に訂正

○ 第２章 p.21, -9 行：右辺の式を = e^bitφX(at)に訂正

○ 第２章 p.22, -6 行：

fX(x) = ¹ 2π

∫ _∞

−∞

e^−itxφX(t)dt に訂正

○ 第２章 p.26, 8 行：問 1 において，正規化定数 c を C に訂正

○ 第２章 p.28, -7 行：問 18 において，Var(X)− を Var(X) = に訂正

○ 第３章 p.34, 10 行：「優しいので」を「易しいので」に訂正

○ 第３章 p.46, 9 行：= e^s+t/e^s= e^tを = e^−λ(s+t)/e^−λs = e^−λtに訂正

○ 第５章 p.88, -3 行：{φX1^(it/n)}ⁿを {φX1^(t/n)}^nに訂正

○ 第５章 p.98, 12 行：= f(h(U))] を = E[f(h(U))] に訂正

○ 第５章 p.100, 3 行：

→d

g^′′(θ) 2 ^σ

2_{Y を →}

d

g^′′(µ) 2 ^σ

2_{Y に訂正}

○ 第５章 p.100, 6 行：g^′′(θ^∗)_→p g^′′(0)を g^′′(θ^∗)_→p g^′′(θ)に訂正

○ 第５章 p.111, 13 行：X +^√n_{− 1S を X +}^√nSに訂正

○ 第６章 p.119, 2 行：θ⁻ⁿI[x₍₁₎>0]I[X_(n)<θ]を θ⁻ⁿI[x₍₁₎>0]I[x_(n)<θ]に訂正

○ 第６章 p.121, -3 行：この式を

∂

∂θi

L(θ1, . . . , θk|X) = 0 もしくは _∂θ^∂

i

ℓ(θ1, . . . , θk|X) = 0, i = 1, . . . , k に変更。尤度方程式はどちらで定義してもいいです。

○ 第６章 p.131, 3 行：{E[f(X)g(X)}^2]を {E[f(X)g(X)]}^2に訂正

○ 第６章 p.142, -13 行：「そのその平均」を「その平均」に訂正

○ 第６章 p.143, 7 行：問 17 については，Xi, Yiの分散を Var(Xi) = σ₁², Var(Yi) = σ²₂, 共分散を Cov(Xi, Yi) = ρσ1σ2に訂正して下さい。

○ 第 7 章 p.146, 12 行：「・・・求めることになる.」の後に次の文章を挿入する。

「X1, . . . , Xnを未知の平均 µ, 既知の分散 σ²を持つ正規母集団からのランダム標本とす る。特定の µ0に対して」

○ 第８章 p.172, 13 行：C(X) = の式を次のように訂正 C(X) ={θ | a ≤ Q(X, θ) ≤ b}

○ 第８章 p.177, -5 行：Beta(a, p) を Beta(a, b) に変更

○ 第９章 p.191, 13∼14 行：(X^⊤X)⁻¹ = OΛ^−2O⊤, (X^⊤X)^−1/2 = OΛ^−1O⊤, U = OP^⊤z_に訂正。

○ 第９章 p.200, 2 行と 4 行：|bI(x)|^p/2を |bI(x)|^{1/2に訂正。}

○ 第９章 p.204, 13 行：「GLS」を「GLM」に訂正

○ 第１０章 p.230, 2 行：「共変推定量は」の後ろに，「Z の適当な関数 ϕ(Z) を用いて」 という文章を挿入して下さい。

○ 第１０章 p.238，-7 行，-８行，p.239, ４行，５行：exp{−θ^2/(2τ2)}を^{√2π exp}{−θ^2/(2τ2)} に変更。p.238, -６行：^√2πτ を τ に変更。

○ 第１１章 p.253, 1 行：「Step 3. X₃^(k) ∼ π(x2|」を「Step 3. X3^(k) ∼ π(x3|」に訂正。 同様に，3 行：「Step m. Xm^(k)_{∼ π(x}2|」を「Step m. X^m(k)∼ π(xm|」に訂正。

○ 第１１章 p.255, 256：「凡関数」を「汎関数」に訂正。

○ 第１ 2 章 p.278, -10 行：最後の「+t = t」を「+t = s」に訂正。

○ 第１２章 p.282, 6 行：真ん中の等式で，Pi^(X1 ^{= k, X}2 ^{= j)}を P (X1 ^{= k, X}2 ⁼

j_|X0 = i)に訂正。

○ 付録 A.1 p.289, 4 行：(a)k= a(a− 1) · · · (a − k + 1) に訂正

○ 付録 A.2 p.303, -7 行：f(x | µ, σ^{2))の右側の ) を削除}

「現代数理統計学の基礎」補足説明

1 第６章∼第１２章：ベイズ推測についてのまとめ

ここでは，ベイズ推測の概観が理解できるように, 教科書の第６章∼第１２章のベイズ 推測に関する項目をまとめてみる。

1.1 ベイズ推測の基本事項

○ ベイズの定理と事前分布・事後分布

ベイズの定理

ベイズの定理とは，逆向きの条件付き確率を計算する公式のことである。この公式の意 味を乳がんの定期検診を例にとって説明してみよう。

(1) 40歳の女性のうち 1% が乳がんに罹患しており，(2) 乳がんに罹患している女性の うち検査で正しく陽性と判定されるのは 80% で，(3) 乳がんでない女性の 10% は間違っ て陽性と判定されることがわかっているとする。このとき，検査で陽性と判定された女性 が本当に乳がんに罹患している確率を求めたい。

以上の内容を確率を用いて書いてみると，A =(40 歳の女性が乳がんに罹患している), B =(検査で陽性と判定される) が事象となり，(1) より A の確率が P (A) = 0.01 であり， (2)より A を与えたときの B の条件付き確率が P (B | A) = 0.8，(3) より A を与えたとき の B の条件付き確率が P (B | A) = 0.1 と表されることがわかる。ここで A は乳がんに罹 患していない事象のことで，A の補事象と呼ばれる。このとき，求めたい確率は，B を与 えたときの A の条件付き確率 P (A | B) であり逆向きの条件付き確率になる。

条件付き確率の定義から，P (A | B) = P (A ∩ B)/P (B), P (A ∩ B) = P (B | A)P (A) で あるから，P (A | B) = P (B | A)P (A)/P (B) と表される。ここで，B = B ∩ (A ∪ A) = (B∩ A) ∪ (B ∩ A), (B ∩ A) ∩ (B ∩ A) = ∅ であるから，P (B) = P (B ∩ A) + P (B ∩ A) = P (B | A)P (A) + P (B | A)P (A) と書ける。従って，

P (A_{| B) =} ^{P (B| A)P (A)}

P (B | A)P (A) + P (B | A)P (A)

と書け，逆向きの条件付き確率が順方向の条件付き確率を用いて表されることになる。 これをベイズの定理という。乳がんの定期検診の例では，P (A) = 0.01, P (A) = 0.99, P (B | A) = 0.8, P (B | A) = 0.1 を代入することにより，P (A | B) ≈ 7.5% となる。

事前分布と事後分布

乳がんの定期検診の例において 40 歳の女性が乳がんに罹患している確率 P (A) は事前 にわかっている確率である。ベイズの定理は，事前の確率を利用して，検査で陽性と判定 された場合に乳がんに罹患しているという事後的な確率を与えていることになる。こうし

た考え方は，統計モデルの母数の推測において医学知識や経験を事前情報として組み入れ るのに役立つ。

例えば，新薬を n = 20 人に処方し，効果があるときには 1，ないときには 0 として集 計をとったところ, x = 15 人に効果が現れたとする。各人に効果が現れる確率を θ とする と，x は 2 項分布 Bin(n, θ) に従う。この分布を

f (x _{| θ) =}^{( n} x

)θ^x(1_{− θ)}^n−x

と表して尤度関数という。尤度関数を最大にする θ の推定値 ˆθを最尤推定値といい，いま の場合 ˆθ = x/n = 15/20 = 0.75 で与えられる。

ここで，これまでの臨床研究により θ は平均 0.4, 分散 0.1²であることが事前情報とし てわかっているものとする。θ は確率であるから 0 ≤ θ ≤ 1 を満たすので，区間 (0, 1) 上 の分布としてベータ分布 Beta(a, b) を想定してみる。この分布の平均と分散は a/(a + b), ab/_{{(a + b)}²(a + b + 1)} なので，それぞれ 0.4, 0.1²とおいて連立方程式を解くと，a = 9.2, b = 13.8になる。すなわち，θ の分布として，a, b にこれらの値を代入したもの Beta(a, b) が得られる。これを θ の事前分布といい，π(θ) という記号で表す。

次に，θ の推測に事前情報を組み入れる方法を述べよう。上述の確率分布は x| θ ∼f(x | θ) = Bin(n, θ)

θ ∼π(θ) = Beta(a, b)

と記述することができる。これは，2 項分布の母数 θ にベータ分布 Beta(a, b) を仮定して いるので，2 項・ベータ・モデルと呼ばれる。θ に関する推測を行うために，x を与えた ときの θ の条件付き分布 π(θ | x) を求める。これは逆向きの条件付き分布なのでベイズの 定理を連続分布の場合に拡張することによって求めることができる。(x, θ) の同時確率分 布は f(x | θ)π(θ) で与えられるので，x を与えたときの θ の条件付き分布は

π(θ_{| x) =} ^{f (x| θ)π(θ)} f_π(x)

と表すことができる。ここで分母は x の周辺分布であり fπ(x) =^∫ f (x| θ)π(θ)dθ で与え られる。これを θ の事後分布という。いまの場合，θ の事後分布は

π(θ| x) = Beta(a + x, b + n − x)

なるベータ分布に従うことが確かめられる。事後分布は尤度関数に事前分布を組み入れた 確率分布であり，事後分布に基づいて θ の推測を行うことをベイズ的推測という。これに 対して尤度関数のみに基づいた推測を頻度論的推測という。

以上では 2 項・ベータ・モデルを例として扱ってきたが，その他の代表的な例として， 正規・正規・モデル，ポアソン・ガンマ・モデルを紹介しよう。

(正規・正規・モデル) これは x の分布も θ の事前分布も正規分布で与えられるモデルで， 例えば x の分布が平均 θ, 分散 σ²の正規分布 f(x | θ) = N (θ, σ²⁾とし, θ に π(θ) = N (µ, τ²⁾ なる正規分布を考えると，事後分布は

π(θ_{| x) = N}^(σ

2_{µ + τ}2_x

σ²+ τ^{2 ,} σ²τ² σ²+ τ²

)

で与えられる。x の代わりに標本サイズ n の標本平均 x を考えるときには，σ²の代わり に σ²/nを代入すればよい。

(ポアソン・ガンマ・モデル) これは x の分布に平均 nλ のポアソン分布 P o(nθ)，θ の事前分布にガンマ分布 Ga(a, 1/b) を想定したモデル，すなわち f(x | θ) = P o(nθ), π(θ) = Ga(a, 1/b)とすると，事後分布は

π(θ| x) = Ga(a + x, 1/(n + b)) で与えられる。

○ ベイズ推定，信用区間，予測分布

ベイズ流点推定

ベイズ推定量は事後分布の平均やモード，メディアンで定義されることが多い。2 項・ ベータ・モデルにおいて θ のベイズ推定量を事後分布の平均 E[θ | x] で与えると，

θˆπ = ^{a + y} a + b + n ⁼

a + b a + b + n ^·

a a + b ⁺

n a + b + n ^·

x n

と書ける。a/(a + b) は事前分布の平均であり，x/n は事前分布を仮定しないときの最尤 推定量になるので，ベイズ推定量は両者の加重平均で表現できることがわかる。標本サイ ズ n が大きければ x/n に近づき，n が小さければ a/(a + b) の方へ近づく。また a + b が大 きくなれば事前分布の分散が小さくなり平均 a/(a + b) への確信が強くなる。このときベ イズ推定量は a/(a + b) の方へ近づくことがわかる。すなわち，加重平均の重みは標本サ イズ n と事前分布の確信 a + b の大小に基づいて調整されている。これに先ほどの数値を 代入すると，ˆθπ = 0.56となり，x/n = 0.75 と a/(a + b) = 0.4 の中間の値をとっている。

正規・正規・モデルにおいては，θ のベイズ推定量は σ²µ + τ²x

σ²+ τ^{2 =}

1/τ²

1/σ²+ 1/τ^{2µ +}

1/σ² 1/σ²+ 1/τ^2x

となり，x と µ を，x の精度 1/σ²と µ への確信度 1/τ²で内分した形をしている。 ポアソン・ガンマ・モデルにおいては，θ のベイズ推定量は

a + x b + n ⁼

b b + n^·

a b ⁺

n b + n^·

x n

で与えられ，事前分布の平均 a/b と最尤推定量 x/n との加重平均で表されることがわかる。 ベイズ信用区間

信用係数 1 − γ のベイズ信用区間は，事後分布に関してその区間が θ を含んでいる確率 が 1 − γ になる区間として与えられる。例えば，L と U を

∫ U L

π(θ| x)dθ = 1 − γ

を満たすようにとると，L, U は x の関数になるので L(x), U(x) と書くことにする。この とき，区間 [L(x), U(x)] は Pπ^(θ ∈ [L(x), U(x)] | x) = 1 − γ を満たすので，信用係数 1 − γ の信用区間になる。

正規・正規・モデルについては，ベイズ推定量 ˆθπ = (σ²µ + τ²x)/(σ²+ τ²)に対して，θ

の事後分布は _√

σ²+ τ²

√σ²τ^{2 (θ− ˆθπ)∼ N (0, 1)}

なる形に変形できる。zγ/2を標準正規分布の上側 100γ/2% 点，すなわち標準正規分布の 分布関数 Φ(·) に対して 1 − Φ(zγ/2) = γ/2を満たす分位点とし，

Iπ(x) = ^[θ^ˆπ ₋

√σ²τ²

√σ²+ τ^{2zγ/2, ˆθπ +}

√σ²τ²

√σ²+ τ^2zγ/2 ]

なる区間を考える。このとき，θ がこの区間に入る事後確率は P_π(θ_{∈ Iπ}(x)_{| x) = Pπ}⁽_{− zγ/2 ≤}

√σ²+ τ²

√σ²τ^{2 (θ− ˆθπ)≤ zγ/2 | x}

)= 1_{− γ}

となるので，Iπ^(x)は信用係数 1 − γ の信用区間になることがわかる。

正規・正規・モデルは事後分布が対称分布になるので信用区間はベイズ推定量を中心に 両側に同じ幅の区間を作ることができる。しかし，2 項・ベータ・モデルやポアソン・ガ ンマ・モデルは事後分布が対称ではないので同様な方法では信用区間を作ることができ ない。この場合は，両側に事後確率が γ/2 になる分位点をとることにする。例えば 2 項・ ベータ・モデルの場合は

∫ L(x) 0

π(θ| x)dθ = γ/2,

∫ 1 U (x)

π(θ| x)dθ = γ/2

を満たすように L(x), U(x) を求めると，得られる区間 [L(x), U(x)] は信用係数 1 − γ の θ の信用区間となる。

信用区間についての注意として，信頼区間とは別の概念であることがあげられる。正規 分布モデル x | θ ∼ N (θ, σ²⁾における，信頼係数 1 − γ の信頼区間は

CI(x) = [x_{− σzγ/2}, x + σz_γ/2]

であり，これは P (θ ∈ CI(x) | θ) = 1 − γ を満たしている。このことは，例えば CI(x) が信頼係数 95% の信頼区間であるとは，100 回 x の乱数を発生させたときに 95 回は区間 CI(x)が θ を含んでいることを意味しており，x の値が与えられたときには CI(x) は θ を 含むか含まないかのどちらかである。これに対して，信用係数 95% の信用区間 Iπ(x)は， xの値が与えられたとき θ の事後分布に関して θ が区間 Iπ(x)に含まれる確率が 95% であ ることを意味する。

ベイズ予測分布

xを観測可能な変量，y を x とは独立な将来の変量とし，それぞれ f(x | θ), f(y | θ) に 従っているとき，観測値 x と事前情報 π(θ) に基づいて予測分布 f(y | θ) を予測する問題 を考える。このとき，ベイズ予測分布は

fˆπ(y _{| x) =}

∫

f (y| θ)π(θ | x)dθ =

∫

f (y | θ)f(x | θ)π(θ)dθ/f^{π(x) (1)} で与えられる。実際，^{∫ ˆ}fπ(y | x)dy = 1 を満たすので確率分布になっている。

○ 仮説検定とベイズ・ファクター及びモデル比較

仮説検定とベイズ・ファクター 母数 θ に関する仮説検定は，一般に

H0 : θ_{∈ Θ}0 vs. H1 : θ_{∈ Θ}1

なる形で表される。ここで Θ0_{∩ Θ}1 =_{∅, Θ}0_{∪ Θ}1 = Θを満たしており，Θ は母数全体の集 合である。頻度論的な仮説検定では，H0を帰無仮説といい有意水準 α を設け H0を間違え て棄却してしまう確率が α 以下になるように検定手法を構成する。これに対して，ベイズ 流仮説検定は，それぞれの仮説に事前確率 P (Hi)を仮定する。これは P (H0) + P (H1) = 1 であるから，確率の比 P (H0)/P (H1)はどちらの仮説が起こりやすかを事前に与えている ことになる。これを事前オッズ比という。事前オッズ比が 1 であることは両方の仮説が同 等に起こりやすいことを事前に与えていることを意味する。各仮説における θ の事前分布 を π(θ | H0), π(θ_{| H}1)とすると，θ の事前分布は

π(θ) = π(θ _{| H0})P (H₀) + π(θ_{| H1})P (H₁) (2) と書ける。x の尤度が x | θ ∼ f(x | θ) で与えられるとき，事後分布は

π(θ_{| x) =} ^{f (x}| θ)π(θ | H0^{)P (H}0^{) + f (x}| θ)π(θ | H1^{)P (H}1⁾

fπ(x_{| H}0)P (H0) + fπ(x_{| H}1)P (H1)

の形で表される。ここで fπ(x_{| H}i) =^∫_θ∈Θif (x| θ)π(θ | H^i)dθである。このことから，各 仮説 Hiの事後確率は

P (Hi | x) =

∫

θ∈Θⁱ

π(θ_{| x)dθ =} ^{fπ(x| Hi)P (Hi)}

fπ(x_{| H}0)P (H0) + fπ(x_{| H}1)P (H1) と書ける。この事後確率の比

P (H0 | x)/P (H1 | x)

を事後オッズ比といい，この値の大小によりどちらの仮説を選択するかを判断する。すな わち，P (H0 | x)/P (H1 | x) < 1 のとき仮説 H0は棄却される。

ベイズ・ファクターは，事後オッズ比と事前オッズ比に基づいて BF01= ^{事後オッズ比}

事前オッズ比 ⁼

P (H0 _{| x)/P (H}1 _{| x)}

P (H0)/P (H1)

で定義される。これに上で与えられている事後確率を代入し整理すると， BF01= fπ(x_{| H}0)/fπ(x_{| H}1)

となり，ベイズ・ファクターはそれぞれの仮説での周辺確率の比として表される。 ベイズ・モデルの比較とモデル平均

ベイズ・ファクターはいくつかのモデルを比較するのに使われる。いま K 個のモデル M₁, . . . , M_Kが候補として考えられ，これらのモデルを比較したいとする。それぞれのモ デルの事前確率を P (Mi)とすると P (M1) +· · · + P (M^{K) = 1}を満たす。各モデル Miに

対して事前分布 π(θ | Mⁱ⁾を設定すると，モデル Miに対するモデル Mjのベイズ・ファク ターは，検定の場合と同様に考えて

BFij = ^{P (Mi | x)/P (Mj | x)} P (Mi)/P (Mj) ⁼

f_π(x_{| Mi}) fπ(x_{| M}j)

で与えられる。ここで fπ(x _{| M}i) = ^∫_Mif (x | θ)π(θ | M^{i)dθである。通常は，M1を最も} 小さいモデルもしくは最も大きなモデルに固定し，k = 2, . . . , K に対してベイズ・ファク ターの値 BF1jを比較して最小になるモデルを選択する。M1を固定すれば BF1kによるモ デルの比較は周辺確率密度関数 fπ(x _{| M}k)に基づいて比較することに等しい。ベイズ情 報量規準 (BIC) は標本サイズ n を大きくとったときの −2 log fπ(x _{| M}k)の近似的な量と して

BICk =− log f(x | ˆθ^{k) + pklog n}

により与えられる。ここで pkはモデル Mkの母数 θkの次元であり，ˆθkは θkの最尤推定値 である。

ベイズファクターや BIC, AIC などの情報量規準に基づいて最適なモデルを選択し，選 択されたモデルの母数の推定を行うことになる。ここで注意すべきことは，最適なモデル の選択には不確実性が伴うため，誤ったモデルの選択が母数推定に影響を与える可能性が ある点である。そこで，すべての候補モデル M1, . . . , M_Kに関して，各モデル Mkの起こ りやすさを事後確率 P (Mk | x) でウェイトづけした推定量

θˆ^{M A} =

∑K k=1

P (Mk| x)E[θ | Mk, x]

が考えられる。これは，各モデルでのベイズ推定量をモデルの事後確率に関して平均を とったもので，ベイズ・モデル平均と呼ばれる。P (Mk)が k に関して均一のとき，n を大 きくとると周辺確率が BIC に基づいて近似できるので

P (M_k | x) ≈ exp{−BICk^/2}/

∑K k=1

exp_{−BICk/2_}

と書ける。ベイズ・モデル平均のウェイトとしてこの近似値を用いることができる。 ベイズモデルの診断

想定したベイズモデルが観測されたデータに当てはまっているか否かを調べるためにはク ロスバリデーションという方法を用いる。これは，観測された n 個のデータ x = (x1, . . . , xn) から xiを除いたもの x_−i = (x₁, . . . , x_i−1, x_i+1, . . . , x_n)を考え，x_−iに基づいたモデルから 将来の値 x^∗_i を予測する。このとき実際の xiとどの程度近いかを調べることにより，想定 したモデルの妥当性を診断する方法である。ベイズモデルにおいては，x_−iを与えたとき の x^∗_i のクロスバリデーション予測分布は

f (x^∗_{i | x−i}) =

∫

f (x^∗_i | θ)π(θ | x−i^)dθ

によって与えられる。このとき，x^∗_i のところへ観測値 xiを代入したもの f(xi _{| x−i})を用 いてモデルのデータへの当てはまりの良さを調べることができ，この値が大きければ当て はまりが良いと判断できる。

○ リスク最適性からのアプローチ

ベイズ推定量，ベイズファクター，ベイズ流予測分布はリスク最適性のアプローチによ る合理的な手法として導くことができる。x の尤度関数を f(x | θ) とし θ の事前分布を π(θ) とする。推定や検定は x に基づいて θ に関するある種の決定を行うので，この関数を一般 に δ(x) で表し決定方式という。決定方式には間違いに対する損失を伴うので，それを損失 関数 L(θ, δ(x)) で評価することを考える。例えば，点推定の場合，L(θ, δ(x)) = (δ(x) − θ)² で δ(x) が θ からどの程度離れているかを測ることができる。損失関数は x に依存するの で，これを x の確率分布で平均化したもの

R(θ, δ) =

∫

L(θ, δ(x))f (x_{| θ)dx}

をリスク関数という。頻度論的にはこのリスク関数に基づいて決定方式 δ(x) の良さを評 価することなるが，ベイズの立場では θ の確率分布に関して更に平均化したもの

r(π, δ) =

∫

R(θ, δ)π(θ)dθ =

∫ ∫

L(θ, δ(x))f (x| θ)π(θ)dxdθ

を評価することになる。事後分布 π(θ | x) と周辺分布 fπ^(x)を用いると f(x | θ)π(θ) = π(θ_{| x)f}π(x)と書き直せるので，

r(π, δ) =^{∫ { ∫} L(θ, δ(x))π(θ_{| x)dθ}^}fπ(x)dx

と変形することができる。{} の中身を事後リスクといい，それを最小にする δ(x) をベイ ズ決定方式という。

点推定の場合，ベイズ決定方式はベイズ推定量と呼ばれ，2 乗損失関数 L(θ, δ(x)) = (θ_{− δ(x))}²に関しては事後リスクは^∫[_{δ(x)}²− 2θδ(x) + θ^2]π(θ | x)dθ と書けるので，こ れを最初にするベイズ推定量は事後平均 δ(x) =^∫ θπ(θ | x)dθ で与えられることがわかる。

仮説検定の場合，ベイズ決定方式はベイズ検定と呼ばれる。帰無仮説を H0 : θ _{∈ Θ}0，対 立仮説を H1 : θ _{∈ Θ}1とし，Θ0∩ Θ1 =∅ とする仮説検定を考えてみよう。検定方式 δ(x) は H0を棄却するとき δ(x) = 1, H0を受容するとき δ(x) = 0 をとるので，仮説検定の損失 関数は，θ ∈ Θ⁰で δ(x) = 1 のとき，もしくは θ ∈ Θ¹で δ(x) = 0 のときに 1 の値をとり， その他の場合に 0 の値をとる。このとき，事前分布 (2) に対する事後リスクは^∫_θ∈Θ0π(θ_| x)dθ_{· δ(x) +}^∫_θ∈Θ1π(θ| x)dθ · (1 − δ(x)) = P (H1 | x) + {P (H0 | x) − P (H1 | x)} · δ(x) と 書けるので，ベイズ検定は

δπ(x) =

{ 1 (P (H1 _{| x)/P (H}0 | x) ≥ 1 のとき) 0 (P (H1 | x)/P (H0 | x) < 1 のとき) で与えられることがわかる。

予測分布の予測問題については，損失関数としてカルバック–ライブラ情報量が使われ るのが一般的である。観測値 x の分布 f(x | θ) と事前分布 π(θ) から予測分布 f(y | θ) を予 測する問題を考える。予測量を ˆf (y | x) とし，これで f(y | θ) を予測するときの損失関数 としてカルバック–ライブラ情報量^∫[log{f(y | θ)/ ˆ^{f (y} | x)}]f(y | θ)dy を用いると，ベイ ズ・リスクは _{∫ ∫ ∫ {}

log ^{f (y| θ)} f (yˆ _{| x)}

}f (y | θ)dyf(x | θ)dxπ(θ)dθ

と書ける。これを最小にするベイズ予測分布を求めると (1) で与えられる。

1.2 事前分布の設定及び階層ベイズと経験ベイズ

○ 事前分布の設定

ベイズ推測は事前分布の設定の仕方に大きく影響を受けるので，どのように設定する かが重要なポイントとなる。事前分布の設定には様々な方法があり，次のような簡単な 設定で概略を説明してみよう。x の尤度関数を x | θ ∼ f(x | θ) とし，θ の事前分布を θ | λ ∼ π(θ | λ) とする。λ は超母数と呼ばれ，θ と λ は多次元でもかまわないとする。こ の場合，θ の事後分布は π(θ | x, λ) = f(x | θ)π(θ | λ)/f^π(x | λ) で与えられ，x の周辺分 布は fπ(x_{| λ) =}^∫ f (x | θ)π(θ | λ)dθ である。

主観的事前分布

事前分布の母数の値 λ を経験や知識から事前に定めておく設定を主観的事前分布とい う。この場合 θ に関するベイズ推測は超母数 λ の事前の設定から影響をうける。例えば， 正規・正規モデルでのベイズ推定量は

σ²µ + τ²x σ²+ τ^{2 =}

1/τ²

1/σ²+ 1/τ^{2µ +}

1/σ² 1/σ²+ 1/τ^2x

であるが，µ と τ²の値に依存して決まり，特に µ の値の取り方に敏感である。τ²の値を 非常に大きくとることで µ の影響を抑えることもできる。

共役事前分布

事前分布 π(θ | λ) とその事後分布 π(θ | x, λ) が同じ分布族に入るような事前分布を共役 事前分布という。共役事前分布の利点はデータの発生による事後分布の更新過程を同じ分 布族の中で構成することができることにある。例えば，正規・正規・モデルを考えてみる と，データが時系列的に観測されており n 時点で構成された事後分布が N (ˆθ^{π(n), ˆτπ2(n))な} る形であるとする。n + 1 時点で xn+1が観測されると，n 時点での事後分布を事前分布と 考えて n + 1 時点での事後分布を求めると，N (ˆθ^{(n+1)π , ˆτπ2(n+1))}と表され，平均と分散は

θˆ⁽ⁿ⁺¹⁾_π =(σ²θ^ˆ(n)_π + ˆτ_π²⁽ⁿ⁾xn+1)/(σ² + ˆτ_π²⁽ⁿ⁾) ˆ

τ_π²⁽ⁿ⁺¹⁾ =σ²τˆ_π²⁽ⁿ⁾/(σ² + ˆτ_π²⁽ⁿ⁾)

で与えられる。このように，共役事前分布については超母数を更新するだけで事後分布が 得られることになり便利である。2 項・ベータ・モデルやポアソン・ガンマ・モデルにお いては，それぞれベータ分布，ガンマ分布が共役事前分布になるが，共役事前分布自体そ れほど多くない。

無情報事前分布

主観的事前分布は超母数 λ の値の取り方に影響を受けると述べたが，このことは解析者が 恣意的に解析結果を操作する可能性があることを示唆する。そこで無情報な事前分布が考 えられる。例えば，x の確率密度関数が位置母数 θ と尺度母数 σ が入った関数 σf((x−θ)/σ) の形をしているときには, θ と σ の代表的な無情報事前分布は

π(θ) = 1, π(σ) = 1/σ

で与えられる。これらは位置変換や尺度変換に関して不変であるという性質をもつ。位 置・尺度母数を持つ確率分布は特別な構造であり，一般にはジェフリーズの事前分布が用 いられる。x の確率関数もしくは確率密度関数 f(x | θ) に対してフィッシャー情報量は,

I(θ) = E^{[{ d}

dθ ^{log f (x| θ)} }2]

で与えられるが，ジェフリーズの事前分布は

π^J(θ) =^√_|I(θ)|

で定義される。θ が多次元のときには I(θ) はフィッシャー情報量行列になり |I(θ)| は行 列式の絶対値になる。例えば，ベルヌーイ分布 Ber(p) の母数 p のジェフリーズの事前分 布は π^J(p) = 1/^√p(1− p) であり Beta(1/2, 1/2) に対応している。2 項分布 Bin(n, p) の ジェフリーズ事前分布も同じ形をする。ポアソン分布 P o(λ) のジェフリーズ事前分布は π^J(λ) = 1/^√λとなる。

こうして得られる無情報事前分布は^∫ π(θ)dθが発散してしまい確率分布にならない場 合が多いことに注意する。^∫ π(θ)dθ =∞ となる事前分布を非正則な事前分布という。こ れに対して主観的事前分布のように^∫ π(θ)dθ <∞ を満たすものを正則な事前分布と呼ん でいる。非正則な事前分布を扱う上で大事な点は事後分布が確率分布になることであり， 事後リスク関数が存在していれば最適解を求めることができる。点推定のときにはこれを 一般化ベイズ推定量という。例えば，正規分布モデル x | θ ∼ N (θ, σ^{2)において無情報事} 前分布 π(θ) = 1 を用いると，事後分布 π(θ | x) は N (x, σ²⁾になり θ の一般化ベイズ推定 量は x で与えられることがわかる。 

○ 階層ベイズと経験ベイズ

ベイズ推測の応用上の有用性は，事前分布に知識や経験に基づいた階層構造を組み入れ ることによりデータを説明する豊かなモデルを作ることができる点である。例えば，正規 分布の分散に逆ガンマ分布を仮定すると t-分布のような裾の厚い分布が得られ，更にその 逆ガンマ分布のパラメータに分布を仮定すると裾の厚さを調整してくれるようになる。

階層的事前分布を考える別の利点は，主観的事前分布において問題となった解析者の恣 意性を緩和することができる点である。このようにベイズ解析に客観性を持たせるアプ ローチを客観的ベイズといい，事前分布の超母数に関してベイズ推測が有界になるときロ バスト（頑健）ベイズと呼んでいる。経験ベイズと呼ばれる手法もこの方向性を指向して いるので，階層ベイズと併せて以下で説明する。

階層的事前分布

階層的事前分布は多段階の階層構造をもつ事前分布で，例えば 2 段階の簡単な階層的事 前分布をもつモデルは次のように表すことができる。

x| θ ∼f(x | θ) θ_{|λ ∼π(θ|λ)}

λ _∼ψ(λ)

(3)

π(θ) = ^∫ π(θ | λ)ψ(λ)dλ と書けるので θ に事前分布を想定することに帰着できる。逆に π(θ)を上のような階層構造に分解することができれば，θ 及び λ の事後分布がよく知られ

ている分布で表されるときには，後述するギブス・サンプリングを用いて容易にベイズ推 測を行うことができる。このように補助変量を加えることにより数値計算を容易にする方 法を拡大法という。また客観的ベイズ推測やロバスト・ベイズ推測の視点からは，1 段階 目の事前分布 π(θ | λ) はより正確な分布を与え，2 段階目の事前分布 ψ(λ) はより曖昧な 分布（例えば無情報事前分布）を与えることが望ましいと考えられている。

例えば，x1, . . . , xpが互いに独立に分布し

xi _{| λ}i _∼f(xi _{| λ}i) = P o(λi),

λi _{| b ∼π(λ}i | b) = Ga(a, b) ⁽⁴⁾ に従っているとし，a は正の既知の値とする。このとき，x = (x1, . . . , xp)とおくと，λiの ベイズ推定量は ˆλi(b) = E[λi | x] = {b/(b + 1)}(a + xi)となり，超母数 b の取り方の影響 を大きく受けることになる。そこで b に

b_{∼ π(b) ∝ b}^α−1/(1 + b)^α+β

となる分布を仮定すると，λiの階層的ベイズ推定量は，x =^∑p_i=1x_i/pに対して ˆλ_i(α, β) = E[λ_{i | x] =} ^{px + α}

px + pa + α + β^{(a + xi)}

と書ける。α, β の取り方に影響を受けるもののベイズ推定量 ˆλi(b)のときよりも緩和され ていることがわかる。さらに，α = β = 0 とおいてみると，ˆλ^HB_i ={x/(x + a)}(a + xi)と なり，超母数の影響を取り除くことができる。この場合，b の事前分布は π(b) = 1/b とな り，無情報事前分布になっている。

経験ベイズ法

主観的事前分布のところで注意したように，事前分布の超母数の設定はベイズ推測に敏 感に反映される。そこで，超母数を未知母数としてこれをデータから推定することによっ て事前分布の設定び客観性をもたせることが考えられる。これを経験ベイズ法と呼んで いる。

具体的には，モデル (3) において超母数 λ に分布を仮定する代わりに λ を未知母数と して扱う。この λ を x の周辺分布 fπ^(x | λ) = ^{∫ f (x} | θ)π(θ | λ)dθ から最尤法などの 方法で推定し推定量 ˆλ を求める。この推定量 ˆλ を主観的ベイズ推測法の中に現れる λ の ところに代入することによって経験ベイズ推測手法が得られる。例えばベイズ推定量が θˆ_π(λ) = E[θ| x, λ] なる形で与えれるとき θ の経験ベイズ推定量は ˆθπ^{(ˆλ) = E[θ}| x, ˆλ] とな る。また事前分布を π(θ | ˆλ) により推定することもできる。

例えば，モデル (4) において b を未知母数としてみる。xi の周辺分布は負の２項分布 になるので，x = (x1, . . . , x_p)の同時周辺分布から b の最尤推定量を求めると，ˆb = x/a となる。これをベイズ推定量 ˆλi(b)に代入すると，得られる経験ベイズ推定量は ˆλ^EB_i = {x/(x + a)}(a + xi)となる。これは b をデータから推定することによって b の取り方の恣 意性を排除していることがわかる。この例では，ˆλ^HB_i と ˆλ^EB_i とが一致しており，このよ うな推定量はベイズ経験ベイズ推定量と呼ばれる。

1.3 マルコフ連鎖モンテカルロ法

階層的事前分布を組み入れてベイズモデルを作り事後分布を求めようとすると，よく 知られている分布以外は容易に求めることができない。また事後分布の平均を求めるに は多重積分を計算する必要があり，モデルが複雑になるにつれて解析的に求めるのは困 難になる。そこで，数値的に事後分布を求めるための方法がマルコフ連鎖モンテカルロ (MCMC)法であり，その代表がメトロポリス・ヘイスティング法とギブス・サンプリン グ法である。まず確率分布から乱数を発生させる方法について説明しよう。

○ 乱数の発生法

区間 [0,1] 上の一様分布に従う一様乱数や正規分布に従う正規乱数についてはソフトウェ アに用意されている。確率分布に従う乱数を一様乱数から構成する原理的な方法が以下で 与えられる。

確率積分変換

連続型確率変数の場合に分布関数 F (x) = P (X ≤ x) の形がわかっていれば，この逆関 数 F⁻¹(·) を用いて分布 F (·) からの乱数を発生させることができる。

Step 1. 一様乱数 U ∼ U(0, 1) を発生させる。 Step 2. X = F⁻¹(U )をおく。

このとき F^′(x) = f (x)とおくと，X ∼ f(x) に従う。例えば, 指数分布 f(x) = e^{−xからの乱} 数を発生させたい場合には，F (x) = 1−e^−xより 1−e^{−X = U}を解いて，X = − log(1−U) とおけばよい。

確率変数 X が離散型で x1 < x2 <· · · < x^kに値をとる場合には，次のようにして離散 分布からの乱数を発生させることができる。

Step 1. U ∼ U(0, 1) を発生させる。

Step 2. F (x_i−1) < U _{≤ F (x}i)ならば，X = xiとおく。

例えば X ∼ Bin(2, 1/2) の場合には，U ∼ U(0, 1) に対して X は次のようになる。

X =







0 (0 < U _{≤ 1/4 のとき)} 1 (1/4 < U _{≤ 3/4 のとき)} 2 (3/4 < U _{≤ 1 のとき)}

受容・棄却法

いま確率分布 π(x) からの乱数を発生させたいとする。π(x) のサポートを含む確率密度 関数 g(x) をとり，M を M = maxx{π(x)/g(x)} で定義し有限であるとする。

Step 1. g(x)から乱数 x^∗を発生させる。また U ∼ U(0, 1) を発生させる。

Step 2. U _{≤ π(x}^∗)/_{Mg(x^∗)_{} ならば x}^∗を π(x) からの標本として受容して X = x^∗と おき，そうでなければ棄却して Step 1 へ戻る。

このとき X ∼ π(x) となる。

π(x)からの乱数発生方法がわからなくても g(x) からの乱数発生法がわかっていれば g(x) からの乱数に基づいて π(x) からの乱数を発生させることができる。ただし M の値が大 きくなると棄却する割合が大きくなり非効率なサンプリング方法になってしまう。特に， M < ∞ という制約は重要で，提案分布の密度 g(x) が目標分布の密度 π(x) より分布の裾 が厚くなる必要がある。例えば π(x) ∼ N (0, 1) の場合には g(x) としてコーシー分布をと ることができるが，π(x) がコーシー分布の場合には候補密度 g(x) を与えることができな い。この場合は次の節で述べるメトロポリス・ヘイスティング法が使われる。

例えば a ≥ 1, b ≥ 1 なるベータ分布 Beta(a, b) から乱数を発生させたい場合を考えよう。 π(x)は Beta(a, b) の確率関数であり，g(x) として一様分布 U(0, 1) の確率関数 g(x) = 1 を とると，M = max_0≤x≤1x^a−1(1_{− x)}^b−1/B(a, b)となり，受容・棄却法は次のようになる。

Step 1. U ∼ U(0, 1), V ∼ U(0, 1) を独立に発生させる。

Step 2. U ≤ π(V )/M ならば V を π(x) からの標本として受容して X = V とおき，そ うでなければ棄却して Step 1 へ戻る。

重点サンプリング

ある関数 h(x) の積分 H = ^∫ h(x)dxを計算する際に確率分布からの乱数が利用できる。 g(x)を乱数発生が可能な確率密度関数とし h(x) のサポートを含むものとする。

H =

∫

h(x)dx =

∫ _h(x)

g(x)^{g(x)dx = Eg}

[h(X) g(X) ]

と書けるので，h(x) の積分は h(x)/g(x) の確率密度関数 g(x) に関する期待値として表さ れることになる。従って次のようにして積分を計算できる。

Step 1. g(x)から n 個の乱数 x1, . . . , x_nを発生させる。

Step 2. bH = n^−1∑n_i=1h(x_i)/g(x_i)として積分^∫ h(x)dxを推定する。

重点サンプリングでは推定精度が g(x) の取り方に依存しており，精度を高める様々な工 夫が提案されている。

○ メトロポリス・ヘイスティングス (MH) 法

いま確率密度関数 π(x) から乱数を発生させたい場合を考える。これを目標分布という。 π(x)から直接乱数を発生させることができないため提案分布の密度 q(x, y) を考えてこの 密度から乱数を発生させることを考える。q(x, y) は^∫ q(x, y)dy = 1を満たしており，本 来ならば条件付き密度 q(y | x) の形で表すべきものであるがマルコフ連鎖との関係から通 常 q(x, y) と表記する。

x0を初期値として与え，以下 x_k−1が与えられたとする。 Step 1. q(x_k−1, y)から ykを発生させ,

α(x_k−1, yk) = min^{1, ^{πu(yk)q(yk, xk−1)} πu(x_k−1)q(x_k−1, yk)

}

を計算する。ただし π(x_k−1)q(x_k−1, yk) = 0のときには α(x_k−1, yk) = 0とする。また πu(x) は π(x) において正規化定数を省いたものである。

Step 2. U ∼ U(0, 1) を発生させ，U ≤ α(xk−1^{, yk)なら ykを受容して xk = ykとし，}

U > α(x_k−1, yk)なら ykを棄却して xk= x_k−1とする。k を k + 1 として Step 1 に戻る。

このとき乱数の系列 {X^k, k = 1, 2, . . .} が構成でき，大きな k に対して x^k ∼ π(x) が成 り立つ。

乱数の最初の部分は初期値に依存するので使用せず，それ以降発生する乱数を用いる。 提案密度として代表的なものは酔歩連鎖と独立連鎖と呼ばれるもので，それぞれ q(x, y) = f (y− x), q(x, y) = f(y) なる形で与えられる。

例えば，確率密度関数 π(x) = Ce^−x4(1 +_|x|³)から乱数を発生させることを考える。こ こで C = ^∫_−∞^∞ e^−x4(1 +_|x|³)dxである。提案密度として y|xk−1 ∼ N (xk−1^{, 1), すなわち}

q(x_k−1, y) = (2π)⁻¹e^{−(y−xk−1)2/2}なる酔歩連鎖を考えると，

α(x, y) = min^{1, e^−y4+x4(1 +_|y|³)/(1 +_|x|³)^} (5) となるので，メトロポリス・ヘイスティングス法は次のようになる。

初期値 x0を与え，以下 x_k−1が与えられているとする。 Step 1. yk ∼ N (xk−1^{, 1), U} ∼ U(0, 1) を発生させる。

Step 2. U _{≤ α(xk−1}, yk)なら xk = ykとし，U > α(x_k−1, yk)なら xk = x_k−1として， Step 1に戻る。

このとき大きな k に対して xkは π(x) からの乱数とみなすことがきる。

○ ギブス・サンプリング法

発生させたい変数が m 個あり，k 回目に発生する乱数を x^(k) = (x^(k)₁ , . . . , x^(k)m )，x^(k) から j 番目の元を除いたものを x^(k)_−j = ((x^(k)₁ , . . . , x_j−1^(k) , x^(k)_j+1, . . . , x^(k)m )とする。また x = (x1, . . . , xm)に対して x_−j を同様に定義する。確率密度関数 π(x) から乱数を発生させる ためのギブス・サンプリング法は次のようなアルゴリズとして与えられる。まず，x_−jを 与えたときの xj の条件付き確率密度関数 π(xj_|x−j)とそれからの乱数の発生法がすべて の j = 1, . . . , m について与えられているとする。

初期値 x⁽⁰⁾を与え，以下 x^(k−1) = (x^(k−1)₁ , . . . , x^(k−1)m )が与えられているとする。 Step 1. x^(k)_{1 ∼ π(x}1_|x^(k−1)₂ , . . . , x^(k−1)m )を発生させる。

Step 2. x^(k)_{2 ∼ π(x}2|x^(k)1 ^{, x(k−1)}3 , . . . , x^(k−1)m )を発生させる。 Step 3. x^(k)_{3 ∼ π(x}2_|x^(k)₁ , x^(k)₂ , x^(k−1)₄ , . . . , x^(k−1)m )を発生させる。 以下同様にして

Step m. x^(k)m _{∼ π(x}2|x^(k)1 , . . . , x^(k)_m−1)を発生させる。以上から x^(k) = (x^(k)₁ , . . . , x^(k)m )が得 られるので，k を k + 1 として Step 1 に戻る。

このとき乱数の系列 {x^(k), k = 1, 2, . . .} が構成でき，大きな k に対して x^(k)j ∼ π(xj), j = 1, . . . , m, が成り立つ。

ベイズ階層モデル (3) について，ギブス・サンプリングを利用した事後分布からのサン プリングを構成してみよう。(x, γ) を与えたときの θ の条件付き分布，(x, θ) を与えたとき の γ の条件付き分布は，

π(θ_{|x, γ) =}∫ ^{f (x|θ)π(θ|γ)} f (x|θ)π(θ|γ)dθ^, π(γ_{|x, θ) =∫} ^{f (x}|θ)π(θ|γ)ψ(γ)

f (x|θ)π(θ|γ)ψ(γ)dγ ^{= π(γ|θ),}

と書ける。いま，このような条件付き分布がわかっていてその分布に従う乱数を発生させ ることができるとする。

初期値として θ0, γ0を決める。k = 1, 2, . . . , M に対して次の要領で乱数を発生させる。 Step 1. θ_{|x, γk−1 ∼ π(θ|x, γk−1})から乱数 θkを発生させる。

Step 2. γ_{|x, θ}k _{∼ π(γ|x, θ}k)から乱数 γkを発生させる。k を k + 1 にして Step 1 へ戻る。 このとき，大きな k に対して，θk ∼ π(θ|x), γ^k ∼ π(γ|x) となる。M を大きくとると， E[h(θ)_{|x] は M}^−1∑M_k=1h(θk)により推定することができる。

例えば，ポアソン・ガンマ階層モデルを考える。x|λ ∼ P o(λ), λ|b ∼ Ga(a, b), b⁻¹ ∼ Ga(k, τ )とし，a, k, τ は既知の値とする。このとき，条件付き分布は

λ|x, b ∼π(λ|x, b) = Ga(a + x, b/(1 + b)) b⁻¹_{|x, λ ∼π(b}⁻¹|x, λ) = Ga(a + k, τ/(1 + λτ)) と書けるので，ギブス・サンプリングは次のようになる。

初期値として λ0^{, b}0を決め，k = 1, 2, . . . , M に対して乱数を発生させる。 Step 1. λ_{|x, bk−1} ∼ Ga(a + x, bk−1^{/(1 + b}k−1^{))から乱数 λkを発生させる。}

Step 2. b⁻¹_{|x, λ}k ∼ Ga(a + k, τ/(1 + λ^{kτ ))から乱数 b−1}k を発生させる。k を k + 1 にし て Step 1 へ戻る。

このとき大きな k に対して λk ∼ π(λ|x), b⁻¹k ∼ π(b⁻¹|x) となり，E[h(λ)|x] はM^−1∑Mi=1^h(λi)

によって推定できることになる。

2 「多次元確率変数の分布」に関して

2.1 ２次元のウィシャート分布について

2次元の確率変数 (X, Y )^⊤ が平均がともに 0 で共分散行列 Σ をもつ２変量正規分布 N^{2(0, Σ)} に従うとき，その同時確率密度関数は，教科書の (4.23), (4.24) で与えられて いるように

f (x, y) = ¹ 2π

1

|Σ|^{1/2 exp} {−¹

2^{(x, y)Σ}

−1_{(x, y)}⊤^}

と書ける。いま，(X1, Y1)^⊤, . . . , (Xn, Yn)^⊤_{が N}2(0, Σ)からのランダム標本とするとき，

V ₌

∑n i=1

(Xi

Yi

)₍

X_i Y_i⁾

の従う分布を自由度 n の Wishart 分布といい，V ∼ W^{2(n, Σ)}と書く。1 次元のときには， この分布は自由度 n のカイ２乗分布になるので，Wishart 分布はカイ２乗分布を多次元へ 拡張したものになっている。ここでは，この同時確率密度関数の導出を行ってみる。

まず，Σ が単位行列 Σ = I の場合を考えよう。すなわち，(X1, Y1)^⊤, . . . , (Xn, Yn)^⊤, i.i.d.

∼ N2(0, I)とする。

X ₌



 X1

... X_n



 , ^{Y =}



 Y1

... Y_n





とおくと，V は

V ₌

(X^⊤X X^⊤Y Y^⊤X Y^⊤Y

)

と表される。ここで，∥X∥ = (X^{⊤X)1/2と定義し，} t11=_{∥X∥, t}22= ^(I ₋ ¹

∥X∥^{2X X}

⊤⁾_{Y , t}

21 = ¹

∥X∥^X

⊤_Y

とおく。このとき，次の性質が成り立つ。

補題. Σ = I のとき，t²₁₁, t²₂₂, t21は独立に分布し，それぞれ t²_{11∼ χ}²_n, t²_{22 ∼ χ}²_n−1, t21 ∼ N (0, 1) に従う。

（証明） t²₁₁ = X^⊤X _{∼ χ}²_nであり，X = x を与えたときの t²₂₂の条件付き分布は， I _{− ∥X∥}−2X X^⊤が巾等行列であることから

t²₂₂= Y^⊤(I ₋ ¹

∥X∥^{2X X}

⊤⁾_{YX = x ∼ χ}²

n−1

となることがわかる。この条件付き分布が条件 X = x に依らないことから，t22は X も しくは t11に独立に分布する。また，X を与えたときの t21の条件付き分布は

t21= ¹

∥X∥^X

⊤_YX = x ∼ N (0, 1)

となり，この条件付き分布も条件 X = x に依らないので，t21は X もしくは t11に独立に 分布することがわかる。t²₂₂と t21が独立になることは，任意の可測集合 A, B に対して

P (t²_{22∈ A, t}21∈ B) = P [P (t²22 ∈ A, t²¹∈ B|X)]

において，X を与えたとき Y^⊤(I_{− ∥X∥}^{−2X X⊤})_{と ∥X∥}^−1X⊤Y とが条件付き独立にな るので，

P (t²_{22∈ A, t}21∈ B|X) = P (t²22∈ A|X) × P (t²¹∈ B|X)

となる。上で注意したように，条件付き分布 t²_{22|X, t}21|X が条件 X に依らないので，P (t²22∈ A_{|X) = P (t}²_{22∈ A), P (t}21∈ B|X) = P (t²¹ ∈ B) と書けるので，結局，

P (t²_{22 ∈ A, t}21 ∈ B) = P (t²22 ∈ A) × P (t21 ∈ B)

となり，独立になることが示される。 _□

この補題から，t²₁₁, t²₂₂, t21の同時確率密度関数は f (t²₁₁, t²₂₂, t21) = ¹

2ⁿΓ2(n/2)^t

n−211 ^tn−322 ^exp{−¹ 2^(t

2

11^{+ t2}21^{+ t2}22⁾}

と書けることがわかる。ここで，Γ2(n/2) = ^√πΓ(n/2)Γ((n− 1)/2) であり，多変量ガン マ関数と呼ばれる。行列 T を下三角行列

T ₌

(t₁₁ 0 t21 t22

)

とおくと，t²₁₁+ t²₂₁+ t²₂₂ = tr (T T^⊤)と書けるので， f (t²₁₁, t²₂₂, t21) = ¹

2ⁿΓ2(n/2)^t

n−211 ^tn−322 ^exp{−¹ 2^{tr (T T}

⊤_)}

となる。ここで，

T T^⊤ ₌

( t²₁₁ t11t21

t11t21 t²₂₁+ t²₂₂ )

と書けており，t²₁₁ = X^⊤X, t11t21= X^⊤Y, t²₂₁+ t²₂₂= Y^⊤Y となることから

V ₌

(v11 v12

v12 v22

)

=

(X^⊤X X^⊤Y Y^⊤X Y^⊤Y

)

=

( t²₁₁ t11t21

t11t21 t²₂₁+ t²₂₂ )

となる。そこで, v11= t²₁₁, v12=^√t₁₁² t21, v22= t²₂₁+ t²₂₂なる変数変換を行うと，ヤコビア ンは

^∂(v_∂(t¹¹2 ^{, v12, v22)} 11^{, t}21^{, t2}22⁾

  = det





1 0 0

t₂₁/(2^√t²₁₁) ^√t²₁₁ 0 0 2t21 1



 = t¹¹

となることがわかる。ヤコビアンを組み込んで，(v11, v22, v21)の同時確率密度関数を書い てみると，V の行列式が |V | = |T T^⊤| = |T |^{2 = t2}11^t222^{と書けることから}

f (v11, v22, v12) = ¹ 2ⁿΓ2(n/2)^t

n−311 ^tn−322 ^exp{−¹

2^{tr V} =} 1

2ⁿΓ2(n/2)^{|V |}

(n−3)/2_exp{−¹

2^{tr V}} と表されることがわかる。従って，Σ = I のときの V ∼ W2^{(n, I)}の確率密度関数は

f (V ) = ¹

2ⁿΓ2(n/2)^{|V |}

(n−3)/2_exp{−¹

2^{tr V}} で与えられる。

最後に，W ∼ W^{2(n, Σ)}の確率密度を導こう。V が下三角行列 T を用いて T T^⊤と表さ れたのと同様にして

A₌

(a11 0 a21 a22

)

なる行列を用いて Σ⁻¹ = AA^⊤と書くことができる。tr [Σ^−1W] = tr [AA^⊤W] = tr [A^{⊤W A}] と書けるので，W を

V _{= A}^⊤W A_, W ₌

(w11 w12

w12 w22

)