2018 年 5 月 8 日 @ 統計モデリング 統計モデリング 第四回配布資料 文献 : a) A. J. Dobson and A. G. Barnett: An Introduction to Generalized Linear Models. 3rd ed., CRC Press. b) H. Dung, et al: Monitoring the Transmission of Schistosoma japonicum in Potential Risk Regions of China, 2008-2012 : Int. J. Environ. Res. Public Health vol. 11(2014), no.2, 2278-2287. c) X. N. Zhou, et al: Epidemiology of Schistosomiasis in the People's Republic of China, 2004. Emerging Infectious Diseases, vol. 13 (2007), no.10, 1470-1476. 配布資料の一部は以下からもDLできます. 短縮 URL http://tinyurl.com/lxb7kb8 担当 : 田中冬彦
第四回ベイズ統計 ( 導入 ) 今後の予定 ガンバ大阪ホームページより http://www2.gamba-osaka.net/stadium/ 第六回数値技法 (1) (2) ( M ) θ, θ,, θ のヒストグラム Frequency 0 5 10 15 sample histogram M distribution 0.0 0.1 0.2 0.3 0.4 π ( θ x) の関数形 posterior distribution 第五回ベイズファクター Google map から転載 -2-1 0 1 2 x -3-2 -1 0 1 2 theta 第七回グループ発表 1 第八回線形モデルのベイズ解析
今日の内容 1. 問題編 ( 中国長江流域架空の都市を題材 ) 2. 統計の復習 1~ ベイズの定理 ( 事象編 ) 3. 統計の復習 2~ ベイズの定理 ( 確率変数編 ) 4. ベイズ統計の基礎 5. 解決編 6. 補足 7. グループワーク ( 計算課題 )
本日の主役 二項モデル + 事前分布 X ~ Bin( n, θ ) θ ~ π ( θ )
2000 年代初頭, 中国都市部の急速な発展とは対照的に長江流域の農村地区では水道さえ整備されていなかった 2004 年の調査では中国内で少なくとも約 70 万人と推定 割合は数 % だが 長江流域では村全体が感染することも普通に起きている風 画像はイメージです
Google map から転載
上流 下流 今回 新たに調査に入った F 村 住血吸虫症が蔓延している村 ( 感染率 > 0.5) 住血吸虫症が蔓延していない村 ( 感染率 < 0.5)
上流 下流 F 村 ( 数百人程度 ) 考えてみよう! 住血吸虫の感染率を早期に判断する必要がある 診断用資材は十分高い精度があるが 10 人分しかない F 村の感染率を推測する上で10 人をどうやって選べばよいか?
上流 下流 感 感 非 感 非 非 感 感 感 感 調査で判明したこと 10 人中, 7 人が住血吸虫に感染! F 村全体での感染率 > 0.5 ( 蔓延 ) と考えてよいか?
θ : F 村の感染率 H 0 ( 帰無仮説 ) : θ 0.5 H ( 対立仮説 ) : θ > 0. 5 1 検定統計量 T = 統計的仮説検定 復習 ˆ θ θ0 θˆ = 0.7 θ 0 = 0.5 ˆ(1 θ ˆ) θ / N N =10 検定方式 T T ( 有意水準 0.05) >1.65 1.65 蔓延していると考えてよい 蔓延しているとまではいえない 0.7 0.5 検定統計量の値を計算すると T = = 1. 38 0.7(1 0.7) /10
上流 下流 F 村 住血吸虫症が蔓延 ( 感染率 > 0.5) 住血吸虫症が蔓延してない ( 感染率 <0.5) 10 人中, 7 人が住血吸虫に感染ただし, 蔓延とはいえない (N=10) しかし, F 村周辺の似たような村の 80% が蔓延 ( 感染率 > 0.5) が既に判明している
考えてみよう 問題提起 他の村の 80% で蔓延, この現地の情報を無視した解析はよいのだろうか. 仮定 : 近隣地域の感染率は, 十分高い精度で判明. 住血吸虫症が蔓延している村 ( 感染率 > 0.5) の割合 80% は十分, 信頼できる情報とする.
2. 統計の復習 1 ~ ベイズの定理 ( 事象編 )
ここでのポイント ベイズ統計の根本 = 条件付き確率を用いた推論
条件付き確率とベイズの定理 条件付き確率事象 Aが起きた時にBも起きている確率 Pr( A B) Pr( B A) Pr( B A) : = = Pr( A) Pr( A) A,Bが独立の場合には Pr( B A) = Pr( B) ベイズの定理 ( ベイズの公式 ) Pr( A B) Pr( B) Pr( B A) = Pr( A) 数学的には下の定義の書き換えにすぎない Pr( A B) = Pr( B A) Pr( A) = Pr( A B) Pr( B)
モデリングを考える上での注意点 以上は 数学的な定義 の仕方 練習してみよう! 実際には, 条件付き確率 Pr( B A) を先に考えることも 練習 犯人が犯行後に 犯行現場にやってくる確率 90% Pr( B A) = 0.90 A, B はどのような事象と解釈できるか. A B
結果を予想してみよう! 例題 : がん診断 * ( 以下は架空のものです ) がんの有無を95% の確率で判別できる診断法があります 検査を受ける人の中でがんである割合は年間 0.5% Aさんの診断結果は陽性でした Aさんの正しい対処方法は? 予想される選択肢 : 1.95% でがんだから, 家族と今後について話し合う 2. 所詮は半分半分 3. 統計的にはがんの人は 0.5% 程度だろ? * 松原望 : 入門ベイズ統計, 東京図書
条件付確率の計算例 (1/3) 1. 診断方法が 95% の精度 これは 条件付き確率 で表現される!! A: 診断で陽性 ; ~A: 診断で陰性 B: がん ; ~B: がんでない Pr( A B) = 0.95, Pr(~ A B) = 0.05, Pr(~ A ~ B) = 0.95, Pr( A ~ B) = 0.05. * 記法 : 補集合 c ~ A は Aの補集合 ( A とかくことが多い )
条件付確率の計算例 (2/3) 2. 検査を受ける人ががんである割合 これは 確率 で表現!! Pr( B) = 0.005, Pr(~ B) = 0.995. 3. 検査を受けて陽性が出る確率 Pr( A) = Pr( A B) + Pr( A ~ B) = Pr( A B) Pr( B) + Pr( A ~ B) Pr(~ B) 0.0545 = ここの計算は時間の都合でとばします.
条件付確率の計算例 (3/3) 4. 陽性が出た時に がんである確率 Pr( B A) = Pr( A B) Pr( B) Pr( A) ここでベイズの公式を用いる! = 0.087
例題 : 迷惑メールフィルタの仕組み モデル化の例 : 以下の条件を条件付き確率で表してみよう. 設定 練習してみよう! 通常メールと迷惑メールの受信比率は 90:10 迷惑メールに特徴的な単語としてアダルトが81% で本文に入る 通常メールでもたまにアダルトが本文に入る(1%) A: 本文にアダルトが含まれる ; ~A: 含まれない B: 迷惑メール ; ~B: 通常メール Pr( A B) = 0.81 Pr( A ~ B) = 0.01.
練習してみよう! 迷惑メールの比率 Pr(B) = Pr(~ B) = 本文にアダルトが含まれる比率 Pr(A) = Pr( A B) Pr( B) + Pr( A ~ B) Pr(~ B) = 本文にアダルトが見つかった場合, 迷惑メールである確率 Pr( B A) =
ここまでのまとめ ベイズの定理 ( ベイズの公式 ) Pr( B A) = Pr( A B) Pr( B) Pr( A) 応用例 : 迷惑メールフィルタ 実際には 1 つの単語のみで判断するのは難しいが 複数の単語 ( や他の条件 ) を組み合わせることで確率は上がっていく 機械的に判断できる条件で確率が高いものを迷惑メールとみなして別のフォルダに振り分ける 機械学習 人工知能 (AI) などの基礎にもなっている!
3. 統計の復習 2 ~ ベイズの定理 ( 確率変数編 )
連続確率変数の独立性 確率変数の独立性 復習 2 つの確率変数 X,Y の同時確率密度 ( 結合確率密度 ) p( x, y) p( x, y) 0, p( x, y)dxdy = 1 2 つの確率変数 X,Y が独立 p ( x, y) = p( x) p( y) p( x) 0, p( x)dx = 1 p( y) 0, p( y)dy = 1 注意 1. 確率変数 X, X, 2, も同様 1 X n p( x,, xn) = p( x1 ) p( x 1 n 2. 離散確率変数の時も同様 )
確率変数の条件付き確率 復習 連続確率変数の条件付き確率 2つの確率変数 X,Y の確率密度 p( x, y) 周辺密度 p ( x) = p( x, y) dy p( y) = p( x, y) dx 確率変数 X の Y=y での条件付き確率密度 p( x y) で定まる p( x, y) p ( x y) = p( y) p( y) X Y (X, Y) の確率分布を条件付き確率などという = p( x, y) dx
確率変数のベイズの定理 復習 連続確率変数のベイズ定理 ) ( ) ( ) ( ) ( x p y p y x p x y p = ここで分母の周辺密度は, 次のように計算できる. = = y y p y x p y y x p x p d ) ( ) ( )d, ( ) ( ( 数理以外の人には ) 抽象的でわかりづらいですが 後でグループで計算します.
4. ベイズ統計の基礎
ここでのポイント ベイズ統計の基本 = 統計モデルのパラメータに 確率分布を設定
統計モデルに基いた分析 第二回 データの統計分析 1. データに応じた統計モデルの設定 ( 母集団分布のモデル化 ) i. i. d. X1,, X n ~ p( x θ ) θ 2. パラメータの推測 点推定 区間推定 ( 信頼区間 ) 仮説検定 ベイズ統計もこの流れは同じ!
パラメータの推測における不確実性 分析者の気持ち prob. 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 θ 0.0 0.2 0.4 0.6 0.8 1.0 theta 0.7 prob. 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 分析者が強い確信を持っている場合 ( データが大量 ) 0.0 0.2 0.4 0.6 0.8 1.0 theta θ 0.7 分析者は左図のような確率分布で表現したい!
パラメータの確率分布の導入 統計モデルは所与とする.( 簡単のため, サンプルサイズ n= 1.) x ~ p( x θ ) ベイズ統計ではパラメータにも確率分布を設定 ( 設定方法は後で ). 事前分布 (*) π ( θ ) θ ~ π ( θ ) dθ = 1, π ( θ ) 0 パラメータに対する不確実性を確率分布で表現. データを得る前 ( 事前 ) の分布で事前分布 (prior) *1 確率分布と確率密度関数 / 確率関数は混同して用いる π p(θ ) *2 データの分布と区別するため, を用いるが, と書いてもよい.
データ x とパラメータ パラメータの条件付き確率分布 θ の同時分布を以下で 定義 p ( x, θ ) = p( x θ ) π ( θ ) データ x のみの確率分布 ( 周辺分布 ) は ( θ ) p( x) = p( x θ ) π dθ データ x が与えられた時のパラメータ θ の条件付き分布は π ( θ x) = p( x θ ) π ( θ ) p( x) この条件付き分布を事後分布という
事前分布と事後分布 未知パラメータの事前分布 ( 分析者が設定 ) π (θ ) 事後分布 ( データ x を代入してベイズ定理から計算 ) π ( θ x) π (θ ) π ( θ x) prob. 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 0.0 0.2 0.4 0.6 0.8 1.0 条件付き分布に変化 θ prob. 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 0.0 0.2 0.4 0.6 0.8 1.0 θ theta theta
従来の統計学との違い やや学術的 従来の統計学 ( 頻度論的統計ともいう ) 統計モデルは数学的に扱いやすいものを設定 各分析 ( 推定 検定 ) に応じて公式を導出 (& 理論上はサンプルサイズ大を暗に仮定 ) ベイズ統計学 π (θ ) パラメータに 初期 の確率分布を設定 パラメータの条件付き分布 π ( θ x) に基いて一貫して考える ( 複雑なモデルでもやり方は変わらない ) ベイズ統計の根幹は条件付き分布 ( 事後分布 )!!
二項モデル 統計モデル ~ 二項モデル (* 学部 1 年統計学 B-2, C-2 の内容 ) θ 毎回, 独立に一定の確率で感染者を選んでいると考えて 10 人中 X 人の感染者が選ばれるとする. X ~ Bin(10, θ ) F 村 非 感感感 感 感 非 非 非 感感感感 感 感 感 非 感 非 非 感感感 感 感感感感 感感 非 非 非 感 感感 非 補足 : 先の統計的仮説検定は上のモデルのみで, 公式を導出していた
事前分布の設定 (1/2) 事前分布 θ 感染率は村ごとに異なるが近隣地域の情報 ( 事前情報 ) から F 村の感染率について以下のような確率分布を設定 Pr( θ > 1/ 2) = Pr( θ 1/ 2) = 0.8 0.2 住血吸虫症が蔓延 ( 感染率 > 0.5) 住血吸虫症が蔓延してない ( 感染率 <0.5)
事前分布の設定 (2/2) 離散化 ( 今回のみ ) 計算の都合上 θ = 0.0, 0.1,... 0.9, 1.0 θ π (θ ) ~ π ( θ ) = θ 0.2 6 は以下の離散値とする, θ = 0.8, θ = 5 0,0.1,...,0.5 0.6,0.7,...,1.0 prob. 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 Prior 0.0 0.2 0.4 0.6 0.8 1.0 theta > theta <- (0:10)/10 > prior <- c( rep(0.2/6, 6), rep(0.8/5, 5)) > plot(theta, prior, pch=3, col=3, lwd=5, xlab="theta", ylab="prob.", main="prior", xlim=c(0,1), ylim=c(0, 0.35)); 以上のモデルで θ の条件付確率 ( 事後分布 ) を数値的に求める
1. 事後分布 図の青点 2. 事後分布を最大にする点 (MAP) ˆ = 0.7 θ π 3. 事後分布による蔓延の確率 Pr( θ > 1/ 2 x = 7) 計算結果 = π ( θ = 0.6 x = 7) + π ( θ = 0.7 x = 7) + π ( θ = 0.8 x = 7) + π ( θ = 0.9 x = 7) + π ( θ = 1.0 x = 7) = 0.95 prob. 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 Prior, Posterior 0.0 0.2 0.4 0.6 0.8 1.0 > X_obs <- 7; N <- 10; like <- dbinom(x_obs, N, theta); > poste <- like*prior/(sum(like*prior)); > plot(theta, prior, pch=3, col=3, lwd=5, xlab="theta", ylab="prob.", main="prior, Posterior", xlim=c(0,1), ylim=c(0, 0.35)); > points(theta, poste, pch=1, col=4, lwd=5); > sum(poste[theta > 0.5]); theta θ 周辺の村を考慮した事前分布により実感とあう蔓延確率が出た!
6. 補足
ベイズ解析 ベイズ解析 / ベイズ分析 (Bayesian Analysis) データ x 所与の下, 事後分布 π ( θ x) を用いて統計解析を行うこと ベイズ解析の強み π ( θ x) が計算できてしまえば, 統一的に統計解析が行える パラメータ推定 ( 点推定 ) 信用区間 仮説検定; 判別 モデル選択 予測分布 概念的な難しさはない ( シンプル!)
ベイズ解析の難点 問題は事後分布 π ( θ x) の計算 統計モデル 事前分布の設定 数値手法
統計モデル p( x θ ) 事前分布 メジャーなもののみ扱う π ( θ ) 本講義のスタンス 事後分布が計算しやすいクラス ( 主に共役事前分布 ) 正当性 階層構造を入れると広いクラス サンプルサイズ大なら, 事前分布の取り方に大きく依存はしない cf) 理論的な課題 事前の情報あり / 情報なし どのように事前分布に反映させるべきか?
次回以降やっていくこと いろいろな統計モデルに事前分布をいれてベイズ解析 基本的な統計モデル ( 第一回 ) 線形モデル ( 第二回 ) Y = α + β + ε i x i Y Po( µ ) n i ~ i log µ = i p j= 1 κ x j ij i i =1,, 2 ε ~ N(0, σ ) ポアソン回帰 ( 一般化線形モデル )( 第三回 ) i i =1,, n たとえば α, β,, σ 2 ;, x i κ j 適切なものに をつけよ に確率分布 ( 事前分布 ) を導入
練習問題 ( 後半にグループで取り組む提出はグループ用シートのみ )
グループで計算してみよう! 問 1: グループメンバーで以下を埋めて, お互いの理解を確認 二項モデルの確率関数 θ 毎回, 独立に一定の確率で感染者を選んでいると考えて 10 人中 x 人の感染者が選ばれる確率は p x 10 = x x 10 x ( θ ) θ (1 θ ) 事前分布の確率関数 θ 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 π (θ )
問 2 グループで計算してみよう! グループメンバーで手分けして x =3,5,7,9 のそれぞれの事後分布を求め以下の表を埋めなさい. 事後分布の確率関数 θ 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 π (θ 3) θ 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 π (θ 5)
θ 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 π (θ 7) θ 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 π (θ 9) 問 3 π ( 0.5) = 1 π ( θ ) = 0, θ 0.5 事前分布のとき 事後分布はどうなるか?