行列変数の積分とウィシャート分布の正規化についてのメモ
基本的の元ネタはMuirheadのAspects of Multivariate Statistical Theoryです。
イントロダクション
統計学で現れる行列を引数にもつ確率密度であるウィシャート分布は$ D \times D $の正定値行列$ X, \Sigma $について
\begin{align}
& W( X | \Sigma , P ) := \frac{ \det(X)^{\frac{ P - D - 1 }{2}} \exp(- \frac{1}{2} \tr( \Sigma^{-1} X ) ) }{ \det( 2 \Sigma )^{P/2} \Gamma_D \qty(\frac{P}{2}) } \\
& \Gamma_D \qty(\frac{P}{2}) := \pi^{D(D-1)/4} \prod_{d=1}^{D} \Gamma \qty( \frac{ P - d - 1 }{2} )
\end{align}
と定義されます。これが多変数のガウス分布の性質を調べる中で登場しますが、その導入は極めて天下りかつ唐突である場合が多いです。ただし、性質としては$P$個の正規分布に従う独立なベクトルのディアド積の分布に等しいため、計算の上では以下のようなベクトル$\vb*x$を導入します。
\begin{align}
& X = \sum_{p=1}^P \vb*x_p \vb*x_p^T \qc \vb*{x}_p \sim N( \vb*x | 0, \Sigma )
\end{align}
これを用いて、例えば期待値は以下のようになることなどが紹介されます。
\begin{align}
& \mathbb{E} \qty[ X ] = \int \dd X ~ X W( X | \Sigma , P ) \to \mathbb{E} \qty[ \sum_{p=1}^P \vb*x_p \vb*x_p^T ] = \sum_{p=1}^P \mathbb{E} \qty[ \vb*x_p \vb*x_p^T ] = P \Sigma
\end{align}
ここで本来の期待値計算には行列の積分$ \dd X $が登場しています。多くの人はこれが具体的にどのような計算をしているのか、そもそも確率密度分布の規格化係数$ \det( 2 \Sigma )^{P/2} \Gamma_D \qty(\frac{P}{2}) $がどのように登場したのかも全くわからないと思います。この資料では、性質云々ではなくその解説を行なうのが目標です。
基本的な流れは以下の通りです、知っている部分については適宜飛ばして読んでもらった方が良いです。
ただし、数学的な構成を行なっていない部分もございますので悪しからず。
外積
外積は内積とは異なり、同じベクトルが$0$となる（分配法則を満たす）二つのベクトルの積です。
\begin{align}
& \vb*a \wedge \vb*a := 0
\end{align}
本来はこの積が棲む空間をどのように名づけるかなどを議論する必要がありますが、ここではそういうものとしてラフに取り扱うことにします。この性質からわかる重要な性質は
\begin{align}
& 0 = ( \vb*a + \vb*b ) \wedge ( \vb*a + \vb*b ) = \vb*a \wedge \vb*a + \vb*b \wedge \vb*a + \vb*a \wedge \vb*b + \vb*b \wedge \vb*b = \vb*b \wedge \vb*a + \vb*a \wedge \vb*b
\end{align}
つまり外積は反交換する二つのベクトルの積であることです。
\begin{align}
& \vb*a \wedge \vb*b = - \vb*b \wedge \vb*a
\end{align}
この計算にさらに結合則を入れてより多くのベクトルの積を考えます。
\begin{align}
& \vb*a \wedge \vb*b \wedge \vb*c = - \vb*a \wedge \vb*c \wedge \vb*b = \vb*c \wedge \vb*a \wedge \vb*b = - \vb*c \wedge \vb*b \wedge \vb*a = - \vb*b \wedge \vb*a \wedge \vb*c = \vb*b \wedge \vb*c \wedge \vb*a
\end{align}
定義から明らかなように、元の$a,b,c$の順番から互換をどれだけ行なったかに対応した符号がついています。これは置換と互換の一般論から、置換を表す関数$\epsilon(\sigma)$とその符号を用いて
\begin{align}
& \vb*v_1 \wedge \vb*v_2 \wedge \vb*v_3 = \epsilon(\sigma) \vb*v_{ \sigma(1) } \wedge \vb*v_{ \sigma(2) } \wedge \vb*v_{ \sigma(3) } \\
& \vb*v_1 = \vb*a \qc \vb*v_2 = \vb*b \qc \vb*v_3 = \vb*c
\end{align}
となることを表しています。この議論は一般化して書けば以下の通りです。
\begin{align}
& \Omega(\vb*v) := \bigwedge_{i=1}^n \vb*v_{i} := \vb*v_1 \wedge \vb*v_2 \cdots \wedge \vb*v_n = \epsilon(\sigma) \vb*v_{ \sigma(1) } \wedge \vb*v_{ \sigma(2) } \cdots \wedge \vb*v_{ \sigma(n) }
\end{align}
外積にまつわる重要な話題は多々ありますが、ここでは面積との関係を紹介します。外積はある$ n $個の$ n $次元ベクトル$ \qty{ \vb*{v}_i } $で作られる（一般化された）面積として導入することができます。例えば三次元では
\begin{align}
& \vb*a \wedge \vb*b \wedge \vb*c = \sum_{ijk=1}^3 a_{i} b_{j} c_{k} \vb*e_{i} \wedge \vb*e_{j} \wedge \vb*e_{k} = \sum_{ijk=1}^3 a_{i} b_{j} c_{k} \epsilon_{ijk} \vb*e_{1} \wedge \vb*e_{2} \wedge \vb*e_{3} = \qty( \vb*a \cdot \qty( \vb*b \times \vb*c ) ) \vb*e_{1} \wedge \vb*e_{2} \wedge \vb*e_{3}
\end{align}
という意味です。
$$ \vb*a \cdot \qty( \vb*b \times \vb*c ) $$
は三重積と呼ばれる三つのベクトルで囲まれる平行六面体の体積を表す式となっています。ただしここではレヴィティビタ記号を導入しその性質を用いました
\begin{align}
& \epsilon_{ \sigma_1 \cdots \sigma_n} =
\begin{cases}
0 & \sigma \notin S_n \\
\epsilon(p) & \sigma \in S_n
\end{cases}
\qc \sigma_i \in \qty{1,2,\cdots,n} \qc S_n : n \mathrm{ ~ th ~ Symmetric ~ group}
\end{align}
詳しくは紹介しません。ここでは一般化された面積と紹介しましたが、その正体は行列式を見ていることとなります。一般に
$$ \vb*v = \sum_{d=1}^{D} v_d \vb*e_d $$
と基底$\vb*e_d$で展開されるベクトルは
\begin{align}
& \Omega( \vb*v ) := \bigwedge_{d=1}^D \vb*v_d = \bigwedge_{d=1}^D \qty( \sum_{d'} v_{dd'} \vb*e_{d'} ) = \sum_{d_1' \cdots d'_D} v_{1d'_1} \cdots v_{Dd'_D} \bigwedge_{d=1}^D \vb*e_{d'_d} \\
& = \sum_{d_1' \cdots d'_D} \epsilon_{d'_1 \cdots d'_D} v_{1d'_1} \cdots v_{Dd'_D} \bigwedge_{d=1}^D \vb*e_{d} = \det(V) \Omega( \vb*e )
\end{align}
という関係を満たします。ここで
$$ (V)_{dd'} := v_{dd'} $$
という行列を定義しました。特に積と和の公式
\begin{align}
& \prod_i \sum_j a_{ij} = \qty( \prod_i \sum_{j_i} ) \prod_{i} a_{ij_i}
\end{align}
を用いることに注意が必要になります。またレヴィチヴィタと行列式の関係も用いました。
\begin{align}
& \sum_{d_1' \cdots d'_D} \epsilon_{d'_1 \cdots d'_D} v_{1d'_1} \cdots v_{Dd'_D} = \sum_{d' \in S_D} \epsilon(d') v_{1d'_1} \cdots v_{Dd'_D} = \det(V)
\end{align}
多次元の場合を考えます。ある$D$次元空間の$D$個のパラメータ$\vb*\xi$によって書かれる点は
\begin{align}
& \vb*r( \vb*\xi ) = \mqty[ x_1(\vb*\xi) & \cdots & x_D(\vb*\xi) ] \in \mathbb{R}^D \qc \vb*\xi = ( \xi_1, \cdots , \xi_D )
\end{align}
です。このパラメータの動くことによる面積要素$ \Omega( \dd{\vb*x} ) $は各パラメータの変化分による体積なので
\begin{align}
& \Omega( \dd{\vb*x} ) := \bigwedge_{ d = 1 }^D \qty( \sum_{ p=1 }^D \pdv{ x_d }{ \xi_p } \dd \xi_p ) = \sum_{p_1 \cdots p_D} \pdv{ x_1 }{ \xi_{p_1} } \cdots \pdv{ x_d }{ \xi_{p_D} } \bigwedge_{ d = 1 }^D \dd \xi_{p_d} = \sum_{p_1 \cdots p_D} \epsilon_{p_1\cdots p_D} \pdv{ x_1 }{ \xi_{p_1} } \cdots \pdv{ x_d }{ \xi_{p_D} } \bigwedge_{ d = 1 }^D \dd \xi_{d}
\end{align}
と式変形をすることができます。これによって変数$\vb*\xi$の場合の面積要素が
\begin{align}
& \Omega( \dd{ \vb*x} ) = \det( \pdv{ x }{ \xi } ) \bigwedge_{ d = 1 }^D \dd \xi_{d} =: J(\vb*\xi) \bigwedge_{ d = 1 }^D \dd \xi_{d} \qc \qty( \pdv{ x }{ \xi } )_{dd'} := \pdv{ x_d }{ \xi_{d'} }
\end{align}
と表せます。ここで
$${J(\vb*\xi)}$$
をヤコビアンと呼びます。座標によらずに
$$ \Omega( \dd {\vb*x} ) $$
は存在することに注意してください。もっとも単純なパラメータとして$ x_d = \xi_d $を採用した場合には
$$ \pdv{ x_d }{ \xi_{d'} } = \delta_{dd'} \to J(\vb*\xi) = 1 $$
であるため
\begin{align}
& \Omega( \dd {\vb*x} ) = \bigwedge_{ d = 1 }^D \dd x_{d} =: \dd[D] x
\end{align}
というもっとも親しみのある面積要素を手に入れることができました。この際に、
$$ \bigwedge_{ d = 1 }^D \dd {\xi_{d}} = \xi_1 \wedge \xi_2 \wedge \cdots \xi_D $$
と添字が右に行くように大きくなるように書いていますが、これは必ずしも一意なことではなくこの了解を空間の向きと呼びます。例えば二次元で以下の面積要素のどちらを表と捉えるかは
\begin{align}
& \xi_1 \wedge \xi_2 \qc \xi_2 \wedge \xi_1
\end{align}
計算する人たちの了解のもとで一意に定まっているとしています。（多くの場合は$ \xi_1 \wedge \xi_2 $を採用します。）
行列の面積要素
一般の行列の体積要素
多変数の積分について定義したので、行列を要素として持つ積分を考えます。ここでは任意の$M \times N$次の行列について考えて、のちにより実用的なものの例をも見ることにします。
\begin{align}
& \dd A := \bigwedge_{ m = 1 }^{N} \bigwedge_{n = 1 }^{M} \dd a_{mn} := \dd a_{11} \wedge \dd{ a_{12} } \cdots \wedge \dd{ a_{1n} } \wedge \dd{ a_{21} } \cdots \wedge \dd a_{MN} \\
& = \dd[M] a_1 \wedge \cdots \wedge \dd[M] a_N = (-1)^{MN} \dd[N] a^1 \wedge \cdots \wedge \dd[N] a^M \\
& A := \mqty[ \vb*a_1 & \cdots & \vb*a_{N} ] = \mqty[ \vb*a^{1} & \cdots & \vb*a^{M} ]^T \qc a_{mn} = (\vb*a_n)_m = (\vb*a^m)_n
\end{align}
例えば
$$ \tr( AA^T ) = \sum_{mn} a_{mn}^2 $$
となるので行列の関数
$$ \exp(- \tr( AA^T )) $$
については
\begin{align}
& \int \dd A ~ \exp(- \tr( AA^T )) = \pi^{ MN/2 }
\end{align}
となります。外積の性質を用いれば次の変数変換に対しての振る舞いを解析できます。
\begin{align}
& A \to A' = L A = L \mqty[ \vb*a_1 & \cdots & \vb*a_{N} ] = \mqty[ L \vb*a_1 & \cdots & L \vb*a_{N} ] \\
& \dd A \to \dd A' = \dd[M]{\qty(La_1)} \wedge \cdots \wedge \dd[M] {\qty(L a_N)} = \det(L)^N \dd A
\end{align}
また同様に行列を左右から作用させた場合にも次のように求められます。
\begin{align}
& A \to A'' = L A R = LA' = L \mqty[ \vb*a_1' & \cdots & \vb*a_{N}' ] = \mqty[ \vb*a_1'' & \cdots & \vb*a_{N}'' ] \\
& \mqty[ \vb*a_1' & \cdots & \vb*a_{N}' ] = \mqty[ \vb*a^1 & \cdots & \vb*a^{N} ]^T R = \qty( \mqty[ R^{T} \vb*a^1 & \cdots & R^T \vb*a^{N} ] )^T \\
& \dd A' = (-1)^{MN} \dd[N] {\qty(R^Ta)^1} \wedge \cdots \wedge \dd[N] {\qty(R^Ta)^M} = \det(R^T)^M \dd A = \det(R)^M \dd A \\
& \dd A \to \dd A'' = \det(L)^N \dd A' = \det(L)^N \det(R)^M \dd A
\end{align}
対称行列の面積要素
さて、ここまでは$MN$個の全てが独立な場合の式です。本題に近づいて$M \times M$の対称行列の場合について考えます。対称行列の独立な成分は、対角要素と非対角要素の自由度の和である$ \frac{M(M+1)}{2} $です。つまり
\begin{align}
& \dd A := \bigwedge_{m \leq m'} \dd a_{mm'} = \dd a_{11} \wedge \dd a_{12} \cdots \wedge \dd a_{MM}
\end{align}
となります。面積要素の定義はおおよそ直感的ですが、問題は変数変換です。
\begin{align}
& A \to A' = S A S^T
\end{align}
と変数変換を行った場合を考えます。変換の形から$ A' $も対称行列です。この場合
\begin{align}
& \dd A' = \det(S)^{M+1} \dd A
\end{align}
が成り立つことが知られています。例えば
\begin{align}
& S = \mqty[ s & 0 & \cdots & 0 \\ 0 & 1 & \cdots & 0 \\ \vdots & & \ddots & \\ 0 & 0 & \cdots & 1 ] \qc A' = \mqty[ s^2a_{11} & sa_{12} & \cdots & sa_{1M} \\ sa_{21} & a_{22} & \cdots & a_{2M} \\ \vdots & & \ddots & \\ sa_{M1} & a_{M2} & \cdots & a_{MM} ]
\end{align}
の場合には特別な場合を確認することができます。次数を忘れた場合の検算に役立ちます。
\begin{align}
& \dd A' = s^{M+1} \dd A
\end{align}
この証明は行列要素の多項式$J(S)$で$ \dd A' = J(S) \dd A $となるとした場合に、任意の正則行列$S$について
\begin{align}
& J(SS') = J(S) J(S') \qc J(I) = 1
\end{align}
が成り立ち、上の性質を満たす行列要素$s_{mm'}$の多項式は必ず行列式のベキで表せること
\begin{align}
& J(S) = \det(S)^{j}
\end{align}
そして先ほどの例から$j=M+1$であることから証明できます。自明でないのは$J(S)$が必ず行列式で表せることですが、その証明は以下の様に考えられます。まず$S$が対角化可能な場合に、$J(S)$は定義から
\begin{align}
& J(S) = J(P \Lambda(S) P^{-1} ) = J( \Lambda(S)) =: f(s_1, \cdots ,s_M )
\end{align}
と固有値$s_m$のみに依存する関数だと考えることができます。この$M$変数関数$f(s_1, \cdots ,s_M )$の関数形は
\begin{align}
& J(SS') = J(S) J(S') \to f(s_1s'_1, \cdots ,s_Ms'_M ) = f(s_1, \cdots ,s_M ) f(s'_1, \cdots ,s'_M )
\end{align}
という規則を満たしますが、これは両辺を$s_m'$で微分を行い$s_m' = 1$とした場合
\begin{align}
& s_m \qty[ \pdv{s''_m} f(s''_1, \cdots ,s''_M ) ]_{s'' = s} = s_m \pdv{s_m} f(s_1, \cdots ,s_M ) \\
& = f(s_1, \cdots ,s_M ) \qty[ \pdv{s'_m} f(s'_1, \cdots ,s'_M ) ]_{s'=1} =: j_m f(s_1, \cdots ,s_M ) \qc j_m := \qty[ \pdv{s'_m} f(s'_1, \cdots ,s'_M ) ]_{s'=1}
\end{align}
という関係を満たします。これを全ての$s_m$について考えると、関数形は以下の形になります。
\begin{align}
& f(s_1, \cdots ,s_M ) = f_0 s_1^{j_1} \cdots s_M^{j_M}
\end{align}
この指数$j_m$については、スワップ行列$W_{ij}$と対角行列に成り立つ次の性質からその特徴がわかります。
\begin{align}
& W_{mm'} := \mqty[ 1 & \cdots & 0 & \cdots & 0 & \cdots & 0 \\
& \ddots & & & & & \\
0 & \cdots & 0 & \cdots & 1 & \cdots & 0 \\
& & & \ddots & & & \\
0 & \cdots & 1 & \cdots & 0 & \cdots & 0 \\
& & & & & \ddots & \\
1 & \cdots & 0 & \cdots & 0 & \cdots & 1 \\ ]
\qc \Lambda(S) := \mathrm{diag}\mqty[ s_1 & \cdots & s_m & \cdots & s_{m'} & \cdots & s_M ] \\
& W_{mm'} \Lambda(S) W_{mm'} = \mathrm{diag}\mqty[ s_1 & \cdots & s_{m'} & \cdots & s_{m} & \cdots & s_M ] \qc W_{mm'}^2 = I
\end{align}
が成り立っているため、関数$f$が持つ以下の性質が明らかになります。
\begin{align}
& f(s_1 , \cdots , s_{m} , \cdots , s_{m'} , \cdots , s_M ) := J( \Lambda(S) ) = J( W_{mm'} \Lambda(S) W_{mm'} ) = f(s_1 , \cdots , s_{m'} , \cdots , s_{m} , \cdots , s_M )
\end{align}
つまり、$f$は$M$個の変数の対称関数ということです。従って$j_m=j_{m'}=j$という等式が任意の$1 \leq m,m' \leq M$について成り立っている必要があります。これを改めて書き下すと
\begin{align}
J(S) = f(s_1 , \cdots , s_{m} , \cdots , s_{m'} , \cdots , s_M ) = f_0 \qty( s_1 \cdots s_M )^{j} = f_0 \det(S)^j
\end{align}
となっている必要があるということです。比例係数$f_0$は定義式から
\begin{align}
& J(SS') = J(S) J(S') \to f_0 \det(SS')^j = f_0^2 \det(S)^j \det(S')^j \to f_0 = f_0^2
\end{align}
が成り立つので、$f_0 = 0,1$です。$f_0 = 0$は全ての変換を$0$とする自明な変換なので除くと
\begin{align}
& J(S) = \det(S)^j
\end{align}
のように、考えられる$J(S)$の関数形は行列式のベキのみであることが証明できたことになります。$~_\Box$
正定値行列の面積要素
次は正定値行列$A$についての行列要素を考えます。正定値行列は対称行列の特別な場合で、面積要素や変数変換は対称行列と同じな一方、その行列要素に制限が存在します。この条件を$ A > 0 $と表します。
\begin{align}
& \dd A = \bigwedge_{m \leq m'} \dd a_{mm'} \qc A' = SAS^T \to \dd A' = \det( S )^{M+1} \dd A
\end{align}
積分範囲の性質を詳しく見るため、最初に正定値行列は
\begin{align}
& A = U^T U
\end{align}
と上三角行列$U$によって（のちに見るように対角成分の符号を除いて）一意に分解できることを紹介します。
\begin{align}
& \mqty[ a_{11} & a_{12} & \cdots & a_{1M} \\ a_{12} & a_{22} & \cdots & a_{2M} \\ \vdots & & \ddots & \\ a_{1M} & a_{2M} & \cdots & a_{MM} ] =
\mqty[ u_{11} & 0 & \cdots & 0 \\ u_{12} & u_{22} & \cdots & 0 \\ \vdots & & \ddots & \\ u_{1M} & u_{2M} & \cdots & u_{MM} ]
\mqty[ u_{11} & u_{12} & \cdots & u_{1M} \\ 0 & u_{22} & \cdots & u_{2M} \\ \vdots & & \ddots & \\ 0 & 0 & \cdots & u_{MM} ]
\end{align}
を具体的に計算して書き下すと、各成分の関係を知ることができます。
\begin{align}
& a_{11} = u_{11}^2 \qc a_{12} = u_{11} u_{12} \qc \cdots \qc a_{1M} = u_{11} u_{1M} \\
& a_{22} = u_{12}^2 + u_{22}^2 \qc \cdots \qc a_{2M} = u_{12} u_{1M} + u_{22} u_{2M} \\
& \vdots \\
& a_{MM} = u_{1M}^2 + \cdots + u_{MM}^2
\end{align}
これは$ a_{mm} > 0 $という正定値行列の性質を用いて、具体的な計算行うことによって証明できます。
\begin{align}
& A = \mqty[ a_{11} & \vb*a_1^T \\ \vb*a_1 & R_1 ] =
\mqty[ \pm \sqrt{a_{11}} & 0 \\ \pm \frac{ \vb*a_1 }{ \sqrt{a_{11}} } & I ]
\mqty[ 1 & 0 \\ 0 & R_{1} - \frac{ \vb*a_{1} \vb*a_{1}^T }{ a_{11} } ]
\mqty[ \pm \sqrt{a_{11}} & \pm \frac{ \vb*a_1^T }{ \sqrt{a_{11}} } \\ 0 & I ] =: U_1^T A_{1} U_{1} \\
& A_1 = U_{2}^T A_2 U_2 \qc \cdots \qc A_M = U_{M}^T A_M U_M \qc A_{M} = I
\end{align}
ここで、上三角行列の積が再度上三角行列になることを用いることによって
\begin{align}
& A = U_{1}^T U_{2}^T \cdots U_{M}^T I U_{M} \cdots U_{2} U_{1} = U^T U \qc U := U_{M} \cdots U_{2} U_{1}
\end{align}
と具体的に知りたい$U$を構成することができます。ただし繰り返すためには
$$ R_1 - \frac{ \vb*a_1 \vb*a_1^T }{ a_{11} } $$
が正定値行列である必要があります。つまり大きさが考える系の大きさより一つ小さい任意のベクトル$ \vb*\psi $について
\begin{align}
& \vb*\psi^T \qty( R - \frac{ \vb*a \vb*a^T }{a} ) \vb*\psi = \vb*\psi^T R \vb*\psi - \frac{ \qty( \vb*a^T \vb*\psi )^2 }{ a } > 0
\end{align}
が成り立って欲しいということです。この証明は$A$がもともと正定値行列であることから
\begin{align}
& \mqty[ \varphi \\ \vb*{\varphi}^T ] := U \mqty[ 0 \\ \vb*\psi^T ] = \mqty[ \pm \sqrt{a} & \pm \frac{ \vb*a^T }{ \sqrt{a} } \\ 0 & I ] \mqty[ 0 \\ \vb*\psi^T ] \\
& 0 \leq \mqty[ \varphi & \vb*\varphi^T ] \mqty[ a & \vb*a^T \\ \vb*a & R ] \mqty[ \varphi \\ \vb*\varphi ] = \mqty[ 0 & \vb*\psi^T ] \mqty[ 1 & 0 \\ 0 & R - \frac{\vb* a \vb*a^T }{a} ] \mqty[ 0 \\ \vb*\psi ] = \vb*\psi^T \qty( R - \frac{ \vb*a \vb*a^T }{a} ) \vb*\psi
\end{align}
と特別なベクトル
$$ \mqty[ \varphi & {\vb*\varphi}^T ] $$
について調べれば実証することが可能です。以下では$U$の対角成分を正とします。
この上三角行列$U$の構成要素は$ A $の持っている自由度を過不足なく上手に表しています。これを用いて正定値行列の面積要素を$u_{mm'}$で書き直すことを考えてみましょう、結果的に以下の式が成り立ちます。
\begin{align}
& \dd A := \bigwedge_{m \leq m'} \dd a_{mm'} = 2^{M} u_{11}^M u_{22}^{M-1} \cdots u_{MM} \bigwedge_{m \leq m'} \dd u_{mm'} = 2^M \prod_{m=1}^M u_{mm}^{M-m+1} \dd U
\end{align}
この証明は、先ほどの成分の関係から以下の式が成り立つことからわかります。
\begin{align}
& \dd a_{11} = 2 u_{11} \dd u_{11} \qc \dd a_{11} = u_{11} \dd u_{12} + \order{ \dd u_{11} } \qc \cdots \qc \dd a_{1M} = u_{11} \dd u_{1M} + \order{ \dd u_{11} } \\
& \dd a_{22} = 2 u_{22} \dd u_{22} + \order{ \dd u_{12} } \qc \cdots \qc \dd a_{2M} = 2 u_{22} \dd u_{22} + \order{ \dd u_{12}, \dd u_{1M} , \dd u_{22} } \\
& \vdots \\
& \dd a_{MM} = 2 u_{MM} \dd u_{MM} + \order{ \dd u_{11}, \dd u_{22}, \cdots , \dd u_{ M-1 M-1} }
\end{align}
ここで$ \order{\dd x} $とは$ \dd x $を含んだ何かしらの一形式を表す記号で
\begin{align}
& \dd x \wedge \order{\dd x} = 0
\end{align}
が成り立つような項です。この関係を導けは、具体的に面積要素を計算すれことで求める式
\begin{align}
\dd A := \bigwedge_{m \leq m'} \dd a_{mm'} = 2^{M} u_{11}^M u_{22}^{M-1} \cdots u_{MM} \bigwedge_{m \leq m'} \dd u_{mm'}
\end{align}
が得られることがわかります。これで正定値行列の面積要素について準備が整いました。
この表示で特に重要なのは、正定値行列$A$は対称な場合の成分の空間を一部分に制限しているため積分区間が複雑であるのに対して、$U$の成分は対角の符号が決まっていれば非対角は任意の値をとれる点です。
\begin{align}
& A > 0 \to A = U^T U \qc u_{mm} \geq 0 \qc - \infty < u_{mm'} < \infty ~~ (m \neq m')
\end{align}
マルチガンマ関数
マルチガンマ関数$\Gamma_M(P)$は$M$次元の正定値行列$A$を用いて以下のように定義することができます。
\begin{align}
& \Gamma_{M}(P) := \int_{A>0} \dd A ~ \det(A)^{P - \frac{ M + 1 }{2}} \exp( -\tr(A) )
\end{align}
特に$M=1$の場合は$ A = a > 0 $なので、よく知るガンマ関数に一致します。
\begin{align}
& \Gamma_{M=1}(P) = \int_{0}^\infty \dd a ~ a^{P - 1 } \exp( - a ) = \Gamma(P)
\end{align}
一方で一般の場合について、以下の関係が成り立つことを示すことができます。
\begin{align}
& \Gamma_{M}(P) = \pi^{ \frac{M(M-1)}{4} } \Gamma \qty( P ) \Gamma \qty( P - \frac{1}{2} ) \cdots \Gamma \qty( P - \frac{ M -1 }{2} ) = \pi^{ \frac{M(M-1)}{4} } \prod_{m=1}^M \Gamma \qty( P - \frac{ m -1 }{2} )
\end{align}
証明は、先ほど示した
$A = U^T U$
を用いて行います。
\begin{align}
& \dd A = 2^{M} u_{11}^M u_{22}^{M-1} \cdots u_{MM} \bigwedge_{m \leq m'} \dd u_{mm'} = 2^M \prod_{m=1}^M u_{mm}^{M-m+1} \dd U \\
& \tr(A) = \tr(U^TU) = \sum_{ m \leq m' } u_{mm'}^2 \qc \det( A ) = \det(U)^2 = \prod_{m=1}^M u_{mm}^2
\end{align}
の関係を用いてマルチガンマ関数を書き直すと
\begin{align}
& \Gamma_{M}(P) = 2^{M} \int \dd U ~ \qty( \prod_{m=1}^M u_{mm}^{M-m+1} ) \qty( \prod_{m=1}^M u_{mm}^{2P-(M+1)} ) \exp( - \sum_{ m \leq m' } u_{mm'}^2 ) \\
& = 2^{M} \prod_{m < m'} \qty( \int_{- \infty}^{\infty} \dd u_{mm'} ~ \exp( - u_{mm'}^2 ) ) \prod_{m = 1}^M \qty( \int_{0}^{\infty} \dd u_{mm} ~ u^{2P-m} \exp( - u_{mm'}^2 ) )
\end{align}
となり、1変数の積分公式でよく知られている
\begin{align}
& 2 \int_{0}^{\infty} \dd x ~ x^{ 2\nu } e^{- x^2} = \int_{0}^{\infty} \dd x^2 ~ \qty(x^2)^{\nu-1/2} e^{- x^2} = \Gamma \qty( \nu + \frac{1}{2} )
\end{align}
を用いることによって目的の式を得ることができます。
\begin{align}
& \Gamma_{M}(P) = \prod_{m < m'} \qty( \sqrt{\pi} ) \prod_{m = 1}^M \Gamma \qty(P - \frac{m}{2} + \frac{1}{2} ) = \pi^{ \frac{M(M-1)}{4} } \prod_{m=1}^M \Gamma \qty( P - \frac{ m -1 }{2} )
\end{align}
ウィシャート分布
最後に本題であるウィシャート分布の正規化について話します。示したいことは正定値行列$\Sigma$について
\begin{align}
& \int_{A>0} \dd A ~ \det(A)^{P - \frac{ M + 1 }{2}} \exp( - \frac{1}{2} \tr( \Sigma^{-1} A) ) = \det( 2 \Sigma )^P \Gamma_{M}(P)
\end{align}
が成り立つことでした。これは変数変換を行ったマルチガンマ関数であることがわかります。実際に
\begin{align}
& A = \sqrt{2\Sigma} B \sqrt{2\Sigma} \qc \dd A = \det( 2 \Sigma )^{\frac{ M(M+1) }{2} } \dd B \\
& \frac{1}{2} \tr( \Sigma^{-1} A ) = \tr( B ) \qc \det( A ) = \det( 2 \Sigma ) \det( B )
\end{align}
が成り立つため、積分は以下の通りに変形できます。
\begin{align}
& \int_{A>0} \dd A ~ \det(A)^{P - \frac{ M + 1 }{2}} \exp( - \frac{1}{2} \tr( \Sigma^{-1} A) ) \\
& = \det( 2 \Sigma )^{P-\frac{M+1}{2}} \det( 2 \Sigma )^{\frac{M+1}{2}} \int_{B>0} \dd B ~ \det(B)^{P - \frac{ M + 1 }{2}} \exp( - \tr(B) ) \\
& = \det( 2 \Sigma )^{P} \Gamma_{M}(P)
\end{align}
つまり、正定値行列$X$を確率変数としてもつウィシャート分布
\begin{align}
& W( X | \Sigma , P ) := \frac{ \det(X)^{\frac{ P - D - 1 }{2}} \exp(- \frac{1}{2} \tr( \Sigma^{-1} X ) ) }{ \det( 2 \Sigma )^{P/2} \Gamma_D \qty(\frac{P}{2}) } \qc X > 0
\end{align}
は、確かに規格化されていることが確かめられたことになります。
終わりに
正規分布の学習の話でいきなり出てきて困惑しますよね、ウィシャート分布。そんな気持ちからの雑なメモでした。
あとMarkdownでインラインのままベクトルを二つ以上使おうとすると数式になってくれないの地味に困るので解決したいです。
$J(S)$の関数形の議論に付き合ってくれた後輩のK君、本当にありがとうございました。
ベーカーキャンベルハウスドルフの公式周りの説明
ベーカーキャンベルハウスドルフの補題
最も有名なのはベーカーキャンベルハウスドルフの補題である。
\begin{align}
& e^{A} B e^{-A} = \sum_{n=0}^\infty \frac{1}{n!} \comm{A}{ \cdots \comm{A}{B} } = \qty( \sum_{n=0}^\infty \frac{\hat{A}^n}{n!} ) B = \exp( \hat{A} ) B
\end{align}
ここで$ \hat{A} B := \comm{A}{B} $という記号を定義した。これをベーカーキャンベルハウスドルフの定理と呼ぶ文献も存在する。この公式は以下の関数を用いて証明することができる。
\begin{align}
& f(\beta) := e^{\beta A} B e^{- \beta A}
\end{align}
この関数$f(\beta)$は次の微分の式を満たす。
\begin{align}
& \dv{\beta} f(\beta) = e^{\beta A} \comm{A}{B} e^{- \beta A} = e^{\beta A} \qty( \hat{A} B ) e^{- \beta A} \\
& \dv[n] {\beta} f(\beta) = e^{\beta A} \qty( \hat{A}^n B ) e^{- \beta A}
\end{align}
これを繰り返し用いて$f(\beta)$のテイラー展開を求めることによって求めたい結果を得る。
\begin{align}
& f(\beta) = \sum_{n=0}^\infty \frac{ \beta^n }{n!} \dv[n] {\beta} f(0) = \sum_{n=0}^\infty \frac{ \beta^n }{n!} \hat{A}^n B = \exp( \beta \hat{A} ) B \\
& \to f(1) = e^{A} B e^{-A} = \exp( \hat{A} ) B
\end{align}
この公式の非自明で興味深い例題は
\begin{align}
& e^{A} e^{B} =: e^{C}
\end{align}
となる$C$による計算をする場合である。一般に$C \neq A + B$である。計算すると
\begin{align}
& e^{C} H e^{-C} = e^{A} e^{B} H \qty(e^{A} e^{B})^{-1} = e^{A} e^{B} H e^{-B} e^{-A} = e^{A} \qty( \exp(\hat{B}) H ) e^{-A} = \exp(\hat{A}) \exp(\hat{B}) H
\end{align}
という式が成り立っている。つまり非自明な関係式$ \exp(\hat{C}) = \exp(\hat{A}) \exp(\hat{B}) $が成り立つ。
ベーカーキャンベルハウスドルフの公式1
非可換な演算子同士の計算が何気ない計算を困難にしてしまう。例えば
\begin{align}
& f(\lambda) := e^{ H(\lambda)} \pdv{\lambda} e^{- H(\lambda)}
\end{align}
という量が単純に$ - \pdv{H(\lambda)}{\lambda} $となるのは特別な場合だけである。なぜならば必ずしも$ \pdv{H(\lambda)}{\lambda} $と$H(\lambda)$が可換であるとは限らないためである。以下の量を考える。
\begin{align}
& F(\lambda, \beta) := e^{ \beta H(\lambda)} \pdv{\lambda} e^{- \beta H(\lambda)} \qc F(\lambda, \beta=1 ) = f(\lambda) \qc F(\lambda, \beta = 0 ) = 0
\end{align}
この量の従う微分方程式は以下の通りである。
\begin{align}
& \pdv{\beta} F(\lambda, \beta) = - e^{ \beta H(\lambda)} \pdv{H(\lambda)}{\lambda} e^{- \beta H(\lambda)}
\end{align}
つまりこの両辺を積分することによって求めたい結果を得ることができる。
\begin{align}
& f(\lambda) = - \int_0^1 \dd \beta ~ e^{ \beta H(\lambda)} \pdv{H(\lambda)}{\lambda} e^{- \beta H(\lambda)} = - \int_0^1 \dd \beta ~ \sum_{n=0}^\infty \frac{ \hat{H}^n(\lambda) }{n!} \beta^n \pdv{H(\lambda)}{\lambda} = - \qty( \sum_{n=0}^{\infty} \frac{ \hat{H}^n(\lambda) }{ (n+1)! } ) \pdv{H(\lambda)}{\lambda}
\end{align}
この無限級数はとある解析関数$\phi$を用いて表すことが可能である。
\begin{align}
& f(\lambda) := e^{ H(\lambda)} \pdv{\lambda} e^{- H(\lambda)} = - \phi\qty( \hat{H}(\lambda) ) \pdv{H(\lambda)}{\lambda} \\
& \phi(z) := \frac{ e^{z} -1 }{z} = \sum_{n=0}^\infty \frac{ z^n }{ (n+1)! }
\end{align}
と表すことも可能である。この形にすると色々な"テクニック"が使える。テクニックとは
\begin{align}
& e^{ H(\lambda)} \pdv{\lambda} e^{- H(\lambda)} = - \phi\qty( \hat{H}(\lambda) ) \pdv{H(\lambda)}{\lambda} = K \\
& \to \pdv{H(\lambda)}{\lambda} = - \phi^{-1} \qty(\hat{H}(\lambda)) K \qc \phi^{-1}(z) \phi(z) = 1 \qc \phi^{-1}(z) := \frac{z}{ e^{z} - 1 }
\end{align}
のような変形である。関数系を用いることによって式変形が容易になる。この有用な例題は
\begin{align}
& e^{A} e^{ \lambda B} =: e^{ C(\lambda)}
\end{align}
についての計算である。次の微分方程式が成り立っていることが確認できる。
\begin{align}
& e^{ - C(\lambda)} \pdv{\lambda} e^{ C(\lambda)} = \phi\qty( - \hat{C}(\lambda) ) \pdv{C(\lambda)}{\lambda} \qc (H(\lambda) = - C(\lambda)) \\
& \phi(- \hat{C}(\lambda) ) \pdv{C(\lambda)}{\lambda} = e^{-C(\lambda)} \pdv{\lambda} e^{ C(\lambda)} = e^{ - \lambda B} e^{ - A } e^{ A } e^{ \lambda B } B = B
\end{align}
つまり、先ほどの公式から次の関係が成り立つ。
\begin{align}
& \pdv{C(\lambda)}{\lambda} = \phi^{-1} \qty( - \hat{C}(\lambda)) B \\
& \to C := C(1) = C(0) + \int_0^1 \dd \lambda ~ \phi^{-1} \qty(- \hat{C}(\lambda)) B = A + \int_0^1 \dd \lambda ~ \phi^{-1} \qty(- \ln(\exp(\hat{C}(\lambda)))) B
\end{align}
ここで$ \psi(z) := \phi^{-1}( - \ln z ) = \frac{ z \ln z }{ z - 1 } $を定義すれば$ \exp(\hat{C}(\lambda)) = \exp(\hat{A}) \exp( \lambda \hat{B}) $という公式を用いて
\begin{align}
& C = A + \int_0^1 \dd \lambda ~ \psi \qty( \exp(\hat{A}) \exp( \lambda \hat{B}) ) B \qc \psi(z) := \frac{ z \ln z }{ z - 1 }
\end{align}
という厳密な公式を得ることができる。これがベーカーキャンベルハウスドルフの公式である。
ベーカーキャンベルハウスドルフの公式2
先ほどの式は整理されている反面、積分が含まれることや対称性の悪さがある。ここでは上記の公式とは別の表現を紹介する。この方法は黒木玄さん（の参考にする長谷川浩司さん）によるものを参考にしました。
\begin{align}
& \pdv{H(\lambda)}{\lambda} = - \phi^{-1} \qty(\hat{H}(\lambda)) \qty( e^{ H(\lambda)} \pdv{\lambda} e^{- H(\lambda)} ) \qc \phi^{-1}(z) \phi(z) = 1 \qc \phi^{-1}(z) = \frac{z}{ e^{z} - 1 }
\end{align}
を便利なように$ \hat{H}(\lambda) = \ln \qty( \exp( \hat{H}(\lambda) ) ) $を用いて次の様に書き直す。
\begin{align}
& \pdv{H(\lambda)}{\lambda} = - \chi \qty( \exp( \hat{H}(\lambda) ) ) \qty( e^{ H(\lambda)} \pdv{\lambda} e^{- H(\lambda)} ) \qc \chi(z) := \phi^{-1}(\ln z) = \frac{\ln z}{ z - 1 }
\end{align}
これと$ H(\lambda) \to - H(\lambda) $としたものに、$ A^{-1}(\lambda) \pdv{\lambda} A(\lambda) = - \qty( \pdv{\lambda} A^{-1}(\lambda)) A(\lambda) $を使うと次の公式を得る。
\begin{align}
& \pdv{H(\lambda)}{\lambda} = \chi \qty( \exp( \hat{H}(\lambda) ) ) \qty( \qty(\pdv{\lambda} e^{ H(\lambda)}) e^{- H(\lambda)} ) \\
& \pdv{H(\lambda)}{\lambda} = \chi \qty( \exp( - \hat{H}(\lambda) ) ) \qty( e^{ - H(\lambda)} \qty( \pdv{\lambda} e^{ H(\lambda)} ) )
\end{align}
この公式を$ H(\lambda) := H(t(\lambda),s(\lambda)) := \ln \qty( e^{t(\lambda)A} e^{s(\lambda)B} ) $について計算すると
\begin{align}
& \qty(\pdv{t} e^{ H(t,s)}) e^{- H(t,s)} = A \qc e^{- H(t,s)} \qty(\pdv{s} e^{ H(t,s)}) = B \\
& \dv{H(\lambda)}{\lambda} = \dv{t(\lambda)}{\lambda} \pdv{H(t,s)}{t} + \dv{s(\lambda)}{\lambda} \pdv{H(t,s)}{s} = \dv{t(\lambda)}{\lambda} \chi \qty( \exp( \hat{H}(\lambda) ) ) A + \dv{s(\lambda)}{\lambda} \chi \qty( \exp( - \hat{H}(\lambda) ) ) B
\end{align}
が成り立つ。また交換子の記法については次の関係が成り立つ。
\begin{align}
& \exp( \hat{H}(\lambda) ) = \exp( t(\lambda) \hat{A} ) \exp( s(\lambda) \hat{B} ) \qc \exp( - \hat{H}(\lambda) ) = \exp( - s(\lambda) \hat{B} ) \exp( - t(\lambda) \hat{A} )
\end{align}
以下では$t(\lambda) = s(\lambda) = \lambda$として、以下の新しい関数$Z(x,y|\lambda)$を定義する。
\begin{align}
& \dv{H(\lambda)}{\lambda} = \chi \qty( \exp( \lambda \hat{A} ) \exp( \lambda \hat{B} ) ) A + \chi \qty( \exp( - \lambda \hat{B} ) \exp( - \lambda \hat{A} ) ) B \\
& = Z(\hat{A},\hat{B}|\lambda) A + Z(\hat{B},\hat{A}|-\lambda) B \qc Z(x,y|\lambda) := \chi \qty( \exp( \lambda x ) \exp( \lambda y ) ) = \sum_{m=0}^\infty \lambda^{m} Z_m(x,y)
\end{align}
すると最終的な$H(\lambda=1) =: C $が以下の様にもとまる。
\begin{align}
& C = \sum_{m=0}^\infty \frac{ 1 }{m+1} \qty( Z_m( \hat{A}, \hat{B} ) A + (-1)^{m} Z_m( \hat{B}, \hat{A} ) B ) \\
& Z_m(\hat{A},\hat{B}) = 1 + \sum_{n=1}^\infty \frac{(-1)^{n}}{n+1} \sum_{ {p_i,q_i} \in \mathcal{Z}_m^{\qty(n)} } \frac{ \hat{A}^{p_1} \hat{B}^{q_1} \cdots \hat{A}^{p_{n}} \hat{B}^{q_{n}} }{ p_1 ! q_1 ! \cdots p_{n} ! q_{n} ! } \\
& \mathcal{Z}_m^{\qty(n)} = \qty{ \qty{ p_i, q_i } ~ s.t. ~ \sum_{i=1}^{n} ( p_i + q_i ) = m , ~ \qty(p_i + q_i) > 0 }
\end{align}
多項分布とディリクレ分布の性質
変分ベイズ学習(MLP,中島伸一)を読んでいていろんな確率分布を初めて見ましたが、重要そうなディルクレ分布や多項分布の性質をせっかく自分で考えたのでメモしておこうかなという次第です。文章中のベクトル文字が（なぜか）太字にならなかったので、ところどころ細いです。あしからず。
多項分布
多項分布は名の通り多項式が元ネタの確率分布です。
\begin{align}
& M_N( \vb*x | \vb*\theta) := N! \prod_{k=1}^K \frac{ \theta_k^{x_k} }{ x_k! } \\
& 0 \leq \theta_k \leq 1 \qc \sum_{k=1}^K \theta_k = 1 \qc 0 \leq x \leq N \qc \sum_{k=1}^K x_k = N
\end{align}
この確率分布が正規化（全ての変数を足し合わせて$1$になる）されていることは多項式展開
\begin{align}
& 1 = \qty( \sum_{k=1}^K \theta_k )^N = \sum_{x_1 \cdots x_K}^{ x_1 + \cdots + x_K = N} \frac{N!}{ x_1 ! \cdots x_K! } \theta_1^{x_1} \cdots \theta_K^{x_K} = \sum_{x_1 \cdots x_K}^{ \sum_{k=1}^K x_k = N} N! \prod_{k=1}^K \frac{ \theta_k^{x_k} }{ x_k! }
\end{align}
によって考えている$x_k$の条件と和の取り方が対応して、分布が正規化されていることが理解できるでしょう。
蛇足かもしれないですが多項式展開についても少し説明します。高校数学で馴染みのある二項展開
\begin{align}
& (\alpha + \beta )^N = \sum_{x=0}^N \frac{N!}{x!(N-x)!} \alpha^x \beta^{N-x}
\end{align}
を繰り返せば証明することができます。まず$ \alpha =\theta_1, ~ \beta = \theta_2 + \cdots \theta_K $として
\begin{align}
& \qty( \sum_{k=1}^K \theta_k )^N = \sum_{x_1=0}^N \frac{N!}{x_1!(N-x_1)!} \theta_1^{x_1} \qty( \theta_2 + \cdots \theta_K )^{N-x_1}
\end{align}
がなりたつので、$\qty( \theta_2 + \cdots \theta_K )^{N-x_1}$の項についてさらに
\begin{align}
& \qty( \theta_2 + \cdots \theta_K )^{N-x_1} = \sum_{x_2=0}^N \frac{(N-x_1)!}{x_2!(N-x_1-x_2)!} \theta_2^{x_2} \qty( \theta_3 + \cdots \theta_K )^{N-x_1-x_2}
\end{align}
と展開して合わせることによって以下の式が得られます。
\begin{align}
& \qty( \sum_{k=1}^K \theta_k )^N = \sum_{x_1,x_2=0}^N \frac{N!}{x_1!(N-x_1)!} \frac{(N-x_1)!}{x_2!(N-x_1-x_2)!} \theta_1^{x_1} \theta_2^{x_2} \qty( \theta_3 + \cdots \theta_K )^{N-x_1-x_2} \\
& = \sum_{x_1,x_2=0}^N \frac{N!}{x_1! x_2! (N-x_1-x_2)!} \theta_1^{x_1} \theta_2^{x_2} \qty( \theta_3 + \cdots \theta_K )^{N-x_1-x_2}
\end{align}
現れた二つの二項係数の一部分が約分されて一般化した二項分布が現れています。これを繰り返せば
\begin{align}
& \qty( \sum_{k=1}^K \theta_k )^N = \sum_{x_1 \cdots x_{K-1} =0}^N \frac{N!}{x_1! \cdots x_{K-1}! (N-x_1 \cdots -x_{K-1})!} \theta_1^{x_1} \cdots \theta_{K-1}^{x_{K-1}} \theta_K^{N-x_1 \cdots -x_{K-1}} \\
& = \sum_{x_1 \cdots x_{K} }^{x_1 + \cdots x_K = N} \frac{N!}{x_1! \cdots x_K!} \theta_1^{x_1} \cdots \theta_{K}^{x_{K}} \qc x_K := N - x_1 - \cdots x_{K-1}
\end{align}
とう多項式展開が得られます。$ N! / x_1! \cdots x_K! $は二項係数の一般化です。
ディリクレ型尤度関数
さてこの多項分布
$$ M_N( \vb*x | \vb*\theta) := N! \prod_{k=1}^K \frac{ \theta_k^{x_k} }{ x_k! } $$
は$\theta_k$をパラメータとして含んでいますが、確率論においては与えられた変数$ x_k $をもとにそれを与える確率を最大化するように$ \theta_k $を選びたくなる場合があります。これは未知のパラメータについて少なくとも観測値$ x_k $が出ているのだから、それを尤もらしくなるような$ \theta_k $にしておけば良いだろうという精神に基づく最適化で、最尤度法と呼ばれます。
そこで多項分布が$ \theta $の関数としてどんな形を持っているのかが重要になるのですが、見てわかるように
\begin{align}
& M_N( \vb*x | \vb*\theta) \propto \prod_{k=1}^K \theta_k^{x_k} \qc 0 \leq \theta_k \leq 1 \qc \sum_{k=1}^K \theta_k = 1
\end{align}
という形をしています。この形を持つ確率分布をディリクレ型尤度関数と言います。
\begin{align}
& D_K ( \vb*\theta | \vb*\phi ) := \frac{ \Gamma \qty( \sum_{k=1}^K \phi_k ) }{ \prod_{k=1}^K \Gamma(\phi_k) } \prod_{k=1}^K \theta_k^{ \phi_k - 1 } \\
& \Gamma(\phi) := \int_0^\infty \dd t ~ t^{\phi-1} e^{-t} \\
& 0 \leq \theta_k \leq 1 \qc \sum_{k=1}^K \theta_k = 1 \qc 0 \leq \phi_k < \infty
\end{align}
合わせて紹介したのはガンマ関数$\Gamma(x)$で、階乗の一般化に対応します。この規格化は
\begin{align}
& \int_D \dd \theta_1 \cdots \dd \theta_{K-1} ~ \theta_1^{\phi_1-1} \cdots \theta_{K-1}^{\phi_{K-1}-1} \qty(
1 - \theta_1 \cdots - \theta_{K-1} )^{\phi_{K}-1} = \frac{ \Gamma( \phi_1 ) \cdots \Gamma( \phi_K ) }{ \Gamma \qty( \phi_1 + \cdots + \phi_K ) } \\
& D := \qty{ \theta_k \geq 0 , 1 - \theta_1 \cdots - \theta_{K-1} \geq 0 }
\end{align}
が成り立つことを示すことができれば良いです。この証明は右辺から行ったほうが楽です。
\begin{align}
& \Gamma( \phi_1 ) \cdots \Gamma( \phi_K ) = \int_0^\infty \dd t_1 \cdots \dd t_{K} ~ t_1^{\phi_1-1} \cdots t_K^{\phi_K-1} e^{ - \qty( t_1 + \cdots + t_K ) }
\end{align}
について$ t_k \geq 0, ~ t_1 + \cdots + t_N = t $となるように以下のような変数変換を行います。
\begin{align}
& t_1 = \theta_1 t \qc \cdots \qc t_{K-1} = \theta_{K-1} t \qc t_K = \qty( 1 - \theta_1 \cdots - \theta_{K-1} )t \\
& \theta_k \in D \qc 0 \leq t \leq \infty
\end{align}
この際の変数変換は以下のヤコビアンで実行されます。
\begin{align}
& J = \det \mqty[
t & 0 & 0 & \cdots & 0 & \theta_1 \\
0 & t & 0 & \cdots & 0 & \theta_2 \\
& & & \ddots & & \\
0 & 0 & 0 & \cdots & t & \theta_{K-1} \\
-t & -t & -t & \cdots & & \qty( 1 - \theta_1 \cdots - \theta_{K-1} ) \\
] \\
& = \det \mqty[
t & 0 & 0 & \cdots & 0 & \theta_1 \\
0 & t & 0 & \cdots & 0 & \theta_2 \\
& & & \ddots & & \\
0 & 0 & 0 & \cdots & t & \theta_{K-1} \\
0 & 0 & 0 & \cdots & & 1 \\
] = t^{K-1}
\end{align}
2行目への変形は行基本変形を行いました。$t$の範囲は$ 0 $から$\infty$であったため
\begin{align}
& \Gamma( \phi_1 ) \cdots \Gamma( \phi_K ) \\
& = \int_D \dd \theta_1 \cdots \dd \theta_{K-1} \int_0^\infty \dd t ~ t^{K-1} \theta_1^{\phi_1-1} \cdots \theta_{K-1}^{\phi_{K-1}-1} \qty( 1 - \theta_1 \cdots - \theta_{K-1} )^{\phi_K -1} t^{ \phi_1 + \cdots + \phi_K - K } e^{ - t } \\
& = \int_D \dd \theta_1 \cdots \dd \theta_{K-1} ~ \theta_1^{\phi_1-1} \cdots \theta_{K-1}^{\phi_{K-1}-1} \qty( 1 - \theta_1 \cdots - \theta_{K-1} )^{\phi_K -1} \int_0^\infty \dd t ~ t^{ \phi_1 + \cdots + \phi_K - 1 } e^{ - t } \\
& = \int_D \dd \theta_1 \cdots \dd \theta_{K-1} ~ \theta_1^{\phi_1-1} \cdots \theta_{K-1}^{\phi_{K-1}-1} \qty( 1 - \theta_1 \cdots - \theta_{K-1} )^{\phi_K -1} \times \Gamma \qty( \phi_1 + \cdots + \phi_{K} )
\end{align}
と変形ができます。これが求めたかった式のため、規格化の証明できました。
ちなみに一般的な場合を示しましたが、に変数の場合には
\begin{align}
& \frac{ \Gamma( \phi_1 ) \Gamma( \phi_2 ) }{ \Gamma( \phi_1 + \phi_2 ) } = \int_0^1 d\theta ~ \theta^{\phi_1-1} (1 - \theta)^{\phi_2-1} =: \beta( \phi_1, \phi_2 )
\end{align}
というベータ関数とガンマ関数の公式になっています。
