岩波データサイエンス-ベイズモデリングの世界 読書メモ

岩波データサイエンス-ベイズモデリングの世界 読書メモ

「岩波データサイエンス-ベイズモデリングの世界」 読書メモ。

はじめに
第I部　 ベイズモデリングの世界
平均値から個性へ　 統計的モデリングのひらく世界像　　 伊庭幸人

概要
できるだけ個性を取り入れる統計へ
たった一人のあなたへ
平均値ではなく、限りなく自分に近いデータを求めると最後は、自分自身になる
粒度を細かくすればするほど、個に行き着く
データのそれぞれは、それ自体では夜空にポツンポツンと光星のような美濃で、その間には無限深い闇がある
何らかの形で似たものをまとめて間を補うものとして「モデリング」がある
モデリングなしに、法則を引き出したり予測はできない
モデリングの例
図2 データの分析例
平均(代表値)を考える際に、2分割、3分割・・と考えていく?
モデリング1
直線的な近似をする
モデリング2
確率モデルで知識を表現する
yi = (場所xiでの平均値) + (ランダムなゆらぎ)
ランダムな揺らぎを「正規分布」とする
平均をm, 分散をσ2
p(Y) = p(y1, y2, …, yN)
pi(yi) = 1/√(2𝛑σ2)exp(-(yi-m)2/2σ2)
ランダムだけで表現したモデル
図2bの2つに分割したモデル
図2dの直線的な近似を表すモデル
一般解は特殊解にはならない
統計科学の仕事
モデルの開発とそれを用いてのモデリング
いろいろな確率分布や関数の形
モデルの推定・評価・利用
AIC基準によるモデル選択
モデルから生成したサンプルをモデルの評価に使える
伝統的な統計学ではこちらに主体
モデル自体が少数の正規分布と線形関数に限られていたため
「滑らかな曲線」とは
複雑な位置依存をどのように表現するか?
一つの方向性
2次、3次の多項式を考えて、AIC最小基準などで最適と思われるモデルを選ぶ
「滑らかな曲線」という概念自体の数式での表現
当てはめる曲線fそのものをある確率分布からのサンプルと考える
pi(yi|f) = 1/√2πσ2exp(-(yi-f(xi))2/2σ2)
p*(f)
fの確率密度関数
本来は無限次元
計算機上で処理するために、 多数の小区間に分割して、 分割した区間のそれぞれの上では一定値とする
例:100分割だと{f1, f2, …,f100}の100次元のベクトルがfを表す
p*(f) = Cexp(-1/2δ2(fi+1 – 2fi + fi-1)2)
2階差分が小さいことを表現
δ2が小さいと、平均は直線的に変化
δ2が大きいと任意の形をとる
あるレベルの「ぐにゃぐにゃ」の曲線がたくさん入っている「壺」からfを選ぶ
空間的、時間的なfでの適用
「個性」を表すための使い方
顧客をK個のグループに分ける
K個のグループの平均値をm={mk} (k=1, …,K)
mkが「各グループの平均値の平均」Mの周りに正規分布する
p*(m)
グループの違いを取り入れつつ、全体の情報を利用した分布
階層ベイズモデル
定数を確率変数と読み替えることで、多段式のデータ生成過程を考える
階層モデルを使った推論をどのようにすれば良いのか
δ2を与えて、まずf、次にyという順にデータが生成される
yとfの同時確率密度関数を与えたことに相当
p(y,f) = p(y|f)p*(f)
p(y|f)
確率の段数が増える
任意の確率変数の組みに関して、同時確率密度関数のp(r,s)が
p(r|s)p(s) = p(r, s) = p(s|r)p(r)
(1)
両辺を積分すると左辺は
∫p(r|s)p(s)dr = p(s)∫p(r|s)dr = p(s)
よって
p(s) = ∫p(r, s)dr = ∫ p(s|r)p(r)dr
(2)
(1) / (2)
p(r|s) = p(r, s) / ∫ p(r, s)dr = p(s|r)p(r) / ∫ p(s|r)p(r)dr
ベイズの法則
確率構造の全体が仮定されている場合、「モデルの推定・評価・利用」の大半が片付く
上記を用いて「段階上の」確率分布は
yのみが与えられた時の、fとσ2の同時分布は
全ての変数の同時確率のモデリングを行って、それにデータを入れれば自動的に答えは出る
「個性」の効用 ミクロとマクロ
階層的なモデルはどのように利用されるか?
「個性」を表すようなパラメータ、 小グループや個体や場所ごとに割り当てられてる 「ミクロなパラメータ」
個体ごと、場所ごとに将来の予測に興味がある場合
個人に対する手術の成功確率
細かく分けた集団ごとの購買行動の予測
各地点での植物の分布密度
σ2やmに相当する系全体の特徴を表すパラメータ 「マクロなパラメータ」
ミクロな構造を無視すると 推定されたマクロなパラメータに 偏りが生じる可能性がある
個体ごとに見るとはっきり存在する相関が、個体を無視するような平均操作を行うことで過小評価される例もある
科学の島々と未知の大陸をつなぐ
諸科学はコアとなる命題により構成
エネルギー保存則
DNAによる遺伝子情報のコード化
原子の存在
それぞれの領域の間に どの領域にもはっきりとは属さない 部分が残されている
問題ごとの判断による、当面の知識の最大限の活用が必要
例
生態系は河川の改修により変化したのか
明日の大売り出しには何を売れば良いのか
赤ちゃんをうつ伏せにするのは危険か
それぞれの体質に合わせた治療法を選ぶにはどうしたら良いか
地球環境の予測
大量のテキストの自動処理
熱力学的な世界観
現象を支配するマクロな変数を直接観測し、その間の法則を解明する
科学や工学の任務
「平均値」の統計学
対抗/補完するものとして「ミクロな」統計がある
memo
再生カーネル法
「平均値から個性へ」の例を試してみる
概要
JAGSやSTANのようなMCMCソフト
KFASのようなカルマンフィルターベースのソフト
人工データ
階層モデルで「個性」をとらえる　　 久保拓弥
観測データと統計モデリング
生態学での例
架空植物の単純化した架空データに見られるパターンは?
ある1個体を選んだときに、それが何個ぐらい種子を作るのか知りたい
植物は胚珠という種子のもとになる器官をどの個体も必ず10個持っている
胚芽が種子になることを結実
ある胚珠が種子になる確率を結実確率と呼ぶ
植物100個体を観察した時の趣旨数についてのデータ
100個体x10=1000個
496個が種子になる
結実確率=496/1000=0.496?
割算推定とその統計モデル
平均値計算(496/1000=0.496)とは何なのか?
架空植物の一つの個体をiとする
全ての個体で結実確率qが同じとする
個体iの10胚珠の中で結実した胚珠数がyi個となる確率は二項分布
二項分布
1,0の確率
100個全体ではf(yi|q)を100個体ぶん掛け合わせたもの
この場合の尤度は結実した全胚珠個数を全葉が個数で割ったものになる
個体差を無視したモデルの予測
上記の例での確率分布
二項分布では現象をうまく表現できていない
「どの個体でも胚珠が結実する確率qは同じ」過程が正しくない
個体iの結実種子数yiのばらつきが二項分布モデルの予測から逸脱する現象
過分散(over dispersion)
二項分布を拡張する必要がある
結実する確率qは植物個体によって異なる
個体差を考慮したモデル
結実する確率をロジスティック関数q(z)=1/(1+exp(-z))で表す
ある個体iについての結実しやすさを表す変数をzi
zi=β + αi
β:全個体共通
α:個体差
2項分布モデルの尤度方程式に常識を入れる
最尤推定は不可能
階層ベイズモデルで表現する個体差
階層ベイズモデルの役割
全100個体の個体差αiをいちいち最尤推定しなくて済ませてしまう手法
αiを確定せずに放置して、確率変数として扱う
観察データをうまく説明できる範囲で、個体たちはできるだけにている(αiがゼロにちかい)となるようにαiを決めようね
{αi}を制約する役目を与えられた確率分布を事前分布と呼ぶ
αiの確率分布には観察データ{yi}と「観察された100個体の結実確率には、全体としてどこか似ている部分がある」というルールを与える
観察データと事前分布で決まるαiやβの確率分布は事後分布
個体差αiの事前分布は簡単のため平均ゼロで標準偏差σの正規分布とする
σはこの植物の個体たちがお互いにどれくらい似ているかを表す
σがゼロに近ければどの個体もお互いに似ている
σが大きければ、αiは各個体の結実数yiに合わせるような値をとる
事前分布に対する事後分布
パラメータσをどうすれば良いか?
先送りして、単にσも何らかの確率分布h(σ)に従うとする
事前分布のパラメータの事前分布なので、超事前分布と呼ぶ
観測データ{yi}のもとでのパラメータの同時分布
分母は全ての場合の席なので定数
事後分布の確率密度は、尤度(観測データのもとでの)と事前分布、長事前分布の確率密度の積
経験ベイズ法による最尤法
階層ベイズのパラメータを推定する方法
経験ベイズ法(empirical Bayesian method)
事後分布p{β, {αi}, σ|{yi})において
全個体共通のパラメータβの事前分布gβ(β)と個体差のばらつきを表すσの(超)事前分布h(σ)を「分散がとても大きな一様分布」とする
βもσも(観察データ似合うように)好き勝手な値をとって良い
各個体の個体差αiは平均がゼロかつ標準偏差σの正規分布である事前分布gα(αi|σ)によって制約される
分母のgβ(β)とh(σ)が定数になるため
事後分布
αiについて積分した量
観測データ{yi}のもとでパラメータβとσの尤度方程式となる
一般化線形モデル(generalized linear mixed model : GLMM)というクラスのモデルと全く同じ形式
統計ソフトウェアR(glmmML)で容易に計算可能
GLMMでの分布
Markov Chain Monte Carlo(MCMC)法
経験ベイズでは不可能なより複雑なモデルに対応
個性の生態学と統計学
個体差の正体は?
観測した個体ごとに体の大きさや年齢が違う、個体ごとの遺伝子が違う
食物が育つ場所の明るい・暗いあるいは土壌中の栄養の多い少ない
「個性」とパラメータの推定　　 伊庭幸人
個体差や個人差、グループ間の差など「個性」を無視してデータをまとめるとパラメータの推定がうまくいかなる例
観測値(x,y)の組からなるデータが4つのグループ(j=1,2,3,4)からなる
それぞれのグループのデータの背後に y = ax + bj + η という直線的な関係がある
ηは期待値ゼロのランダムな雑音
傾きaは全部共通
切片bjはグループによって違う
データを観測する人は各観測値がそれぞれどのグループに族しているか知っている
データ
グループごとに見るとxとyに相関があることがわかる
グループの情報を捨てて全て混ぜてしまうと全体が団子になって、相関らしいものは見えなくなる
傾きaのような対極的なパラメータを 誤差を含めて正しく効率的に推定するためには
「グループごとに異なるbjを取り入れて、 それらを緩く関係づけて推定するモデル」 が必要
グループごとに適当な基準点 (例えばグループに所属するデータx,yの平均値)を考えて、 そこからの差をとってから、全体を混ぜて解析すれば良い
個体差・個人差・グループ差 がある時の古典的な手法
差や比などをうまくとって個体差や個人差を消してしまう
医療で同じ人の処置前と処置後の差を取ることで、個人差を消す
問題が複雑になると困難
無理に消去しようとせずに、 「個体差・個人差・グループ差を表す直接観測されない変数」 を使ってモデリング
この例ではb
導入した多数のパラメータを事前分布で縛っておいて、 後からMCMCなどで積分(周辺化)して消す
階層ベイズモデリングの考え

個人差・地域差をとりこむ統計科学　 医学分野の事例　　 丹後俊郎

はじめに
「薬が効く」≠「病気が治る」
全ての患者に一様に効くわけでない
「予測不能な個人差」が存在
個人によって変化する変量あるいは確率変数(random variable)
混合モデル(mixed-effects model)
性、年齢などのように個人によらず 効果が一定と考えることが自然な母数効果(fixed-effects)
薬の効果のように個人によって変わる変量効果(randomーeffects)
ベイズモデル(Bayes model)
全ての要因効果に確率変数(事前分布)を仮定
個人差を取り込んだシステム
検査の基準値は健常者の約95%が含まれる集団
個人の生理的変動幅は集団のそれに比較して著しく狭い
ある検診センターでま赤血球の5回の測定値をプロット
個人の変動よりも、集団内の変動が大きい
個人差が大きい
任意の個人の検査データの分布が適当な変数変換をすることで正規分布N(μi, σi2)に従う
個体差があるということ
H0: μi = μj, σi2 = σj2 の帰無仮説が否定される
線形変量効果モデル (linear random-effects model)
xijを個人iのj回目の測定値
xij = μi + εij =(μ +βi) +εij
i = 1, …,n (個人): j = 1,2, …,r(反復)
βi:個人差を示す個人iの変量効果
βi ~ N(0, σβ2) (σβ2は個人間分散)
εij:反復誤差
εij ~ N(0, σε2) (σε2は個人内分散)
集団の分散σ2
σ2 = σβ2 + σε2
個人差指数
検査項目の個人差の大きさを評価
η = σβ / σε
個人間、孤児内の平方和Vβ,Vε
個人差を取り込んだ治療効果の評価
地域差を取り込んで疾病リスクの推定
階層ベイズモデル
治療薬Progabideの有効性補評価するモデルの階層ベイズモデル
yji ~ Poisson (μji) log μji = log(T) + b0i + b1ixji2 + γxji1xji2 b0i ~ N(μλ, σλ2), b1i ~ N(μθ, σθ2)
確率変数であるパラメータの事前分布は「無情報ぶり」を示すた
γ, μλ, μθ ~ N(0, 1002) 1/σλ2, 1/σθ2 ~ Gamma(0.001, 0.001)
MCMC法で解く

全体モデルから局所モデルへ　 状態空間モデルとシミュレーション　　 樋口知之

はじめに
時間的に局在化した「局所モデル」
局所モデルをうまくつなぎ合わせたモデル
東京の気温データのモデル化
全体モデルを用いたデータ解析
少数の固定したパラメータを持つ統計モデル
局所線形モデル
時間的に局在化した情報を取り出す
局所非線形モデル
非線形への拡張
確率差分方程式
等式で表現される制約条件からの確率的なズレを許容する”柔らかな”モデル
確立的な揺らぎを生み出すノイズ項が従う分布形を、非ガウス系に一般化
ジャンプや異常値といったまれに生成する確立事象にも対応
非ガウス性
一般状態空間モデル
粒子フィルタ
一般状態空間モデルの数値的解放に提案
一般空間モデルの持つ非線形性・非ガウス性を忠実に取り扱う
計算実装の著しい簡便さ
ロボティクス、ITS(高度交通システム)、ファイナンス、マーケティングとうで実用化が進んでいる
粒子フィルタを用いたシミュレーション
シミュレーションモデルに含まれる変量で構成される位相空間
シミュレーションの数値解をプロットすると、初期値・境界条件を与えた時点では天であったものが、計算時間とともに”紐”になる
初期条件を変えると違った軌跡を位相空間に描く
シミュレーションを続けると複数のパス(解の多様性)が生じる
観測データを用いることでシミュレーションの選択肢の中から絞り込みを行える
全体モデルから局所モデルへ
東京の平均気温の経年変化
データ数をN個
Y={y1, y2, …,yN}
気温は上昇傾向
直線a・n + bを当てはめる
直線と観測値の差分項をwn
平均0、分散σ2のガウス分布
wn ~ N(0, σ2)
データセットYの確率分布
データynを平均μn = a・n+β、分散σ2のガウス分布からのN個の独立したサンプル
未知のパラメータをベクトルθ=(a,b,σ2)T
Θの関数
尤度関数
尤度関数あるいは対数尤度関数の最大化によりθの値を定める
直線の傾きが一定ではない
1980年以降は急峻
データ全体に対して一直線のモデルではなく、μnに対して”局所的にほぼ直線”を表現
μn = 2μn-1 – μn-2 + vn, vn ~ N(0, τ2)
vnは連続する3点の一直線からのズレを表すノイズ
λ = τ2 / σ2
λ=0だと直線
λ=±∞(σ2→0)では全く直線ではなくなる
状態空間モデル
システムモデル
Xn = Fn(xn-1, vn)
xn=(μn, μn-1)T, vn=(vn)
Xn:状態ベクトル
Vn:システムノイズベクトル
q(v|θsys)に従う白色ノイズ
Θsys:分布を記述するパラメータベクトル
Fn:非線形関数
観測モデル
yn = Hn(xn, wn)
yn=(yn), wn=(wn)
Wn:観測ノイズベクトル
r(w|θobs)に従う白色ノイズ
Θobs:分布を記述するパラメータベクトル
Hn:非線形関数
局所非線形モデル
局所的増減値が、都市のヒートアイランド化効果で増幅される
観測できない量ρnを導入
“ランダムウォークの概要とアルゴリズム及び実装例“でも述べている離散時間ランダムウォーク
Xn =(μn, μn-1, ρn)T
Vn =(vμ,n, vρ,n)T
逐次ベイズフィルタ
Z1:n : 最初の時刻から時刻nまでのベクトルzを全て並べた量
条件付き分布の解説
予測分布
去年までのデータを使って今年を予測
フィルタ分布
今年までのデータに基づいた今年の状態ベクトルの分布
平滑化分布
全て手元にデータがある元での今年の状態ベクトルの分布
条件推定のための漸化式の模式図
(1) 一期先予測
手元に去年のフィルタ分布p(xn-1|y1:n-1)がある
p(xn|y1:n-1) = ∫ p(xn|xn-1)p(xn-1|y1:n-1)dx
p(xn|xn-1)は一般状態空間モデルのシステムモデル
(2) フィルタリング
今年の予測分布が得られると、今年のデータが入ってきて、ベイズの定理によるフィルタリング計算で、フィルタ分布が得られる
モンテカルロ近似と粒子フィルタ
一般状態空間モデルでは条件付き分布p(xj|y1:k)はあらゆる形状を示す可能性があるため、解析関数を利用した表現は無理
状態ベクトルが高次元の場合にどのように表現するか
逐次式で出てくる状態ベクトルの次元の積分への対応も課題
超高次元のp(xj)の表現をコンピューター上で可能にしつつ、逐次更新式の実現をシンプルにするには?
条件付き分布を、 そこから得られたとみなす独立な多数の実現値 (例えば、数百〜100万個)で持って近似すれば良い
モンテカルロ近似
一つ一つの実現値を「粒子」と呼ぶ
例
予測分布 p(xn|y1:n-1)
Xn|n-1 = {xn|n-1(1), xn|n-1(2),…,xn|n-1(m)}
予測粒子
フィルタ分布 p(xn|y1:n)
Xn|n = {xn|n(1), xn|n(2),….,xn|n(m)}
フィルタ粒子
xj|k(i)
時刻jの状態ベクトル
状態ベクトルの推定に利用した観測データの最後の時刻がk
I番目の粒子
粒子近似された系での逐次更新式
規格化した尤度の確立で予測粒子をリサンプリング(復元抽出)し、得られた粒子を時刻nのフィルタ粒子とする
適合度が低いと粒子は死滅
適合度が高いと分裂して仲間を増やす
シミュレーションとデータ同化
時間的発展形式を取るシミュレーションモデルの数理モデル形式
連続時間・空間の偏微分方程式
空間上では経度緯度格子系(グリッド)
科学・物理でも格子系上で様々な変数を定義して計算

生きた言葉をモデル化する　 自然言語処理と数学の接点　　 持橋大地

はじめに
言語学者の経験と主観による仮説と検証
言葉を統計的に考える分野
計算言語学
工学的な立場からは自然言語処理
言語を統計的に捉えることによって、複雑で膨大な言語現象を計算機で自動的にモデル化
規則では捉えきれない曖昧性や例外、文脈依存性を数学てに適切に扱う
言語の統計モデル
言語は記号列
まずは言語は単語からなるとした仮定
単語の頻度には大きな偏りがある
順位と頻度は反比例
Zpifの法則
自然界の多くの離散現象に共通する累乗法則
全体でN語の文章の中で、単語iがni回現れたとすると確率は
Pi = ni / N
各単語の出現率をV次元(Vは語彙数)のベクトルで表したp=(p1,p2,…pV)
例:語彙が(w1,w2,w3)=(“花束”,”航海”,”バイアグラ”)の3つしかない場合
p1=(0.3, 0.7, 0)
p2=(0.4, 0.3, 0.2)
p3=(0.1,0.1,0.8)
広告メールはp3, 通常のメールはp1,p2のような確率分布から生成
Pの確率分布:ディリクレ分布
単語の確率分布がディリクリ分布から生まれたと考える
手順
p ~ Dir(p|α)を生成
~ :「~の確率分布に従って」
単語wi ~ p(i=1,…N)を生成
wの確率は様々なpの可能性について 積分を行なって期待値を計算する
niは単語がwの中に現れた関数
多数の文章w1,…,wDに関する確率の積
αに関して凸であり、ニュートンほうを用いて、データの確率を最大にする事前分布のパラメータαを求めることができる
Wが与えられれば、それを生んだpの確率分布はねベイズの定理からディリクレ事後分ぷとして推定
期待値は
E[pi|w] = (ni + αi ) / (N + α)
ディリクレ過程による無限離散確率分布Gの生成
ディリクレ過程を使うと未知の語の生成確率も仮定される
Nグラムモデルと無限nグラムモデル
言葉は独立に生起することはない(前述の仮定は違う)
Nグラムモデル
言葉がその前の(n-1)語の言葉に依存する(隣接したn語の間の関係を捉える)モデル
1グラムのモデル
言葉が互いに独立に生成
Nグラムモデルの、文w=w1w2・・・wTの確率は
条件付き確立の積
例:3グラムの場合
p(彼女 が 見る 夢)
= p(彼女) x p(が|彼女) x p(見る|彼女 が) x p(夢|が 見る)
言葉の(n-1)次のマルコフ過程
非常に単純なモデルではあるが、音声認識や統計的機械翻訳などで言語的に不適格な分の確率を小さくするために極めて有効
Nグラムモデルではnを増やすほど、言葉の間の関係をより緻密に捉えることができる
nを増やしすぎると、データがないために単純な推定では条件付き確率が0になってしまう
例「お魚くわえた三毛猫」
p(三毛猫 | お魚 くわえた) = n(お魚 くわえた 三毛猫)/n(お魚 くわえた) = 0
階層ディリクレ過程
1グラム分布p(*)
これを基底測度G0としたディリクレ過程による2グラム分布p(*|w1)
さらにミレを基底測度G0とした3グラム分布p(*|w2w1)が生成される
概念図
言語の場合はディリクレ過程の当てはまりは完全ではない
その拡張
2パラメータ・ポアソン-ディリクレ過程(Pitman-Yo過程)
階層Pitman-Yo過程
文脈長を隠れ変数として確率モデル化することでnグラムのn(木の階層の深さ)モ可変とする
ライダムウォーク生成文
言葉の意味の統計モデル
言葉の持つ”意味”を統計的に扱うためには?
文章によって異なる確率分布θ=(θ1,…,θk)によって表され、 ディリクレ分布から各文章ごとに生成される
K:潜在的な話題の総数(通常100程度)
LDAモデル
話題分布θを選び
文章の持つ話題分布θ~Dir(θ|α)を生成
Θに従って話題をランダムに選び
For n=1,…,N
ある話題kn~θを選択
その話題から文が生成される
ある話題knから、単語wn~p(w|kn)を生成
それぞれの文章w=w1w2…wNの確率は
ベイズにより事後確率推定
数百万〜数億語のテキストに対してサンプリングを行うことで、各単語の生成された正しい話題と、文書の話題分布を全て推定できる
川端康成の雪国でのデータ処理結果
K=100
話題58:トンネル-ガラス-信号-窓
毎日新聞テキストのLDA処理結果
参考文献
トピックモデルによる統計的潜在意味解析
自然言語処理シリーズ
トピックモデル
機械学習プロフェッショナルシリーズ
岩波データサイエンスvol2
言葉の意味とword2vec
Word2vecの単語ベクトル
テキスト中で各単語の前後の数単語の単語ベクトルから自分自身の単語ベクトルが予測される
2016年のプリンストン大学のAroraによる、 「潜在的な意味空間での乱ダムウォークによってテキストが生成された」 数統計モデルにより説明されることの示唆
A Latent Variable Model Approach to PM 1-based Word Embeddings
テキストの時刻tでの単語wtが、 その座標wtの「文脈ベクトル」Ctの近さに従って 生成されたと仮定する
p(wt|ct) = exp(wtTct) / ∑exp(vTct)
動的な線形対数モデル(ロジスティック回帰)
終わりに
言語には宇文構造があり、これを線形な単語列から(セイカイナシデ)獲得することが、研究の最先端

ポスト近代科学としての統計科学　　 田邉國士

統計学はわからない?統計学は「推論の方法」
「合理的な推論は何か?」
仮説(モデル)設定と推論の仕方に関する多様な接近法に対するパースペクティブを持つ
ニュートン・デカルトパラダイム
近代の科学的推論の手法は「仮説演繹法」に因る
仮説定立→演繹→実験という手続き
演繹主導型の推論法
(要素)還元主義
複雑な対象の性質はその構成要素の性質の総和である
各要素を調べれば全体が把握できる
対象を構成要素に分解し、その個々の要素に関するモデルを構築し検証すれば十分である
演繹推論の限定性
数学的知識の証明においては、演繹的推論のみが許される
演繹的推論が数学の絶対的確かさを保証する
数学においては「∞」や「0に収束する」といった概念を
「任意の(すべての)」といった形容詞や
Ε-δ論法などのような有限の表記(表現)と有限の操作によって
有限の世界と無限の世界
特定の数学の証明が正しく行われたことの証明は、当該の数学の外にある
数学は閉じた形式的体系ではないので、内添い的に心理であることを保証されているわけではない
帰納という原罪をてことして
演繹的手法以外の知識を獲得する方法
帰納的推論
限られた経験データから経験外の事態を予測するための汎化された知を導く
機能による推論は論理的には誤り
黒いカラスを何羽見た経験があろうと「カラスは黒い」と結論することはできない
結論の誤りは「カラスは黒い」という命題を確率により表現しても変わらない
我々はカラスの総数を知らないのみならず、見たカラスがランダムだったかを知らない
過去から「帰納推論」に「客観性」をもたらす試みがなされてきた
アリストテレスの「単純牧拳法」
ドゥンス・スコトスの「一致法」
ウィリアム・オッカムの「差異法」
R.グロステストとR.ベーコンによる「試験(テスト)」という手続き
仮説検定法
有限個のデータから、無限の事象について言及する一般命題を導出することは誤り
人間が獲得する知識
数学の場合のように定義に基づいて記号表現された対象に対する演繹的推論の操作によって得られるもの
古典物理学のの場合のように対象の先験的文節化とその記号表現を前提とする仮説演繹法によるもの
意識化され記号化できる知識
意識上に上がってこない暗黙の直感知
人間も「演繹的」ではなく「帰納的」に生きている
数理科学の中の「帰納的推論」
経験科学の「原因から結果を演繹的に推論する問題」→「純問題」
「結果から原因を訴求する問題」→「逆問題」
例:現在の温度分布を知って過去の温度分布を求めるために「熱方程式」を逆にとく
因から果への推移の間にエントロピーの増大が起こる
因果の連鎖を逆にたどる逆問題の多くは「非適切逆問題」
J.S.アダマール
因から果を導く作用素が悪条件あるいは非可逆となる場合
因果関係を方程式に表現して、因を道として、果の部分に観測データを代入して方程式を解こうとしても、誤差が生じて解が求められなくなる
実験科学やエンジニアリングでの問題にも非適切逆問題は多くある
電場の情報からその場を作る電荷分布を求める
重力場の情報からその場を作る質量の分布を求める
統計学は帰納推論である
擬似演繹
フィッシャー、ピアソン、ネイマンによる「仮説検定論」
数理統計学
統計学に確率論を導入し、推論手続きの一部を数式化
推論過程の一部を数式化したのみ
データの有権性に由来する統計学の帰納という本質は変わらない
優位性検定
観測できる事象の構成要素に関する特定の確率分布(モデル)を措定
措定した確率分布は人間が勝手に想定(モデリング)したもの
観測できる生の事象の確率分布をこれから導き出すことはできない
仮説が間違う可能性もある
帰納仮説
このモデルが真であると仮定したときにまれにしか起こらない事象群を確率の計算に基づいて求める
観測データがこの事象群に入ったならばこの仮説を「棄却」する
どのくらい小さければ稀と判断すべきかも根拠が判然としない
拡大する知と統計科学
現代科学の対象は異なる時間スケールで変化する多様な要素が複雑かつ階層的に結合した多自由度を持つ系
統計学は、データの収集と整理のための技術に関する記述的な学問だけではない
現代の統計学は、 単一の事象の数量的把握にのみ関わるものではなく、 相互に絡み合った複雑な事象間に伏在する関係の構造をモデリングし、 線形的知識および有限の経験データを統合し 事象の認識・予測・制御を行う方法を提供する
統計科学の帰納推論の3つの要素
統計・確率モデル
対象の構造や法則性を柔軟に表現することのできる確率・統計モデル
データ
アルゴリズム
非適切問題の解法の発展
「方程式モデル」→「最適化(変分)モデル」→「統計・確率モデル」
複数の部分系を持つシステムのモデリングの場合
この部分はよくわかっているが、そこの部分はかなり怪しい、あそこは全くわからない
方程式アプローチだと「わからない」部分があると頓挫する
最適化モデルだと、 よくわかっていないまたはわからない部分を 大雑把にモデリングして重み付けして 最適化すべき関数に加えこめば処理できるが 重みは値をどう取るべきか指針が得られない
最適化(最小化)すべき関数と指数関数を統合し、正規化するとギブス分布という統計モデルを作ることができる
経験ベイズほう
統計の目的はパラメータの特定ではなく、確率分布自体
罰金付きロジスティック回帰モデル
[補論]ベイズ統計と機械学習
ニュートン・デカルト・パラダイムからの離脱としての機械学習
機械学習における推論には精密な測定データは必ずしも必要としない
機械学習では変数の間に措定されるべき機序の先験的な発見を行わない
認知と推論の道具としてのモデル
モデルは世界の事象を意識的に選択・解釈・表現しその帰結を推論するために必須な道具
明示的なモデルだけではなく、無意識のモデルもある
モデルと現実世界は乖離している
機械学習は、 本質的な要因の文節やデータの次元圧縮等もなしで、 無関係なものを含む現象データを単に並べ上げた超多次元のデータセットから 意味ある共起構造を見えない形で掴み、それに基づいて推論するのみ
機械学習におけるモデルとベイズ統計
機械学習におけるモデルでは、内部パラメータを調節すればどんなデータにも適合しうる「ぐにゃぐにゃ」の可塑性を持ったモデルを用いるので過剰適合する可能性がある
これを回避して予測性(汎化性)を損なわない仕組みが機械学習には備わっている
ベイズモデルでは内部パラメータを間接的に柔らかく制限することで汎化性を担保する
階層化することでハイパーパラメータを作りより自動的に処理することも可能
機械学習の適用領域と限界
学習データから大きく外れた道データに固いする予測力はない
外れたデータであるか否かも予めわからない
機械学習による推論では機序は人間には解釈不能

第II部　 階層ベイズ講義　　 伊庭幸人

はじめに
はじめに
階層ベイズモデリング
3つの特徴
データの生成過程を確率分布で表現する
直接観測されない隠れた要素(潜在変数)を多数考える
線形モデルや正規分布に拘らず、離散と連続を同等に扱う
構成と文献の引用
岩波DS(岩波データサイエンス)
統フロ(統計科学のフロンティア)
久保緑本(データ解析のための統計モデリング)
PRML(パターン認識と機械学習)
Bayesian Data Analysis
StanとRでベイズ統計モデリング
記号と用語

講義0　 ベイズ・階層ベイズ・経験ベイズ

はじめに
岩波DS1のベイズ超速習コース
ベイズ
データyが確率分布p(y|x)から生成される
確率分布を決めるパラメータxも別の確率分布p(x)からのサンプル
p(x)をxの事前分布
→x→y
p(x|y) = p(x,y) / ∫ p(x,y)dx = p(y|x)p(x) /∫ p(y|x)p(x)dx
積分の範囲はxの定義されている範囲全体
xが多変量のベクトルであれば積分は多重積分になる
離散変数なら和に読み替えられる
事後分布
xのみに興味があるのであれば、分母を1/Cとして
p(x|y) =Cp(y|x)p(x)と表記できる
Cは事後分布の正規化定数
事後分布から情報を抽出する方法
興味のある統計量A(x)の事後分布のもとでの期待値や中央値、四分位数を計算
事後分布の確率密度を最大にするxの値をxの推定量とする(MAP推定量)
最尤推定量
パラメータxそのものではなく、 推定したxに基づく未来のデータzの予測が目的
p(z|y) = ∫ p(z|x)p(x|y)dx = ∫ p(z|x)p(y|x)p(x)dx / ∫ p(y|x)p(x)dx
事前分布の影響
無情報事前分布(散漫事前分布)
「パラメータxについての事前知識がない」ことを表現する大きな広がりを持つ事前分布
適切な設定は難しい
事前分布に積極的に知識を取り込むことも可能
パラメータさの成分の個数を固定して、データyのサンプルサイズを大きくすると、事前分布の効果は少なくなって、ベイズ推定でも最尤推定でも同じような結果に近づく
サンプルサイズが小さい時は事前分布p(x)とデータを含む部分p(y|x)の釣り合いで推定結果が決まる
階層ベイズ
事前分布p(x)にパラメータγを入れてp(x|γ)とし
このγにまた事前分布p(γ)を仮定する
→γ→x→y
xの事前分布を単なる主観ではなく、データに適応的に決める
xとγの同時事後分布
p(x,γ|y) = p(y|x)p(x|γ)p(γ) / ∫ p(y|x)p(x|γ)p(γ)dxdγ
経験ベイズ
xとγの同時事後分布p(x,γ|y)を直接使う推定法
フルベイズ法
周辺尤度・エビデンス
階層ベイズの目的
局所的なパラメータxに興味があり、そのための事前分布p(x|γ)をデータから求める
大局的なパラメータγやそれから定まるデータ生成プロセスの全体(混合分布)に興味がある
p(x|y)の部分に「回帰曲線の傾き」のようなパラメータが含まれていて、主に知りたいのはその値

講義1　 階層ベイズの2つの顔

1 スタイン推定量から階層ベイズへ
はじめに
階層ベイズモデリングの一つの起源
「縮小推定」の問題
推定に関係ない特徴に対応するパラメータの影響を 縮小するように推定する方法を パラメータ縮小推定 ( parameter shrinkage) や， 単に 縮小推定(shrinkage method) という
問題設定
N個の観測値{yi}
各yiは別々の期待値θiを持つ 正規分布から得られたと仮定
p(yi|θi) = 1/√(2πσ2)exp(-(yi-θi)2/2σ2)
例:yiが「重さ」とするとy1,y2..は 「あるクラスの生徒の体重」、 「いろいろな種類のペットの体重」、 「家の中にある家具の重さ」
各正規分布の分散(測定誤差)は全て同じσ2で既知と仮定
サンプルサイズnは4以上
推定量θi*({yi})の良さの尺度
パラメータの真値θiと推定値の2乗誤差の期待値
スタイン推定量
上記の設定では、iごとに一つの測定値yiしかないので、θiの推定量はそれぞれの測定値そのままのθに=yi
どのような{θi}についても
スタイン推定量
直感的意味
それぞれのyiを、全部の{yi}の平均値の方向にaだけ引っ張ってやる
{yi}の平均値
引っ張る程度のaを人間が決めずに、 2行目の式でデータから適応的に求める
スタイン推定は「縮小推定量」の一種
「平均との差」を「縮小」する
スタイン推定量の仕組み
「関係のないものを一緒に図ると結果が良くなる」?
期待値θiがバラバラ(測定したもの同士が関係ない)場合はS2がσ2と比べて大きくなりaはゼロにになる
期待値θiが近いもの(測定したものが「関係がある」)場合は、s2が小さくなり、aが無視できなくなり全体の平均に引き寄せられる
逆説から応用へ
「バイアスとばらつきのバランスをとる」
例
大リーグの打者の打率をハーズンはじめの少ない打席から推定する
シーズン末の打率(推定に使った打席の分は除く)が「真の打率」
地域ごとに分けるとサンプルサイズが不足して誤差が大きくなるが、 全部平均してしまうと地域差が全く見えなくなって困るケースへの応用
「小地域推定」の問題
ベイズ統計による解釈
上述では「θiのばらつき」は考慮しているが{θi}自体は正規分布のパラメータで、確率変数ではない
{θi}を確率変数と見る
Θiの事前分布を、平均θ0、分散δ2の正規分布と仮定
p(θi|δ2, θ0) = 1/ √(2πδ2)exp(-(θi-θ0)2/2δ2)
{yi}はp(yi|θi)からiごとに独立に生成
ベイズの公式から、θiの事後分布は
p(θi|yi, δ2, θ0) = Cp(yi|θi)p(θi|δ2, θ0) = C’ exp[-(yi-θi)2/2σ2-(θi-θ0)2/2δ2] C,C’はδ2を含む正規化定数
[ ]の中を最大にするθiが事後密度を最大にする推定量(MAP推定量)
Θiで微分してゼロと置くと θiMAP = (1-b)yi + bθ0, b=σ2/(δ2 + σ2)
スタイン推定量と同様の式
経験ベイズ
事前分布のθ0と分散δ2が未知
Θ0とδ2をデータから求める
ベイズモデリングの観点ではΘ0とδ2にさらに事前分布を仮定
経験ベイズ法による近似
Θ0とδ2の事後分布を直接扱う代わりに、周辺尤度を最大化する
2つの考え方の違い
その後の展開
ベイズモデルでは
「少なくとも悪くはならない」証明を諦めることで、適用範囲が転がる
モデルが不適切だと「何も考えない単純な推定方式」より悪い推定結果になる
スタインの補題に基づく議論は、 訓練データのサンプルサイズが有限の場合に 厳格な結果が得られる論理的手法として貴重
グループ構造と複数の要因を含む階層ベイズモデルの例
2 過分散から階層ベイズへ
はじめに
階層ベイズモデルの別の道筋
「過分散」の問題から始まるランダム効果モデル、混合分布モデル
過分散とは
2項分布やポアソン分布では「平均と分散の間に関係があり、母平均を与えると母分散も決まってしまう」という性質がある
ポアソン分布では、母平均がθから母分散もそれと等しくθになる
サンプルサイズが大きければ、サンプルの平均と分散はほぼ等しくなる
例:ガイガーカウンターで30秒のカウント数を10回ずつ計ったもの
データAの平均は21.2、不偏分散は24.2でほぼ同じ
データAは固定点での測定
データBは平均28なのに普遍分散は54.4
ポアソン分布の過程と乖離
データBは家の中のあちこちに移動して測定
ポアソン分布のばらつきのほかに、測定場所のばらつきが加わる
過分散
病気の発生数で個人差や地域差が効く場合
植物の開花数で個体差や土壌の差が問題になる場合
隠れた非統一性がある場合
混合分布
「直接観測されていない非統一性」に対して確率分布を考えて、 はじめに仮定した分布を「混合」したものを考える
非負の整数のデータ{yi}について観測値yiがポアソン分布に従うとする
p(yi|θi) =θiyi / yi! exp(-θi)
強度を表すパラメータθiに対してパラメータγを含んだ分布p(θi|γ)を考える
γ→θi→yi
p(yi|γ) = ∫ p(yi|θi)p(θi|γ)dθ
共役事前分布の利用
事前分布p(θi|γ)の選び方
典型的な2通りの方法の一つ
混合分布の式の形が簡単になるように選ぶ
共役事前分布
ポアソン分布に対してはガンマ分布を選ぶ
2項分布に対してはベータ分布
多項分布に対してはディリクレ分布
自然言語処理ではディリレクレ分布の無限次元への拡張である「ディリクレ過程」が事前分布として活用される
リンク関数の利用
もう一つの方法
リンク関数を用いて、強度θiを「μi=logθi」と変換してから、
μiを「μi =β + γi」のように定数βとその周りのばらつきγiの和で表現
γiの分布として正規表現を仮定
γiは直接観測されない量
ランダム効果
特徴:回帰(例:ポアソン回帰)への拡張が容易
共役事前分布でもポアソン回帰への拡張は可能(負の二項回帰)だが、モデルが複雑になると扱いにくくなる
リンク関数を使う方法は融通が効く
一般化線形混合モデル(Generalized Linear Mixed Model, GLMM)
ポアソン回帰などの一般化線形モデル(GLM)ランダム効果を取り入れたモデル
渋滞空間モデルやCARモデルに組み込む場合にも便利
階層ベイズモデルとみなす意味
複雑な分布を表現するのに{θi}を導入
{θi}をこのしてMCMCで数値計算していく方が現実的
実際の問題はγやαやβ等のパラメータの混合分布の最尤推定
高次元数値積分と最適化の組み合わせでモデルを当てはめるのは困難
その後の展開
過小分散と過大分散
過小分散の原因
人間や機械の働きで意図的に均一にされている
要素の間に反発力が働いている
例:川岸に並んだカップルの一定区間ごとの人数
講義1 のまとめ
ともに階層ベイズで扱うことが可能
個体差や非一様性を表現する局所的なパラメータ
γ,α,βのような大域的なパラメータ

講義2　 相関を表現する事前分布

1 状態空間モデル
はじめに
状態空間モデルは、時系列解析の手法を統計モデリングの視点から統一的に見直すもの
状態空間モデルとは
状態空間モデル
時系列データ{yt},t=1,…,nの背後に、 直接観測されない状態の列{xt},t=1,…,nを仮定したもの
xtの時間発展を記述する式(システム方程式)
xt+1 = F(xt) + ηt
データytの発生(観測プロセス)を表現する式(観測方程式)
yt = H(xt) + εt
H,Fは任意の関数、システム雑音ηtと観測雑音εtは各時点で独立の確率変数
システムモデルの条件付き確率分布
p(xt+1|xt) =1/√(2πδ2)exp(-(xt+1-F(xt))2/2δ2), t=0,…,n-1
観測モデルの条件先確率分布
p(yt|xt) =1/√(2πδ2)exp(-(yt-H(xt))2/2δ2), t=0,…,n
例:xtがスカラー(1成分)でF(xt)=xt
システム方程式
xt+1=xt + ηt
ランダムウォークの式
ローカルレベルモデル
隣接する時刻tとt+1でのxの値に差が少ない
一般化とパラメータ推定
雑音の分布としてコーシー分布などを採用しても良い
観測方程式の部分をポアソン分布や2項分布で置き換えることもできる
一般化線形モデル(GLM)の時系列モデル版
yt ~ p(yt|g(xt))
複数の時刻の状態の組みを新しい状態と定義することで 「xt+1の条件付き確率がxtとxt-1に依存する」 というモデルを作ることができる
遅延座標
時間遅れ座標
「xが滑らかに変化する」ことの表現
2階の差分方程式
xt+1 = 2xt – xt-1 +ηt
一階の差分の形式に書ける
ローカルトレンドモデル
任意のARモデルやATMAモデルを状態空間モデルの枠内に取り込める
σ2やδ2のような雑音の大きさを決めるパラメータや H,Fに含まれるパラメータをデータから推定する
これらをまとめてα,βとする
p(yt|xt:α),p(xt+1|xt:β)のような 条件付き確率にパラメータが含まれている
最尤法アプローチ
式を最大化するα,β
ベイズ的な解釈
状態空間モデルのベイズ的な解釈
状態x={xt}を局所的なパラメータとする
xからのデータ生成を表す部分
ベイズの式
カルマンフィルタでの「線形最小分散推定量の計算」
カルマンフィルタは線形モデル(システム方程式と観測方程式がどちらも線形で、雑音が全て独立で正規分布に従う)でのベイズ
1次元での平滑化、関数や曲線の推定
状態空間モデルは直線上に観測点が配置されているような空間データの解析にも使える
2 空間への拡張
はじめに
相関を表現する事前分布という考え方は空間構造を表現するもでるにも拡張可能
CARモデル
空間一次元では時間のモデルと同様に組めるが2次元ではうまくいかない
アプローチとして、 各点iに対して、近傍N(i)のxの値を固定した時の 条件付き確率p(xi|{xj}j∈N(i))を求めるモデルで表現する
うまくいかない場合もある
「ミクロな条件付き確率密度関数の積による表現」という発想をやめて 「事前確率密度の対数をミクロな項の和」として与える
マルコフ確率場 (Markov Random Field, MRF)
定義された分布が多変量正規分布になる場合
ガウス型マルコフ確率場 (GMRF)
正方格子上で定義されたCARモデル
ガウス過程とカーネル回帰
ガウス過程を事前分布として使う
計算量がサンプルサイズで決まり、次元によらないので、高次元で有利になる
一般のマルコフ確率場モデル
Geman達の論文
離散確率変数のマルコフ確率場(統計物理でのボッツモデル) およびそれに先過程と呼ばれる補助的な離散変数を 導入したモデルを事前分布として用いる
MCMCによって事後分布を最大化してMAP推定を求める
問題点
自然画像へのアプローチは難しい
3 非適切逆問題
逆問題の例いろいろ
非適切逆問題の重要な問題には時間や空間に関連したものが多い
代表的な例
コンピュータ断層撮影(CT)
地震波を使って地球の内部を探る
影から月の山の高さを測る
両眼の視差から距離を求める
物体の表面の再構成や物体の運動の認識
「錯覚」の解釈
非適切逆問題と正則化
逆問題の多くでは、観測値岳からは一位に定まらなかったり、 微小な雑音の影響で解が大きく変動したりする
非適切(ill-used)
偏微分方程式の初期値問題に端を発する
推定する対象についての事前知識を何らかの形で組み込むことが必要
対象のモデリング
Yをデータ、xを推定したい対象とした時、 対数尤度logp(y|x)に罰則項(罰金項、正則化項)f(x)を加える
l(x) = logp(y|x) – λf(x)
λ:罰則の強さ(正則化の強さ)
罰則付き最尤推定
非適切逆問題と階層ベイズ
ベイズで解釈するために、xの事前分布を
p(x|λ) = exp(-λf(x)) / Z(λ)
罰則項f(x)がxの成分の2次式の場合は、xの事前分布はガウス型になる

講義3　 外れ値・クラスター分け・欠測

はじめに
階層ベイズモデルは、外れ値、確率モデルに基づくクラスター分け、欠陥などを組織的に吸収可能
1 離散値を取るパラメータの利用
はじめに
離散値をとるパラメータ(潜在変数、ラベル)を導入することで、 外れ値やクラスターーの分類のような問題を 階層ベイズモデルの枠組みで扱うことができる
ベイズの枠組みでは離散と連続の垣根が低い
外れ値のモデル
外れ値例
xi∈{0,1}が観測iの結果が「外れ値」であるかどうかを決める 「ラベル」(あるいは「離散パラメータ」)
データ{yi}i=1,…,nの条件付き確率(尤度関数)
xi=0であれば「外れ値でない時のモデルp(y|γ)からyiが発生
xi=1であれば「外れ値の時のモデル」pi(y)(例えば幅の大きな一様分布)からyiが発生
γの事前分布をp(γ)とすると事後確率は
p(γ,{xi}|{yi}) = Cp({yi}|γ,{xi})p({xi}|Q)p(γ)
MCMCでサンプリングすることで求まる
I番目が外れ値である事後確率
外れ値の影響を考慮したγの推定値
有限混合分布モデル
「外れ値か否か」だけではなく 2個のクラスターに分けるモデリング一般に使える
観測値iが一方のクラスターに属せばxi=0、 他方のクラスターに属せばxi=1
その値によって、 観測値yiが分布p(yi|γ0)と p(yi|γ1)のいずれから得られたか決まる
K個のクラスターに分ける場合
p(yi|γk)は各クラスターに対応する分布
事後分布
有限混合分布モデル
隠れマルコフモデル
データ{yi}が時系列
背後にある潜在変数{xi}も時系列
{xi}の事前分布
隠れマルコフモデル
直接観測されない状態に従う事前分布が、遷移確率p(xt+1|xt)で決まる「マルコフ連鎖」
応用例
音声認識
病気の進行や治療の効果を表す
経済事例列で「景気拡張期」「景気後退期」のような 「レジーム」を隠れた状態と考えるモデル
マルコフ・スイッチングモデル
データの従うARモデルの一部が隠れた状態によって切り替わる
「系列事象」のモデリングが可能
自然言語処理
Nグラムモデル
DNAの塩基配列(配列のアライメント)
この先の発展と問題点
隠れマルコフモデルと状態空間モデルのハイブリッド
状態xiとして離散値を取る変数と連続変数の両者を含むような拡張
課題
複雑なモデリングは多数の極大値を持つ事後分布(多峰性分布)を持つことが多い
MCMCの計算に困難を生じる
多峰性分布に強いMCMCのアルゴリズムもあるがJAGS等のMCMCツールに組み込まれていない
Stanはそもそも離散変数のサンプリングを直接サポートしていない
ラベル・スイッチングの問題
2 欠測
はじめに
「欠測」をアドホックな扱いではなく、「観測できない状態」とみなして統計的に扱う
階層ベイズモデルでの取り扱い
ランダムな欠測
欠測の例
{(yi,zi)}というデータに対してのモデル
パラメータγを推定する問題
ただしziの値がランダムに欠測
{zi}の背後に{xi}という「隠れた状態」がある
ziはxiから以下の確率で生成
確率q NA
NA:欠測
確率1-q xi
欠測しているiの集合をNA={i|zi=NA}
NAの要素の個数#NAをm
ディラックのδ関数
超関数
鋭いピークを持つ確率関数
正規分布の密度関数で分散が弧状に小さいもの
拘束条件s=s’の表現
Yi,zi,xiの同時関数はγ→(yi,xi)→(yi,zi)という生成プロセスを考えて
p(yi,zi,xi|γ) = p(yi,xi|γ)[(1-q)xδ(zi-xi) + qxδ(zi=NA)] [ ]の中は[確率1-qでzi=xi][確率qでzi=NA]を式で表現
Γの事前分布をp(γ)とすると
{yi},{zi}, {xi}, γの同時分布は
完全にランダムな欠測(MCAR)
Ziの欠測の有無は独立
「打ち切り」を含むモデル
結束の条件が明確に分かっている場合
例:ziの値が既知の閾値ξを超えると必ず欠測し、 そうでない場合は結束しない
測定時間の上限が一時間でそこで打ち切られる
zi={NA, xi>ξ {xi, xi≤ξ
I(条件式) = {1, (条件が成立) {0, (条件が不成立)
p(yi, zi, xi|γ) = p(yi, xi|γ) [I (xi ≤ ξ)δ(zi-xi) + I ( xi > ξ)δ(zi=NA)] xi ≤ ξなら zi=xi
xi > ξなら zi=NA
同時分布
欠測とベイズモデリング
欠陥値の背後にある状態を局所的なパラメータ{xi}と考える
欠陥のモデリングは階層ベイズモデリングとほぼ同じ
「パラメータの事前分布」とみなすか「欠陥が生じる前の完全データの分布」とみなすか
融合して考える利点
欠陥のモデリングと階層ベイズモデリングを一体化して考える利点
欠陥の処理とモデリングにもともと含まれている「隠れた状態」の処理を合わせて 「直接観測されていない変数について積分する」という統一したアプローチができる
数学的手法の共通化ができる
ギブス・サンプラーの適用
データ拡大法 (data augmentation)
考えるべき点
ベイズモデリングの背後には「バイアスとばらつきをバランスさせる」という考え方がある
欠陥に関しては「バイアスの補正を優先する」
因果推論
ルービンの枠組みでは因果推論における「反事実」を欠測の一種として捉える
交路によるバイアスの補正

付録A　 階層ベイズモデルの予測分布

2種類の予測分布
「興味のあるパラメータの違い」が 階層ベイズの予測分布の違いにどのように反映されるか
局所的なパラメータxに興味がある場合
「xに興味がある」ということは 「予測」に焼き直すと 「xがデータから推定した値をとるとして、そこから未来のデータzが独立に発生する」」
xを事後分布からサンプリングして、それを用いて新しいデータを生成した時の分布
対極的なパラメータγに興味がある場合
現在のxの値は偶発的なものとして無視
Γを事後分布からサンプルして、 そのγの値を入れた条件付き分布から新たにx’を生成し、 そのx’を用いて未来のデータzを生成した時の分布
「回帰モデルの直線の傾き」のようなパラメータに興味がある場合
予測分布の評価
性能評価によく使われる手法
サンプルサイズの大きい極限を考える方法
漸近理論
階層ベイズモデルの場合に、 サンプル祭がを大きくした時漸近理論の結果が 良い近似になるかどうか
例:グループ構造を持つモデルでの2つの極限
グループ数を増やして、各グループのメンバーは一定にする
グループ数を一定にして、各グループのメンバー数を増やす

付録B　 スタイン推定量が2乗誤差の期待値を改良することの証明

一般的な考察
交差項の評価
与えられた値に縮小する場合
データから求めた平均値に縮小する場合

付録C　 事前分布が指数型分布族の場合の経験ベイズ推定

ここで考える問題
式(2)の直感的意味
生成モデルが良くない場合の経験ベイズ法の問題点
式の導出