はじめに
こちら統計検定準1級のテキストでおなじみ、「統計検定準1級対応 統計学実践ワークブック」の例・例題解答をまとめたものです。
テキストはこちら
今回の記事ではワークブックの21〜32章について解答をまとめています。
個人的なメモを記事にまとめたものなので、記述ミスなどがあるかもしれません。
その場合は優しく教えていただけるとありがたいです笑
その他注意事項
- 著作権の兼ね合いで問題は記述していません
- 問題の解答方法は1つでは無いので、ワークブックの解法とは違う解き方をしているものも存在します
- 一部理解が追いついてないためスキップしている問題があります(試験に出ないわけではありません)
関連記事
また、以下の記事内容も一緒に学習していたので併せて参考にしてみてください
1.例題解説 1~10章 統計学実践ワークブック
2.例題解説 11~20章 統計学実践ワークブック
3.統計用語の解説がわかりやすかったサイトチートシート
今後受験される方の参考になれば幸いです。
21章 標本調査法
例
-
例1
- 標本平均の分散
- 非復元抽出:$\frac{N-n}{N-1} \cdot \frac{1}{n} \sigma^2$
- 復元抽出:$\frac{1}{n} \sigma^2$
標本平均の分散
\displaylines{ \frac{N - n}{N-1} \cdot \frac{1}{n} \sigma^2 \\ = \frac{1600 - 100}{1600-1} \cdot \frac{1}{100} 120\\ = 1.125 }
標本平均が0.5以上になるには
\displaylines{ 0.5= \frac{1600 - n}{1600-1} \cdot \frac{1}{n} 120\\ n = \frac{1600 \cdot 120}{120 + (1600 - 1)\frac{1}{2}}\\ n = 208.809 }
- 標本平均の分散
-
例2
-
ネイマン分配法:推定量の分散を最小にする
\displaylines{ n_h = \frac{N_h \cdot \sigma_h \cdot \sqrt{\frac{N_h}{N_h-1}}}{\sum_{h=1}^L N_h \cdot \sigma_h \cdot \sqrt{\frac{N_h}{N_h-1}}} \cdot n }
$h〜L$:層の数
有限修正項を無視すると
\displaylines{ n_h = \frac{N_h \cdot \sigma_h }{\sum_{h=1}^L N_h \cdot \sigma_h } \cdot n }
\displaylines{ \frac{500 \cdot 6 \cdot \sqrt{\frac{500}{500 - 1}}}{(500\cdot 6+500\cdot 18)\sqrt{\frac{500}{500 - 1}}}\cdot 20\\ = \frac{3000\cdot 1.001}{12000\cdot 1.001}\cdot 20\\ = 5 }
-
例題
-
問 21.1
- 無作為抽出
- 集落抽出
- 層化抽出
- 優位抽出
- 二段抽出
より、2が正解
-
問 21.2
-
(1)
\displaylines{ \hat{Y} = \left(\frac{20}{2}\cdot 2\cdot 15\right)+\left(\frac{10}{2}\cdot 2\cdot 150\right)+\left(\frac{5}{2}\cdot 2\cdot 510\right)+\left(\frac{5}{2}\cdot 2\cdot 1010\right)\\ = 9,400\\ \frac{\hat{Y}}{40} = 235 }
-
(2)
ネイマン分配法は
\displaylines{ n_h = \frac{N_h \cdot \sigma_h \cdot \sqrt{\frac{N_h}{N_h-1}}}{\sum_{h=1}^L N_h \cdot \sigma_h \cdot \sqrt{\frac{N_h}{N_h-1}}} \cdot n }
より以下が分散の少ない分配比率となり、Cが近い
層No 分配比率 1 1.02 2 1.83 3 4.02 4 1.11
-
-
問 21.3
-
(1)
- $n=\frac{N}{2}$として計算すると$\frac{\sigma^2}{N-1}$となるので✗
- 無限の場合は$\frac{\sigma^2}{n}$、有限の場合は$\frac{N-n}{N-1}\frac{\sigma^2}{n}$より、常に無限の場合が大きくなるので◯
- どちらも$\frac{\sigma^2}{n}$となるので◯
-
(2)
-
非復元抽出
\displaylines{ \frac{N-n}{N-1}\frac{1}{n}\sigma^2\\ = \frac{9585-600}{9585-1}\frac{1}{600}\cdot 420\\ = 0.656 }
-
復元抽出
\displaylines{ \frac{\sigma^2}{n}= \frac{420}{600}\\ = 0.7 }
-
-
22章 主成分分析
例
-
例1
分散:$\frac{1}{n-1}\sum^n_i (x - \bar{x}) = 9.1$共分散:$\frac{1}{n-1}\sum_i^n (x_{1i} - \bar{x_1})(x_{2i} - \bar{x_2})=\frac{-3.5}{6-1} = -0.7$
共分散行列:対角は分散、それ以外はそれぞれの変数の共分散
相関行列:
対角は1それ以外は相関係数となる行列:$\frac{s_{i,j}}{\sqrt{s_{i,i} s_{j,j}}} = \frac{共分散}{\sqrt{分散_i 分散_j}}$
-
例2
-
用語
-
寄与率:主成分全体に対する割合$\frac{\lambda_i}{\sum_i^n \lambda_i}$
-
主成分得点:i番目の実数(or標準化した数)にそれぞれの固有ベクトルを掛けて足したもの
$z_{ik} = \sum^p_{j=1} x_{ji} \cdot \phi_{jk}$
$z_{ik}$: 第k主成分のi番目の得点(主成分得点)
$x_{ij}$: j変数のi番目の標準化データ(標準化前の場合もある)
$\phi_{jk}$: j変数に対する第 k主成分の固有ベクトル
-
主成分負荷量:もとの変数と主成分の相関係数
$r_{y_j, x_k} = \frac{Cov[x_k, y_j]}{\sqrt{V[y_j] s_P{k,k}}} = \frac{\lambda_ju_{k,j}} {\sqrt{s_{k,k}}}$
-
-
(1)
$\frac{20.2+19.4}{20.2+19.4+0.85+0.18} = \frac{39.6}{40.63} = 0.974$
97%以上説明できているので問題ない
-
(2)
$\frac{\sqrt{\lambda} \cdot u_i(固有ベクトル)}{\sqrt{s_{k,k}(その列の分散)}} = \frac{\sqrt{19.4} \cdot 0.564}{\sqrt{9.1}} = 0.823$
-
例題
- 問 22.1
-
(1)
主成分負荷量:$r_{jk} = \frac{\sqrt{\lambda} u_k} {s_{k,k}}$、相関行列の場合は$\sqrt{\lambda} \cdot u_k$となるので
\displaylines{ u_1 = \sqrt{2.26} \cdot -0.497 = -0.74\\ u_2 = \sqrt{1.5} \cdot 0.217 = 0.265 }
-
(2)
ワークブックの回答の通りなので省略
-
(3)
第一主成分:足の遅さ
第二主成分:肩の強さ
なので左上のほうが良い結果となり、A
-
23章 判別分析
例題
-
問 23.1
-
(1)
明らかに線形で判別できないのでSVMの精度が高くなる
-
(2)
適切にガウシアンカーネルのパラメータやソフトマージンパラメータが調整されておらず過学習している状態と想定される
クロスバリデーションなどを利用して適切に調整する必要がある
-
-
問 23.2
-
E, F, G の重心はそれぞれ
E:$(W^T x_E, W^T z_E) = (6.34, 4.24)$
F:$(W^T x_F, W^T z_F) = (0.95, 0.22)$
G:$(W^T x_G, W^T z_G) = (9.93, 0.22)$
$X$ を計算
$X_A$:$(W^T x_L, W^T z_L) = (0.44, -1.24)$Eとの距離を計算すると
\displaylines{ (0.44 - 6.34, -1.24 - 4.24 )\\ = (-5.9, -5.48)\\ (-5.9, -5.48)^T(-5.9, -5.48)\\ = (-5.9)^2 + (-5.48)^2 \\ = 34.81 + 30.0304 = 64.8404\\ }
よって
$\sqrt{64.8404} \approx 8.05$Fとの距離:$\sqrt{4.4117} \approx 2.1$
Gとの距離:$\sqrt{92.1917} \approx 9.6$
一番距離が小さいFに分類される
-
-
問 23.3
ベイズの定理を用いて分母と分子をそれぞれ表す
$P(y=1| X) = \frac{P(X|y=1)P(y=1)}{P(X|y=1)P(y=1) + P(X|y=-1)P(y=-1)}$
$P(y=-1|X) = \frac{P(X|y=-1)P(y=-1)}{P(X|y=1)P(y=1) + P(X|y=-1)P(y=-1)}$
上記より
\displaylines{ f(x) = \log \frac{P(y=1|X)}{P(y=-1|X)}\\ = \log \frac{P(X|y=1)P(y=1)}{P(X|y=-1)P(y=-1)} }
問題より
事前分布:$P(y=1), P(y=-1)$はそれぞれ$\pi_1, \pi_2$
\displaylines{ \log \frac{P(X|y=1)\pi_1}{P(X|y=-1)\pi_2}\\ = \log \frac{P(X|y=1)}{P(X|y=-1)} + \log \frac{\pi_1}{\pi_2} }
問題より
$P(X|y=1), P(X|y=-1)$はそれぞれ$N(\mu_1, \Sigma_1), N(\mu_2, \Sigma_2)$に従うので、確率は
\displaylines{ P(X|y=1) = \frac{1}{\sqrt{(2\pi)^k |\Sigma_1|}} e^{-\frac{1}{2} (x-\mu_1)^T \Sigma_1^{-1} (x-\mu_1)}\\ P(X|y=-1) = \frac{1}{\sqrt{(2\pi)^k |\Sigma_2|}} e^{-\frac{1}{2} (x-\mu_2)^T \Sigma_2^{-1} (x-\mu_2)} }
$e^{-\frac{1}{2} (x-\mu_1)^T \Sigma_1^{-1} (x-\mu_1)}$の$\mu_1$はベクトルなので$(x-\mu_1)^T(x-\mu_1)$となります。
\displaylines{ \log \frac{\frac{1}{\sqrt{(2\pi)^k |\Sigma_1|}} e^{-\frac{1}{2} (x-\mu_1)^T \Sigma_1^{-1} (x-\mu_1)}}{ \frac{1}{\sqrt{(2\pi)^k |\Sigma_2|}} e^{-\frac{1}{2} (x-\mu_2)^T \Sigma_2^{-1} (x-\mu_2)}} + \log \frac{\pi_1}{\pi_2}\\ }
\displaylines{ \log \frac{|\Sigma_2|^{\frac{1}{2}}}{|\Sigma_1|^{\frac{1}{2}}} + \log e^{{-\frac{(x-\mu_1)^T \Sigma_1^{-1} (x-\mu_1)}{2}} - (- \frac{(x-\mu_2)^T \Sigma_2^{-1} (x-\mu_2)}{2})} + \log \frac{\pi_1}{\pi_2} }
$log\ e$の部分を計算すると
\displaylines{ \frac{1}{2}\Big( (x-\mu_2)^T \Sigma_2^{-1} (x-\mu_2) - (x-\mu_1)^T \Sigma_1^{-1} (x-\mu_1)\Big)\\ =\frac{1}{2} \Big( (x^T\Sigma_2^{-1}x - x^T\Sigma_2^{-1}\mu_2 - \mu_2^T\Sigma_2^{-1}x + \mu_2^T\Sigma_2^{-1}\mu_2)\\ - (x^T\Sigma_1^{-1}x - x^T\Sigma_1^{-1}\mu_1 - \mu_1^T\Sigma_1^{-1}x + \mu_1^T\Sigma_1^{-1}\mu_1) \Big) }
計算すると
\displaylines{ f(x) = -\frac{1}{2}x^T (\Sigma_1^{-1} - \Sigma_2^{-1})x \\ + (\Sigma_1^{-1}\mu_1 - \Sigma_2^{-1}\mu_2)^Tx\\ + \frac{1}{2}(\mu_2^T \Sigma_2^{-1} \mu_2 - \mu_1^T \Sigma_1^{-1} \mu_1) + \log \frac{|\Sigma_2|^{\frac{1}{2}}}{|\Sigma_1|^{\frac{1}{2}}} + \log \frac{\pi_1}{\pi_2} }
-
(2)
計算は省略
-
24章 クラスター分析
例題
-
問 24.1
-
(1)
(b)最近隣法:近いものから結合していく
(a)再遠燐法:距離が小さいもの同士でクラスタを作成し、作成した中から大きい数値をクラスタ距離
-
最遠燐法の具体例
例:6×6 の距離行列
A B C D E F A 0 2 6 10 9 12 B 2 0 4 8 7 10 C 6 4 0 7 6 9 D 10 8 7 0 5 6 E 9 7 6 5 0 3 F 12 10 9 6 3 0
-
最も近いクラスタを探す
- $d(E,F)=3d(E, F) = 3$ が最小なので、クラスタ EとF を統合。
-
新しいクラスタ CEFC_{EF} の距離を計算
-
各クラスタとの距離は 最大値 をとる:
$$
d(C_{EF}, A) = max(9, 12) = 12\
d(C_{EF}, B) = max(7, 10) = 10\
d(C_{EF}, C) = max(6, 9) = 9\
d(C_{EF}, D) = max(5, 6) = 6
$$
-
-
距離行列を更新
A B C D EF A 0 2 6 10 12 B 2 0 4 8 10 C 6 4 0 7 9 D 10 8 7 0 6
┌──(E, F) │ ├──(D, (E, F)) │ ├──(C, (D, (E, F))) │ ├──(B, (C, (D, (E, F)))) │ └──(A, (B, (C, (D, (E, F)))))
-
-
(2)
ユークリッド距離:$d_m (x,y)= \sqrt{\sum_{i=1}^p |x_i - y_i|^2 }$
k-means
初期の中心点と各データのユークリッド距離を求める
$x_1^{(0)}$
データ 平均$\bar{x}_1^{(0)}$との差 $(x_i - \bar{x}_1^{(0)})^T(x_i - \bar{x}_1^{(0)})$ $\sqrt{(x_i - \bar{x}_1^{(0)})^T(x_i - \bar{x}_1^{(0)})}$ 7.7, 5.8, 7.7 -1.3, -0.2, 0.7 2.22 1.489 10.1, 5.9, 7.3 1.1, -0.1, 0.3 1.31 1.144 4.9, 4.1, 7.9 -4.1, -1.9, 0.9 21.23 4.607 4.8, 5.0, 9.0 -4.2, -1.0, 2.0 22.64 4.758 6.4, 3.7, 8.0 -2.6, -2.3, 1.0 13.05 3.612 4.0, 5.2, 11.6 -5.0, -0.8, 4.6 46.8 6.841 $x_2^{(0)}$
5,5,10
データ 平均$\bar{x}_1^{(0)}$との差 $(x_i - \bar{x}_1^{(0)})^T(x_i - \bar{x}_1^{(0)})$ $\sqrt{(x_i - \bar{x}_1^{(0)})^T(x_i - \bar{x}_1^{(0)})}$ 7.7, 5.8, 7.7 -2.7, 0.8, -2.3 13.22 3.635 10.1, 5.9, 7.3 5.1, 0.9, -2.7 34.11 5.84 4.9, 4.1, 7.9 -0.1, -0.9, -2.1 5.23 2.28 4.8, 5.0, 9.0 -0.2, 0, -1 1.04 1.019 6.4, 3.7, 8.0 1.4, -1.3, -2 7.65 2.76 4.0, 5.2, 11.6 -1, 0.2, 1.6 3.6 1.897 上記より
ユーザー クラスター 1 1 2 1 3 2 4 2 5 2 6 2
-
25章 因子分析・グラフィカルモデル
例
-
例1
-
用語
- 共通因子:数学の能力のような観測できない共通の要因
- 因子スコア($f_i$):数学力のような各情報を集約したもの
- 因子負荷量:モデルを $x_{ij}= a_j f_i +d_j u_{ij}$ とした時、$a_j$ の部分で共有因子と各設問の関係を持つ係数
- $d_j$:独自係数
- $u_{ij}$:独自性因子
-
(1)
方程式、微積分、複素数が高いので、計算力と推定できる
-
(2)
共通性と独自性はそれぞれ、$V[x_{ij}] = a_j^2 + d_j^2 = 1$ より
因子負荷量$a_j$ が0.92より
共通性:$0.92^2 = 0.846..$
$0.92^2 + b_j^2 = 1$ より
独自性:$1- 0.84 = 0.1536$
-
(3)
数学の計算力と捉えられるので
得点 因子スコア A 55 -1.613 B 80 1.462 C 72 0.193 平均 70
-
-
例2
両辺に$X_1$を掛けて期待値を取る\displaylines{ X_3 = b \cdot X_2 + v\\ E[X_3 \cdot X_1] = E[b\cdot X_2 \cdot X_1] +E[ v \cdot X_1]\\ p_{13} = b \cdot p_{12} + 0\\ b = \frac{p_{13}}{p_{12}} }
例題
-
問 25.1
\displaylines{ a^2 = \text{第一因子}^2 + \text{第二因子}^2 }
より、
\displaylines{ a_2 = \sqrt{0.8545 -(-0.92^2)}\\ = 0.09 }
第一因子の酸味が弱く、第二因子の味の濃さが高いJが近い
-
問 25.2
式1:$Y= aX + bW + u$式2:$Z = cX + dY +v$
-
式1
両辺に$X$を掛ける
\displaylines{ XY = aX^2 + b XW + uX }
期待値を取る
\displaylines{ E[XY] = aE[X^2]+bE[XW]+uE[X] }
ここで、
-
$Cov[X,Y](=p_{xy}) = E[XY] - E[X]E[Y]$($E[X],E[Y]$ は問題文より標準化されて期待値0)
-
$V[X] = E[X^2] - (E[X]^2)$→$E[X]$は$N(0,1)$より期待値0なので$V[X]=1$
-
$E[XW]$:問題文より独立なので0
より
$p_{xy} = a$
両辺に$W$を掛けて期待値を取ると
\displaylines{ YW = aXW + bW^2 + uW\\ E[YW] = aE[XW] +bE[W^2]+uE[W] }
$p_{yw} = b$
-
-
式2
両辺に$X$を掛けて期待値を取ると
\displaylines{ E[XZ] = cE[X^2] + dE[XY]+vE[X] }
$p_{xz} = c + d\ p_{xy}$
両辺に$Y$を掛けて期待値を取ると
\displaylines{ E[YZ] = c E[XY] + d E[Y^2] + v E[Y] }
$p_{yz} = c\ p_{xy} + d \cdot 1 + 0$
2つの式を計算すると
\displaylines{ d = \frac{p_{xz}\cdot p_{xy} - p_{yz}}{(p_{xy})^2 -1} }
\displaylines{ c = \frac{p_{xz} - p_{yz}\cdot p_{xy}}{1-(p_{xy})^2} }
-
(2)
\displaylines{ p_{xz}= \frac{p_{xz} - p_{yz}\cdot p_{xy}}{1-(p_{xy})^2} + p_{xy} \cdot \frac{p_{xz}\cdot p_{xy} - p_{yz}}{(p_{xy})^2 -1}\\ = \frac{p_{xz} - p_{yz}\cdot p_{xy}}{1-(p_{xy})^2} + \frac{p_{xy} (p_{xz} \cdot p_{xy} - p_{yz})}{1-(p_{xy})^2}\\ = \frac{(p_{xz} - p_{yz} \cdot p_{xy}) + (p_{xy}\cdot p_{xz} \cdot p_{xy} - p_{xy} \cdot p_{yz})}{1-(p_{xy})^2 }\\ = \frac{p_{xz} - p_{xz} \cdot (p_{xy})^2}{1-(p_{xy})^2}\\ = \frac{p_{xz}(1 - (p_{xy})^2)}{1-(p_{xy})^2}\\ = p_{xz} }
-
26章 その他の多変量解析手法
例題
-
問 26.1
-
計量MDSの手順
-
距離行列Dを準備
-
二重中心化して行列Bを作成
\displaylines{ B = -\frac{1}{2} (I_n - \frac{1}{n} J_n)D (I_n - \frac{1}{n} J_n) }
$I_n$:対角が1の単位行列
$J_n$:全てが1の行列
-
Bの固有値、固有ベクトルを求める
$|B - \lambda I| = 0$
-
次元を選択し、固有値の大きい順に並べる
固有値 ( A = \begin{Vmatrix} \lambda_1 & 0\ 0 & \lambda_2 \end{Vmatrix}),
固有ベクトル ( U=(u_1, u_2) )
-
座標軸Xを求めて可視化する
$X = UA^{\frac{1}{2}}$
-
-
(1)
\displaylines{ B = -\frac{1}{2}\Big( \begin{Vmatrix} 1&0&0&0\\ 0&1&0&0\\ 0&0&1&0\\ 0&0&0&1\\ \end{Vmatrix} - \frac{1}{4} \begin{Vmatrix} 1&1&1&1\\ 1&1&1&1\\ 1&1&1&1\\ 1&1&1&1\\ \end{Vmatrix} \Big) \begin{Vmatrix} 0&5&5&16\\ 5&0&4&5\\ 5&4&0&5\\ 16&5&5&0\\ \end{Vmatrix} \Big( \begin{Vmatrix} 1&0&0&0\\ 0&1&0&0\\ 0&0&1&0\\ 0&0&0&1\\ \end{Vmatrix} - \frac{1}{4} \begin{Vmatrix} 1&1&1&1\\ 1&1&1&1\\ 1&1&1&1\\ 1&1&1&1\\ \end{Vmatrix} \Big)\\ }
計算すると
\displaylines{ \begin{Vmatrix} 4&0&0&-4\\ 0&1&-1&0\\ 0&-1&1&0\\ -4&0&0&4\\ \end{Vmatrix} }
Bの固有値を解くと
\displaylines{ |B - \lambda I|= 0 }
\displaylines{ \Big| \begin{Vmatrix} 4&0&0&-4\\ 0&1&-1&0\\ 0&-1&1&0\\ -4&0&0&4\\ \end{Vmatrix} - \begin{Vmatrix} \lambda &0&0&0\\ 0&\lambda &0&0\\ 0&0&\lambda &0\\ 0&0&0&\lambda \\ \end{Vmatrix}\Big|\\ = \begin{Vmatrix} 4-\lambda &0&0&-4\\ 0&1-\lambda &-1&0\\ 0&-1&1-\lambda &0\\ -4&0&0&4-\lambda \\ \end{Vmatrix} }
余因子展開を行う
\displaylines{ (4-\lambda ) \begin{Vmatrix} 1-\lambda &-1&0\\ -1&1-\lambda &0\\ 0&0&4-\lambda \\ \end{Vmatrix} - 0 \begin{Vmatrix} 0&-1&0\\ 0&1-\lambda &0\\ -4&0&4-\lambda \\ \end{Vmatrix}\\ + 0 \begin{Vmatrix} 0&1-\lambda &0\\ 0&-1&0\\ -4&0&4-\lambda \\ \end{Vmatrix} - (-4) \begin{Vmatrix} 0&1-\lambda &-1\\ 0&-1&1-\lambda \\ -4&0&0 \\ \end{Vmatrix}\\ = (4-\lambda ) \begin{Vmatrix} 1-\lambda &-1&0\\ -1&1-\lambda &0\\ 0&0&4-\lambda \\ \end{Vmatrix} +4 \begin{Vmatrix} 0&1-\lambda &-1\\ 0&-1&1-\lambda \\ -4&0&0 \\ \end{Vmatrix} }
行列式の計算$\begin{bmatrix}a & b\ c& d \end{bmatrix} = ad- bc$ より
\displaylines{ (4-\lambda )\big[(1-\lambda)(1-\lambda) -(4-\lambda) \big]\\ - (-4)\big[(1-\lambda)(1-\lambda)(-4) - (-1)(-1)(-4) \big]\\ = (4-\lambda)^2[ (1-\lambda)^2-1] - 16[(1-\lambda)^2 -1)]\\ = [(1-\lambda)^2-1][(4-\lambda)^2-16] }
$[(1-\lambda)^2-1][(4-\lambda)^2-16]=0$ を解くとよって固有値は
8,2,0,0
$B=\begin{Vmatrix}
4&0&0&-4\
0&1&-1&0\
0&-1&1&0\
-4&0&0&4\
\end{Vmatrix}$より
\displaylines{ \begin{Vmatrix} 4-\lambda &0&0&-4\\ 0&1-\lambda &-1&0\\ 0&-1&1-\lambda &0\\ -4&0&0&4-\lambda \\ \end{Vmatrix} }
$\lambda = 8$:
\displaylines{ \begin{bmatrix} -4 & 0 & 0 & -4 \\ 0 & -7 & -1 & 0 \\ 0 & -1 & -7 & 0 \\ -4 & 0 & 0 & -4 \end{bmatrix} \begin{bmatrix} x \\ y \\ z \\ w \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\ 0 \\ 0 \end{bmatrix} }
- $−4x−4w=0$
- $−7y−z=0$
- $-y - 7z = 0$
- $−4x−4w=0$
よって$\begin{bmatrix} 1 \ 0 \ 0 \ 1 \end{bmatrix}$
標準化(=1)するには
\displaylines{ \sqrt{a^2 + 0^2 + 0^2+ a^2} = 1\\ \sqrt{2 a^2} =1 \\ a = \sqrt{\frac{1}{2}} }
\displaylines{ \begin{bmatrix} \sqrt{\frac{1}{2}} \\ 0 \\ 0 \\ -\sqrt{\frac{1}{2}} \end{bmatrix} }
-
(2)
次元を1とした場合の座標X
第1固有ベクトル8の場合のXは
$X = UA^{\frac{1}{2}}$ より
\displaylines{ \begin{bmatrix} \sqrt{\frac{1}{2}} \\ 0 \\ 0 \\ -\sqrt{\frac{1}{2}} \end{bmatrix} \cdot \sqrt{8} = \begin{bmatrix} \sqrt{\frac{1}{2}} \\ 0 \\ 0 \\ -\sqrt{\frac{1}{2}} \end{bmatrix} \cdot 2\sqrt{2} = \begin{bmatrix} 2 \\ 0 \\ 0 \\ -2\end{bmatrix} }
この距離行列を求めると
\displaylines{ \begin{bmatrix} (2-2)^2 & (2-0)^2 & (2-0)^2 & (2-(-2))^2\\ (0-2)^2 & (0-0)^2 & (0-0)^2 & (0-(-2))^2\\ (0-2)^2 & (0-0)^2 & (0-0)^2 & (0-(-2))^2 \\ (-2-2)^2 & (-2-0)^2 & (-2-0)^2 & (-2 - (-2))^2 \end{bmatrix}\\ \begin{bmatrix} 0 & 4 & 4 & 16\\ 4 & 0 & 0 & 4\\ 4 & 0 & 0 & 4\\ 16 & 4 & 4 & 0 \end{bmatrix}\\ }
Dと一致して無いので再現できていない
-
(3)
次元を2とした場合の$X = UA^{\frac{1}{2}}$は
\displaylines{ \begin{bmatrix} \sqrt{\frac{1}{2}} & 0 \\ 0 & \sqrt{\frac{1}{2}}\\ 0 & -\sqrt{\frac{1}{2}}\\ -\sqrt{\frac{1}{2}} & 0 \end{bmatrix} \begin{bmatrix} \sqrt{8} & 0\\ 0 & \sqrt{2} \end{bmatrix}\\ = \begin{bmatrix} 2 & 0\\ 0 & 1\\ 0 & -1\\ -2 & 0 \end{bmatrix}\\ }
距離行列を計算すると
\displaylines{ \begin{bmatrix} (2-2)^2 + (0-0)^2 & (2-0)^2+ (0-1)^2& (2-0)^2+ (0-(-1))^2 & (2-(-2))^2 + (0-0)^2\\ (0-2)^2+ (1-0)^2 & (0-0)^2 + (1-1)^2 & (0-(0)^2 + (1-(-1))^2 & (0-(-2))^2 + (-1-0)^2\\ (0-2)^2 + (-1-0)^2 & (0-0)^2 + (-1-1)^2 & (0-0)^2 + (-1- (-1))^2 & (0-(-2) )^2 + (-1-0)^2\\ (-2-2)^2 + (0-0)^2 & (-2-2)^2 + (0-0)^2 & (-2- 0)^2 + (0-(-1))^2 & (-2 -(-2))^2 + (0-0)^2 \end{bmatrix}\\ = \begin{bmatrix} 0 & 5 & 5 & 16\\ 5 & 0 & 4 & 5 \\ 5 & 4 & 0 & 5\\ 16 & 5 & 5 & 0 \end{bmatrix}\\ }
Dと一致するので再現できている
-
27章 時系列解析
例題
-
問 27.1
-
用語
-
定常性
強定常:
$f(y_t, ... y_{t+l}) = f(y_{t+k}, ..., y_{t+k+l})$:時点の数が同じならどの集まりの同時分布も同一
弱定常:
$E(y_t) = \mu$:どの時点においても$y$の期待値は一定
$Cov(y_t, y_{t-j})$:時刻$t$によらず$j$次の自己共分散が一定
-
自己回帰仮定(AR仮定):
$Y_t = \phi_1 Y_{t-1} +U_t$
$\phi= \pm 1$:過去は影響せずに一定(or符号が変えるだけ)
$|\phi|>1$:過去に遡るほど影響が多きい
$|\phi|<0$:弱定常となる
自己共分散 $\gamma_h = \phi_1^h \frac{\sigma^2}{1-\phi_1^@}$
自己相関係数 $p_h = \phi_1^h$
$\phi<0$:符号が交互に変わる
$\phi > 0$:数値が大きくなる事に滑らかになる
-
-
回答
-
(1)
負の場合は符号が交互になり、正の場合は数が大きいほどなめからになるので、
a 0.7 b 0 c 1 d -0.8 -
(2)
最小二乗推定量は普遍推定量だが、最小分散を取ることは無いので、自己相関構造を考慮した最尤法や一般化最小二乗法などを用いる必要がある
DW比は2より小さくなる
-
-
-
問 27.2
-
用語
-
移動平均過程(MA過程)
1次の式:$Y_t = \mu + U_t + \theta_1 U_{t-1}$
$q$次の移動平均過程
分散 $\gamma_0 = (1+\theta_1^2+ ...+\theta_q^2 )\sigma^2$
自己共分散:$(\theta_h + \theta_1 \theta_{h=1},...+\theta_{q-h}\theta_q)\sigma^2$
-
自己回帰移動平均過程(ARMA過程)
階差($\Delta Y=Y_t - Y_{y-1}$)を用いた自己回帰移動平均過程は「自己回帰和分移動平均過程(ARIMA)」
-
ディッキーフラー検定(AR過程は拡張ディッキーフラー検定)
自己回帰過程において階差を取る必要があるのかを検定する
帰無仮説:$\phi = 1$
対立仮説:$|\phi| < 1$
-
次数の選択
自己相関係数(AR過程でゆっくり減衰) 編自己相関係数(MA過程でゆっくり減衰) 選択するモデル 2次以降0 ゆっくり減衰 MA(1) 3次以降0 ゆっくり減衰 MA(2) ゆっくり減衰 2次以降0 AR(1) ゆっくり減衰 3次以降0 AR(2) ゆっくり減衰 ゆっくり減衰 ARMA(1,1
-
-
回答
自己共分散
-
ラグ0の場合
\displaylines{ Cov(Y_t, Y_t) = V[Y_t]\\ = E[Y_t^2] - (E[Y_t])^2 }
$E[Y_t]$ の計算
\displaylines{ E[Y_t] = E[U_t + \theta_1 U_{t-1} + \theta_2 U_{t-2}]\\ = E[U_t]+ E[\theta_1 U_{t-1}] +E[ \theta_2 U_{t-2}]\\ = 0+0+0 }
$U_t$はホワイトノイズなので期待値0になる
よって共分散は
\displaylines{ Cov(Y_t, Y_t) = E[Y_t^2]\\ =E[(U_t + \theta_1 U_{t-1} + \theta_2 U_{t-2})^2]\\ = E[U_t^2 + \theta_1^2 U_{t-1}^2 +\theta_2^2 U_{t-2}^2 + 2\theta_1 U_t U_{t-1} + 2\theta_2U_t U_{t-2} + 2\theta_1 U_{t-1} U_{t-2}] }
ホワイトノイズは以下の性質がある
$E[U_t] = 0$
$E[U_t^2] = \sigma^2$
$E[U_t U_{t-k}]=0$
上記に適応すると
\displaylines{ \sigma^2 + \theta_1^2 \sigma^2 + \theta_2^2 \sigma^2\\ = \sigma^2(1+\theta_1^2 + \theta_2^2) }
-
ラグ1の場合
\displaylines{ \gamma(1) = Cov(Y_t, Y_{t-1})\\ = E[(Y_t - E[Y_t])(Y_{t-1} - E[Y_{t-1}])] }
$U_t$はホワイトノイズのため、$E[Y_t]=0$となるので
\displaylines{ E[Y_t\ Y_{t-1}] }
$Y_t = U_t + \theta_1 U_{t-1} + \theta_2 U_{t-2}$
$Y_{t-1} = U_{t-1} + \theta_1 U_{t-2} + \theta_2 U_{t-3}$
\displaylines{ Y_t Y_{t-1}\\ = (U_t + \theta_1 U_{t-1} + \theta_2 U_{t-2})(U_{t-1} + \theta_1 U_{t-2} + \theta_2 U_{t-3})\\ = U_t U_{t-1} + \theta_1 U_t U_{t-2} + \theta_2 U_t U_{t-3}\\ + \theta_1 U_{t-1}^2 + \theta_1^2 U_{t-1} U_{t-2} + \theta_1 \theta_2 U_{t-1} U_{t-3}\\ + \theta_2 U_{t-2} U_{t-1} + \theta_1 \theta_2 U_{t-2}^2 + \theta_2^2 U_{t-2} U_{t-3} }
ホワイトノイズの性質
$E[U_t] = 0$
$E[U_t^2] = \sigma^2$
$E[U_t\ U_{t-k}]= 0$
\displaylines{ \theta_1 E[U^2_{t-1}] + \theta_1 \theta_2 E[U^2_{t-2}]\\ = \theta_1 \sigma^2 + \theta_1 \theta_2 \sigma^2\\ = \sigma^2 (\theta_1 + \theta_1 \theta_2) }
-
ラグ2の場合
ラグ1と同様に$E[Y_t\ Y_{t-2}]$
\displaylines{ E[ (U_t + \theta_1 U_{t-1} + \theta_2 U_{t-2}) (U_{t-2} + \theta_1 U_{t-3} + \theta_2 U_{t-4}) ]\\ E[ U_t U_{t-2} + \theta_1 U_t U_{t-3} + \theta_2 U_t U_{t-4}\\ \theta_1 U_{t-1}U_{t-2}+ \theta_1^2 U_{t-1}U_{t-3} + \theta_1 \theta_2 U_{t-1}U_{t-4}\\ \theta_2 U_{t-2}^2 +\theta_1 \theta_2 U_{t-2} U_{t-3} + \theta^2 U_{t-2}U_{t-4} ] }
ホワイトノイズの性質
$E[U_t] = 0$
$E[U_t^2] = \sigma^2$
$E[U_t\ U_{t-k}]= 0$
より
\displaylines{ \theta_2 \sigma^2 }
ラグ3以降は0
-
-
-
問 27.3
-
DW検定
\displaylines{ DW = \frac{\sum^T_{t=2} (\hat{U_t} - \hat{U}_{t-1})^2}{\sum^T_{t=1}\hat{U_t}^2} }
自己相関$p$ 自己共分散$\hat{\gamma}$ DW $0$ $0$ 2 $p>0$ $\hat{\gamma}>0$ $DW<2$ $p<0$ $\hat{\gamma}<0$ $DW>2$ DW比の近似値:$2(1-\hat{\gamma}_1)$
-
次数の選択
自己相関係数(AR過程でゆっくり減衰) 偏自己相関係数(MA過程でゆっくり減衰) 選択するモデル 2次以降0 ゆっくり減衰 MA(1) 3次以降0 ゆっくり減衰 MA(2) ゆっくり減衰 2次以降0 AR(1) ゆっくり減衰 3次以降0 AR(2) ゆっくり減衰 ゆっくり減衰 ARMA(1,1) -
回答
-
(1)
自己回帰(AR)過程の場合は標本編自己相関係数を用いる
0になっている手前の次数を選択するので、AR(2)を選択する
-
(2)
DW比は1次自己相関を検出するものなので
$2(1-\gamma_1) = 2(1-0.691) = 0.618$
-
-
-
問 27.4
AICが一番小さいものを選択するのでAR(3)
28章 分割表
例
-
例1
帰無仮説:新薬と従来薬の効果は同等対立仮説:新薬のほうが効果が高い
効果あり 効果なし 計 新薬 6 1 7 従来薬 3 4 7 計 9 5 14
効果あり 効果なし 計 新薬 $X_{11}$ $7-X_{11}$ 7 従来薬 $9-X_{11}$ $7-(9-X_{11})=X_{11}-2$ 7 計 9 5 14 Xの取りうる範囲 $2〜7$
$X_{11}$が6以上の確率を超幾何分布にて計算する
超幾何分布の確率
\displaylines{ \frac{{}_{効果あり}C_{新薬の効果あり} \cdot {}_{効果無し}C_{新薬の効果なし}}{{}_{全体の数}C_{新薬グループの数}} }
\displaylines{ P(P=6) = \frac{{}_{9}C_{6} \cdot {}_{5}C_{1}}{{}_{14}C_{7}}\\ = \frac{84\cdot 5}{3432}= 0.12237 }
\displaylines{ P(P=7) = \frac{{}_{9}C_{7} \cdot {}_{5}C_{0}}{{}_{14}C_{7}}\\ = \frac{36}{3432}= 0.01049 }
0.1328…
よって5%で有意では無い
例題
-
問 28.1
- コントロールによって操作可能なので✗
- こちらもコントロールのサンプル数によって変動するので✗
- ◯
-
問 28.2
-
当てはめ値:最尤推定値より
$x_{ij} = 48\frac{10}{48}\frac{12}{48}= 2.5$ 10-2.5=7.5 12-2.5=9.5 38- (12-2.5) = 28.5 -
適合度$X^2$:
\displaylines{ X^2 = \sum \frac{(O_i - E_i)^2}{E_i} }
$O_i$:観測データ
$E_i$:期待度数
\displaylines{ \frac{(4-2.5)^2}{2.5}+\frac{(6-7.5)^2}{7.5}+\frac{(8-9.5)^2}{9.5}+\frac{(30-28.5)^2}{28.5}\\ = 0.9 + 0.3 + 0.2368 + 0.078 \\ = 1.516... }
自由度:カテゴリ数-推定するパラメータ数(今回は$x_{ij}$のみ)→$4-1=3$
-
逸脱度
\displaylines{ G^2 = 2 \log\ \Lambda = 2 \sum_i \sum_j x_{ij} \log \frac{x_{ij}}{E_{ij}} }
より
\displaylines{ 2( 4 \log \frac{4}{2.5} + 6 \log\frac{6}{7.5} + 8 \log\frac{8}{9.5} + 30 \log\frac{30}{28.5} )\\ = 2( 4 \cdot 0.470 + 6 \cdot -0.223 + 8 \cdot -0.172 + 30 \cdot 0.051 )\\ = 2 \cdot 0.696 = 1.392 }
自由度:同様に3
どちらも自由度3のΧ二乗により棄却できない
-
信頼区間
標本対数オッズ比
\displaylines{ \frac{\frac{4}{8}}{\frac{6}{30}} = \frac{1}{2}\frac{5}{1}= 2.5\\ \log\ 2.5 = 0.916 }
対数オッズ比の標準誤差
\displaylines{ \sqrt{\frac{1}{4} + \frac{1}{6} + \frac{1}{8} + \frac{1}{30}}\\ = \sqrt{0.57499} = 0.758 }
よって対数オッズ比の信頼区間
\displaylines{ 0.916 \pm 1.96 \cdot 0.758\\ (-0.569 , 2.401) }
対数を戻すと
$(0.5655, 11.04)$
信頼区間が1を含まれていることは検定が有意では無いことに対応する
-
-
問 28.3
-
グラフィカルモデルの因子の独立性
-
性質1:$A, B$ が連結でない場合、因子$A,B$は独立である
-
性質2:$A$と$B$の間に直接エッジがない $(A, B) \notin E$ 場合、他の全てのノードを条件にすると$A$と$B$は条件付き独立
$A \perp B \mid V \setminus {A, B}$
-
性質3:グラフ上で$A$と$B$を分離するノード集合$S$が存在するとき、$A$と$B$は$S$を条件にすると条件付き独立。
$A \perp B \mid S$
-
-
(1) 4つの頂点は連結しているので独立な因子は存在しない
-
(2)
- 性質2より$A\perp C | B,D$
- 性質3より$A\perp C | B$
-
29章 不完全データの統計処理
例
-
例1
ID 全てのデータ 観測されたデータ 1 35 35 2 38 38 3 42 42 4 45 45 5 49 49 6 52 52 7 53 53 8 56 56 9 63 60+ 10 65 60+ 全てのデータが観測されたとすると
$\bar{x} = 49.8$
これが最尤推定値となる
観測データの平均は$46.25$
打ち切りの更新アルゴリズム
\displaylines{ \mu^{t+1}= \frac{1}{n} \Big( \sum_{i=1}^m x_i + (n-m) \big(\mu^t + \frac{\phi(a^t)}{1- \Phi(a^t)}\sigma \big)\Big) }
ここで
$n$:データ数
$m$:観測できた数
$a$:打ち切りポイントのZ値($\frac{c- \mu^t}{\sigma}$)
$c$:打ち切りポイント
問題文より標準偏差10として
$a^0 = \frac{60-46.25}{10}=1.375$
$1-\Phi(1.375)= 0.084$ (標準正規分布表より)
$\phi(1.375) = \frac{1}{\sqrt{2 \cdot 3.14 }\cdot 10}\exp(- \frac{1.375^2}{2}) = 0.39904 \cdot 0.3885 = 0.155$
1回目の更新
\displaylines{ \frac{1}{10} \Big(370+ (10-8) \big(46.25 + \frac{0.155}{0.084}10 \big)\Big)\\ = \frac{1}{10}\Big(395.5 + 142.628\Big)\\ = 49.94 }
$a^1 = \frac{60-49.94}{10} = 1.005$
$1-\Phi(1.005) = 0.1587$
$\phi(1.005) = \frac{1}{\sqrt{2 \cdot 3.14 }\cdot 10} \exp(- \frac{1.005^2}{2}) = 0.39904 \cdot 0.6035 = 0.2408$
\displaylines{ \frac{1}{10} \Big(370+ (10-8) \big(49.94 + \frac{0.2408}{0.1587}10 \big)\Big)\\ = 50.02 }
-
例2
x y 35 40 37 38 40 48 45 41 47 60 54 46 57 57 60 62 63 $x$ の全てを利用した統計量
$\bar{x} = 50$
$\bar{\sigma_x^2} = 100.6$
欠損が無い7組のデータのそれぞれの統計量は
$\bar{x} = 45$
$\bar{s_x^2} = 59.71$
$\bar{y} = 47.14$
$\bar{s_y^2} = 62.41$
$r = \frac{39.14}{\sqrt{59.71 \cdot 62.41}}=0.6412$
回帰パラメータの推定値
誤差分散 $\hat{t} = s_Y^2 - \frac{(s_{xy})^2}{s_x^2} = 62.41 - \frac{(39.14)^2}{59.71} = 36.75$
$\hat{\beta} = \frac{s_{xy}}{s_x^2} =\frac{39.14}{59.71}= 0.655$
$\hat{a} = \bar{y} - \hat{\beta} \bar{x} = 47.14 - 0.655 \cdot 45 = 17.655$
上記より$X$の周辺分布以外のパラメータの最尤推定値は以下のようになる
$\hat{\mu_Y} = \hat{\alpha} + \hat{\beta} \hat{x} = 17.65 + 0.655 \cdot 50 = 50.65$
$\hat{\sigma}_{XY} = \hat{\beta} \cdot \hat{\sigma^2}_X =0.655 \cdot 100.6 = 65.94$
$\hat{\sigma}^2_Y = \hat{t}^2 + \hat{\sigma}_{XY} = 36.75 + 0.655 \cdot 100.6 = 80.57$
$\hat{p} = \frac{\hat{\sigma}_{XY}}{\sqrt{\hat{\sigma_X}\hat{\sigma_Y}}}= \frac{65.94}{\sqrt{100.6\cdot 80.57}} = 0.74$
例題
- 問 29.1
-
(1)⑤
高得点のデータが存在しないので$\bar{\mu}$は低くなる
同様に高得点のデータが存在しないので一部データのみとなり、$r_A$も過少評価となる
-
(2)②
回帰のデータにより$\mu$は偏り無く推定するが、分散を考慮できないので、$r$は過大評価する
-
30章 モデル選択
例題
-
問 30.1
最大対数尤度:
$$
L = \big(\frac{2\pi S_e}{n} \big)^{-\frac{n}{2}} e^{-\frac{n}{2}}
$$$S_e$:残差平方和
-
問 30.2
最大尤度 $L = \big(\frac{2\pi S_e}{n} \big)^{-\frac{n}{2}} e^{-\frac{n}{2}}$
AIC:$-2 \log L + 2k$
BIC:$-2 \log L + k \log n$
より
説明変数無しの場合のAIC\displaylines{ \log\ L = -\frac{20}{2} \log \big(\frac{2 \pi S_e}{20} \big) - \frac{20}{2}\\ = -10 (\log 2\pi + \log S_e - \log 20)- 10\\ = -10(1.838 + 10.881 - 2.996) - 10\\ = -10 \cdot 9.723 -10\\ =-107.23 }
\displaylines{ -2 \log\ L\ + 2k\\ =-2 \cdot -107.23 + 2\cdot 2\\ = 214.46 + 4 = 218.46 }
説明変数無しの場合のBIC
\displaylines{ -2 \log L + k \log n\\ = -2 \cdot -107.23 + 2 \log 20\\ = 214.46 + 2\cdot 2.996\\ = 220.452 }
その他それぞれを計算すると
説明変数 AIC BIC なし 218.46 220.45 x1 215.50 218.49 x2 209.24 212.23 x3 214.98 217.97 x1-x2 210.84 214.82 x1-x3 208.36 212.34 x2-x3 208.86 212.84 x1-x2-x3 208.46 213.44 よって選択されるのは
AIC:$x_1, x_3$
BIC:$x_2$
31章 ベイズ法
例
-
例1
- ベータ分布
期待値:$\frac{a}{a+b}$
モード:$\frac{a-1}{a+b-2}$
事後分布はサンプル数を$n$有効数を$x$として
$a$:$a+x$
$b$:$b+n-x$
上記から、
$(2+x, 20+20-x) = (2+x, 40-x)$
- ベータ分布
-
例2
例1の設定より$(2+x, 40-x)$-
(1)
ベイズ推定量は期待値なので、ベータ分布の期待値
\displaylines{ \frac{a}{a+b} = \frac{2+x}{42} }
-
(2)
\displaylines{ \frac{x}{20}\cdot \frac{20}{42} + \frac{1}{11}\cdot \frac{22}{42}\\ = \frac{2+x}{42} }
-
(3)
MAP推定量はモードなので
\displaylines{ \frac{a-1}{a+b-2} = \frac{1+x}{40} }
-
-
例3
ガンマ・ポアソンモデルの推定-
(1)
$\sum x=20$
$n= 5$
ガンマ分布のパラメータ$G(a, 1/b)$
事後分布の推定方法は
$a$:$a+\sum x$
$b$:$b + n$
より
$G(a+\sum x, \frac{1}{b+n}) = G(22, 1/6)$
-
(2)
ガンマポアソンモデルのベイズ推定量はガンマ分布の期待値
$$
a/b = 22/6 = 11/3
$$ -
(3)
ガンマ分布のモードとなるので
$$
(a-1)/ b = 21/6 = 7/2
$$
-
例題
-
問 31.1
-
(1)
Aくんのベイズ推定量
事前分布は一様分布$U(0,1)$→ベータ分布にすると$Be(1,1)$
10回くじを引いて3回あたりが出たので
$a$:$a+ x =1+3=4$
$b$:$b+n-x = 1+10-3 =8$
事後分布は$(4, 8)$となる
ベイズ推定量は
\displaylines{ \frac{a}{a+b} = \frac{4}{12} = \frac{1}{3} }
-
別解
$U(0,1)$の確率密度関数$f(x) = \frac{1}{b-a} = 1$
\displaylines{ 事後分布 = \frac{尤度 \cdot 事前分布 }{周辺尤度(省略してOK)}\\ =f(x|\theta)\cdot \pi }
$f(x|\theta)$において、xは$Bin(n(=10), \theta)$に従う
3回あたりが出ているので
{}_{10}C_{3} \theta^3 (1-\theta)^7
ベータ分布の$x^{a-1} (1-x)^{b-1}$と一致しているため
$Be(4,8)$
-
(2)
Bさんの事前分布は$Be(1,4)$
$a$:$a+x = 4$
$b$ $b+n-x = 4 + 10 -3 =11$
よって$(4, 11)$となる
モードは
\displaylines{ \frac{a-1}{a+b-2} = \frac{4-1}{4+11-2}\\ = \frac{3}{13} }
-
-
問 31.2
-
(1)
事前分布 $G(2,1/3)$
$\sum x = 15$
$n =5$
より
事後分布 $G(17, 1/8)$
確率密度関数は
$f(x) = \frac{b^a}{\Gamma(a)}x^{a-1} e^{-xb}$ より
\displaylines{ \frac{8^{17}}{\Gamma(17)}x^{17-1} e^{-8x}\\ = \frac{8^{17} \lambda^{16} e^{-8\lambda}}{\Gamma(17) } }
-
(2)
\displaylines{ \int^{\infty}_0 f(x_6|\lambda) \pi(\lambda | x) d\lambda\\ }
$f(x_6|\lambda)$ :6日目の確率密度関数PDF→$\frac{\lambda^{x_6} e^{-\lambda}}{x_6!}$
$\pi(\lambda | x)$ :5日までのPDF→ $\frac{8^{17} \lambda^{16} e^{-8\lambda}}{\Gamma(17) }$
\displaylines{ \int^{\infty}_0 \frac{\lambda^{x_6} e^{-\lambda}}{x_6!}\frac{8^{17} \lambda^{16} e^{-8\lambda}}{\Gamma(17) }d\lambda\\ = \int^{\infty}_0 \frac{8^{17}\cdot \lambda^{x_6+16}\cdot e^{-9\lambda}}{x_6!\ \Gamma(17)}\\ = 8^{17} \frac{\int^{\infty}_0\lambda^{x_6+16}\cdot e^{-9\lambda}}{x_6!\ \Gamma(17)} }
$\int^{\infty}_0\lambda^{x_6+16}\cdot e^{-9\lambda}$の計算
以下2パターンのどちらかに変換できる
- $\Gamma(a) = \int^{\infty}_0x^{a-1}\cdot e^{-x}dx$
- $\frac{\Gamma(a)}{b^a} = \int^{\infty}_0x^{a-1}\cdot e^{-bx}dx$
よって
\displaylines{ \int^{\infty}_0\lambda^{x_6+16}\cdot e^{-9\lambda}\\ =\frac{\Gamma(x_6 + 17)}{9^{x_6 + 17}} }
\displaylines{ \frac{8^{17} }{x_6!\ \Gamma(17)} \frac{\Gamma(x_6 + 17)}{9^{x_6 + 17}}\\ = \frac{\Gamma(x_6 + 17)}{x_6!\ \Gamma(17)}(\frac{8}{9})^{17} (\frac{1}{9})^{x_6} }
この形は負の二項分布
-
-
問 31.3
-
メトロポリス法
「提案して、受け入れるか決める」方式のサンプリング
メトロポリス法は、ある確率分布 $p(x)$ からサンプルを得る方法
-
手順
-
初期値 $x_0$ を決める。
適当な値からスタートする。
-
新しい値 $x'$を提案する。
- 提案分布$q(x' | x)$を用いて、新しいサンプル候補 $x'$を生成する。
- 例えば、現在の $x$に正規分布の乱数を加えるなど。
-
受け入れるか決める。
- 受け入れ率(受理確率) $A$を計算:
$A = \min\left(1, \frac{p(x')}{p(x)}\right)$ - 確率$A$に基づいて、次のステップで$x'$に移るか、現在の$x$を維持する
- 受け入れ率(受理確率) $A$を計算:
-
繰り返す。
- ステップ2から3を繰り返し、十分な回数サンプリングを行う
-
-
-
(1)
メトロポリス・ヘイスティング法の採択率$C$の導出
採択率:$\alpha(x^{(t)},y) = \min \big[ 1, \frac{\pi(y)}{\pi(x^{(t)})}\big]$
$N(\mu,1)$ は確率密度関数$\phi(y-\mu)$と書けるので
\displaylines{ \frac{\pi(y)}{\pi(x^{(t)})}\\ = \frac { 1/8 \phi(y-(-5)) + 3/4\phi(y-0)+ 1/8\phi(y-5) } { 1/8 \phi(x^{(t)}-(-5)) + 3/4\phi(x^{(t)}-0)+ 1/8\phi(x^{(t)}-5) } }
よって
\displaylines{ \min\Big[ 1, \frac { 1/8 \phi(y+5) + 3/4\phi(y)+ 1/8\phi(y-5) } { 1/8 \phi(x^{(t)}+5) + 3/4\phi(x^{(t)})+ 1/8\phi(x^{(t)}-5) } \Big] }
-
(2)
ステップ幅が小さいと他の山に移動しにくく、幅が大きすぎると採択率が低くとどまる
今回は3つの分布の期待値幅が10なので0.1〜10に連れて安定すると考えられる
- ア:その次に安定しているので1
- イ:正規分布の期待値幅と同様に安定して分散しているので10
- ウ:採択率が低くとどまっていることが多いので$100$
- エ:分散が最も低いので$0.1$
-
(3)
初期値の影響が大きく乱数が安定していない状態のため
-
-
問 31.4
-
ギブスサンプリング
「変数を1つずつ更新する」方式のサンプリング
ギブスサンプリングは、複数の変数 $x_1, x_2, ..., x_n$ を持つ確率分布からサンプリングする方法です。特に、条件付き分布が簡単に扱える場合に便利。
-
手順(2変数の場合)
- 初期値$(x_1^{(0)}, x_2^{(0)})$を決める。
-
1つの変数を条件付き分布からサンプリングする。
- 例えば、$x_1$を$p(x_1 | x_2)$ からサンプリング。
-
もう1つの変数を更新する。
- $x_2$ を $p(x_2 | x_1)$からサンプリング。
-
繰り返す。
- $(x_1, x_2)$を何度も更新し、十分にサンプルを集める。
-
ポイント
- 1つの変数ずつ確実に更新できるため、条件付き分布が簡単に計算できる場合に便利
- 例えば、ベイズ統計のように「条件付き分布が解析的に求まる」場合に向いている
- ただし、変数間の依存が強い場合、収束が遅くなることがある
-
交互にサンプリングするので条件付き確率を求める
$\mu = (0, 2)^T$, $\Sigma=
\begin{bmatrix}
1 & 0.5 \
0.5 & 1
\end{bmatrix}$ より$E[X]=0,\ V[X]= 1,\ E[Y] = 2,\ V[Y] = 1$
相関係数 $p=\frac{Cov(x,y)}{\sqrt{1\cdot 1}}=0.5/1 =0.5$
-
(1)
条件付き期待値の公式:$E[x|Y=y] = E[X] + p(\sigma_x/\sigma_y)(y^{(t)} - E[Y])$ より
\displaylines{ E[X|Y=y]= 0 + 0.5(1/1)(y^{(t)}-2)\\ = 0.5(y^{(t)}-2)\\ = \frac{y^{(t)}}{2}-1 }
条件付き分散の公式:$V[X|Y=y] = \sigma_x^2 (1-p^2)$ より
\displaylines{ V[X|Y=y]= 1(1-(0.5)^2)\\ = 1-0.25 = 0.75 }
よって$N(\frac{y^{(t)}}{2}-1, 0.75)$
-
(2) $x^{t+1}$を所与とした時、
\displaylines{ E[Y|X=x] = 2+0.5(1/1)(x^{(t+1)} - 0)\\ = 2+\frac{1}{2}x^{(t+1)} }
\displaylines{ V[Y|X=x] = 1(1-(0.5)^2)\\ = 0.75 }
よって$N(2+\frac{1}{2}x^{(t+1)}, 0.75)$
-
32章 シミュレーション
例
-
例1
半径1の円の面積は$\pi r^2 = \pi 1^2 = \pi$$U^2+V^2 \leq 1$ となる確率は$\frac{\pi/4(円の面積の1/4)}{1(正方形の面積)}=\frac{\pi}{4}$ より
$U^2+V^2 \leq 1$となる個数$M$は二項分布$Bin(N, \pi/4)$に従う
二項分布の分散は
$V[X] = p(1-p)/N = \frac{\pi/4(1-\pi/4)}{N}$
標準偏差が0.01となるのは
\displaylines{ 4\sqrt{\frac{\pi/4(1-\pi/4)}{N}} = 0.01\\ 4\sqrt{\frac{\pi(4-\pi)}{4^2 N}}\\ = \frac{1}{\sqrt{N}}\sqrt{\pi (4-\pi)}=0.01\\ \sqrt{N} = 100\sqrt{\pi(4-\pi)}\\ N=10000\cdot 3.14(4-3.14)\\ = 27004... }
よって27000以上のNが必要
-
例2
確率変数$X$の累積分布関数(CDF)を$F(x)$とする
$U \sim \text{Uniform}(0,1)$の乱数を使って、$X$の乱数を生成する方法は以下の手順
-
一様乱数の生成:$U \sim \text{Uniform}(0,1)$
-
$X$ の累積分布関数を求める
指数分布:$2e^{-2t}$から累積分布関数を求める
\displaylines{ F(x) = \int^{x}_0 2e^{-2t}dt\\ = [-e^{-2t}]^x_0\\ = -e^{-2x} + 1\\ = 1-e^{-2x} }
-
逆関数を求める
$x=$の形にして$X$を求める
\displaylines{ F(x) = 1-e^{-2x}\\ e^{-2x} = 1-F(x)\\ \text{対数を取る}\\ -2x = \log(1-F(x))\\ x = -\frac{1}{2}\log(1-F(x)) }
-
$F(x)$に一様乱数の$U \sim \text{Uniform}(0,1)$を代入
$X = -\frac{1}{2}\log(1-U)$
$1-U=U$なので
$X = -\frac{1}{2}\log(U)$ となる
-
-
例3
ベータ分布$Be(2,2)$に従う乱数を採択棄却法で生成する-
確率密度関数
\displaylines{ f(x) = 6x(1-x), \quad 0<x<1 }
-
正規化されていない関数$l(x)$
$l(x) = x(1-x)$
-
正規化定数
$c=\int^1_0 x(1-x) dx = 1/6$
-
提案分布$g(x)$とスケール係数$M_1, M_2$
提案分布:一様分布$g(x) = 1$
-
スケール係数
$M_1=1$
$M_2 =1/4$
-
-
手順
- 一様乱数$U(0,1)$から生成した$x$を取得
- $r= x(1-x)/M_i \cdot g(x)$として$x$とは独立な一様乱数$U(0,1)$から生成した$u$と比較
- $u\leq r$なら$x$を出力(採択)、小さければ何も出さない(棄却)
- これを繰り返してサンプリングすると$x$は確率密度$x(1-x)$に収束する
-
$M_i$効率
$x(1-x)$の極値から最大値を取得する
微分する
\displaylines{ \frac{\partial}{\partial x} x(1-x)\\ = 1(1-x) + x(-1)\\ 1-2x }
これを0として
\displaylines{ 1-2x = 0\\ x = 1/2 }
$x=1/2$を代入して
\displaylines{ 1/2(1-1/2) = 1/4 }
最大値は$1/4$ となる
採択率が高いほうが効率が良い=閾値が近いほうが良く、それは面積比で決まる
$l(x)$の面積は問題より$1/6$なので
$M_1$の面積比:$\frac{1}{6\cdot 1} =1/6$
$M_2$の面積比:$\frac{1}{6}\cdot \frac{4}{1} = 2/3$
$M_2$が4倍効率が良い
-
-
例4
モンテカルロ積分:乱数を用いて積分を近似する方法$I=\int^5_2 e^{-x} dx$
ただし、一様乱数はU(2,5)を使うので密度は1/3となる
そのため $I$ は$I = 3\int^5_2 e^{-x}/3 dx$ としておく
$X_i$ の確率密度関数は
$f(x) = 1/3$ , $2\leq x \leq 5$
積分 $I$ を期待値の形で表すと
$$
I = 3\int^5_2 e^{-x}/3 dx= 3E[g(X)]
$$ここで $g(X) = e^{-x}/3$
$\hat{g} = \frac{1}{m}\sum^m_{i=1}g(x_i)$ とする
$3\hat{g}$ を推定量として出力する
-
例5
理解が追いついてないので詳細なし
-
例6
理解が追いついてないので詳細なし
-
例7
-
ブーストラップ法:元のデータから復元抽出したサンプルを多数作成し、統計量の分布や信頼区間を推定する方法
$x_1, …x_n$ のデータがあるとき、その中から復元抽出でランダムに取得した標本をブートストラップ標本といい、 $i$ 番目のブートストラップ標本を $x^{*(b)}_i$ とする
$x^{*(b)}_i$を用いて計算した統計量の $\hat{\theta}^*(b)$ を用いて標本誤差の推定は以下にて求められる
\displaylines{ \hat{se}_B = \sqrt{\frac{1}{B-1}\sum^B_{b=1}(\hat{\theta}^*(b)-\bar{\hat{\theta^*}})} }
-
ブートストラップ標本:$x^{*(b)}_i$
ブートストラップ標本を元にした推定量:$\hat{\theta}^*(b)$
-
ブートストラップ推定量:ブートストラップ法によって得られた推定値
$\bar{\hat{\theta^*}} = \frac{1}{B}\sum^B_{b=1} \hat{\theta}^*(b)$
-
ブートストラップバイアス:ブートストラップで推定した統計量期待値と、元のサンプルの統計量との差
推定量 $\hat{\theta}$ のバイアス: $Bias = E[\hat{\theta}^*] - \hat{\theta}$
-
バイアス修正後の推定値:バイアスを補正して得られる推定値
元の推定量 $\hat{\theta}$ から、ブートストラップバイアスを引いて補正する
$\hat{\theta}_{corrected} = \hat{\theta}- Bias$
標準誤差の推定でBは50〜200程度で問題無い
-
-
回答
-
(1)
標本平均のバイアス補正後の推定値は
$11.083 - 0.021 = 11.062$
-
(2)
標本分散:13.538
標本平均の標本誤差 $\hat{se}(\bar{x})=\sigma/\sqrt{x}$ より
\displaylines{ \sqrt{\frac{\sigma^2}{n}}\\ = \sqrt{13.538/12} = 1.0621... }
$\hat{se}_B/\hat{\bar{x}} = 1.010/1.0621 = 0.951$
-
-
-
例8
-
ジャックナイフ法:データの一部を順番に除外しながら統計量を再計算し、その分散やバイアスを推定する手法
$x_1, …x_n$ のデータがある時、$i$ を除いたパラメータ推定量 $\theta$ を$\hat{\theta}_i$ とする
-
ジャックナイフでの推定量の平均を $\hat{\theta}_.$ :
\displaylines{ \hat{\theta}_. = \frac{1}{n}\sum_{i=1}^n \hat{\theta}_i }
-
推定量のバイアス:
\displaylines{ \hat{bias} = (n-1)(\hat{\theta}_. - \hat{\theta}) }
-
ジャックナイフ推定量 $\hat{\theta}_{jack}$ :
\displaylines{ \hat{\theta}_{jack} = \hat{\theta}-\hat{bias} }
- 推定量の標準誤差のジャックナイフ推定量:
\displaylines{ \hat{se}_{jack} = \sqrt{\frac{n-1}{n}\sum^n_{i=1}(\hat{\theta}_i - \hat{\theta}_. )^2} }
-
-
回答
-
(1)
それぞれを計算すると
番号 死者数 iを除く11個の平均 $\bar{\hat{\theta}_.}$ 1 11 11.09090909 11.08333333 2 9 11.27272727 3 9 11.27272727 4 12 11 5 5 11.63636364 6 8 11.36363636 7 7 11.45454545 8 12 11 9 13 10.90909091 10 14 10.81818182 11 18 10.45454545 12 15 10.72727273 -
(2)
\displaylines{ \hat{se}_{jack} = \sqrt{\frac{n-1}{n}\sum^n_{i=1}(\hat{\theta}_i - \hat{\theta}_. )^2} }
より
\displaylines{ \sqrt{(12-1) s^2} = \sqrt{11 \cdot 0.10256} = 1.062148.. }
-
-
例題
-
問 32.1
-
(1)
$V[x] = E[x^2] - (E[x])^2$ より
\displaylines{ E[x^2]= \int^1_0 1-u^2 du \\ = \int^1_0 1 - \int^1_0 u^2\\ =[u]^1_0 - [\frac{1}{3}u^3 ]^1_0\\ = 1- \frac{1}{3} = \frac{2}{3} }
\displaylines{ V[\sqrt{1-u^2}]= \frac{2}{3} - (\frac{\pi}{4})^2\\ = 0.666... - 0.6262..\\ = 0.05.. }
$\hat{\pi} = 4\cdot \frac{1}{n} \sum^n_{i=1} \sqrt{1-U^2}$ より標準偏差 $SD[\hat{\pi}]$ は
\displaylines{ V[\hat{\pi}] = V[4\cdot \frac{1}{n} \sum^n_{i=1} \sqrt{1-U^2}]\\ = 4^2 \frac{1}{n^2}\sum^n_{i=1} V[\sqrt{1-U^2}]\\ = 4^2 \frac{1}{n}\cdot 0.05 }
\displaylines{ SD[\hat{\pi}] = \sqrt{4^2\cdot \frac{1}{n} 0.05} \\ = \frac{4\sqrt{0.05}}{\sqrt{n}} }
こちらが0.01の場合を考える
\displaylines{ \frac{4\sqrt{0.05}}{\sqrt{n}} =0.01\\ \frac{16\cdot 0.05}{n} =0.0001\\ n = 8000 }
-
(2)
$V[\hat{\pi}] = V[4\cdot \frac{1}{4}\sum^{n}_{i=1} \sqrt{1-U^2_i} - (1/2 - U_i)]$
\displaylines{ V[\hat{\pi}] = V[4\cdot \frac{1}{n}\sum^{n}_{i=1} \sqrt{1-U^2_i} - (1/2 - U_i)]\\ \frac{4^2}{n^2}\sum^{n}_{i=1} V[\sqrt{1-U^2_i} - (1/2 - U_i)] }
ワークブックp313モンテカルロ積分の公式より
$V[x]= \sigma^2 = \int^b_a (x-E[X])^2 f(x) dx$
を用いて\displaylines{ V[\sqrt{1-U^2_i} - (1/2 - U_i)]\\ = \int^1_0\Big[(\sqrt{1-U^2_i}-(1/2-U_i)) - E[\sqrt{1-U^2_i}-(1/2-U_i)] \Big]^2 \cdot \frac{1}{1-0}dx\\ = \int^1_0\Big[(\sqrt{1-U^2_i}-(1/2-U_i)) - E[\sqrt{1-U^2_i}]+E[(1/2-U_i)] \Big]^2 \cdot \frac{1}{1-0}dx }
$E[\sqrt{1-U^2_i}]$ :(1)より$4/\pi$
$E[1/2 - U_i]$ :問題より0
こちらを用いて
\displaylines{ \int^1_0\Big[(\sqrt{1-U^2_i}-(1/2-U_i)) - E[\sqrt{1-U^2_i}]+E[(1/2-U_i)] \Big]^2 \cdot \frac{1}{1-0}dx\\ = \int^1_0\Big[(\sqrt{1-U^2_i} +U - \frac{1}{2}) - \frac{\pi}{4} - 0 \Big]^2\\ = \int^1_0\big(\sqrt{1-U^2_i} + u -\frac{1}{2} - \frac{4}{\pi}\big)^2\\ \text{問題より}\\ = 0.0144 }
\displaylines{ V[\hat{\pi}] = \frac{4^2}{n^2}\sum^n_{i=1} 0.0144\\ = \frac{16}{n} 0.0144 }
標準偏差が0.01になるには
\displaylines{ \sqrt{\frac{16}{n} 0.0144} =0.01\\ \frac{16}{n} 0.0144 = 0.0001\\ n = 16\cdot 0.0144 \cdot 10000\\ = 2304 }
-
おわりに
いかがでしたでしょうか。
こちらで一部を除くワークブックの全章について例題解答をまとめました。
ワークブックの解答に抜けがちだった途中計算の公式などもなるべく記述するようにしていますので皆さんの参考になれば幸いです。