２変量正規分布の累積分布関数
2変量正規分布の確率密度関数は平面上の局面であることが考えられる。確率密度関数の外形は地形図の等高線と同じであり、山の形をしている。\(\rho>0\)のとき\(y_2=y_1\)状にそって山ができる。山のほとんどは第１、3象限にある。\(x_i=\sigma_i y_i+\mu_i\)の逆変換をするとき、楕円をそれぞれの軸に対して\(\sigma_i\)倍させ、中心を\((\mu_1, \mu_2)\)にずらしたものと一致する。
単変量正規分布の累積分布関数の数値は数値表からわかる。次に示される
\begin{align}F(x_1, x_2)&=\mathrm{Pr}\{X_1\geq x_1, X_2\geq x_2\}\\&=\mathrm{Pr}\left\{\cfrac{X_1-\mu_1}{\sigma_1}\geq\cfrac{x_1-\mu_1}{\sigma_1}, \cfrac{X_2-\mu_2}{\sigma_2}\geq\cfrac{x_2-\mu_2}{\sigma_2}\right\}\\&=\mathrm{Pr}\left\{\cfrac{X_1-\mu_1}{\sigma_1}\geq y_1, \cfrac{X_2-\mu_2}{\sigma_2}\geq y_2\right\}\label{eq1}\tag{1}\end{align}
の数値はPearson(1931)により得られる。ここに上の\(y_1=(x_1-\mu_1)/\sigma_1\)、\(y_1=(x_2-\mu_2)/\sigma_2\)とする。また数地表についてはNasional Bureanof Standards(1959)に記載されている。これらの数地表に関する文献はGupta(1963)を参照されたい。またPearsonは
\begin{align}F(x_1,x_2)=\sum_{j=0}^{\infty}\rho^j \tau_j(y_1)\tau_j(y_2)\label{eq2}\tag{2}\end{align}
を示した。ここにtetrachoric function \(\tau_j(y)\)はPearson(1930)に記載されている。またHarris and Sons(1980)は\eqref{eq2}の一般化を行った。\(\tau_j(y)\)は次で定義される。
\begin{align}\tau_j(y)=\cfrac{(-1)^{j-1}}{\sqrt{j!}}\cfrac{d^{j-1}}{dx^{j-1}}\left(\cfrac{1}{\sqrt{2\pi}}e^{-\frac{y^2}{2}}\right),\ \ \ \ j=0,1,2,\ \ldots .\end{align}
Pearson(1930)により、\(\tau_0(y)\)から\(\tau_{19}(y)\)まで得られる。\eqref{eq2}より\eqref{eq1}の2変量正規分布の累積分布関数の数値を近似的に計算することができる。以下、\eqref{eq2}の証明をする。
累積分布関数の近似
確率ベクトル\(\boldsymbol{Y}\)の平均ベクトル、共分散行列をそれぞれ次で与える。
\begin{align}\boldsymbol{\mu}=\begin{pmatrix}0\\0\end{pmatrix},\ \ \ \ \begin{pmatrix}1&\rho\\\rho&1\end{pmatrix}.\end{align}
このとき特性関数\(\psi(t_1,t_2)\)は
\begin{align}\psi(t_1, t_2)&=\exp\left\{i(t_1\cdot0Lt_2\cdot0)-\cfrac{1}{2}(t_1^2+2\rho t_1t_2+t_2^2)\right\}\\&=\exp\left\{-\cfrac{t_1^2}{2}-\cfrac{t_2^2}{2}-\rho t_1t_2\right\}.\end{align}
である。ここで\(f(y_1, y_2, \rho)\)を次の確率密度関数とする。
\begin{align}f(y_1, y_2, \rho) =\cfrac{1}{2\pi\sqrt{1-\rho^2}}\exp\left\{-\cfrac{1}{2(1-\rho^2)}(y_1^2-2\rho y_1y_2+y_2^2)\right\}.\end{align}
次に特性関数を無限級数に展開する。\(\rho=0\)周りで\(\rho\)のべき級数に展開する。まず特性関数が
\begin{align}f(y_1, y_2, \rho)&=\exp\left\{-\cfrac{t_1^2}{2}-\cfrac{t_2^2}{2}-\rho t_1t_2\right\}\\&=\exp\left\{-\cfrac{t_1^2}{2}-\cfrac{t_2^2}{2}\right\}\sum_{k=0}^{\infty}\cfrac{\rho^k}{k!}(it_1)^k(it_2)^k \end{align}
と展開できることに注意する。補足として\(f(y_i)=(1/\sqrt{2\pi})e^{-\frac{y_i^2}{2}}\)の特性関数は\(\exp\{-\frac{t_i^2}{2}\}, i=1, 2\)であることから、上記の式の指数部分は標準正規分布の確率密度関数の積であることがわかる。ここでエルミーと多項式の性質より、確率密度関数が次のように展開される。
\begin{align}\label{eq3} f(y_1, y_2, \rho)&=\phi(y_1)\phi(y_2)\left\{\sum_{k=0}^{\infty}\cfrac{\rho^k}{k!}H_k(y_1)H_k(y_2)\right\},\tag{3}\end{align}
ここでエルミーと多項式は\begin{align}H_k(x)=(-1)^k\cfrac{\partial^k}{\partial x^k}\phi(x)/\phi(x)\end{align}で定義される。ここに\(\phi\)は標準正規分布の確率密度関数である。上記を示すのに、まずエルミーと多項式に対応する特性関数を求める。まず\(\int e^{it_jy_j}(-1)^k\tfrac{\partial^k}{\partial y_j^k}\phi(y_j)dy_j, j=1, 2\)は\(k=0, 1, 2, \ldots\)に対して次を満たすことを帰納法により示す。
\begin{align}\int e^{it_j y_j}(-1)^k\cfrac{\partial^k}{\partial y_j^k}\phi(y_j)dy_j&= (it_j)^k\exp\left\{-\cfrac{t_j^2}{2}\right\}.\end{align}
\(k=0\)のとき
\begin{align}\int e^{it_jy_j}(-1)^0\cfrac{\partial^0}{\partial y_j^0}\phi(y_j)dy_j&=\int e^{it_jy_j}\phi(y_j)dy_j\\&=e^{-\frac{t_j^2}{2}}\\&=(it_j)^0e^{-\frac{t_j^2}{2}}.\end{align}
である。\(k=l-1\)が成り立つと仮定し、\(k=l\)のとき
\begin{align}&\int e^{it_jy_j}(-1)^l\cfrac{\partial^l}{\partial y_j^l}\phi(y_j)dy_j\\&=\int e^{it_jy_j}(-1)^l\cfrac{\partial}{\partial y_j}\left(\cfrac{\partial^{l-1}}{\partial y_j^{l-1}}\phi(y_j)\right)dy_j\\&=\left[e^{it_jy_j}(-1)^l\cfrac{\partial^{l-1}}{\partial y_j^{l-1}}\phi(y_j)\right]_{-\infty}^{\infty}-\int(it_j)e^{it_jy_j}(-1)^l\cfrac{\partial^{l-1}}{\partial y_j^{l-1}}\phi(y_j)dy_j\\&=[T_{l-1}(y_j)e^{-\frac{1}{2}y_j^2+it_jy_j}]_{-\infty}^{\infty}-(-1)(it_j)\int e^{it_jy_j}(-1)^{l-1}\cfrac{\partial^{l-1}}{\partial y_j^{l-1}}\phi(y_j)dy_j\\&=[T_{l-1}(y_j)e^{-\frac{1}{2}y_j^2+it_jy_j}]_{-\infty}^{\infty}+(it_j)(it_j)^{l-1}\exp\left\{-\cfrac{y_j^2}{2}\right\}\\&=[T_{l-1}(y_j)e^{-\frac{1}{2}y_j^2+it_jy_j}]_{-\infty}^{\infty}+(it_j)^l\exp\left\{-\cfrac{t_j^2}{2}\right\}.\end{align}
このとき\(T_m(y_j), m=0, 1, 2 , \ldots\)は\(y_j\)のべき乗の和から成る多項式である。\(y_j\to\infty\)、または\(-\infty\)のとき\(T_m(y_j)\)と\(e^{-\frac{1}{2}y_j^2+it_jy_j}\)の積は\(0\)に収束する。具体例として、\(x>1\)に対して\(\lim_{a\to\infty}a^r/x^a=0\)であることと、\(y_j\to\infty\)のとき\(e^{-\frac{1}{2}y_j^2+it_jy_j}\to0\)であることから、\(\lim_{y_j\to\infty}y_j^r/e^{\frac{1}{2}y_j^2-it_jy_j}=\lim_{y_j\to-\infty}y_j^r/e^{\frac{1}{2}y_j^2-it_jy_j}=0\)であることがいえる。ここに\(r\)は整数値である。例として
\begin{align}T_0(y_j)&=\cfrac{1}{\sqrt{2\pi}}(-1),\\T_1(y_j)&=\cfrac{1}{\sqrt{2\pi}}(-y_j),\\T_2(y_j)&=\cfrac{1}{\sqrt{2\pi}}(-y_j^2+1),\\T_3(y_j)&=\cfrac{1}{\sqrt{2\pi}}(-y_j^3+3y_j)\end{align}である。よって\begin{align}&[T_{l-1}(y_j)e^{-\frac{1}{2}y_j^2}]_{-\infty}^{\infty}+(it_j)^l\exp\left\{-\cfrac{t_j^2}{2}\right\}\\&=0+(it_j)^l\exp\left\{-\cfrac{t_j^2}{2}\right\}\\&=(it_j)^l\exp\left\{-\cfrac{t_j^2}{2}\right\}\end{align}
である。したがって
\begin{align}\int e^{it_jy_j}(-1)^k\cfrac{\partial^k}{\partial y_j^k}\phi(y_j) dy_j=(it_j)^k\exp\left\{-\cfrac{t_j^2}{2}\right\}\end{align}
が示せた。このことから
\begin{align}&\int\int e^{it_1y_1+it_2y_2}H_k(y_1)H_k(y_2)\phi(y_1)\phi(y_2)dy_1dy_2\\&=\int\int e^{it_1y_1+it_2y_2}(-1)^k\cfrac{\partial^k}{\partial y_1^k}\phi(y_1)\times (-1)^k\cfrac{\partial^k}{\partial y_2^k}\phi(y_2)dy_1dy_2\\&=\int e^{it_1y_1}(-1)^k\cfrac{\partial^k}{\partial y_1^k}\phi(y_1)dy_1\int e^{it_2y_2}(-1)^k\cfrac{\partial^k}{\partial y_2^k}\phi(y_2)dy_2\\&=(it_1)^k(it_2)^k\exp\left\{-\cfrac{t_1^2}{2}-\cfrac{t_2^2}{2}\right\}.\end{align}
したがって\(\rho^k/k!\)を掛けて、\(k=0,1,2,\ldots\)について、和をとると次を得る。
\begin{align}\exp\left\{-\cfrac{t_1^2}{2}-\cfrac{t_2^2}{2}\right\}\sum_{k=0}^{\infty}\cfrac{\rho^k}{k!}(it_1)^k(it_2)^k&=\int\int e^{it_1y_1+it_2y_2}\left\{\sum_{k=0}^{\infty}\cfrac{\rho^k}{k!}H_k(y_1)H_k(y_2)\right\}\\&\ \ \ \ \times\phi(y_1)\phi(y_2)dy_1dy_2\\&=\int\int e^{it_1y_1+it_2y_2}\phi(y_1)\phi(y_2)\\&\ \ \ \ \times\left\{\sum_{k=0}^{\infty}\cfrac{\rho^k}{k!}H_k(y_1)H_k(y_2)\right\}dy_1dy_2.\end{align}
よって\eqref{eq3}の確率密度関数\(f(y_1, y_2, \rho)\)の展開が示された。
次に\(Y_1, Y_2\)が独立であり、同時分布が周辺分布の積で表せることに注意し、級数展開した確率密度関数を\(y_1, y_2\)についてこう別積分をすると次の\(Y_1, Y_2\)の累積分布関数を得る。
\begin{align}&\int\int f(y_1, y_2 , \rho)dy_1dy_2 \\&= \int\int \phi(y_1)\phi(y_2)\left\{\sum_{k=0}^{\infty}\cfrac{\rho^k}{k!}H_k(y_1)H_k(y_2)\right\}dy_1dy_2\\&=\int\int \phi(y_1)\phi(y_2)\Biggl\{1+\cfrac{\rho}{1!}(-1)\cfrac{\partial \phi(y_1)}{\partial y_1}\cfrac{1}{\phi(y_1)}(-1)\cfrac{\partial \phi(y_2)}{\partial y_2}\cfrac{1}{\phi(y_2)}\\&\ \ \ \ + \cfrac{\rho^2}{2!}(-1)^2\cfrac{\partial^2 \phi(y_1)}{\partial y_2^2}\cfrac{1}{\phi(y_1)}(-1)^2\cfrac{\partial^2\phi(y_2)}{\partial y_2^2}\cfrac{1}{\phi(y_2)}+\cdots\Biggr\}dy_1dy_2\\&=\int\int \phi(y_1)\phi(y_2)+\left\{\cfrac{\rho}{1!}\phi(y_1)\phi(y_2)+\cfrac{\rho^2}{2!}\cfrac{\partial^2\phi(y_1)}{\partial y_1^2}\cfrac{\partial^2\phi(y_2)}{\partial y_2^2}+\cdots\right\}dy_1dy_2\\&=\Phi(y_1)\Phi(y_2)+\left\{\cfrac{\rho}{1!}\phi(y_1)\phi(y_2)+\cfrac{\rho^2}{2!}\cfrac{\partial\phi(y_1)}{\partial y_1}\cfrac{\partial\phi(y_2)}{\partial y_2}+\cdots\right\}\\&=\Phi(y_1)\Phi(y_2)+\phi(y_1)\phi(y_2)\sum_{k=0}^{\infty}\cfrac{\rho^{k+1}}{(k+1)!}H_k(y_1)H_k(y_2).\end{align}
ここに\(\Phi\)は標準正規分布の累積分布関数である。これは\eqref{eq2}であり、近似的に類正規分布関数の値を計算することができる。
