今日はもう更新しないかもしれないと思ったけれども、一階線形微分方程式の一般解について考えてみたので記録する。一階微分方程式の一般形を
\[ y(x)^{\prime} + p(x) y(x) = g(x) \]
と書くと、斉次のときの一般解は
\[ y(x) = C e^{-\int p(x) dx}\]
である。ここまでの議論はスキップして、非斉次な場合を考える。斉次なときは、\( z(x) = y(x) e^{\int p(x) dx}\)は任意定数が一般解になるから\(z^{\prime} = 0\)が\(z\)の微分方程式。では、非斉次の場合も微分方程式がそのまま\(y\)で考えるより単純になりそうだと思って\(z\)の方程式を計算する。結果は
\[ z(x)^{\prime} = g(x) e^{\int p(x) dx}\]
になるので\(z\)はすぐに次のようになると分かる。
\[ z(x) = \int g(x) e^{\int p(x) dx} + C^{\prime} \]
ところで、\(z(x) = y(x) e^{\int p(x) dx}\)だから一般解は
\[ y(x) = \left( \int g(x) e^{\int p(x) dx} + C^{\prime} \right) e^{-\int p(x) dx}\]
である。これは定数変化法と同等だけれども、このように変形すると、定数を変化しているようにはちょっと見えない。
2013年6月2日日曜日
2013年6月1日土曜日
2013年5月31日金曜日
2013年5月30日木曜日
2013年5月29日水曜日
2013年5月28日火曜日
2013年5月27日月曜日
2013年5月26日日曜日
2013年5月23日木曜日
2013年5月22日水曜日
2013年5月21日火曜日
2013年5月20日月曜日
2013年5月19日日曜日
2013年5月18日土曜日
無相関検定に使う推定量の確率分布について(3)
問題
\[ \begin{align*} \boldsymbol{u}_{1} &= \frac{1}{\sqrt{n}}\sum_{i=1}^{n} \boldsymbol{e}_{i}\\ \boldsymbol{u}_{2} &= \sum_{i=1}^{n} c_{i}\boldsymbol{e}_{i} \end{align*} \]
を含む規格直交基底を構成できることを示せ。
解答
\(n\)個の\(\boldsymbol{u}_{1}, \boldsymbol{u}_{2}\)を含む線形独立なベクトルの組を見い出せば、グラムシュミット直交化法により、規格直交基底を構成することができる。そこでまず、\(\{\boldsymbol{e}_{i}\}\)の\(\boldsymbol{e}_{1}\)を\(\boldsymbol{u}_{1}\)に置き換えたものが線形独立であることを示してから、\(c_{k} \neq c_{1}\)になる添字を\(k\)として、\(\boldsymbol{e}_{k}\)と\(\boldsymbol{u}_{2}\)を置き換えたものが線形独立になることを示す。
\(\{\boldsymbol{u}_{1}, \boldsymbol{e}_{2}, \boldsymbol{e}_{3}, \ldots , \boldsymbol{e}_{n}\}\)が線形従属と仮定すると
\[ \sum_{i=2}^{n} \alpha_{i} \boldsymbol{e}_{i} + \alpha_{1} \boldsymbol{u}_{1} = 0\]
となる\(\boldsymbol{\alpha}\)が0ベクトル以外に存在する。\(\alpha_{1} = 0\)のときは、\(\{\boldsymbol{e}_{2}, \boldsymbol{e}_{3}, \ldots , \boldsymbol{e}_{n}\}\)が線形独立であることに矛盾し、\(\alpha_{1} \neq 0\)のときは、\(\boldsymbol{u}_{1}\)からくる\(\boldsymbol{e}_{1}\)成分を打ち消すベクトルが存在しない。従ってこの仮定は矛盾するので、\(\{\boldsymbol{u}_{1}, \boldsymbol{e}_{1}, \ldots , \boldsymbol{e}_{n}\}\)は線形独立である。
次に\(\boldsymbol{u}_{2}\)と\(\boldsymbol{e}_{k}\)を置き換え、添字の\(k\)と\(2\)を交換する。すると先ほどと同様にn個のベクトルの組は線形従属だと仮定すると
\[ \boldsymbol{V} = \sum_{i=3}^{n} \alpha_{i} \boldsymbol{e}_{i} + \alpha_{1} \boldsymbol{u}_{1} + \alpha_{2} \boldsymbol{u}_{2} = 0\]
が成り立つ。\(\alpha_{2} = 0\)とすると、残りのベクトルの組が線形従属になり矛盾する。\(\alpha_{2} \neq 0\)のときは、\(\boldsymbol{u}_{i} \cdot \boldsymbol{V} = 0\)から
\[
\begin{align*}
\alpha_{1} &= - \frac{n-2}{\sqrt{n}}\\
\alpha_{2} &= c_{1} + c_{2} \end{align*}
\]
が得られる。1成分と2成分以外は\(\{\boldsymbol{e}_{i}\}\)で打ち消せるので、1, 2成分が0になるように係数を調整する。\( \boldsymbol{e}_{1} \cdot \boldsymbol{V} = \boldsymbol{e}_{2} \cdot \boldsymbol{V} = 0\)から条件はもとまり
\[ \begin{align*} c_{1} ( c_{1} + c_{2} ) = \frac{n - 2}{n}\\ c_{2} ( c_{1} + c_{2} ) = \frac{n - 2}{n} \end{align*} \]
となる。従って\(c_{1} = c_{2}\)ならば、線形従属になる。もし置き換えとして選んだ\(\boldsymbol{e}_{k}\)がたまたまこの条件を満たすとしても、\(k = 2, 3, \ldots n\)全てについて成り立つことはない。なぜならば、これが成り立つと\(c_{1} = c_{2} = c_{3} = \ldots = c_{n}\)となり、\(\boldsymbol{u}_{2}\)が\(\boldsymbol{u}_{1}\)のスカラー倍になってしまい、互いに直交しなくなる。以上より、必ず\(c_{1} \neq c_{k}\)にできるので、\(\{\boldsymbol{u}_{1}, \boldsymbol{u}_{2}, \boldsymbol{e}_{3}, \boldsymbol{e}_{4} , \ldots , \boldsymbol{e}_{n} \}\)は線形独立である。
が得られる。1成分と2成分以外は\(\{\boldsymbol{e}_{i}\}\)で打ち消せるので、1, 2成分が0になるように係数を調整する。\( \boldsymbol{e}_{1} \cdot \boldsymbol{V} = \boldsymbol{e}_{2} \cdot \boldsymbol{V} = 0\)から条件はもとまり
\[ \begin{align*} c_{1} ( c_{1} + c_{2} ) = \frac{n - 2}{n}\\ c_{2} ( c_{1} + c_{2} ) = \frac{n - 2}{n} \end{align*} \]
となる。従って\(c_{1} = c_{2}\)ならば、線形従属になる。もし置き換えとして選んだ\(\boldsymbol{e}_{k}\)がたまたまこの条件を満たすとしても、\(k = 2, 3, \ldots n\)全てについて成り立つことはない。なぜならば、これが成り立つと\(c_{1} = c_{2} = c_{3} = \ldots = c_{n}\)となり、\(\boldsymbol{u}_{2}\)が\(\boldsymbol{u}_{1}\)のスカラー倍になってしまい、互いに直交しなくなる。以上より、必ず\(c_{1} \neq c_{k}\)にできるので、\(\{\boldsymbol{u}_{1}, \boldsymbol{u}_{2}, \boldsymbol{e}_{3}, \boldsymbol{e}_{4} , \ldots , \boldsymbol{e}_{n} \}\)は線形独立である。
2013年5月17日金曜日
2013年5月16日木曜日
無相関検定に使う推定量の確率分布について(2)
問題
\[
\begin{align*}
\boldsymbol{u}_{1} &= \frac{1}{\sqrt{n}}\sum_{i=1}^{n} \boldsymbol{e}_{i}\\
\boldsymbol{u}_{2} &= \sum_{i=1}^{n} c_{i}\boldsymbol{e}_{i}
\end{align*}
\]
を含む規格直交基底を構成できることを\(n=3\)について示せ。
解答
まず線形独立なベクトルの組を選ぶ。そうすれば、グラムシュミットの直交化によって、正規直交基底を構成することができる。ベクトルの組が線形独立かどうかの判定は、グラム行列の行列式が0で無いことにより確認する。
ベクトルの組として\(\boldsymbol{u}_{1}, \boldsymbol{u}_{2}, \boldsymbol{e}_{i}\)を選んだとする(\(i=1, 2, 3\))。するとグラム行列\(G_{i}\)は次のようになる。
\[ G_{i} = \left(
\begin{array}{ccc}
1 & 0 & 1/\sqrt{n}\\
0 & 1 & c_{i}\\
1/\sqrt{n}& c_{i} & 1
\end{array}
\right)
\]
行列式の計算は、3行目に\(-1/\sqrt{n} \text{1行目} - c_{i} \text{2行目}\)をたすと、行列式の値を変えることなく三角行列に変形できるので、対角成分の積から行列式は\(1 - 1/n - c_{i}^{2}\)だと分かる。仮に\(|G_{i}| = |G_{j}| = 0~~(i \neq j)\)とすると、\(c_{i}^{2} + c_{j}^{2} > 1\)になるので、\(\sum_{i=1}^{3} c_{i}^{2} = 1\)に矛盾する。従って、2つ以上のグラム行列の行列式が0になることは無いので、線形独立になる組は\(i = 1, 2, 3\)のいずれかだと言える。
\(\boldsymbol{u}_{1}, \boldsymbol{u}_{2}, \boldsymbol{e}_{3}\)が線形独立とする。この3つのベクトルから、正規直交基底をグラムシュミットの直交化の方法によって、正規直交基底を構成する。
\[ \boldsymbol{v}^{(0)} = \boldsymbol{e}_{3} - (\boldsymbol{u}_{1} \cdot \boldsymbol{e}_{3}) \boldsymbol{u}_{1} \]
は\(\boldsymbol{u}_{1}\)と直交する。\(\boldsymbol{u}_{1}^{(0)} = \boldsymbol{v}^{(0)}/|\boldsymbol{v}^{(0)}|\)として
\[ \boldsymbol{v}^{(1)} = \boldsymbol{u}^{(0)} - (\boldsymbol{u}^{(0)}\cdot \boldsymbol{u}_{2}) \boldsymbol{u}_{2}\]
を定義する。これは\(\boldsymbol{u}_{1}, \boldsymbol{u}_{2}\)と直交するので、これを大きさ1に規格化したものと\(\boldsymbol{u}_{3}\)とする。
\[ \boldsymbol{u}_{3} = \frac{\boldsymbol{v}^{(2)}}{|\boldsymbol{v}^{(2)}|}\]
こうして構成したベクトル\(\{\boldsymbol{u}_{i}\}\)は規格直交基底をなす。
行列式の計算は、3行目に\(-1/\sqrt{n} \text{1行目} - c_{i} \text{2行目}\)をたすと、行列式の値を変えることなく三角行列に変形できるので、対角成分の積から行列式は\(1 - 1/n - c_{i}^{2}\)だと分かる。仮に\(|G_{i}| = |G_{j}| = 0~~(i \neq j)\)とすると、\(c_{i}^{2} + c_{j}^{2} > 1\)になるので、\(\sum_{i=1}^{3} c_{i}^{2} = 1\)に矛盾する。従って、2つ以上のグラム行列の行列式が0になることは無いので、線形独立になる組は\(i = 1, 2, 3\)のいずれかだと言える。
\(\boldsymbol{u}_{1}, \boldsymbol{u}_{2}, \boldsymbol{e}_{3}\)が線形独立とする。この3つのベクトルから、正規直交基底をグラムシュミットの直交化の方法によって、正規直交基底を構成する。
\[ \boldsymbol{v}^{(0)} = \boldsymbol{e}_{3} - (\boldsymbol{u}_{1} \cdot \boldsymbol{e}_{3}) \boldsymbol{u}_{1} \]
は\(\boldsymbol{u}_{1}\)と直交する。\(\boldsymbol{u}_{1}^{(0)} = \boldsymbol{v}^{(0)}/|\boldsymbol{v}^{(0)}|\)として
\[ \boldsymbol{v}^{(1)} = \boldsymbol{u}^{(0)} - (\boldsymbol{u}^{(0)}\cdot \boldsymbol{u}_{2}) \boldsymbol{u}_{2}\]
を定義する。これは\(\boldsymbol{u}_{1}, \boldsymbol{u}_{2}\)と直交するので、これを大きさ1に規格化したものと\(\boldsymbol{u}_{3}\)とする。
\[ \boldsymbol{u}_{3} = \frac{\boldsymbol{v}^{(2)}}{|\boldsymbol{v}^{(2)}|}\]
こうして構成したベクトル\(\{\boldsymbol{u}_{i}\}\)は規格直交基底をなす。
2013年5月15日水曜日
無相関検定に使う推定量の確率分布について(1)
問題
\[ \begin{align*} r &= \frac{S_{XY}}{S_{XX}S_{YY}}\\ S_{AB} &= \frac{\sum_{i=1}^{n}(A_{i}-\bar{A})(B_{i}-\bar{B})}{n} \end{align*} \]で定義される相関係数(ピアソンの積率相関係数)\(r\)使って定義した量
\[t = \sqrt{n-2}\frac{r}{\sqrt{1-r^{2}}}\]
は、\(X,Y\)に相関が無い場合、自由度n-2のt分布に従うことを示せ。
考えたこと
\[
\begin{align*}
&c_{i} = \frac{x_{i}-\bar{x}}{\sqrt{nS_{XX}}}\\
&\sum_{i=1}^{n}c_{i} = 0, ~~~\sum_{i=1}^{n}(c_{i})^{2} = 1
\end{align*}
\]
で定義した\(c_{i}\)を使って
\[ \begin{align*} Z_{1} &= \sqrt{n}\bar{Y} = Z_{1}\\ Z_{2} &= \sum_{i=1}^{n}c_{i}Y_{i} \end{align*} \]
を定義する。直交行列\(\mathcal{O}\)を
\[ \mathcal{O}^{t} = (\boldsymbol{u}_{1}, \boldsymbol{u}_{2}, \ldots ,\boldsymbol{u}_{3})\]
として、\(\mathcal{O}\boldsymbol{Y} = \boldsymbol{Z}\)とすることができる。実際、
\[ \begin{align*} \boldsymbol{u}_{1} &= \frac{1}{\sqrt{n}}\sum_{i=1}^{n} \boldsymbol{e}_{i}\\ \boldsymbol{u}_{2} &= \sum_{i=1}^{n} c_{i}\boldsymbol{e}_{i} \end{align*} \]
となり、これらは互いに直交し、大きさも1に規格化されている。残りのベクトルについても、適当に基底を選んで、全体として規格直交ベクトルになるようにできれば、3成分以降の\(Z_{i}\)も定義されるので、 \[\sum_{i=1}^{n} Y_{i}^{2} = \sum_{i=1}^{n}Z_{i}^{2}\] となる。これを使うと
\[ t = \frac{Z_{2}}{\sqrt{\sum_{i=3}^{n}(Z_{i})^{2}/(n-2)}}\]
と書くことができる。\(Y_{i} \sim \mathrm{N}(0, 1)\)だから、正規分布の再生性から\(Z_{i} \sim \mathrm{N}(0, 1)\)なので、\(t\)は自由度\(n-2\)のt分布に従う。
で定義した\(c_{i}\)を使って
\[ \begin{align*} Z_{1} &= \sqrt{n}\bar{Y} = Z_{1}\\ Z_{2} &= \sum_{i=1}^{n}c_{i}Y_{i} \end{align*} \]
を定義する。直交行列\(\mathcal{O}\)を
\[ \mathcal{O}^{t} = (\boldsymbol{u}_{1}, \boldsymbol{u}_{2}, \ldots ,\boldsymbol{u}_{3})\]
として、\(\mathcal{O}\boldsymbol{Y} = \boldsymbol{Z}\)とすることができる。実際、
\[ \begin{align*} \boldsymbol{u}_{1} &= \frac{1}{\sqrt{n}}\sum_{i=1}^{n} \boldsymbol{e}_{i}\\ \boldsymbol{u}_{2} &= \sum_{i=1}^{n} c_{i}\boldsymbol{e}_{i} \end{align*} \]
となり、これらは互いに直交し、大きさも1に規格化されている。残りのベクトルについても、適当に基底を選んで、全体として規格直交ベクトルになるようにできれば、3成分以降の\(Z_{i}\)も定義されるので、 \[\sum_{i=1}^{n} Y_{i}^{2} = \sum_{i=1}^{n}Z_{i}^{2}\] となる。これを使うと
\[ t = \frac{Z_{2}}{\sqrt{\sum_{i=3}^{n}(Z_{i})^{2}/(n-2)}}\]
と書くことができる。\(Y_{i} \sim \mathrm{N}(0, 1)\)だから、正規分布の再生性から\(Z_{i} \sim \mathrm{N}(0, 1)\)なので、\(t\)は自由度\(n-2\)のt分布に従う。
次回
適当に\(\boldsymbol{u}_{i}\)をとれば\(t\)分布に従うことは分かったので、次回は \(n=3\)の場合に、そのようなベクトルを構成できることを確かめる。2013年5月14日火曜日
火曜日
多項分布が多次元正規分布に近似する計算は、ある1つの属性のカテゴリの生起確率が与えられているときの、適合度検定に使う推定量が従う確率分布の自由度を求めるためだった。2つの属性が独立で、生起確率が不明なとき、最尤推定値に置き換えて計算した確率分布も、2つの属性の独立性の検定に使う推定量が従う確率分布の自由度を求めるためだった。
最初は、検定の説明から始めて、適合度検定と独立性の検定について、説明を書こうとも思っていた。あまり、それを書く意味について考えていなかったけれども、あらためて考えてみると、今の自分の理解では、書く価値は無さそうだと感じた。
かわりに、新しい計算として、無相関検定に使う推定量が従う確率分布の計算を始めた。上の2つの計算もそうだったけれども、今回の計算も、手持ちの教科書には、しっかりとした議論、または計算が載っていないので、またちびちびと計算していこうと思う。
最初は、検定の説明から始めて、適合度検定と独立性の検定について、説明を書こうとも思っていた。あまり、それを書く意味について考えていなかったけれども、あらためて考えてみると、今の自分の理解では、書く価値は無さそうだと感じた。
かわりに、新しい計算として、無相関検定に使う推定量が従う確率分布の計算を始めた。上の2つの計算もそうだったけれども、今回の計算も、手持ちの教科書には、しっかりとした議論、または計算が載っていないので、またちびちびと計算していこうと思う。
2013年5月13日月曜日
2種の属性が独立であるときの確率分布
問題
独立な2種の属性\(A, B\)のカテゴリー数をそれぞれ\(a, b (\ge 2)\)とする。\(n\)個のサンプルを取った時、\(A = A_{i},~~B = B_{j}\)をとる度数はそれぞれ\(n_{i j}\)とする。\(A = A_{i},~~B = B_{j}\)の生起確率は、属性は独立なので\(p_{ij} = p_{i\circ} p_{\circ j}\)である。\(p_{i \circ}, p_{\circ j}\)の真の値は不明として、最尤推定値
\[ p_{i \circ} = \frac{n_{i\circ}}{n},~~p_{\circ j} = \frac{n_{\circ j}}{n}\]
を使う。このとき
\[ z_{ij} = \frac{n_{ij} - np_{ij}}{\sqrt{np_{ij}}}\]
を実現値とする、確率変数\(Z_{ij}\)の2乗の和
\[Z^{2} = \sum_{i=1}^{a}\sum_{j=1}^{b} (Z_{ij})^{2}\]
の従う確率分布は自由度\((a-1)(b-1)\)の\(\chi^{2}\)分布であることを示せ。
\[ p_{i \circ} = \frac{n_{i\circ}}{n},~~p_{\circ j} = \frac{n_{\circ j}}{n}\]
を使う。このとき
\[ z_{ij} = \frac{n_{ij} - np_{ij}}{\sqrt{np_{ij}}}\]
を実現値とする、確率変数\(Z_{ij}\)の2乗の和
\[Z^{2} = \sum_{i=1}^{a}\sum_{j=1}^{b} (Z_{ij})^{2}\]
の従う確率分布は自由度\((a-1)(b-1)\)の\(\chi^{2}\)分布であることを示せ。
解答
\(n\)回の測定によって得られた度数分布が得られる確率は、多項分布に従っている。\(n\)が大きいことを使うと更に、これを正規分布に近似することができる。指数関数の肩は\(-z^{2}/2\)になる。
今回\(p_{ij}\)を推定値に置き換えたことによって、\((a-1) + (b-1)\)個の条件式が課されている。この条件を\(z_{ij}\)で書き直すと
\[
\begin{align*}
&\sum_{j=1}^{a} z_{ij} \sqrt{p_{\circ j}} = 0\\
&\sum_{i=1}^{b} z_{ij} \sqrt{p_{i\circ}} = 0 (*)\end{align*}
\]
である。これは、\(np_{ij} = n_{i\circ}p_{\circ j} = n_{\circ j}p_{i \circ}\)を使うと確率変数が
\[ z_{ij} = \frac{n_{ij} - n_{i\circ}p_{\circ j}}{\sqrt{n_{i\circ}p_{\circ j}}} = \frac{n_{ij} - n_{\circ j}p_{i \circ}}{\sqrt{n_{\circ j}p_{i \circ}}}\]
と書けるので、あとは\(\sum_{i} p_{i\circ} = \sum_{j} p_{\circ j} = 1\)と\(\sum_{i}n_{ij} = n_{\circ j},~~\sum_{j} n_{ij} = n_{i \circ}\)を使うと得られる。
この得られた式を使うと、\(z_{ia} (i=1, \ldots , b)\)は\(z_{ij} (i=1, \ldots , b),~~(j = 1, \ldots , a-1)\)の和で書ける。これを\(z_{i\circ}^{2} := \sum_{j=1}^{a} z_{ij}^{2}\)に使うと
\[ z_{i\circ}^{2} = \sum_{j, k=1}^{a-1} c_{jk} z_{ij}z_{ik}\]
である。特に\(z_{b\circ}^{2}\)は\(z_{bj}~~(j = 1, \ldots , a-1)\)の二次式の和である。式(*)のまだ使っていない、もう1つの条件式を使うと\(z_{bj}\)は\(z_{ij}~~(i=1, \ldots, b-1)(j=1, \ldots , a-1)\)の和で書けるから、
\[ z_{b\circ}^{2} = \sum_{i,k=1}^{b-1}\sum_{j, l=1}^{a-1} c_{ijkl}z_{ij}z_{kl}\]
と書ける。\(z_{i\circ}^{2}~~(i=1, \ldots , b-1)\)はすでに\(z_{ij}~~(i=1,\ldots , b-1)(j=1, \ldots , a-1)\)の2次式の和であったから\(z^{2}\)も、その\((a-1)(b-1)\)個の変数の二次式になる。ここで2つの添字をまとめて1つの\((a-1)(b-1)\)次元ベクトル\(\boldsymbol{y}\)を導入すると、対称行列\(\mathrm{A}\)を使って、次のようにまとめることができる。
\[ z^{2} = \boldsymbol{y}^{t} \mathrm{A} \boldsymbol{y}\]
対称行列なので、直交行列\(\mathrm{O}\)で対角化できる。\(\mathrm{A}\)の固有値を\(\lambda_{i}\)とすると、\(x_{i} = \frac{1}{\sqrt{\lambda_{i}}}\sum_{j}\mathrm{O}_{ij}y_{j}\)となる変数で\(z^{2}\)を書き直すと
\[ z^{2} = \sum_{i=1}^{(a-1)(b-1)}x_{i}^{2}\]
になる。すると確率密度は\(\prod_{i=1}^{(a-1)(b-1)}e^{-x_{i}^{2}/2}\)に比例することが分かるが、規格化条件より
\[ p(\boldsymbol{x}) = \prod_{i=1}^{(a-1)(b-1)}\frac{1}{\sqrt{2\pi}} e^{-x_{i}^{2}/2}\]
となることが分かる。以上より、実現値\(x_{i}\)に対応する確率変数\(X_{i}\)が\(\mathrm{N}(0, 1)\)に従うことが分かった。\(\mathrm{N}(0, 1)\)に従う、独立変数の二乗和は、その変数の数を自由度とする\(\chi^{2}\)分布に従うので、\(Z^{2}\)は自由度\((a-1)(b-1)\)の\(\chi^{2}\)分布に従う。
である。これは、\(np_{ij} = n_{i\circ}p_{\circ j} = n_{\circ j}p_{i \circ}\)を使うと確率変数が
\[ z_{ij} = \frac{n_{ij} - n_{i\circ}p_{\circ j}}{\sqrt{n_{i\circ}p_{\circ j}}} = \frac{n_{ij} - n_{\circ j}p_{i \circ}}{\sqrt{n_{\circ j}p_{i \circ}}}\]
と書けるので、あとは\(\sum_{i} p_{i\circ} = \sum_{j} p_{\circ j} = 1\)と\(\sum_{i}n_{ij} = n_{\circ j},~~\sum_{j} n_{ij} = n_{i \circ}\)を使うと得られる。
この得られた式を使うと、\(z_{ia} (i=1, \ldots , b)\)は\(z_{ij} (i=1, \ldots , b),~~(j = 1, \ldots , a-1)\)の和で書ける。これを\(z_{i\circ}^{2} := \sum_{j=1}^{a} z_{ij}^{2}\)に使うと
\[ z_{i\circ}^{2} = \sum_{j, k=1}^{a-1} c_{jk} z_{ij}z_{ik}\]
である。特に\(z_{b\circ}^{2}\)は\(z_{bj}~~(j = 1, \ldots , a-1)\)の二次式の和である。式(*)のまだ使っていない、もう1つの条件式を使うと\(z_{bj}\)は\(z_{ij}~~(i=1, \ldots, b-1)(j=1, \ldots , a-1)\)の和で書けるから、
\[ z_{b\circ}^{2} = \sum_{i,k=1}^{b-1}\sum_{j, l=1}^{a-1} c_{ijkl}z_{ij}z_{kl}\]
と書ける。\(z_{i\circ}^{2}~~(i=1, \ldots , b-1)\)はすでに\(z_{ij}~~(i=1,\ldots , b-1)(j=1, \ldots , a-1)\)の2次式の和であったから\(z^{2}\)も、その\((a-1)(b-1)\)個の変数の二次式になる。ここで2つの添字をまとめて1つの\((a-1)(b-1)\)次元ベクトル\(\boldsymbol{y}\)を導入すると、対称行列\(\mathrm{A}\)を使って、次のようにまとめることができる。
\[ z^{2} = \boldsymbol{y}^{t} \mathrm{A} \boldsymbol{y}\]
対称行列なので、直交行列\(\mathrm{O}\)で対角化できる。\(\mathrm{A}\)の固有値を\(\lambda_{i}\)とすると、\(x_{i} = \frac{1}{\sqrt{\lambda_{i}}}\sum_{j}\mathrm{O}_{ij}y_{j}\)となる変数で\(z^{2}\)を書き直すと
\[ z^{2} = \sum_{i=1}^{(a-1)(b-1)}x_{i}^{2}\]
になる。すると確率密度は\(\prod_{i=1}^{(a-1)(b-1)}e^{-x_{i}^{2}/2}\)に比例することが分かるが、規格化条件より
\[ p(\boldsymbol{x}) = \prod_{i=1}^{(a-1)(b-1)}\frac{1}{\sqrt{2\pi}} e^{-x_{i}^{2}/2}\]
となることが分かる。以上より、実現値\(x_{i}\)に対応する確率変数\(X_{i}\)が\(\mathrm{N}(0, 1)\)に従うことが分かった。\(\mathrm{N}(0, 1)\)に従う、独立変数の二乗和は、その変数の数を自由度とする\(\chi^{2}\)分布に従うので、\(Z^{2}\)は自由度\((a-1)(b-1)\)の\(\chi^{2}\)分布に従う。
2013年5月12日日曜日
2013年5月11日土曜日
4次元正規分布の条件付き確率が1次元正規分布になる話
問題
2つの特徴、例えば、背の高さと足の速さのようなものを考える。これらはそれぞれ、{高い, 低い}と{速い、 遅い}を取るとする。小学4年生男子を、この特徴で分類すると4つに別れる。これらの数をそれぞれ\(n_{11}, n_{12}, n_{21}, n_{22}\)と書くことにする。ここで左のインデックスを背の高さ、右のインデックスを足の速さとし、
\[
\begin{align*}
&(\text{高い}, \text{低い}) \to (1, 2)\\
&(\text{速い}, \text{遅い}) \to (1, 2)
\end{align*}
\]
という対応をとるとする。このようにすると、特徴の詳細に全く依らない書き方になる。
\(n_{ij}\)の合計は\(n\)に固定しておく。これらの分類がどれだけの割合になるかは分からないとする。
この2つの特徴は独立であるとすると、特徴1が\(i\)、特徴2が\(j\)に分類される確率\(p_{ij}\)は次のように2つの確率の積で書ける。
\[ p_{ij} = p_{i\circ} p_{\circ j}\]
最尤推定から、これらの確率を推定すると結果は次のようになる。
\[ \begin{align*} &p_{i\circ} = \frac{n_{i \circ }}{n},~~~~p_{\circ j} = \frac{n_{\circ j}}{n}\\ &n_{i \circ} := n_{i 1} + n_{i 2}\\ &n_{\circ j} := n_{1 j} + n_{2 j} \end{align*} \]
このときの確率分布を求めよ。
最尤推定によって加わった条件式を除けば、4項分布の場合にあたる。\(n\)は大きいとして、正規分布に近似する。\(z_{ij} = \frac{n_{ij} - n p_{ij}}{\sqrt{np_{ij}}}\)を確率変数とすると、確率密度関数は、次のように書ける。
\[ p(z_{11}, z_{21}, z_{12}, z_{22}) = \prod_{i=1, j=1}^{2} \frac{1}{\sqrt{2\pi}} e^{-z_{ij}^{2}/2}\]
ここで\(p_{ij}\)の定義に注意すると\(z_{ij}\)は、次のように書ける。
\[ \begin{align*} &z_{11} = \frac{n_{11} - n_{1 \circ} p_{\circ 1}}{\sqrt{n_{1\circ} p_{\circ 1}}}\\ &z_{12} = \frac{n_{12} - n_{1 \circ} p_{\circ 2}}{\sqrt{n_{1\circ} p_{\circ 2}}}\\ &z_{21} = \frac{n_{21} - n_{2 \circ} p_{\circ 1}}{\sqrt{n_{2\circ} p_{\circ 1}}}\\ &z_{22} = \frac{n_{22} - n_{2 \circ} p_{\circ 2}}{\sqrt{n_{2\circ} p_{\circ 2}}} \end{align*} \]
上2つと下2つのそれぞれのペアは2次元正規分布の条件付き確率を求めたときの計算を\(n \to n_{\circ i},~~ p_{i} \to p_{i \circ}\)に置き換えたものに一致するから、以前\(\delta\)で\(z_{2}\)を書きなおしたように\(z_{i 2}\)を\(\delta_{i}\)で書き直せる。
\[p(z_{11}, z_{21}, \delta_{1}, \delta_{2}) = \left(\frac{1}{\sqrt{2\pi}}\right)^{4}\prod_{i=1}^{2} \exp \left[-\frac{1}{2p_{\circ 2}}(z_{i1} - \mu(\delta_{i}))^{2} + \phi(\delta_{i}) \right] \]
更に\(z_{21} = -\sqrt{\frac{p_{1\circ}}{p_{2\circ}}} z_{11} + \delta_{3}\)と置くと
\[p(z_{11}, \boldsymbol{\delta}) = \left(\frac{1}{\sqrt{2\pi}}\right)^{4} \exp \left[-\frac{1}{2p_{2\circ}p_{\circ 2}}(z_{11} - \mu(\boldsymbol{\delta}))^{2} + \phi(\boldsymbol{\delta}) \right]\]
と書ける。条件確率は定義から、次のように計算できる。
\[ \begin{align*} p(Z_{11} = z_{11}, \boldsymbol{\Delta} = 0) &= \frac{p(Z_{11} = z_{11}, \boldsymbol{\Delta} = 0)}{p(\boldsymbol{\Delta} = 0)}\\ &= \frac{1}{\sqrt{2\pi p_{2\circ} p_{\circ 2}}} e^{-\frac{1}{2p_{2\circ} p_{\circ 2}}z_{11}^{2}} \end{align*} \]
この式から、\(Z := Z_{11}/\sqrt{p_{2\circ} p_{\circ 2}}\)は1次元の標準正規分布に従うことが分かる。
\[ p_{ij} = p_{i\circ} p_{\circ j}\]
最尤推定から、これらの確率を推定すると結果は次のようになる。
\[ \begin{align*} &p_{i\circ} = \frac{n_{i \circ }}{n},~~~~p_{\circ j} = \frac{n_{\circ j}}{n}\\ &n_{i \circ} := n_{i 1} + n_{i 2}\\ &n_{\circ j} := n_{1 j} + n_{2 j} \end{align*} \]
このときの確率分布を求めよ。
解答
\(p_{i \circ}, p_{\circ j}, n\)は定数なので、下の2式が条件になる。\(\sum_{i=1}^{2} n_{i \circ} = \sum_{i=1}^{2} n_{\circ i} = n\)だから独立なのは、2つである。今考えている問題は\(\sum_{i, j} n_{i, j} = n\)と合わせれば、条件式は計3つになる。以後の計算では、\(n_{1\circ}, n_{2\circ}, n_{\circ 1}\)が定数になるのが条件だと思って計算する。最尤推定によって加わった条件式を除けば、4項分布の場合にあたる。\(n\)は大きいとして、正規分布に近似する。\(z_{ij} = \frac{n_{ij} - n p_{ij}}{\sqrt{np_{ij}}}\)を確率変数とすると、確率密度関数は、次のように書ける。
\[ p(z_{11}, z_{21}, z_{12}, z_{22}) = \prod_{i=1, j=1}^{2} \frac{1}{\sqrt{2\pi}} e^{-z_{ij}^{2}/2}\]
ここで\(p_{ij}\)の定義に注意すると\(z_{ij}\)は、次のように書ける。
\[ \begin{align*} &z_{11} = \frac{n_{11} - n_{1 \circ} p_{\circ 1}}{\sqrt{n_{1\circ} p_{\circ 1}}}\\ &z_{12} = \frac{n_{12} - n_{1 \circ} p_{\circ 2}}{\sqrt{n_{1\circ} p_{\circ 2}}}\\ &z_{21} = \frac{n_{21} - n_{2 \circ} p_{\circ 1}}{\sqrt{n_{2\circ} p_{\circ 1}}}\\ &z_{22} = \frac{n_{22} - n_{2 \circ} p_{\circ 2}}{\sqrt{n_{2\circ} p_{\circ 2}}} \end{align*} \]
上2つと下2つのそれぞれのペアは2次元正規分布の条件付き確率を求めたときの計算を\(n \to n_{\circ i},~~ p_{i} \to p_{i \circ}\)に置き換えたものに一致するから、以前\(\delta\)で\(z_{2}\)を書きなおしたように\(z_{i 2}\)を\(\delta_{i}\)で書き直せる。
\[p(z_{11}, z_{21}, \delta_{1}, \delta_{2}) = \left(\frac{1}{\sqrt{2\pi}}\right)^{4}\prod_{i=1}^{2} \exp \left[-\frac{1}{2p_{\circ 2}}(z_{i1} - \mu(\delta_{i}))^{2} + \phi(\delta_{i}) \right] \]
更に\(z_{21} = -\sqrt{\frac{p_{1\circ}}{p_{2\circ}}} z_{11} + \delta_{3}\)と置くと
\[p(z_{11}, \boldsymbol{\delta}) = \left(\frac{1}{\sqrt{2\pi}}\right)^{4} \exp \left[-\frac{1}{2p_{2\circ}p_{\circ 2}}(z_{11} - \mu(\boldsymbol{\delta}))^{2} + \phi(\boldsymbol{\delta}) \right]\]
と書ける。条件確率は定義から、次のように計算できる。
\[ \begin{align*} p(Z_{11} = z_{11}, \boldsymbol{\Delta} = 0) &= \frac{p(Z_{11} = z_{11}, \boldsymbol{\Delta} = 0)}{p(\boldsymbol{\Delta} = 0)}\\ &= \frac{1}{\sqrt{2\pi p_{2\circ} p_{\circ 2}}} e^{-\frac{1}{2p_{2\circ} p_{\circ 2}}z_{11}^{2}} \end{align*} \]
この式から、\(Z := Z_{11}/\sqrt{p_{2\circ} p_{\circ 2}}\)は1次元の標準正規分布に従うことが分かる。
関連する日記
2013年5月10日金曜日
金曜日
再び、独立性検定のための計算。昨日、やった計算と全く同じようにやると、答えが1次元正規分布になる場合は期待していた答えが出てきた。ただし、説明するときは、どれが与えられている定数なのかをはっきりさせておかないと混乱のもとになりそう。
では、まず考えていた問題をもう少し正確に書いておこう。
この問題は明日、計算しよう。これがうまくいけば次は、最終的な確率変数が2つの場合をやる。
では、まず考えていた問題をもう少し正確に書いておこう。
問題
2つの特徴、例えば、背の高さと足の速さのようなものを考える。これらはそれぞれ、{高い, 低い}と{速い、 遅い}を取るとする。小学4年生男子を、この特徴で分類すると4つに別れる。これらの数をそれぞれ\(n_{11}, n_{12}, n_{21}, n_{22}\)と書くことにする。ここで左のインデックスを背の高さ、右のインデックスを足の速さとし、
\[ \begin{align*} &(\text{高い}, \text{低い}) \to (1, 2)\\ &(\text{速い}, \text{遅い}) \to (1, 2) \end{align*} \]
という対応をとるとする。このようにすると、特徴の詳細に全く依らない書き方になる。
\(n_{ij}\)の合計は\(n\)に固定しておく。これらの分類がどれだけの割合になるかは分からないとする。
この2つの特徴は独立であるとすると、特徴1が\(i\)、特徴2が\(j\)に分類される確率\(p_{ij}\)は次のように2つの確率の積で書ける。
\[ p_{ij} = p_{i\circ} p_{\circ j}\]
最尤推定から、これらの確率を推定すると結果は次のようになる。
\[ \begin{align*} &p_{i\circ} = \frac{n_{i \circ }}{n},~~~~p_{\circ j} = \frac{n_{\circ j}}{n}\\ &n_{i \circ} := n_{i 1} + n_{i 2}\\ &n_{\circ j} := n_{1 j} + n_{2 j} \end{align*} \]
このときの確率分布を求めよ。
\[ p_{ij} = p_{i\circ} p_{\circ j}\]
最尤推定から、これらの確率を推定すると結果は次のようになる。
\[ \begin{align*} &p_{i\circ} = \frac{n_{i \circ }}{n},~~~~p_{\circ j} = \frac{n_{\circ j}}{n}\\ &n_{i \circ} := n_{i 1} + n_{i 2}\\ &n_{\circ j} := n_{1 j} + n_{2 j} \end{align*} \]
このときの確率分布を求めよ。
この問題は明日、計算しよう。これがうまくいけば次は、最終的な確率変数が2つの場合をやる。
2013年5月9日木曜日
正規分布の条件付き確率。2次元の場合
問題
2次元正規分布の確率密度関数を次のように定義する。
\[
\begin{align*}
&p(z_{1}, z_{2})dz_{1}dz_{2} = \prod_{i=1}\frac{dz_{i}}{\sqrt{2\pi}} e^{-z_{i}^{2}/2}\\
&p_{1} + p_{2} = 1~~(p_{i} \ge 0)
\end{align*}
\]
このとき、\(z_{1} + z_{2} = 0\)としたときの条件付き確率\(p(Z_{1} = z_{1}; \Delta = 0)\)を求めよ。
ただし、\(\Delta = \sqrt{\frac{p_{1}}{p_{2}}}Z_{1} + Z_{2}\)とする。
このとき、\(z_{1} + z_{2} = 0\)としたときの条件付き確率\(p(Z_{1} = z_{1}; \Delta = 0)\)を求めよ。
ただし、\(\Delta = \sqrt{\frac{p_{1}}{p_{2}}}Z_{1} + Z_{2}\)とする。
解答
条件付き確率の定義から
\[ p(Z_{1} = z_{1}; \Delta = 0) = \frac{p(Z_{1} = z_{1}, \Delta = 0)}{p(\Delta = 0)}~~(1)\]
と書ける。まず、\(z_{2}\)を\(- \sqrt{\frac{p_{1}}{p_{2}}}z_{1} + \delta\)で\(z^{2} = z_{1}^{2} + z_{2}^{2}\)を書き直すと、
\[z^{2} = \frac{1}{p_{2}}\left(z_{1} - \sqrt{p_{1}p_{2}}\delta \right)^{2} + p_{2}\delta^{2}\]
となる。これを確率密度関数に代入して、\(dz_{2}\)を\(d\delta\)に書き換えると、確率密度関数を\(z_{1}, \delta\)の関数として書き直すことができる。
\[p(Z_{1} = z_{1}, \Delta = \delta) = \left(\frac{1}{\sqrt{2\pi}}\right)^{2}\exp\left[- \frac{1}{2p_{2}}\left(z_{1} - \sqrt{p_{1}p_{2}}\delta\right)^{2} - \frac{p_{2}}{2} \delta^{2}\right]~~(2)\]
これを\(z_{1}\)について積分して周辺確率が計算できる。
\[ p(\Delta = \delta) = \sqrt{\frac{p_{2}}{2\pi}} e^{-\frac{1}{2p_{2}}\delta^{2}}~~(3)\]
(2),(3)を(1)の右辺に代入すると、条件付き確率が得られる。
\[ p(Z_{1} = z_{1}, \Delta = 0) = \frac{1}{\sqrt{2\pi p_{2}}} e^{-\frac{1}{2p_{2}}z_{1}^{2}}\]
関連する日記
2013年5月8日水曜日
水曜日
条件付き確率でうまくいくはずなんだけれども、なかなかうまくいかない。よくよく考えると、二項分布から1次元の正規分布に近似されるのも、2次元正規分布に条件\(z_{1} + z_{2} = 0\)を付けた場合にあたるなぁと気付いた。これについて、ちょっとまじめに考えてみよう。問題をもう少し、しっかり書くと次のようになる。
問題
2次元正規分布の確率密度関数を次のように定義する。
\[
\begin{align*}
&p(z_{1}, z_{2})dz_{1}dz_{2} = \prod_{i=1}\frac{dz_{i}}{\sqrt{2\pi}} e^{-z_{i}^{2}/2}\\
&p_{1} + p_{2} = 1~~(p_{i} \ge 0)
\end{align*}
\]
このとき、\(z_{1} + z_{2} = 0\)としたときの条件付き確率\(p(Z_{1} = z_{1}; \Delta = 0)\)を求めよ。
ただし、\(\Delta = \sqrt{\frac{p_{1}}{p_{2}}}Z_{1} + Z_{2}\)とする。
問題は書けたので、次回はこの問題を解こう。
このとき、\(z_{1} + z_{2} = 0\)としたときの条件付き確率\(p(Z_{1} = z_{1}; \Delta = 0)\)を求めよ。
ただし、\(\Delta = \sqrt{\frac{p_{1}}{p_{2}}}Z_{1} + Z_{2}\)とする。
問題は書けたので、次回はこの問題を解こう。
2013年5月7日火曜日
2013年5月6日月曜日
2013年5月5日日曜日
日曜日
昨日までで多項分布がサンプルサイズが大きいときに正規分布に近づくことの説明を書きました。適合度検定と独立性検定で、そこで示したことを使うのですが、まだ自分の言葉でまとめることができていません。
適合度検定については、だいたい理解するところまで来ているので、今日は独立性の検定の方を調べていました。本を斜め読みしていたせいか、事前に事象が発生する確率が分かっているときと、そうでないときの2通りあることに気づかなくて、しばらく考えこんでしまいました。
そこに気付いてからは、\(\chi^{2}\)分布の自由度の勘定の辻褄が合うかなぁ程度の理解まで進んだところで、今日はおしまい。
適合度検定については、だいたい理解するところまで来ているので、今日は独立性の検定の方を調べていました。本を斜め読みしていたせいか、事前に事象が発生する確率が分かっているときと、そうでないときの2通りあることに気づかなくて、しばらく考えこんでしまいました。
そこに気付いてからは、\(\chi^{2}\)分布の自由度の勘定の辻褄が合うかなぁ程度の理解まで進んだところで、今日はおしまい。
更新スタイルについて
確かな理解に至るまでは、日記という形でほぼ文章だけで途中経過を記録していこうと思います。理解が進んで、下書きもだいたいできた段階から、記事公開していこうかなと思っています。そうすれば、根本的に考えが間違っていて、後戻りということも少なくなるんじゃないかなぁと期待しています。
2013年5月4日土曜日
多項分布から多変量正規分布への近似について(9)
問題
\[\begin{align*}
&|\mathrm{A}| =
\left| \begin{array}{cccc}
1 + \frac{p_{1}}{p_{s+1}} & \frac{\sqrt{p_{1}p_{2}}}{p_{s+1}} & \ldots & \frac{\sqrt{p_{1}p_{k-1}}}{p_{s+1}} \\
\frac{\sqrt{p_{1}p_{2}}}{p_{s+1}} & 1 + \frac{p_{2}}{p_{s+1}} & \ldots & \frac{\sqrt{p_{2}p_{s}}}{p_{s+1}} \\
\vdots & \vdots & \ddots & \ldots \\
\frac{\sqrt{p_{1}p_{s}}}{p_{s+1}} & \frac{\sqrt{p_{2}p_{s}}}{p_{s+1}} & \ldots & 1 + \frac{p_{s}}{p_{s+1}}\\
\end{array}
\right| = \left|(\mathrm{C}_{s}, 1 + \frac{\sqrt{p_{s}}}{p_{s+1}}\boldsymbol{v}^{\prime}) \right|\\
&\mathrm{C}_{i} = \left(
\begin{array}{c}
\mathrm{I}_{i\times i} \\
\mathrm{O}_{(s-i) \times (i)}
\end{array}
\right)\\
&\boldsymbol{v}^{\prime} = \frac{1}{p_{s+1}}
\left(
\frac{p_{1}}{\sqrt{p_{1}}},
\frac{p_{1} + p_{2}}{\sqrt{p_{2}}},
\ldots ,
\frac{\sum_{l=1}^{i}p_{l}}{\sqrt{p_{i}}},
\ldots ,
\frac{\sum_{l=1}^{s}p_{l}}{\sqrt{p_{s}}}
\right)
\end{align*}
\]
を示せ。これが示されれば、三角行列になっているので、対角成分の積が行列式に等しい。従って
\[|A| = 1 + \frac{\sqrt{p_{s}}}{p_{s+1}}\frac{\sum_{l=1}^{s}p_{l}}{\sqrt{p_{s}}} = \frac{\sum_{l=1}^{s+1}p_{l}}{p_{s+1}}\]
となる。特に\(\sum_{l=1}^{s+1} p_{l} = 1\)のとき\(|A| = 1/p_{s+1}\)になる。
を示せ。これが示されれば、三角行列になっているので、対角成分の積が行列式に等しい。従って
\[|A| = 1 + \frac{\sqrt{p_{s}}}{p_{s+1}}\frac{\sum_{l=1}^{s}p_{l}}{\sqrt{p_{s}}} = \frac{\sum_{l=1}^{s+1}p_{l}}{p_{s+1}}\]
となる。特に\(\sum_{l=1}^{s+1} p_{l} = 1\)のとき\(|A| = 1/p_{s+1}\)になる。
考えたこと
\(s = 3\)のときと同様の操作を\(i\)回繰り返すと、\(|\mathrm{A}|\)は次のように変形できることをまず示す。
\[
\begin{align*}
&|\mathrm{A}| = |(\mathrm{C}_{i}, \boldsymbol{e}_{i+1} + \sqrt{p_{i+1}}\boldsymbol{v}^{(i+1)}, \boldsymbol{e}_{i+2} + \sqrt{p_{i+2}}\boldsymbol{v}^{(i+1)}, \ldots , \boldsymbol{e}_{s} + \sqrt{p_{s}}\boldsymbol{v}^{(i+1)})|\\
&\boldsymbol{v}^{(i+1)} =\boldsymbol{v}^{(i)} + \sqrt{\frac{p_{i}}{p_{i+1}}}(v^{(i)})_{i} \boldsymbol{e}_{i+1}
\end{align*}
\]
ただし\(\boldsymbol{v}^{(1)} = \boldsymbol{v}\)である。 \(i = 1\)のときは、直接計算することで確認できるので、\(i = k\)のとき成り立つと仮定すれば、\(i = k+1\)のときも成り立つことを示すことによって、これを証明する。
まず\(\boldsymbol{v}^{(i)}\)は漸化式から
\[ \left(\boldsymbol{v}^{(i)}\right)^{t} =\frac{1}{p_{s+1}} \left(\frac{p_{1}}{\sqrt{p_{1}}}, \frac{p_{1} + p_{2}}{\sqrt{p_{2}}}, \ldots , \frac{\sum_{l = 1}^{i-1} p_{l}}{\sqrt{p_{i}}}, v_{i}, v_{i+1}, \ldots , v_{s}\right)\]
となることが確認できるので、これを以降の証明で使う。
次の操作で、\(k\)列以前は変化しないので、\(k+1\)列以降の成分がどう変化するかを確認する。
\[ \begin{align*} \text{(k+1)列} &= \boldsymbol{e}_{k+1} + \sqrt{p_{k+1}} \boldsymbol{v}^{(k+1)} - \sqrt{\frac{p_{k+1}}{p_{k+2}}}\left( \boldsymbol{e}_{k+2} + \sqrt{p_{k+2}} \boldsymbol{v}^{(k+1)} \right)\\ &= \boldsymbol{e}_{k+1} - \frac{p_{k+1}}{p_{k+2}}\boldsymbol{e}_{k+2} \end{align*} \]
次に、操作2では(k+1)列以降の全ての列が変化するが、それ以前は\(k+1\)行の成分が全て0なので、変化しない。\(k+1\)行目はk+1行目のk+1成分は1なので、次のようになる。
\[ \begin{align*} (\text{k+1行}) &= \boldsymbol{e}_{k+1} - \sqrt{\frac{p_{k+1}}{p_{k+2}}}\boldsymbol{e}_{k+2} + \frac{p_{k+1}}{p_{k+2}}(\text{k+1行})_{k+1} \boldsymbol{e}_{k+2}\\ &= \boldsymbol{e}_{k+1} \end{align*} \]
\(k + 2 \le i \le s\)列目は次のようになる。
\[ \begin{align*} (\text{i行}) &= \boldsymbol{e}_{i} + \sqrt{p_{i}} \boldsymbol{v^{(k+1)}} + \sqrt{\frac{p_{k+1}}{p_{k+2}}}\sqrt{p_{i}}(v^{(k+1)})_{k+1} \boldsymbol{e}_{k+2}\\ &= \boldsymbol{e}_{i} + \sqrt{p_{i}} \boldsymbol{v^{(k+2)}} \end{align*} \]
最後の等式には\(\boldsymbol{v}^{(i)}\)の漸化式を使った。この結果を使うと、次のようになり、\(i = k+1\)も成り立つことが分かる。
\[ \begin{align*} |\mathrm{A}| &= |(\mathrm{C}_{k}, \boldsymbol{e}_{k+1}, \boldsymbol{e}_{k+2} + \sqrt{p_{k+2}}\boldsymbol{v}^{(k+2)}, \boldsymbol{e}_{k+3} + \sqrt{p_{k+3}}\boldsymbol{v}^{(k+2)}, \ldots , \boldsymbol{e}_{s} + \sqrt{p_{s}}\boldsymbol{v}^{(k+2)})|\\ &= |(\mathrm{C}_{k+1},\boldsymbol{e}_{k+2} + \sqrt{p_{k+2}}\boldsymbol{v}^{(k+2)}, \boldsymbol{e}_{k+3} + \sqrt{p_{k+3}}\boldsymbol{v}^{(k+2)}, \ldots , \boldsymbol{e}_{s} + \sqrt{p_{s}}\boldsymbol{v}^{(k+2)})| \end{align*} \]
\(i\)回操作を繰り返した後の式が得られたので、これを\(s\)回繰り返すと、欲しかった式が得られる。これにより、一般に\(|A| = 1/p_{s}\)だと分かった。
ただし\(\boldsymbol{v}^{(1)} = \boldsymbol{v}\)である。 \(i = 1\)のときは、直接計算することで確認できるので、\(i = k\)のとき成り立つと仮定すれば、\(i = k+1\)のときも成り立つことを示すことによって、これを証明する。
まず\(\boldsymbol{v}^{(i)}\)は漸化式から
\[ \left(\boldsymbol{v}^{(i)}\right)^{t} =\frac{1}{p_{s+1}} \left(\frac{p_{1}}{\sqrt{p_{1}}}, \frac{p_{1} + p_{2}}{\sqrt{p_{2}}}, \ldots , \frac{\sum_{l = 1}^{i-1} p_{l}}{\sqrt{p_{i}}}, v_{i}, v_{i+1}, \ldots , v_{s}\right)\]
となることが確認できるので、これを以降の証明で使う。
次の操作で、\(k\)列以前は変化しないので、\(k+1\)列以降の成分がどう変化するかを確認する。
- \(\sqrt{\frac{p_{k+1}}{p_{k+2}}} \times \text{k+2列}\)をk+1列から引く
- \(\sqrt{\frac{p_{k+1}}{p_{k+2}}} \times \text{k+1行}\)をk+2行に足す
\[ \begin{align*} \text{(k+1)列} &= \boldsymbol{e}_{k+1} + \sqrt{p_{k+1}} \boldsymbol{v}^{(k+1)} - \sqrt{\frac{p_{k+1}}{p_{k+2}}}\left( \boldsymbol{e}_{k+2} + \sqrt{p_{k+2}} \boldsymbol{v}^{(k+1)} \right)\\ &= \boldsymbol{e}_{k+1} - \frac{p_{k+1}}{p_{k+2}}\boldsymbol{e}_{k+2} \end{align*} \]
次に、操作2では(k+1)列以降の全ての列が変化するが、それ以前は\(k+1\)行の成分が全て0なので、変化しない。\(k+1\)行目はk+1行目のk+1成分は1なので、次のようになる。
\[ \begin{align*} (\text{k+1行}) &= \boldsymbol{e}_{k+1} - \sqrt{\frac{p_{k+1}}{p_{k+2}}}\boldsymbol{e}_{k+2} + \frac{p_{k+1}}{p_{k+2}}(\text{k+1行})_{k+1} \boldsymbol{e}_{k+2}\\ &= \boldsymbol{e}_{k+1} \end{align*} \]
\(k + 2 \le i \le s\)列目は次のようになる。
\[ \begin{align*} (\text{i行}) &= \boldsymbol{e}_{i} + \sqrt{p_{i}} \boldsymbol{v^{(k+1)}} + \sqrt{\frac{p_{k+1}}{p_{k+2}}}\sqrt{p_{i}}(v^{(k+1)})_{k+1} \boldsymbol{e}_{k+2}\\ &= \boldsymbol{e}_{i} + \sqrt{p_{i}} \boldsymbol{v^{(k+2)}} \end{align*} \]
最後の等式には\(\boldsymbol{v}^{(i)}\)の漸化式を使った。この結果を使うと、次のようになり、\(i = k+1\)も成り立つことが分かる。
\[ \begin{align*} |\mathrm{A}| &= |(\mathrm{C}_{k}, \boldsymbol{e}_{k+1}, \boldsymbol{e}_{k+2} + \sqrt{p_{k+2}}\boldsymbol{v}^{(k+2)}, \boldsymbol{e}_{k+3} + \sqrt{p_{k+3}}\boldsymbol{v}^{(k+2)}, \ldots , \boldsymbol{e}_{s} + \sqrt{p_{s}}\boldsymbol{v}^{(k+2)})|\\ &= |(\mathrm{C}_{k+1},\boldsymbol{e}_{k+2} + \sqrt{p_{k+2}}\boldsymbol{v}^{(k+2)}, \boldsymbol{e}_{k+3} + \sqrt{p_{k+3}}\boldsymbol{v}^{(k+2)}, \ldots , \boldsymbol{e}_{s} + \sqrt{p_{s}}\boldsymbol{v}^{(k+2)})| \end{align*} \]
\(i\)回操作を繰り返した後の式が得られたので、これを\(s\)回繰り返すと、欲しかった式が得られる。これにより、一般に\(|A| = 1/p_{s}\)だと分かった。
次回
今回の話題は、これで一旦終了です。
2013年5月3日金曜日
多項分布から多変量正規分布への近似について(8)
問題
\[
|\mathrm{A}| = \left| \begin{array}{cccc}
1 + \frac{p_{1}}{p_{5}} & \frac{\sqrt{p_{1}p_{2}}}{p_{5}} & \frac{\sqrt{p_{1}p_{3}}}{p_{5}} & \frac{\sqrt{p_{1}p_{4}}}{p_{5}} \\
\frac{\sqrt{p_{1}p_{2}}}{p_{5}} & 1 + \frac{p_{2}}{p_{5}} & \frac{\sqrt{p_{2}p_{3}}}{p_{5}} & \frac{\sqrt{p_{2}p_{4}}}{p_{5}} \\
\frac{\sqrt{p_{1}p_{3}}}{p_{5}} & \frac{\sqrt{p_{2}p_{3}}}{p_{5}} & 1 + \frac{p_{3}}{p_{5}} & \frac{\sqrt{p_{3}p_{4}}}{p_{5}} \\
\frac{\sqrt{p_{1}p_{4}}}{p_{5}} & \frac{\sqrt{p_{2}p_{4}}}{p_{5}} & \frac{\sqrt{p_{1}p_{3}}}{p_{5}} & 1 + \frac{p_{4}}{p_{5}} \\
\end{array}
\right| = \frac{1}{p_{5}}
\]
になることを示せ。
になることを示せ。
考えたこと
大きな行列式を機械的に計算する1つの方法に、行列式の値を変えない操作を繰り返して、三角行列を作るという方法がある。このときの行列式の値は、三角行列の対角成分の積になる。
今回の場合は2つの操作を合わせて1回と数えたとき、\(i\)回目の操作を、次のように定義する。
- \(i+1\)列の\(\sqrt{\frac{p_{i}}{p_{i+1}}}\)倍を\(i\)列から引く
- \(i\)行の\(\sqrt{\frac{p_{i}}{p_{i+1}}}\)倍を\(i+1\)行に足す
\[
\begin{align*}
|\mathrm{A}|&=
\left|
\begin{array}{cccc}
1 + \frac{\sqrt{p_{1}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right) & \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right) & \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right) & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right)\\
\frac{\sqrt{p_{1}}}{p_{5}}\left(\frac{p_{2}}{\sqrt{p_{2}}}\right) & 1 + \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{2}}{\sqrt{p_{2}}}\right) & \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{2}}{\sqrt{p_{2}}}\right) & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{2}}{\sqrt{p_{2}}}\right)\\
\frac{\sqrt{p_{1}}}{p_{5}}\left(\frac{p_{3}}{\sqrt{p_{3}}}\right) & \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{3}}{\sqrt{p_{3}}}\right) &
1 + \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{3}}{\sqrt{p_{3}}}\right) & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{3}}{\sqrt{p_{3}}}\right)\\
\frac{\sqrt{p_{1}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right) & \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right) & \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right) & 1 + \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right)\\
\end{array}
\right|\\
&=
\left|
\begin{array}{cccc}
1 & \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right) & \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right) & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right)\\
-\sqrt{\frac{p_{1}}{p_{2}}} & 1 + \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{2}}{\sqrt{p_{2}}}\right) & \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{2}}{\sqrt{p_{2}}}\right) & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{2}}{\sqrt{p_{2}}}\right)\\
0 & \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{3}}{\sqrt{p_{3}}}\right) &
1 + \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{3}}{\sqrt{p_{3}}}\right) & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{3}}{\sqrt{p_{3}}}\right)\\
0 & \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right) & \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right) & 1 + \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right)\\
\end{array}
\right|\\
&=
\left|
\begin{array}{cccc}
1 & \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right) &
\frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right) & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right)\\
0 &
1 + \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{\sum_{l=1}^{2}p_{l}}{\sqrt{p_{2}}}\right) &
\frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{\sum_{l=1}^{2}p_{l}}{\sqrt{p_{2}}}\right) &
\frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{\sum_{l=1}^{2}p_{l}}{\sqrt{p_{2}}}\right)\\
0 & \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{3}}{\sqrt{p_{3}}}\right) &
1 + \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{3}}{\sqrt{p_{3}}}\right) & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{3}}{\sqrt{p_{3}}}\right)\\
0 & \frac{\sqrt{p_{2}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right) & \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right) & 1 + \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right)\\
\end{array}
\right|\\
&=
\left|
\begin{array}{cccc}
1 & 0 & \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right) & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right)\\
0 & 1 &
\frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{\sum_{l=1}^{2}p_{l}}{\sqrt{p_{2}}}\right) &
\frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{\sum_{l=1}^{2}p_{l}}{\sqrt{p_{2}}}\right)\\
0 & 0 &
1 + \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{\sum_{l=1}^{3}p_{l}}{\sqrt{p_{3}}}\right) &
\frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{\sum_{l=1}^{3}p_{l}}{\sqrt{p_{3}}}\right)\\
0 & 0 & \frac{\sqrt{p_{3}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right) & 1 + \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{4}}{\sqrt{p_{4}}}\right)\\
\end{array}
\right|\\
&=
\left|
\begin{array}{cccc}
1 & 0 & 0 & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{p_{1}}{\sqrt{p_{1}}}\right)\\
0 & 1 & 0 & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{\sum_{l=1}^{2}p_{l}}{\sqrt{p_{2}}}\right)\\
0 & 0 & 1 & \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{\sum_{l=1}^{3}p_{l}}{\sqrt{p_{3}}}\right)\\
0 & 0 & 0 & 1 + \frac{\sqrt{p_{4}}}{p_{5}}\left(\frac{\sum_{l=1}^{4}p_{l}}{\sqrt{p_{4}}}\right)\\
\end{array}
\right|
\end{align*}
\]
三角行列の対角成分の積から行列式が\(\sum_{l=1}^{5}p_{l}/p_{5}\)だと分かる。更に\(\sum_{i=1}^{5} p_{i} = 1\)を使うと\(1/p_{5}\)が得られる。これが示したいことであった。
次回
次回、最終回は今回の計算を参考に一般の場合の計算を扱います。三角行列を得る手続きは、次元に依らず変わらないことから一般化できます。
2013年5月2日木曜日
多項分布から多変量正規分布への近似について(7)
問題
\[
\left| \begin{array}{cccc}
1 + \frac{p_{1}}{p_{s+1}} & \frac{\sqrt{p_{1}p_{2}}}{p_{s+1}} & \ldots & \frac{\sqrt{p_{1}p_{k-1}}}{p_{s+1}} \\
\frac{\sqrt{p_{1}p_{2}}}{p_{s+1}} & 1 + \frac{p_{2}}{p_{s+1}} & \ldots & \frac{\sqrt{p_{2}p_{s}}}{p_{s+1}} \\
\vdots & \vdots & \ddots & \ldots \\
\frac{\sqrt{p_{1}p_{s}}}{p_{s+1}} & \frac{\sqrt{p_{2}p_{s}}}{p_{s+1}} & \ldots & 1 + \frac{p_{s}}{p_{s+1}}\\
\end{array}
\right| = \frac{1}{p_{s+1}} \]
を\(\sum_{l=1}^{s+1}p_{l} = 1\)として示せ。
を\(\sum_{l=1}^{s+1}p_{l} = 1\)として示せ。
考えたこと
まず\(s=2\)のときから考えるのが一番単純になる。2次の行列式の公式を使うと\(|\mathrm{A}| = \frac{p_{1} + p_{2} + p_{3}}{p_{3}}\)になるが、\( p_{1} + p_{2} + p_{3} = 1\)であるから\(|\mathrm{A}| = 1/p_{s+1}\)が成り立つ。
ここで、行列式は固有値の積であることを使うと、\(\mathrm{A}\)の固有値は\(1,~~1/p_{s+1}\)だと分かる。実は\(\mathrm{A} = 1 + \mathrm{B}\)であるから、\(|\mathrm{B}| = 0\)である(そうでなければ固有値$1$を取ることと矛盾する)。このことは、次のように書けば明らかである。
\[
\begin{align*}
&\mathrm{B} = (\sqrt{p_{1}}\boldsymbol{v}, \sqrt{p_{2}} \boldsymbol{v})\\
&\boldsymbol{v}^{t} = \frac{1}{p_{3}}(\sqrt{p_{1}}, \sqrt{p_{2}})
\end{align*}
\]
このことは、一般に成り立つ。
\[ \begin{align*} &\mathrm{B} = (\sqrt{p_{1}}\boldsymbol{v}, \sqrt{p_{2}} \boldsymbol{v}, \ldots , \sqrt{p_{s}}\boldsymbol{v})\\ &\boldsymbol{v}^{t} = \frac{1}{p_{s+1}}(\sqrt{p_{1}}, \sqrt{p_{2}}, \ldots , \sqrt{p_{s}}) \end{align*} \]
このような特殊な形をしていることに気づいた。このことを利用して、一般の場合の行列式が計算できないだろうか。
\[ \begin{align*} &\mathrm{B} = (\sqrt{p_{1}}\boldsymbol{v}, \sqrt{p_{2}} \boldsymbol{v}, \ldots , \sqrt{p_{s}}\boldsymbol{v})\\ &\boldsymbol{v}^{t} = \frac{1}{p_{s+1}}(\sqrt{p_{1}}, \sqrt{p_{2}}, \ldots , \sqrt{p_{s}}) \end{align*} \]
このような特殊な形をしていることに気づいた。このことを利用して、一般の場合の行列式が計算できないだろうか。
次回
2次、3次の正方行列の場合は、公式を使えばとりあえず\(1/p_{s+1}\)が成り立つことは示せます。次回は4次の場合の計算を扱って、一般的な場合の計算方法を模索します。
2013年5月1日水曜日
多項分布から多変量正規分布への近似について(6)
問題
\(k\)項分布は、サンプルサイズが大きいときに\(k-1\)次元正規分布に近似できる。
考えたこと
\(k = 2, 3\)のときと同様に、スターリングの公式を使うと次のように書ける。
\[ W \simeq \prod_{l=1}^{k-1} \frac{1}{\sqrt{2\pi n p_{l}}} \frac{1}{\sqrt{p_{k}}}\left[\prod_{l=1}^{k}\left(1 + \frac{m_{l}}{np_{l}}\right)\right]^{-n}\prod_{l=1}^{k}\left(1 + \frac{z_{l}^{2}}{m_{l}}\right)^{-m_{l}}~~(1)\]
\([\ldots]\)の部分は、
\[ \left(1 + \frac{m_{l}}{np_{l}}\right)^{p_{l}} \simeq 1 + \frac{m_{l}}{n} + \frac{m_{l}^{2}}{2n^{2}} - \frac{m_{l}^{2}}{2p_{l}n^{2}}~~(*)\]
なので、\(s_{k} = \sum_{i=1}^{k} m_{i}\)とすると一般に
\[ [\ldots] \simeq 1 + \frac{s_{k}}{n} + \frac{s_{k}^{2}}{2n^{2}} - \sum_{l=1}^{k}\frac{m_{l}^{2}}{2n^{2}p_{l}}~~(2)\]
になる。これを数学的帰納法を使って、証明できる。\(k=1\)のときは、\((*)\)から成り立つ。\(k = t\)のときに成り立つと仮定して、\(k = t + 1\)も成り立つことを示す。
\[
\begin{align*}
[\ldots] &\simeq \left(1 + \frac{s_{t}}{n} + \frac{s_{t}^{2}}{2n^{2}} - \sum_{l=1}^{t}\frac{m_{l}^{2}}{2n^{2}p_{l}}\right)
\left(1 + \frac{m_{t+1}}{n} + \frac{m_{t+1}^{2}}{2n^{2}} - \frac{m_{t+1}^{2}}{2n^{2}p_{t+1}} \right)\\
&\simeq 1 + \frac{s_{t} + m_{t+1}}{n} + \frac{s_{k}^{2} + m_{k+1}^{2} + 2s_{t}m_{t+1}}{2n^{2}} - \sum_{l=1}^{t+1}\frac{m_{l}^{2}}{2n^{2}p_{l}}\\
&= 1 + \frac{s_{k+1}}{n} + \frac{s_{k+1}^{2}}{2n^{2}} - \sum_{l=1}^{t+1}\frac{m_{l}^{2}}{2n^{2}p_{l}}
\end{align*}
\]
以上より\((2)\)が一般に成り立つことが示せた。\(s_{k} = 0\)と\(z_{l}^{2} = m_{l}^{2}/(np_{l})\)を使うと
\[ W \simeq \frac{1}{\sqrt{p_{k}}}\prod_{l=1}^{k-1} \frac{1}{\sqrt{2\pi np_{l}}} e^{-x_{k}^{2}/2}\]
が得られる。\(x_{i}\)の定義は三項分布のときと同じである。この式から確率密度関数を求めると、次のようになる。
\[ p(x_{1}, x_{2}, \ldots x_{k-1}) = \frac{1}{\sqrt{\det [\mathrm{A}] p_{k}}}\prod_{l=1}^{k-1} \frac{1}{\sqrt{2\pi}} e^{-x_{i}^{2}/2}\]
式中の行列\(\mathrm{A}\)は以前、三項分布のときに定義したように
\[
\left(
\begin{array}{cccc}
1 + \frac{p_{1}}{p_{k}} & \frac{\sqrt{p_{1}p_{2}}}{p_{k}} & \ldots & \frac{\sqrt{p_{1}p_{k-1}}}{p_{k}} \\
\frac{\sqrt{p_{1}p_{2}}}{p_{k}} & 1 + \frac{p_{2}}{p_{k}} & \ldots & \frac{\sqrt{p_{2}p_{k-1}}}{p_{k}} \\
\vdots & \vdots & \ddots & \ldots \\
\frac{\sqrt{p_{1}p_{k-1}}}{p_{k}} & \frac{\sqrt{p_{2}p_{k-1}}}{p_{k}} & \ldots & 1 + \frac{p_{k-1}}{p_{k}}\\
\end{array}
\right)
\]
である。確率が規格化されているためには、\(\det[\mathrm{A}] = 1/p_{k}\)になるべきだが、実際にそうなっているだろうか。これが次の疑問である。
以上より\((2)\)が一般に成り立つことが示せた。\(s_{k} = 0\)と\(z_{l}^{2} = m_{l}^{2}/(np_{l})\)を使うと
\[ W \simeq \frac{1}{\sqrt{p_{k}}}\prod_{l=1}^{k-1} \frac{1}{\sqrt{2\pi np_{l}}} e^{-x_{k}^{2}/2}\]
が得られる。\(x_{i}\)の定義は三項分布のときと同じである。この式から確率密度関数を求めると、次のようになる。
\[ p(x_{1}, x_{2}, \ldots x_{k-1}) = \frac{1}{\sqrt{\det [\mathrm{A}] p_{k}}}\prod_{l=1}^{k-1} \frac{1}{\sqrt{2\pi}} e^{-x_{i}^{2}/2}\]
式中の行列\(\mathrm{A}\)は以前、三項分布のときに定義したように
\[
\left(
\begin{array}{cccc}
1 + \frac{p_{1}}{p_{k}} & \frac{\sqrt{p_{1}p_{2}}}{p_{k}} & \ldots & \frac{\sqrt{p_{1}p_{k-1}}}{p_{k}} \\
\frac{\sqrt{p_{1}p_{2}}}{p_{k}} & 1 + \frac{p_{2}}{p_{k}} & \ldots & \frac{\sqrt{p_{2}p_{k-1}}}{p_{k}} \\
\vdots & \vdots & \ddots & \ldots \\
\frac{\sqrt{p_{1}p_{k-1}}}{p_{k}} & \frac{\sqrt{p_{2}p_{k-1}}}{p_{k}} & \ldots & 1 + \frac{p_{k-1}}{p_{k}}\\
\end{array}
\right)
\]
である。確率が規格化されているためには、\(\det[\mathrm{A}] = 1/p_{k}\)になるべきだが、実際にそうなっているだろうか。これが次の疑問である。
次回
行列\(\mathrm{A}\)の行列式を計算しようと思います。
2013年4月30日火曜日
多項分布から多変量正規分布への近似について(5)
問題
三項分布がサンプルサイズ\(n\)が大きいとき、2次元正規分布に近似できることを示せ
考えたこと
\[ z_{3} = - \sqrt{\frac{p_{1}}{p_{3}}} z_{1}- \sqrt{\frac{p_{2}}{p_{3}}} z_{2}\]
を使って\(z^{2}\)を\(z_{1}, z_{2}\)だけで書き直すと次のようになる。
\[ \begin{align*} &z^{2} = \boldsymbol{z}_{2}^{t}\mathrm{A}\boldsymbol{z}_{2}\\ &\mathrm{A} = \left( \begin{array}{cc} 1 + \frac{p_{1}}{p_{3}} & \frac{\sqrt{p_{1}p_{2}}}{p_{3}}\\ \frac{\sqrt{p_{1}p_{2}}}{p_{3}} & 1 + \frac{p_{2}}{p_{3}} \end{array} \right)\\ &\boldsymbol{z}_{2}^{t} = (z_{1}, z_{2}) \end{align*} \]
\(\mathrm{A}\)は対称行列なので、必ず対角化できる直交行列\(\mathrm{O}\)が存在する。対角行列を\(\mathrm{D}\)とすると
\[\mathrm{A} = \mathrm{O}\mathrm{D}\mathrm{O}^{t}\]
と書ける。\(\boldsymbol{x}_{2} = \mathrm{O}\boldsymbol{z}_{2}\)とすると
\[z^{2} = \sum_{i=1}^{2}\lambda_{i} (x_{2})_{i}\]
\[ \begin{align*} &z^{2} = \boldsymbol{z}_{2}^{t}\mathrm{A}\boldsymbol{z}_{2}\\ &\mathrm{A} = \left( \begin{array}{cc} 1 + \frac{p_{1}}{p_{3}} & \frac{\sqrt{p_{1}p_{2}}}{p_{3}}\\ \frac{\sqrt{p_{1}p_{2}}}{p_{3}} & 1 + \frac{p_{2}}{p_{3}} \end{array} \right)\\ &\boldsymbol{z}_{2}^{t} = (z_{1}, z_{2}) \end{align*} \]
\(\mathrm{A}\)は対称行列なので、必ず対角化できる直交行列\(\mathrm{O}\)が存在する。対角行列を\(\mathrm{D}\)とすると
\[\mathrm{A} = \mathrm{O}\mathrm{D}\mathrm{O}^{t}\]
と書ける。\(\boldsymbol{x}_{2} = \mathrm{O}\boldsymbol{z}_{2}\)とすると
\[z^{2} = \sum_{i=1}^{2}\lambda_{i} (x_{2})_{i}\]
\[ p(z_{1}, z_{2})\Delta z_{1} \Delta z_{2} = \frac{1}{2\pi\sqrt{p_{3}}}e^{-z^{2}/2}\Delta z_{1} \Delta z_{2}\]
になるが、ここから\(\Delta (x_{2})_{1} \Delta (x_{2})_{2}\)に書きなおそうとするとヤコビアンがかかるけれども、これは\(|\det \mathrm{O}|\)なので因子として寄与しない。最後に\(\Delta x_{1} \Delta x_{2}\)に書き換えると、\(\prod_{i=1}^{2}\lambda_{i}\)が出てくるが、これは\(\det \mathrm{A}\)に等しい。まとめると、次のようになる。
\[ \begin{align*} p(x_{1}, x_{2}) \Delta (x_{1})_{1} \Delta (x_{2})_{2} &= \frac{1}{2\pi}\frac{1}{\sqrt{\det [\mathrm{A}]p_{3}}}\exp \left(-\sum_{i=1}^{2}x_{i}^{2}/2\right)\Delta x_{1} \Delta x_{2}\\ &= \prod_{i=1}^{2} \frac{1}{\sqrt{2\pi}}e^{-x_{i}^{2}/2}\Delta x_{i} \end{align*} \]
これは、2次元正規分布である。
次回
ここまでで、「三項定理が2次元正規分布に近似できる」ことが分かりました。次回は一旦、\(k\)項分布の場合に挑戦してみます。そこでまた、分からないところが出れば、立ち返って考えてみます。
2013年4月29日月曜日
多項分布から多変量正規分布への近似について(4)
問題
三項分布がサンプルサイズ\(n\)が大きいとき、2次元正規分布に近似できることを示せ
考えたこと
問題に答えるには、\(n\)が大きいときの確率関数\(W\)の式を求めて、
\[ W(n_{1}, n_{2}) \Delta n_{1} \Delta n_{2} \simeq p(z_{1}, z_{2}) \Delta z_{1} \Delta z_{2}\]
を使って、確率密度関数を求めれば良い。この式は
\[
\begin{align*}
&\mathrm{Prob}[z_{1} \le Z_{1} \le z_{1} + \Delta z_{1}, z_{2} \le Z_{2} \le z_{2} + \Delta z_{2}] \\
&= \mathrm{Prob}[n_{1} \le N_{1} \le n_{1} + \Delta n_{1}, n_{2} \le N_{2} \le n_{2} + \Delta n_{2}]
\end{align*}
\]
として、1次元の時と同様に考えると得られる。
ここでは、\(W\)の近似式を得るところまでをやる。確率関数は
\[
\begin{align*}
&W = \frac{n!}{n_{1}! n_{2} !n_{3}!} p_{1}^{n_{1}} p_{2}^{n_{3}} p_{3}^{n_{3}}\\
&n_{1} + n_{2} + n_{3} = n\\
&p_{1} + p_{2} + p_{3} = 1
\end{align*}
\]
である。二項分布と同様に変数
\[z_{i} = \frac{n_{i} - np_{i}}{\sqrt{np_{i}}} = \frac{m_{i}}{\sqrt{np_{i}}}\]
を導入して、スターリングの公式を使って整理すると、次のようになる。
\[ \begin{align*} W &\simeq \frac{1}{2\pi}\sqrt{\frac{n}{n_{1}n_{2} n_{3}}} \prod_{i=1}^{3} \left(\frac{n_{i}}{np_{i}}\right)^{-n_{i}}\\ &= \frac{1}{2\pi}\sqrt{\frac{n}{n_{1}n_{2} n_{3}}} \left[\prod_{i=1}^{3} \left(1 + \frac{m_{i}}{np_{i}}\right)^{p_{i}}\right]^{-n} \prod_{i=1}^{3} \left(1 + \frac{z_{i}^{2}}{m_{i}}\right)^{-m_{i}} \end{align*} \]
\([\ldots ]\)の部分は、\(m_{i}/(np_{i})\)でテイラー展開して\(\sum_{i=1}^{3} m_{i} = 0\)を使うと
\[ [\ldots ] \simeq 1 + \frac{1}{n} \frac{z^{2}}{2}\]
なので、これを\(W\)に代入して\(\sqrt{\frac{n}{n_{1}n_{2} n_{3}}} \simeq \frac{1}{\sqrt{np_{1} np_{2} p_{3}}}\)を使うと、次のようにまとめられる。
\[ \begin{align*} &W \simeq \prod_{i=1}^{2}\frac{1}{\sqrt{2\pi np_{i}}} \frac{1}{\sqrt{p_{3}}} e^{-z^{2}/2}\\ &z^{2} = z_{1}^{2} + z_{2}^{2} + z_{3}^{2} \end{align*} \]
\[z_{i} = \frac{n_{i} - np_{i}}{\sqrt{np_{i}}} = \frac{m_{i}}{\sqrt{np_{i}}}\]
を導入して、スターリングの公式を使って整理すると、次のようになる。
\[ \begin{align*} W &\simeq \frac{1}{2\pi}\sqrt{\frac{n}{n_{1}n_{2} n_{3}}} \prod_{i=1}^{3} \left(\frac{n_{i}}{np_{i}}\right)^{-n_{i}}\\ &= \frac{1}{2\pi}\sqrt{\frac{n}{n_{1}n_{2} n_{3}}} \left[\prod_{i=1}^{3} \left(1 + \frac{m_{i}}{np_{i}}\right)^{p_{i}}\right]^{-n} \prod_{i=1}^{3} \left(1 + \frac{z_{i}^{2}}{m_{i}}\right)^{-m_{i}} \end{align*} \]
\([\ldots ]\)の部分は、\(m_{i}/(np_{i})\)でテイラー展開して\(\sum_{i=1}^{3} m_{i} = 0\)を使うと
\[ [\ldots ] \simeq 1 + \frac{1}{n} \frac{z^{2}}{2}\]
なので、これを\(W\)に代入して\(\sqrt{\frac{n}{n_{1}n_{2} n_{3}}} \simeq \frac{1}{\sqrt{np_{1} np_{2} p_{3}}}\)を使うと、次のようにまとめられる。
\[ \begin{align*} &W \simeq \prod_{i=1}^{2}\frac{1}{\sqrt{2\pi np_{i}}} \frac{1}{\sqrt{p_{3}}} e^{-z^{2}/2}\\ &z^{2} = z_{1}^{2} + z_{2}^{2} + z_{3}^{2} \end{align*} \]
次回
得られた\(W\)は一見すると3次元正規分布のようだが、実際は2次元である。これはよく、
\[ z_{3} = - \sqrt{\frac{p_{1}}{p_{3}}} z_{1}- \sqrt{\frac{p_{2}}{p_{3}}} z_{2}\]
だからだとか言われるが、次回はこれについて具体的に確かめる。
2013年4月28日日曜日
多項分布から多変量正規分布への近似について(3)
問題
二項分布は、サンプルサイズ\(n\)が大きいときに正規分布に近似できることを示せ。
である。
使う定理や公式
- スターリングの公式
- 指数関数
考えたこと
前回の確率関数\(W\)を\(C = \sqrt{2\pi}\)として、\(p_{i},~~m\)で書くと次のようにまとめられる。
\[ W \simeq \sqrt{\frac{n}{2\pi n_{1}n_{2}}} \left[ \left(1 + \frac{m}{np_{1}}\right)^{p_{1}}\left(1 - \frac{m}{np_{2}}\right)^{p_{2}} \right]^{-n} \left(1 + \frac{x_{1}^{2}}{m}\right)^{-m}\left(1 + \frac{-x_{2}^{2}}{m}\right)^{m}\]
\([\ldots]\)以外は前回と同じ。\(x^{2} = x_{1}^{2} + x_{2}^{2}\)として、
\[ [\ldots] \simeq \exp \left(\frac{x_{1}^{2}}{2p_{1}p_{2}}\right) = e^{x^{2}/2} ~~(*)\]
になることを認めれば、
\[ \begin{align*} &W \simeq \frac{1}{\sqrt{2\pi n p_{1}p_{2}}} e^{-x^{2}/2}\\ &x = \frac{n_{1} - np_{1}}{\sqrt{np_{1}p_{2}}}\\ \end{align*} \]
となる。従って、二項分布は\(n >> 1\)のとき正規分布に近づく。
\([\ldots]\)以外は前回と同じ。\(x^{2} = x_{1}^{2} + x_{2}^{2}\)として、
\[ [\ldots] \simeq \exp \left(\frac{x_{1}^{2}}{2p_{1}p_{2}}\right) = e^{x^{2}/2} ~~(*)\]
になることを認めれば、
\[ \begin{align*} &W \simeq \frac{1}{\sqrt{2\pi n p_{1}p_{2}}} e^{-x^{2}/2}\\ &x = \frac{n_{1} - np_{1}}{\sqrt{np_{1}p_{2}}}\\ \end{align*} \]
となる。従って、二項分布は\(n >> 1\)のとき正規分布に近づく。
(*)の評価
(*)を示すには
\[ \lim_{n\to \infty} \left| e^{x^{2}/2} - \left[\left(1 + \frac{m}{np_{1}}\right)^{p_{1}}\left(1 - \frac{m}{np_{2}}\right)^{p_{2}} \right]^{n}\right| = 0\]
になることを確認すれば良い。
\[
\begin{align*}
a &= e^{x^{2}/(2n)}\\
b &= \left(1 + \frac{m}{np_{1}}\right)^{p_{1}} \left(1 - \frac{m}{np_{2}}\right)^{p_{2}}
\end{align*}
\]
としたとき、\(a^{n} - b^{n} = (a-b)\sum_{i=0}^{n-1}a^{i}b^{n-1-i}\)なので\(a \le c,~~b \le c\)な\(c\)をとると
\[ |a^{n} - b^{n}| \le n c^{n}|a - b|\]
である。\(b\)を\(m/n\)で展開して、\(x_{1} = m/\sqrt{n} = \mathcal{O}(1)\)であることを使うと
\[ \begin{align*} a - b &= e^{y/n} - \left( 1 + \frac{y}{n} + \mathcal{O}(N^{-3/2})\right)\\ &= \left(\frac{y}{n}\right)^{2}\sum_{k=2}^{\infty}\frac{1}{k!}\left(\frac{y}{n}\right)^{k} + \mathcal{O}(N^{-3/2})\\ \end{align*} \]
としたとき、\(a^{n} - b^{n} = (a-b)\sum_{i=0}^{n-1}a^{i}b^{n-1-i}\)なので\(a \le c,~~b \le c\)な\(c\)をとると
\[ |a^{n} - b^{n}| \le n c^{n}|a - b|\]
である。\(b\)を\(m/n\)で展開して、\(x_{1} = m/\sqrt{n} = \mathcal{O}(1)\)であることを使うと
\[ \begin{align*} a - b &= e^{y/n} - \left( 1 + \frac{y}{n} + \mathcal{O}(N^{-3/2})\right)\\ &= \left(\frac{y}{n}\right)^{2}\sum_{k=2}^{\infty}\frac{1}{k!}\left(\frac{y}{n}\right)^{k} + \mathcal{O}(N^{-3/2})\\ \end{align*} \]
となる。\(y/n \le 1\)なので右辺第一項は
\[ \text{第一項} \le \left(\frac{y}{n}\right)^{2}\sum_{i=2}^{\infty} \frac{1}{k!} = \left(\frac{y}{n}\right)^{2}(e - 2)\]
であるから全体として\(a - b = \mathcal{O}(N^{-3/2})\)である。これを使うと\( |a^{n}-b^{n}| = \mathcal{O}(N^{-1/2})\)になることが分かるので\(n \to \infty\)で\(0\)に収束する。
\[ \text{第一項} \le \left(\frac{y}{n}\right)^{2}\sum_{i=2}^{\infty} \frac{1}{k!} = \left(\frac{y}{n}\right)^{2}(e - 2)\]
であるから全体として\(a - b = \mathcal{O}(N^{-3/2})\)である。これを使うと\( |a^{n}-b^{n}| = \mathcal{O}(N^{-1/2})\)になることが分かるので\(n \to \infty\)で\(0\)に収束する。
次回
次回は三項分布について同じ解析をする予定です。$k$項分布の場合の式を予想するためです。
2013年4月27日土曜日
多項分布から多変量正規分布への近似について(2)
問題
二項分布は、サンプルサイズ\(n\)が大きいときに正規分布に近似できることを示せ。二項分布の確率関数\(W\)は
\[
\begin{align*}
&W = \frac{n!}{n_{1}! n_{2}!} p_{1}^{n_{1}} p_{2}^{n_{2}}\\
&n_{1} + n_{2} = n\\
&p_{1} + p_{2} = 1
\end{align*}
\]
である。
使う定理や公式
- 不完全なスターリングの公式
- 指数関数
考えたこと
不完全なスターリングの公式を使うと
\[
\frac{n!}{n_{1}!n_{2}!} \simeq \frac{1}{C} \sqrt{\frac{n}{n_{1}n_{2}}}\left(\frac{n_{1}}{n}\right)^{- n_{1}}\left(\frac{n_{2}}{n}\right)^{-n_{2}}
\]
である。これを使うと、次のようにまとめることができる。
\[ W \simeq \frac{1}{C} \sqrt{\frac{n}{n_{1}n_{2}}} \left(\frac{n_{1}}{np_{1}}\right)^{- n_{1}}\left(\frac{n_{2}}{np_{2}}\right)^{-n_{2}} \]
連続分布に移行したとき、\(X_{i} \sim \mathrm{N}(0, 1)\)になるように定義する。
\[ E[n_{i}] = V[n_{i}] = np_{i}\]
なので、
\[ x_{i} = \frac{n_{i} - np_{i}}{\sqrt{np_{i}}} = \frac{m_{i}}{\sqrt{np_{i}}}\]
とする。以下では、計算を容易にするために\(m = m_{1} = - m_{2}\)として、\(p_{1} = p_{2} = p (= 1/2)\)の場合を考える。すると次のように書ける。
\[ W \simeq \frac{1}{C} \sqrt{\frac{n}{n_{1}n_{2}}} \left(1 + \frac{-2x_{1}^{2}}{n}\right)^{-n/2} \left(1 + \frac{x_{1}^{2}}{m}\right)^{-m}\left(1 + \frac{-x_{1}^{2}}{m}\right)^{m}\]
\(\sqrt{\frac{n}{n_{1}n_{2}}} \simeq 2 \frac{1}{\sqrt{n}}\)で、指数関数の定義を使えば
\[ \begin{align*} W &\simeq \frac{2}{C\sqrt{n}} e^{-x^{2}/2}\\ x^{2} &= x_{1}^{2} + x_{2}^{2} \end{align*} \]
となる。前回得た式と\(\Delta n = \frac{\sqrt{n}}{2} \Delta x\)を使うと
\[ p(x) = \frac{1}{C}e^{-x^{2}/2}\]
を得る。確率密度関数の規格化条件から\(C = \sqrt{2\pi}\)も得られるので、完全なスターリングの公式も得られる。
次回
\(p_{1} = p_{2}(=1/2)\)の場合のみ得られたので、次回はこの制限を取り払った場合を考える。
2013年4月26日金曜日
多項分布から多変量正規分布への近似について(1)
理解したいこと
多項分布はサンプルサイズが大きいとき、多変量正規分布に近づく
考えたこと
多項分布は離散分布で、多変量正規分布は連続分布なので、どう対応付くのかまず知る必要がある。ここでは1変数になる場合を考える。離散確率変数が\(n\)から\(n + \Delta n\)に変わるのに対応して、近似によって得られた確率分布の連続確率変数が\(x\)から\(x + \Delta x\)に変わるとすると
\[ \mathrm{Prob}[x \le X \le x + \Delta x] = \mathrm{Prob}[n \le N \le n + \Delta n]\]
が成り立つ。左辺を確率密度関数\(p(x)\)で書くと、
\[ \mathrm{Prob}[x \le X \le x + \Delta x] = \int_{x}^{x + \Delta x} p(x^{\prime}) dx^{\prime}\]
となる。\(F(x) = \int p(x) dx\)とすると右辺は\(F(x + \Delta x) - F(x)\)になる。\(\Delta x\)が小さいとして、\(F(x+\Delta x)\)を\(x\)まわりでテイラー展開すると、次のようになる。
\[ \begin{align*} F (x + \Delta x) &= F (x) + \frac{dF(x)}{dx} \Delta x + \mathcal{O}(\Delta x^{2})\\ &= F(x) + p(x)\Delta x + \mathcal{O}(\Delta x^{2}) \end{align*} \]
従って
\[ \mathrm{Prob}[x \le X \le x + \Delta x] = p(x) \Delta x + \mathcal{O}(\Delta x^{2})\]
だと分かる。次に離散確率の場合は、確率\(W(n)\)を使って書くと
\[ \mathrm{Prob}[n \le N \le n + \Delta n] = \sum_{k=n}^{n + \Delta n} W(k)\]
である。\(k\)が\(n\)から\(n + \Delta n\)に変化する間、\(\Delta x\)しか変化しないので、\(W(k)\)の値もあまり変化しない。従って
\[ \mathrm{Prob}[n \le N \le n + \Delta n] \simeq W(n)\Delta n\]
である。以上より、離散分布を連続分布に近似できるとき、次の関係が成り立つことが分かった。
\[ W(n)\Delta n \simeq p (x) \Delta x \]
\[ \mathrm{Prob}[x \le X \le x + \Delta x] = \mathrm{Prob}[n \le N \le n + \Delta n]\]
が成り立つ。左辺を確率密度関数\(p(x)\)で書くと、
\[ \mathrm{Prob}[x \le X \le x + \Delta x] = \int_{x}^{x + \Delta x} p(x^{\prime}) dx^{\prime}\]
となる。\(F(x) = \int p(x) dx\)とすると右辺は\(F(x + \Delta x) - F(x)\)になる。\(\Delta x\)が小さいとして、\(F(x+\Delta x)\)を\(x\)まわりでテイラー展開すると、次のようになる。
\[ \begin{align*} F (x + \Delta x) &= F (x) + \frac{dF(x)}{dx} \Delta x + \mathcal{O}(\Delta x^{2})\\ &= F(x) + p(x)\Delta x + \mathcal{O}(\Delta x^{2}) \end{align*} \]
従って
\[ \mathrm{Prob}[x \le X \le x + \Delta x] = p(x) \Delta x + \mathcal{O}(\Delta x^{2})\]
だと分かる。次に離散確率の場合は、確率\(W(n)\)を使って書くと
\[ \mathrm{Prob}[n \le N \le n + \Delta n] = \sum_{k=n}^{n + \Delta n} W(k)\]
である。\(k\)が\(n\)から\(n + \Delta n\)に変化する間、\(\Delta x\)しか変化しないので、\(W(k)\)の値もあまり変化しない。従って
\[ \mathrm{Prob}[n \le N \le n + \Delta n] \simeq W(n)\Delta n\]
である。以上より、離散分布を連続分布に近似できるとき、次の関係が成り立つことが分かった。
\[ W(n)\Delta n \simeq p (x) \Delta x \]
次回
今回得られた関係を使って、二項分布を正規分布に近似できることを示す。
登録:
投稿 (Atom)