平均・分散が未知の正規分布のベイズ推定に引き続き，報酬が平均・分散が未知の正規分布に従うと仮定したThompson Samplingの解説第３回目です．

第３回目となる今回は，多腕バンディット問題に対する腕選択アルゴリズムの１つである，Thompson Samplingについて，具体的な式を交えて解説したいと思います．

TL;DR

• Thompson Samplingは，腕選択におけるexplorationとexploitationのバランスを取ることができる多腕バンディット問題に対する方策の１つ
• Thompson Samplingは，各腕が最適である事後確率に従って腕を選択
• そのためには，報酬の期待値に対する事後分布からの乱数が必要になる
  • 報酬が平均・分散が未知の正規分布に従っていると仮定し，事後分布を導出
  • 平均に対する事前分布は正規分布，分散に対する事前分布はscaled inverse chi-squared分布をチョイス
  • これによって事前分布が共役事前分布になってくれる
• 結局は，Thompson Samplingのアルゴリズムでは平均に対する事後分布のパラメータと，分散に対する事後分布のパラメータを更新するだけでよくなる

Thompson Samplingとは

Thompson Samplingは，多腕バンディット問題に対する方策（どの腕をプレイするかを選択するアルゴリズム）の１つです（スロットマシンの台のことを腕とも呼びます）．

多腕バンディット問題の解説記事で触れましたが，多腕バンディット問題では，探索（exploration）と活用（exploitation）のバランスを取りながら腕を選択できるか，ということが重要になります．

Thompson Samplingは，比較的簡単なアルゴリズムでありながら，問題に合わせて探索と活用のバランスを取ってくれる，という強力なアルゴリズムです． UCB戦略と同様に，多腕バンディット問題の方策の代表例として挙げられることが多いです．

さて，以降では，$K$本の腕があるとし，各腕$i(1\leq i\leq K)$からの報酬が平均（期待値）$\mu_i$と分散$\sigma^2_i$が未知の正規分布$N(\mu_i,\sigma^2_i)$に従って生成されると仮定します． この仮定の元で，Thompson Samplingのアルゴリズムについて数式を交えて解説していきます．

Thompson Samplingの腕選択方法

Thompson Samplingは，各腕$i$が最適である事後確率に従って腕を選択します． すなわち，次の式のように，「ある腕の報酬の期待値が他のどの腕の報酬の期待値よりも高い」確率に従って，その腕を選択します．

$$\pi(\mu_i \geq \max_{j\neq i}\mu_j \vert \{X_i\}_{i=1}^K) = \int_0^1 \pi(\theta_i \vert X_i) \left( \prod_{j\neq i} \int_0^{\theta_i} \pi(\theta_j \vert X_j) d\theta_j \right) d\theta_i$$

ここで，$\mu_i$は腕$i$の報酬の期待値，$X_i$は腕$i$をプレイしたことによって観測した$T_i$個の報酬の履歴$X_i=\{r_t\}_{t=1}^{T_i}$です．

上の式は積分計算が入っているため，計算を行うことは一般的には困難です． しかし，各腕$i$から乱数$\theta_i$を事後分布$\pi(\mu_i \vert X_i)$に従って生成し，$\theta_i$が最大となった腕をプレイすれば，上の式に従って腕を選択したことと同等の動作となります．

つまり，事後分布$\pi(\mu_i \vert X_i)$からの乱数を生成することができれば，「ある腕の報酬の期待値が他のどの腕の報酬の期待値よりも高い確率に従って腕を選択する」ということが実現できるのです．

以上をまとめると，Thompson Samplingでは，以下のような流れを繰り返すことで，各ステップごとに腕を選択します．

1. 各腕$i$の報酬の期待値に関する事後分布$\pi(\mu_i \vert X_i)$から乱数$\theta_i$を生成
2. 最大の乱数を生成した腕$i_{max}={\rm arg~max}_i\theta_i$をプレイし，報酬$r$を観測
3. 観測した報酬の履歴$X_{i_{max}}$に$r$を追加

事後分布$\pi(\mu_i \vert X_i)$からの乱数生成

これまでに記述したように，Thompson Samplingのアルゴリズムを実現するためには，事後分布$\pi(\theta_i \vert X_i)$からの乱数$\theta_i$を生成することが必要になります．

報酬の分布である正規分布$N(\mu_i,\sigma^2_i)$の分散$\sigma^2_i$が既知であるならば，単純に$\mu_i$に対する事後分布$\pi(\mu_i \vert X_i)$を導出し，その分布から乱数を生成すれば乱数$\theta_i$を得ることができます．

しかし，今回の場合，平均$\mu_i$だけでなく分散$\sigma^2_i$も未知としているため，$\mu_i$のみに対する事後分布$\pi(\mu_i \vert X_i)$を導出するのは困難です．

そこで，平均$\mu_i$と分散$\sigma^2_i$の組$(\mu_i,\sigma^2_i)$に対する事後分布$\pi(\mu_i,\sigma^2_i \vert X_i)=\pi(\mu_i\vert \sigma^2_i ; X_i)\pi(\sigma^2_i \vert X_i)$を導出し，この分布から乱数$\theta_i$を以下のようにして生成するようにします．

1. $\pi(\sigma^2_i \vert X_i)$に従って，乱数$\phi_i$を生成
2. $\pi(\mu_i\vert \sigma^2_i ; X_i)=\pi(\mu_i\vert \phi_i ; X_i)$に従って，乱数$\theta_i$を生成

事後分布$\pi(\mu_i,\sigma^2_i \vert X_i)$の計算

乱数$\theta_i$を生成するのに事後分布$\pi(\mu_i,\sigma^2_i \vert X_i)=\pi(\mu_i\vert \sigma^2_i ; X_i)\pi(\sigma^2_i \vert X_i)$を用いるため，$\pi(\mu_i\vert \sigma^2_i ; X_i)$と$\pi(\sigma^2_i \vert X_i)$を導出する必要があります．

これは，平均・分散が未知の正規分布のベイズ推定をしていることにほかなりません．

こちらの記事で記述しているように，平均$\mu_i$に対する事前分布$\pi(\mu_i \vert \sigma^2_i)$を正規分布$N(\hat{\mu}_{i,0}, \frac{\sigma^2_i}{k_{i,0}})$とし，分散$\sigma_i^2$に対する事前分布$\pi(\sigma^2_i)$をscaled inverse chi-squared分布$Scale-inv-\chi^2(v_{i,0},\hat{\sigma}^2_{i,0})$とすると，事後分布$\pi(\mu_i\vert \sigma^2_i ; X_i)$と$\pi(\sigma^2_i \vert X_i)$は以下のように導出することができます．

$$\begin{aligned} \pi(\mu_i \vert \sigma^2_i; X_i)&=\frac{1}{\sqrt{2\pi \sigma^2_i/k_{i,T_i}}}\exp{\left\{-\frac{(\mu_i-\hat{\mu}_{i,T_i})^2}{2\sigma^2_i/k_{i,T_i}}\right\}}, \\ \pi(\sigma^2_i \vert X_i)&=\frac{(\hat{\sigma}^2_{i,T_i}v_{i,T_i}/2)^{v_{i,T_i}/2}}{\Gamma(v_{i,T_i}/2)}\frac{\exp{\left(\frac{-v_{i,T_i}\hat{\sigma}^2_{i,T_i}}{2\sigma^2_i}\right)}}{(\sigma^2_i)^{1+v_{i,T_i}/2}} \end{aligned}$$

ここで，

$$\begin{aligned} v_{i,T_i}&=v_{i,0}+T_i, \\ v_{i,T_i}\hat{\sigma}^2_{i,T_i}&=v_{i,0}\hat{\sigma}^2_{i,0}+\sum_{t=1}^{T_i}(r_t-\frac{\sum_{t^{\prime}=1}^{T_i} r_{t^{\prime}}}{T_i})^2+\frac{k_{i,0}T_i}{k_{i,0}+T_i}(\frac{\sum_{t=1}^{T_i} r_t}{T_i}-\hat{\mu}^2_{i,0})^2, \\ k_{i,T_i}&=k_{i,0}+T_i, \\ \hat{\mu}_{i,T_i}&=\frac{k_{i,0}}{k_{i,0}+T_i}\hat{\mu}_{i,0}+\frac{\sum_{t=1}^{T_i} r_t}{k_{i,0}+T_i} \end{aligned}$$

です． また，$\hat{\mu}_{i,0}, k_{i,0}, v_{i,0}, \hat{\sigma}^2_{i,0}$は事前分布のパラメータとします．

導出方法の詳細は，記事を参照してください．

アルゴリズム

これまでの解説をまとめると，Thompson Samplingのアルゴリズムは以下のような流れになります．

これまでの解説で長ったらしい数式が出てきていたわりには，スッキリとしたアルゴリズムになっているのではないでしょうか．

これがThompson Samplingの良いところで，事後分布の導出は面倒ですが，最終的なアルゴリズムはシンプルになり，かつそこそこのリグレットで抑えられることが多いです．

報酬の分布には，今回は正規分布を仮定しましたが，それ以外にもいろいろな分布を仮定したThompson Samplingが存在します． 特に，ベルヌーイ分布を仮定したThompson Samplingは有名でコードも様々な方が公開しているため，初めて試すにはちょうどいいかも知れません．