東京大学・株式会社Nospareの菅澤です．
今回はブートストラップ法を用いた事後分布と一般化ベイズ法との関連について紹介します．

はじめに

$X=(x_1,\dots,x_n)$を未知の分布$F_0$からの独立なサンプルとします．
以前に紹介したweighted likelihood boostrap (WLB)での設定とは異なり，今回は興味のあるパラメータ$\theta$の真値$\theta_0$は，以下のように一般的な損失関数$\ell(\theta, x)$によって定まるとします．

\theta_0\equiv \theta(F_0)={\rm argmin}_{\theta}\int \ell(\theta, x)dF_0(x)

データに対して統計モデル$p(x;\theta)$を考え，損失関数として負の対数尤度$\ell(\theta, x)=-\log p(x;\theta)$を考えると，上記の$\theta_0$はモデル$p(x;\theta)$とデータ生成分布$F_0$の間のカルバックライブラー距離を最小化する値として解釈することができます．その意味で，上記のような$\theta_0$の定義は，統計モデルが誤特定された状況や，モデルを特定せずに損失関数のみによってパラメータ推定を実行するなど，様々な状況を含んだ一般的な設定になっています．

WLBは$p(x;\theta)$が生成分布$F_0$を含んでいる状況 (モデルが正しく特定化されている状況) で事後分布の近似を与えたのに対し，Lyddon et al. (2019)では，上記で定義された一般的な$\theta_0$の事後分布をターゲットにしたLoss-likelihood bootstrap (LLB) 法を提案しています．

Loss-likelihood bootstrap (LLB)

もし未知の分布$F_0$を知ることができれば，$\theta_0=\theta(F_0)$から真の値を「計算」することができます．しかし，実際には$F_0$は未知です．LLBは，$\theta$の代わりに分布$F$を未知のものとして扱い，ターゲットパラメータ$\theta_0$はあくまで分布$F_0$によって自動的に決まる量という考え方を採用します．

このように考えることで，$\theta$の事後分布を考える問題が$F$の事後分布を求める問題に話がすり替わります．分布$F$をベイズ推定するために，ディリクレ過程を事前分布として使うことを考えます．するとデータ$X$が与えられたもとで，$F$の事後分布は

\pi(F \mid X) \sim {\rm DP}\left(\alpha + \sum_{i=1}^n \delta_{x_i}\right)

となります．ここで，$\delta_x$は$x$上の一点分布で，$\alpha$はディリクレ過程を定める測度（$F$に関する事前情報を反映したもの）です．ここでは無情報な状況を想定し，任意の集合$\Omega$に対して$\alpha(\Omega)\to 0$のケースを考えます．この場合，$F$の事後分布から生成される分布はデータ点${x_1,\ldots,x_n}$上にサポートを持ち，各点の重みが${\rm Dir} (1,\ldots,1)$で決まる離散分布になります．

LLBは$F$の事後分布から$F$を生成し，$\theta=\theta(F)$として$\theta$を生成することで$\theta$の事後分布を与えるアルゴリズムです．ちなみにこのような事後分布の計算方法はnonparametric Bayesian learningとも呼ばれています．

LLBの具体的な手順としては，以下のステップを$B$回繰り返します．

${\rm Dir} (1,\ldots,1)$からランダムな重み$(w_1,\ldots,w_n)$を生成する．
重み付き損失関数$\sum_{i=1}^n w_i\ell(\theta,x_i)$を最小化する$\theta$を計算する．

上記のアルゴリズムはWLBと非常に似ていることがわかります．実際にアルゴリズム上では，WLBとLLBは重み付き尤度関数を最大化するか重み付き損失関数を最小化するかの違いしかありません．

一般化事後分布とWLBの違い

損失関数$\sum_{i=1}^n \ell(\theta,x_i)$に基づく一般化事後分布は

\pi(\theta\mid X)\propto \pi(\theta)\exp\left\{-\omega \sum_{i=1}^n \ell(\theta,x_i)\right\}

で与えられます．この事後分布の形の妥当性については一般化ベイズ法の記事で紹介しております．ここで$\omega>0$はlearning rateと呼ばれ，事後分布の広がりを調整するチューニングパラメータです．

ここで，損失関数として負の対数尤度を用いた場合，上記で$\omega=1$とした事後分布は通常の尤度による事後分布に一致し，LLBはWLBに一致します．WLBの記事で紹介したように，WLBによって得られる分布と事後分布は漸近的に（$n\to\infty$の状況で）一致します．このように，尤度が正しく特定化されている状況では，ブートストラップによって与えられる分布と事後分布に近似的な対応関係があるのですが，LLBのように一般的な損失関数のもとではその関係性が変わってきます．

Lyddon et al. (2019)ではLLBによって与えられる$\theta$の分布の漸近分布が与えられています．一般化事後分布の漸近分布と並べて書くと以下のようになります．

LLBの漸近分布: $N(\widehat{\theta}, J(\theta_0)^{-1}I(\theta_0)J(\theta_0)^{-1}/n)$
一般化事後分布の漸近分布: $N(\widehat{\theta}, \omega^{-1}J(\theta_0)^{-1}/n)$

ここで$\widehat{\theta}$は経験損失$\sum_{i=1}^n \ell(\theta, x_i)$の最小化によって与えられる推定量であり

I(\theta)=\int \left(\frac{\partial}{\partial\theta}\ell(\theta, x)\right)\left(\frac{\partial}{\partial\theta}\ell(\theta, x)\right)^\top dF_0(x), \ \ \ \ \ \ 
J(\theta)=\int \frac{\partial^2}{\partial\theta\partial\theta^\top}\ell(\theta, x) dF_0(x)

です．両者の漸近分布は共分散行列の部分が異なっていることがわかります．LLBの漸近分布では，いわゆるサンドイッチ分散の形が現れており，経験損失に基づくM推定量の漸近的な性質と同じような結果になっています．一方で，一般化事後分布の漸近分散はサンドイッチ型になっていないため，一般の損失関数では不確実性を誤って見積もってしまう可能性があることを示唆しています．

一般化事後分布のカリブレーション

一般化事後分布を使う場合はlearning rateである$\omega$を決める必要がありました．前述のLLBの漸近的な結果の副産物として，$\omega$のカリブレーションの方法を与えることができます．

$\theta$が多次元の場合にスカラー$\omega$によって2つの行列$J(\theta_0)^{-1}I(\theta_0)J(\theta_0)^{-1}$と$\omega^{-1}J(\theta_0)^{-1}$を一致させることはできませんが，代替案として両者のトレースを一致させるように$\omega$をカリブレーションする方法が考えられます．具体的には，$\omega$を以下のように設定します．

\widehat{\omega}=\frac{{\rm tr}\{J_n(\widehat{\theta})^{-1}I_n(\widehat{\theta})J_n(\widehat{\theta})^{-1}\}}{{\rm tr}\{J_n(\widehat{\theta})^{-1}\}} \\
I_n(\theta)=\frac{1}{n}\sum_{i=1}^n\left(\frac{\partial}{\partial\theta}\ell(\theta, x_i)\right)\left(\frac{\partial}{\partial\theta}\ell(\theta, x_i)\right)^\top \ \ \ \ \ \ \ 
J_n(\theta)=\frac{1}{n}\sum_{i=1}^n\frac{\partial^2}{\partial\theta\partial\theta^\top}\ell(\theta, x_i)

このように定めた$\omega$を使った一般化事後分布を考えることで，LLBと同等程度の不確実性を持った一般化事後分布を得ることができます．

おわりに

今回は一般的な損失に対するブートストラップ法と一般化ベイズ法との関連について紹介しました．

株式会社Nospareでは統計学の様々な分野を専門とする研究者が所属しております．統計アドバイザリーやビジネスデータの分析につきましては株式会社Nospareまでお問い合わせください．