楕円型分布の特性関数とモーメントについて解説する。
楕円型分布が直交変換に関し不変であるように、楕円型分布の特性関数も不変であることをみていく。
また、楕円型分布のモーメントの導出をおこなう。多変量正規分布と似た性質を持つことが分かる。
楕円型分布の特性関数
楕円型分布を持つ確率ベクトル\(\boldsymbol{Y}\)の特性関数\(\mathrm{E}[e^{i\boldsymbol{t}^T\boldsymbol{Y}}]\)は直交変換に関して不変である性質を持つ。\(\boldsymbol{OO}^T = \boldsymbol{I}\)を満たす直交行列\(\boldsymbol{O}\)に対して、\(\boldsymbol{OY}\)の特性関数が不変であることを示す。\(|\partial(y_1, \ldots, y_p)/ \partial(z_1, \ldots, z_p)| = 1\)より、\(\boldsymbol{Z}\)の特性関数は次のように導出される。
\begin{align}\mathrm{E}[e^{i\boldsymbol{t}^T\boldsymbol{OY}}] &= \int_{-\infty}^{\infty}\cdots\int_{-\infty}^{\infty}e^{i\boldsymbol{t}^T\boldsymbol{OY}}g(\boldsymbol{y}^T\boldsymbol{y}) dy_1\cdots dy_p\\&=\int_{-\infty}^{\infty}\cdots\int_{-\infty}^{\infty} g\bigl[(\boldsymbol{O}^T\boldsymbol{z})^T\boldsymbol{O}^T\boldsymbol{z}\bigr]dz_1\cdots dz_p\\&=\int_{-\infty}^{\infty}\cdots\int_{-\infty}^{\infty}e^{i\boldsymbol{t}^T\boldsymbol{OO}^T\boldsymbol{z}}g(\boldsymbol{z}^T\boldsymbol{OO}^T\boldsymbol{z}) dz_1 \cdots dz_p\\&=\int_{-\infty}^{\infty}\cdots\int_{-\infty}^{\infty}e^{i\boldsymbol{t}^T\boldsymbol{z}}g(\boldsymbol{z}^T\boldsymbol{z}) dz_1 \cdots dz_p \\\label{eq1}&=\mathrm{E}[e^{i\boldsymbol{t}^T\boldsymbol{Z}}],\tag{1}\end{align}
ここに\(\boldsymbol{Z} = \boldsymbol{OY}\)は確率密度関数\(g(\boldsymbol{y}^T\boldsymbol{y})\)をもつ。したがって、直交変換に関して特性関数は不変である。すべての直交行列\(\boldsymbol{O}\)に対して\eqref{eq1}が成り立つことは、\(\mathrm{E}[e^{i\boldsymbol{t}^T\boldsymbol{Z}}]\)は\(\boldsymbol{t}^T\boldsymbol{t}\)の関数であることを意味する。例として\(\phi(\boldsymbol{t}^T\boldsymbol{t}) = \mathrm{E}[e^{i\boldsymbol{t}^T\boldsymbol{Z}}]\)とおき、直交変換\(\boldsymbol{O}^T\boldsymbol{t}\)を行うと、\(\phi[(\boldsymbol{O}^T\boldsymbol{t}))^T\boldsymbol{O}^T\boldsymbol{t}] = \phi(\boldsymbol{t}^T\boldsymbol{t})\)となることが分かる。次に、確率ベクトル\(\boldsymbol{Y}\)の特性関数を次のように書くことにする。
\begin{align}\label{eq2}\mathrm{E}[e^{i\boldsymbol{t}^T\boldsymbol{Y}}] = \phi(\boldsymbol{t}^T\boldsymbol{t})\tag{2}.\end{align}
このとき、\(\boldsymbol{\Lambda} = \boldsymbol{CC}^T\)であるような正則行列\(\boldsymbol{C}\)を用いて、\(\boldsymbol{X} = \boldsymbol{\mu} + \boldsymbol{CY}\)の変換を考える。\(\boldsymbol{X}\)の特性関数は次となる。
\begin{align}\mathrm{E}[e^{i\boldsymbol{t}^T\boldsymbol{X}}] &= \mathrm{E}[e^{i\boldsymbol{t}^T(\boldsymbol{\mu} + \boldsymbol{CY})}]\\&=e^{i\boldsymbol{t}^T\boldsymbol{\mu}} \mathrm{E}[e^{i\boldsymbol{t}^T\boldsymbol{CY}}]\\&= e^{i\boldsymbol{t}^T\boldsymbol{\mu}} \phi\bigl(\boldsymbol{t}^T\boldsymbol{C}(\boldsymbol{t}^T\boldsymbol{C})^T\bigr)\\&= e^{i\boldsymbol{t}^T\boldsymbol{\mu}} \phi(\boldsymbol{t}^T\boldsymbol{C}\boldsymbol{C}^T\boldsymbol{t})\\\label{eq3}&= e^{i\boldsymbol{t}^T\boldsymbol{\mu}} \phi(\boldsymbol{t}^T\boldsymbol{\Lambda t}).\tag{3}\end{align}
このことから、\(e^{i\boldsymbol{t}^T\boldsymbol{\mu}}\phi(\boldsymbol{t}^T\boldsymbol{\Lambda t})\)の形の任意の特性関数は、楕円型分布の(2)式の確率密度関数を持つ確率ベクトル\(\boldsymbol{X}\)に対応する。
楕円型分布のモーメント
楕円型分布を持つ\(\boldsymbol{X}\)のモーメントは特性関数\(e^{i\boldsymbol{t}^T\boldsymbol{\mu}}\phi(\boldsymbol{t}^T\boldsymbol{\Lambda t})\)または\(\boldsymbol{X} = \boldsymbol{\mu} + R\boldsymbol{CU}\)から求めることができる。ここに、\(\boldsymbol{C}^T\boldsymbol{\Lambda}^{-1} \boldsymbol{C} = \boldsymbol{I}\)である。このとき、次が成り立つ。
定理1 楕円型分布のモーメント
証明 \(R\)の\(2m\)次モーメントを求めるために、\(\boldsymbol{Y}\)の特性関数を次のように書き換える。
\begin{align}\mathrm{E}[\exp(i\boldsymbol{t}^T\boldsymbol{Y})] &= \pi(\boldsymbol{t}^T\boldsymbol{t})\\ &= \sum_{k=0}^{\infty}\cfrac{\pi^{(k)}(0)}{k!}(\boldsymbol{t}^T\boldsymbol{t})^k\\&= \mathrm{E}\left[_0F_1\left(\tfrac{1}{2}o; -\tfrac{1}{4}\boldsymbol{t}^T\boldsymbol{t}\boldsymbol{Y}^T\boldsymbol{Y}\right)^k\right]\\&=\mathrm{E}\left[\sum_{k=0}^{\infty}\cfrac{1}{k!(p/2)_k} \left(-\tfrac{1}{4}\boldsymbol{t}^T\boldsymbol{t}\boldsymbol{Y}^T\boldsymbol{Y}\right)^k\right]\\&=\sum_{k=0}^{\infty}\cfrac{1}{k!(p/2)_k}\left(-\tfrac{1}{4}\right)^k\mathrm{E}\bigl[(\boldsymbol{Y}^T\boldsymbol{Y})^k\bigr](\boldsymbol{t}^T\boldsymbol{t})^k,\end{align}
ここに、\(_0F_1\)はベッセル関数（Bessel function）である。両辺の\((\boldsymbol{t}^T\boldsymbol{t})^k\)の係数を比べると
\begin{align}&\sum_{k = 0}^{\infty}\cfrac{\phi^{(k)}(0)}{k!}(\boldsymbol{t}^T\boldsymbol{t})^k = \sum_{k=0}^{\infty} \cfrac{1}{k!(p/2)_k} \left(-\tfrac{1}{4}\right)^k\mathrm{E}\bigl[(\boldsymbol{Y}^T\boldsymbol{Y})^k\bigr](\boldsymbol{t}^T\boldsymbol{t})^k\\ &\Leftrightarrow \sum_{k=0}^{\infty}\cfrac{\phi^{(k)}(0)}{k!} = \sum_{k=0}^{\infty}\cfrac{1}{k!(p/2)_k}\left(-\tfrac{1}{4}\right)^k\mathrm{E}\bigl[(\boldsymbol{Y}^T\boldsymbol{Y})^k\bigr]\\&\Leftrightarrow \phi^{(k)}(0) = \cfrac{1}{(p/2)_k}\left(-\tfrac{1}{4}\right)^k\mathrm{E}\bigl[(\boldsymbol{Y}^T\boldsymbol{Y})^k\bigr]\\&\Leftrightarrow \mathrm{E}\bigl[(\boldsymbol{Y}^T\boldsymbol{Y})^k\bigr] = (-4)^k\left(p/2\right)_k\phi^{(k)}(0).\end{align}
よって\(\mathrm{E}[R^{2m}]\)は
\begin{align}\mathrm{E}[R^{2m}] &= \mathrm{E}\bigl[(\boldsymbol{y}^T\boldsymbol{y})^m\bigr]\\&=(-4)^m\left(p/2\right)_m\pi^{(m)}(0)\end{align}
である。このことから、\(\mathrm{E}[R^2]\)と\(\mathrm{R^4}\)はそれぞれ
\begin{align}\mathrm{E}[R^2] &= (-4)^1 (p/2)_1 \phi'(0)\\&=(-4)(p/2)\phi'(0)\\&=-2p\phi'(0),\end{align}\begin{align}\mathrm{E}[R^4] &= (-4)^2(p/2)_2\phi''(0)\\&=16\bigl\{(p/2)(p/2 + 1)\bigr\}\phi''(0)\\&=4p(p+2)\phi''(0)\end{align}
となる。□
\(\boldsymbol{Y}= R\boldsymbol{U}\)の高次モーメントを考える。\(R = \boldsymbol{Y}^T\boldsymbol{Y}\)の二次形式の奇数次モーメントを計算することで、\(R\)の奇数次モーメントは\(0\)であることが分かる。故に、\(\boldsymbol{Y}\)の奇数次モーメントは\(0\)である。次がいえる。
\begin{align}\label{eq6}\mathrm{E}\bigl[(X_i-\mu_i)(X_j-\mu_j)(X_k-\mu_k)\bigr] = 0.\tag{6}\end{align}
実際に、\(\boldsymbol{X} - \boldsymbol{\mu}\)の全ての奇数次モーメントは\(0\)である。
\(\mathrm{E}[U_iU_jU_kU_l]\)を考える。\(\boldsymbol{U}^T\boldsymbol{U} =1\)より
\begin{align} 1 &= \mathrm{E}\bigl([\boldsymbol{U}^T\boldsymbol{U})^2\bigr]\\&= \sum_{i,j = 1}^p\mathrm{E}[U_i^2U_j^2]\\&=\sum_{i=1}^p\mathrm{E}[U_i^4] + \sum_{\substack{i,j=1\\i\neq j}}^p\mathrm{E}[U_i^2U_j^2]\\\label{eq7}&= p\mathrm{E}[U_1^4] + p(p-1)\mathrm{E}[U_1^2U_2^2]\tag{7}\end{align}
\(\mathrm{E}[U_1^4]\)を求める。\(\mathrm{E}[U_1^4] =\mathrm{E}[\sin^4\Theta_1]\)であることから、このモーメントを計算する。
\begin{align}\mathrm{E}[\sin^4\Theta_1^4] &= \int_{-\pi/2}^{\pi/2} \sin^4\theta_1 \cfrac{\Gamma(\frac{1}{2}p)\cos^{p-2}\theta_1}{\Gamma(\frac{1}{2})\Gamma\bigl[\frac{1}{2}(p-1)\bigr]}d\theta_1\\&= \cfrac{\Gamma(\frac{1}{2}p)}{\Gamma(\frac{1}{2})\Gamma\bigl[\frac{1}{2}(p-1)\bigr]}\int_{-\pi/2}^{\pi/2}\sin^4\theta_1\cos^{p-2}\theta_1 d\theta_1\\&=\cfrac{\Gamma(\frac{1}{2}p)}{\Gamma(\frac{1}{2})\Gamma\bigl[\frac{1}{2}(p-1)\bigr]} \int_0^1u^{\frac{1}{2}(p-1) - 1}(1-u)^{\frac{5}{2}-1}du\\&=\cfrac{\Gamma(\frac{1}{2}p)}{\Gamma(\frac{1}{2})\Gamma\bigl[\frac{1}{2}(p-1)\bigr]} \cfrac{\Gamma\bigl[\frac{1}{2}(p-1)\bigr]\Gamma(\frac{5}{2})}{\Gamma\bigl[\frac{1}{2}(p+4)\bigr]} \\&= \frac{3}{p(p+2)}.\end{align}
したがって、\eqref{eq7}について次が成り立つ
\begin{align}&1 = p\mathrm{E}[U_1^4] + p(p-1)\mathrm{E}[U_1^2U_2^2]\\ & \Leftrightarrow 1 = \cfrac{3p}{p(p+2)} + p(p-1)\mathrm{E}[U_1^2U_2^2] \\&\Leftrightarrow \cfrac{p(p+2) - 3p}{p(p+2)}= p(p-1)\mathrm{E}[U_1^2U_2^2]\\ &\Leftrightarrow \mathrm{E}[U_1^2U_2^2] = \cfrac{1}{p(p+2)}\end{align}
よって、\(Y_i^4\)と\(Y_i^2Y_j^2\)のモーメントが次のように導出される。
\begin{align}\mathrm{E}[Y_i^4] &= \mathrm{E}\bigl[(RU_i)^4\bigr] \\&= \mathrm{E}[R^4]\mathrm{U_i^4}\\&= \cfrac{3\mathrm{E}[R^4]}{p(p+2)},\end{align}\begin{align}\mathrm{E}[Y_i^2Y_j^2] &= \mathrm{E}[R^4]\mathrm{E}[U_i^2U_j^2]\\&=\cfrac{\mathrm{E}[R^4]}{p(p+2)}.\end{align}
注意として、\(i=j=k=l\)、\(i=j\neq k=l\)、\(i=k\neq l=j\)、\(i=l\neq j=k\)出ない限り、\(\mathrm{E}[U_iU_jU_kU_l]=0\)である。これは黒ネッカーのデルタによって次のように表現できる。
\begin{align}\mathrm{E}[U_iU_jU_kU_l]= \cfrac{\delta_{ij}\delta_{kl} + \delta_{ik}\delta_{jl}+\delta_{il}\delta_{jk}}{p(p+2)}\end{align}
一般の場合の4次モーメント\(\mathrm{E}[(X_i - \mu_i )(X_j - \mu_j )(X_k- \mu_k )(X_l- \mu_l )]\)は次の楕円型分布のモーメントとキュムラントで導出する。
