はじめに
https://yuhi-sa.github.io/posts/20260223_bayesian_optimization/1/ では、ガウス過程を代理モデルとして、EI・UCB・PIという3種類の獲得関数で「次にどこを評価すべきか」を決める仕組みを解説しました。これらの獲得関数はいずれも、GPの予測平均と予測分散を明示的な数式に代入して次の点を決定します。本記事では、これとは異なるアプローチ——事後分布からサンプリングした値をそのまま使って意思決定する「トンプソンサンプリング(Thompson Sampling)」を、最も基本的な設定である多腕バンディット問題から出発して解説します。トンプソンサンプリングは1933年に提案された歴史上最も古い探索戦略の一つでありながら、現代の推薦システム・広告配信・A/Bテストの逐次最適化で標準的に使われており、後半ではガウス過程と組み合わせることでベイズ最適化の獲得関数の1つ(GP-TS)としても機能することを実装で確認します。
多腕バンディット問題
多腕バンディット(Multi-Armed Bandit, MAB)問題は、\(K\) 本の「腕(arm)」があり、腕 \(i\) を引くと未知の確率分布に従う報酬が得られる、という逐次意思決定問題です。本記事では最も基本的なBernoulliバンディットを扱います。腕 \(i\) を引くと確率 \(p_i\) で報酬1、確率 \(1-p_i\) で報酬0が得られます(\(p_i\) はエージェントには未知)。
\(T\) ステップの試行を通じて累積報酬を最大化したい、というのがこの問題の目標です。これは以下の**累積リグレット(cumulative regret)**の最小化と等価です。
\[ R(T) = \sum_{t=1}^{T} \left( p^* - p_{a_t} \right) \tag{1} \]ここで \(p^* = \max_i p_i\) は最良の腕の成功確率、\(a_t\) は時刻 \(t\) でエージェントが選んだ腕です。式(1)は「もし最初から最良の腕だけを引いていたら得られたはずの報酬」と「実際に得た報酬」の差を表しており、優れたアルゴリズムほど \(R(T)\) の増加が緩やか(理想的には \(O(\log T)\) )になります。
多腕バンディットは、探索(exploration)——まだよく知らない腕を試して情報を得る——と活用(exploitation)——これまでの情報から最良と思われる腕を選び続ける——のトレードオフという、強化学習・逐次実験計画に共通する中心的な課題を、最もシンプルな形で凝縮した問題です。
3つの戦略:ε-greedy・UCB1・トンプソンサンプリング
ε-greedy
最も単純な戦略です。確率 \(1-\varepsilon\) でこれまでの経験平均が最大の腕を選び(活用)、確率 \(\varepsilon\) でランダムな腕を選びます(探索)。
\[ a_t = \begin{cases} \arg\max_i \hat{p}_i(t) & \text{確率 } 1-\varepsilon \\ \text{一様ランダムな腕} & \text{確率 } \varepsilon \end{cases} \tag{2} \]ここで \(\hat{p}_i(t)\) は時刻 \(t\) までの腕 \(i\) の経験成功率です。\(\varepsilon\) を固定値にすると、収束後も一定確率で「無駄な」ランダム探索を続けるため、累積リグレットは漸近的に \(\Theta(T)\) (線形)で増加し続けます。
UCB1(Upper Confidence Bound)
「不確実性の中では楽観的であれ(optimism in the face of uncertainty)」という原則に基づき、経験平均に信頼区間の幅を加えたスコアで腕を選びます。
\[ a_t = \arg\max_i \left( \hat{p}_i(t) + \sqrt{\frac{2\ln t}{n_i(t)}} \right) \tag{3} \]ここで \(n_i(t)\) は腕 \(i\) をそれまでに引いた回数です。引いた回数が少ない腕ほどボーナス項 \(\sqrt{2\ln t / n_i(t)}\) が大きくなり、自動的に探索が促されます。理論的には \(O(\log T)\) の累積リグレットが保証されています(Auer, Cesa-Bianchi & Fischer, 2002)。
トンプソンサンプリング
トンプソンサンプリングは、各腕の成功確率 \(p_i\) に対するベイズ事後分布を保持し、その事後分布から1つサンプルを引いて、サンプル値が最大だった腕を選ぶという手続きです。Bernoulli報酬の共役事前分布はBeta分布なので、腕 \(i\) の事後分布は
\[ p_i \mid \mathcal{D} \sim \text{Beta}(\alpha_i, \beta_i), \qquad \alpha_i = 1 + s_i, \quad \beta_i = 1 + f_i \tag{4} \]と閉形式で更新できます(\(s_i, f_i\) はそれぞれ腕 \(i\) の成功・失敗回数、事前分布は無情報事前分布 \(\text{Beta}(1,1)\) =一様分布)。各ステップで
\[ \theta_i \sim \text{Beta}(\alpha_i, \beta_i) \quad (i=1,\ldots,K), \qquad a_t = \arg\max_i \theta_i \tag{5} \]としてサンプル値最大の腕を引き、観測結果に応じて式(4)の \(\alpha_{a_t}\) または \(\beta_{a_t}\) を1増やします。この手続きは事後分布の広がりそのものが探索量を自動調整するという美しい性質を持ちます。データが少なく事後分布が広い(不確実性が高い)腕ほどサンプル値が大きく振れて選ばれやすくなり、データが蓄積して事後分布が真の値の周りに集中すると自然に選ばれにくくなる仕組みです。UCB1の信頼区間ボーナスが「決め打ちの数式」で探索量を制御するのに対し、トンプソンサンプリングはベイズ事後分布のサンプリングという確率的な仕組みで同じ効果を実現しているという対比が重要です。
Python実装
3手法をnumpyでフルスクラッチ実装します。
import numpy as np
def eps_greedy(rng, K, T, true_probs, eps=0.1):
counts, sums = np.zeros(K), np.zeros(K)
chosen = np.zeros(T, dtype=int)
for t in range(T):
if t < K:
a = t # 全腕を1回ずつ引いてから開始
elif rng.random() < eps:
a = rng.integers(K)
else:
a = np.argmax(sums / np.maximum(counts, 1))
r = float(rng.random() < true_probs[a])
counts[a] += 1
sums[a] += r
chosen[t] = a
return chosen
def ucb1(rng, K, T, true_probs):
counts, sums = np.zeros(K), np.zeros(K)
chosen = np.zeros(T, dtype=int)
for t in range(T):
if t < K:
a = t
else:
bonus = np.sqrt(2 * np.log(t + 1) / counts)
a = np.argmax(sums / counts + bonus)
r = float(rng.random() < true_probs[a])
counts[a] += 1
sums[a] += r
chosen[t] = a
return chosen
def thompson_sampling(rng, K, T, true_probs):
alpha, beta = np.ones(K), np.ones(K)
chosen = np.zeros(T, dtype=int)
for t in range(T):
theta = rng.beta(alpha, beta) # 式(5):各腕の事後分布からサンプリング
a = np.argmax(theta)
r = float(rng.random() < true_probs[a])
if r > 0:
alpha[a] += 1 # 式(4):成功なら alpha を更新
else:
beta[a] += 1 # 式(4):失敗なら beta を更新
chosen[t] = a
return chosen
eps_greedy・ucb1・thompson_sampling はいずれも式(2)・式(3)・式(4)(5)にそのまま対応しています。counts[a] += 1 などの更新は各手法の統計量(経験平均、Beta分布のパラメータ)をオンラインで維持する部分です。
数値実験:累積リグレットの比較
\(K=8\)
本の腕、成功確率 \(p_i \in \{0.10, 0.20, 0.30, 0.35, 0.40, 0.45, 0.50, 0.55\}\)
(最良腕 \(p^*=0.55\)
)のBernoulliバンディットで、\(T=5000\)
ステップ・200シード平均で累積リグレット(式1)を比較しました(seed=1000〜1199)。
K = 8
true_probs = np.array([0.10, 0.20, 0.30, 0.35, 0.40, 0.45, 0.50, 0.55])
T, N_SEEDS = 5000, 200
best_p = true_probs.max()
regret_curves = {name: np.zeros((N_SEEDS, T)) for name in
["eps_greedy", "ucb1", "thompson"]}
for seed in range(N_SEEDS):
for name, fn in [("eps_greedy", eps_greedy), ("ucb1", ucb1),
("thompson", thompson_sampling)]:
rng = np.random.default_rng(1000 + seed)
chosen = fn(rng, K, T, true_probs) if name != "eps_greedy" \
else fn(rng, K, T, true_probs)
regret_curves[name][seed] = np.cumsum(best_p - true_probs[chosen])
実行結果
K=8 arms, true_probs=[0.1, 0.2, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55], best_p=0.55
T=5000 steps, N_SEEDS=200
eps_greedy : cumulative regret at T=5000: mean=167.53 std=95.31 %pulls on best arm=70.0%
ucb1 : cumulative regret at T=5000: mean=297.35 std=26.22 %pulls on best arm=51.7%
thompson : cumulative regret at T=5000: mean=81.25 std=31.75 %pulls on best arm=84.2%
トンプソンサンプリングの累積リグレットは81.25で、ε-greedy(167.53)の約半分、UCB1(297.35)の約4分の1という結果になりました。最良腕を選んだ割合も、トンプソンサンプリングが84.2%と最も高く、ε-greedy(70.0%)・UCB1(51.7%)を上回っています。
途中経過を見ると、この差がどう開いていくかが分かります。
t= 100: eps_greedy= 10.14 ucb1= 15.56 thompson= 12.39
t= 500: eps_greedy= 30.68 ucb1= 61.23 thompson= 35.52
t= 1000: eps_greedy= 50.41 ucb1= 104.26 thompson= 49.75
t= 2000: eps_greedy= 83.82 ucb1= 169.98 thompson= 64.22
t= 5000: eps_greedy= 167.53 ucb1= 297.35 thompson= 81.25
\(t=100\) 時点ではトンプソンサンプリングとε-greedyはほぼ互角ですが、\(t=1000\) を境にトンプソンサンプリングが優位に転じ、差が開き続けます。実際、各手法の「終盤の傾き」(1ステップあたりの追加リグレット)を \(t\in[500,1000]\) と \(t\in[4000,5000]\) で比較すると、この傾向が定量的に確認できます。
eps_greedy slope[500:1000]=0.0395/step slope[4000:5000]=0.0259/step
ucb1 slope[500:1000]=0.0860/step slope[4000:5000]=0.0351/step
thompson slope[500:1000]=0.0285/step slope[4000:5000]=0.0041/step
トンプソンサンプリングの傾きは0.0285から0.0041へと約7倍縮小しており、リグレットの増加が対数的に鈍化している(累積リグレットの理論的下限 \(O(\log T)\) に近い挙動)ことが伺えます。一方ε-greedyは0.0395から0.0259と1.5倍程度しか縮小せず、固定の探索率 \(\varepsilon\) による一定確率のランダム探索が終盤まで続くため、傾きが完全にはゼロに近づきません。このデータからだけでは真に線形かどうかは断定できませんが、傾きの減少がトンプソンサンプリングよりはるかに緩やかであることは明確です。
なお本実験のUCB1がε-greedyより悪化している点は一見意外に見えますが、原因は腕同士の成功確率の差(ギャップ)が最大でも0.05刻みと小さいことにあります。UCB1の信頼区間ボーナス \(\sqrt{2\ln t/n_i(t)}\) は理論保証のために全腕に対して十分保守的に設計されており、ギャップが小さい設定では「まだ十分自信が持てない」と判断して次善の腕への引き直しを続けやすく、5000ステップ程度では収束しきりません。ギャップが大きい設定(例えば \(p_i\) の差を0.2刻みにする)ではUCB1の性能は大きく改善します。

Beta事後分布の収束過程
トンプソンサンプリングの核心である「事後分布の更新」を可視化します。真の成功確率 \(p=0.35\)
の腕を固定して繰り返し引き、試行回数 \(n \in \{0, 1, 5, 20, 100, 500\}\)
の各時点でのBeta事後分布 \(\text{Beta}(1+s, 1+f)\)
をプロットしました(seed=20260729)。
from scipy.stats import beta as beta_dist
true_p = 0.35
rng = np.random.default_rng(20260729)
outcomes = (rng.random(500) < true_p).astype(int)
cum_successes = np.cumsum(outcomes)
for n in [0, 1, 5, 20, 100, 500]:
s = cum_successes[n - 1] if n > 0 else 0
f = n - s
alpha, beta_param = 1 + s, 1 + f
mean = alpha / (alpha + beta_param)
print(f"n={n}: Beta({alpha},{beta_param}) mean={mean:.4f}")
n= 0 successes= 0 Beta(1,1) mean=0.5000 std=0.2887
n= 1 successes= 1 Beta(2,1) mean=0.6667 std=0.2357
n= 5 successes= 3 Beta(4,3) mean=0.5714 std=0.1750
n= 20 successes= 7 Beta(8,14) mean=0.3636 std=0.1003
n= 100 successes= 35 Beta(36,66) mean=0.3529 std=0.0471
n= 500 successes= 188 Beta(189,313) mean=0.3765 std=0.0216

\(n=0\) (一様事前分布)では標準偏差0.289と非常に広い分布ですが、\(n=20\) で標準偏差0.100、\(n=500\) で標準偏差0.022まで収縮し、事後平均も真値0.35の近くに収束しています(\(n=20\) 時点ではたまたま7/20=0.35とちょうど一致する経験率が観測されたため事後平均0.364が真値に近く出ていますが、\(n=1\) の時点では1/1=1.0という極端な観測から事後平均0.667とかなりずれていることにも注意してください——これは事前分布 \(\text{Beta}(1,1)\) の効果で極端な値には抑制がかかっているものの、少数の観測ではまだ大きなブレが残るという、ベイズ推論の基本的な性質を表しています)。この「データが少ないうちは事後分布が広く、多くのデータを見ると自然に収束する」という性質こそが、トンプソンサンプリングが探索と活用を自動でバランスさせる仕組みの直接的な可視化になっています。
発展:ガウス過程 + トンプソンサンプリング(GP-TS)
ここまでのBernoulliバンディットは腕が離散的で有限個でしたが、https://yuhi-sa.github.io/posts/20260223_bayesian_optimization/1/ で扱ったベイズ最適化は連続空間 \(\mathcal{X} \subseteq \mathbb{R}^d\) 上の最適化です。トンプソンサンプリングの考え方は、離散的な腕の代わりにhttps://yuhi-sa.github.io/posts/20260502_gaussian_process/1/を代理モデルとして使うことで、連続空間にそのまま拡張できます。これが**GP-Thompson Sampling(GP-TS)**です。
GP-TSの手続きは次の通りです。
- これまでの観測 \(\mathcal{D} = \{(\mathbf{x}_i, y_i)\}\) でガウス過程を学習する
- 候補点集合 \(X_{\text{cand}}\) 上でGP事後分布の関数を1つサンプリングする(多変量正規分布 \(\mathcal{N}(\boldsymbol{\mu}, \boldsymbol{\Sigma})\) からのサンプリングであり、Bernoulliバンディットの式(5)でBeta分布から1点サンプリングしたことの連続空間版に相当します)
- サンプルされた関数を最小化(最大化)する点を次の評価点として選ぶ
EI・UCB・PIが予測平均・予測分散を明示的な数式(式7-10、https://yuhi-sa.github.io/posts/20260223_bayesian_optimization/1/参照)に代入して次の点を決めるのに対し、GP-TSは事後分布から実現値をサンプリングし、そのサンプル自体を最適化するという点が本質的に異なります。ここで重要なのは、式(6)のサンプリングには予測分散だけでなく、候補点間の共分散(GP posteriorの非対角成分)が必要になることです。EIなどの獲得関数は各点の周辺分散だけで計算できましたが、GP-TSは「候補点全体で1つの一貫した関数」をサンプリングする必要があるため、共分散行列全体を使います。
Python実装
class GaussianProcess:
def __init__(self, length_scale=1.0, signal_var=1.0, noise_var=1e-6):
self.l, self.sf2, self.sn2 = length_scale, signal_var, noise_var
def kernel(self, X1, X2):
dist_sq = (np.sum(X1**2, axis=1, keepdims=True)
- 2.0 * X1 @ X2.T + np.sum(X2**2, axis=1))
return self.sf2 * np.exp(-0.5 * dist_sq / self.l**2)
def fit(self, X, y):
self.X_train, self.y_train = X.copy(), y.copy()
K = self.kernel(X, X) + self.sn2 * np.eye(len(X))
self.K_inv = np.linalg.inv(K)
def predict_full_cov(self, X_new):
"""予測平均と*full*共分散行列(GP-TSの関数サンプリングに必要)"""
k_star = self.kernel(self.X_train, X_new)
mu = k_star.T @ self.K_inv @ self.y_train
cov = self.kernel(X_new, X_new) - k_star.T @ self.K_inv @ k_star
return mu, cov + 1e-8 * np.eye(len(X_new)) # ジッタで数値安定化
def gp_thompson_step(gp, X_cand, rng):
"""式(6):GP事後分布から1関数サンプリングし、その最小点を返す"""
mu, cov = gp.predict_full_cov(X_cand)
f_sample = rng.multivariate_normal(mu, cov)
return X_cand[np.argmin(f_sample)]
GaussianProcess クラスはhttps://yuhi-sa.github.io/posts/20260223_bayesian_optimization/1/のクラスと同じkernel・fitを使い、predict(周辺分散のみ)の代わりにpredict_full_cov(共分散行列全体)を追加しています。gp_thompson_step が式(6)の実装で、rng.multivariate_normal が「事後分布から1つの関数を引く」操作に対応します。
数値実験:EIとの比較
https://yuhi-sa.github.io/posts/20260223_bayesian_optimization/1/と同じテスト関数 \(f(x) = \sin(3x) + 0.1x^2 - 0.5\cos(7x)\)
、探索範囲 \(x\in[-3,3]\)
、初期点3点+逐次評価12回(計15回評価)、seed=42で、EIとGP-TSを同一の初期点から実行し比較しました。
EI : best_x=1.6894 best_f=-1.0210 regret=0.0172 n_evals=15
GP-TS : best_x=-0.7876 best_f=-0.9993 regret=0.0388 n_evals=15
真の大域最小値は \(f(\mathbf{x}^*) = -1.038189\) (https://yuhi-sa.github.io/posts/20260223_bayesian_optimization/1/参照)なので、この1条件ではEIのregretが0.0172、GP-TSのregretが0.0388と、EIの方が真の最適値に近い解を得ました。

上図の上段は、8個の観測点を持つGP事後分布から3本の関数をサンプリングした様子です。観測点付近(不確実性が低い領域)では3本のサンプルがGP平均に近い値に収束する一方、観測がまばらな領域(\(x\in[-0.5, 1]\) 付近)ではサンプルごとに大きくばらついています。この「不確実性が高い領域でサンプル値が大きく散らばる」という性質が、GP-TSの探索メカニズムそのものです。下段は、同一の初期点からEIとGP-TSがどう収束するかを示しており、この1回の実行ではGP-TSが評価8回目・10回目で改善を続けEIより早く \(f\approx-0.999\) に到達しているものの、最終的にはEIが評価12回目でより深い最小値 \(-1.0210\) を発見しています。
複数シードでの追加検証
単一シードの結果だけから「EIの方が優れている」と結論づけるのは早計です(https://yuhi-sa.github.io/posts/20260223_bayesian_optimization/1/の実践上の注意点でも述べた通り)。そこで同じ設定(初期点3・逐次評価12回)をseed=100〜129の30シードで実行し、regretの分布を確認しました。
n_seeds=30
EI regret: mean=0.0146 median=0.0022 std=0.0174
GP-TS regret: mean=0.0640 median=0.0388 std=0.0918
GP-TS better in 26.7% of seeds
30シード中GP-TSがEIを上回った(regretが小さかった)のはわずか26.7%で、平均regretもEI(0.0146)がGP-TS(0.0640)よりおよそ4.4倍小さいという結果になりました。この特定の低評価予算(15回)・1次元・単峰性の低いテスト関数の設定では、EIの方が頑健に良い解を見つける傾向があると言えます。これは、GP-TSが式(6)で「関数を1本サンプリングしてその最小点を選ぶ」ため、たまたま外れ値的なサンプルを引いてしまうと明後日の方向を評価してしまうのに対し、EIは予測平均・分散という2つの統計量に基づく決定論的な基準で常に「改善期待値が最大の点」を選ぶため、評価回数が少ない設定ではばらつきが小さくなりやすいことに起因すると考えられます。一方でGP-TSは、多腕バンディットの節で見た通り逐次的な意思決定を多数回繰り返す設定(評価予算が大きい、あるいは非同期・並列評価が必要な設定)ではEIやUCBに劣らない性能を示すことが知られており(Chowdhury & Gopalan, 2017; Kandasamy et al., 2018)、本実験の「評価回数15回のみ」という設定はGP-TSにとって不利な条件である可能性が高い点には注意してください。
トンプソンサンプリングとUCB系獲得関数の使い分け
最後に、本記事で見た内容を踏まえて、トンプソンサンプリングが実務でどう位置づけられるかを整理します。
| ε-greedy | UCB1 / UCB | トンプソンサンプリング | |
|---|---|---|---|
| 探索の仕組み | 固定確率でランダム探索 | 決め打ちの信頼区間ボーナス | 事後分布からのサンプリング |
| 実装の複雑さ | 最も単純 | 単純(式が閉形式) | 中程度(事後分布の管理が必要) |
| 理論的リグレット | \(\Theta(T)\) (線形) | \(O(\log T)\) | \(O(\log T)\) (Bernoulliの場合) |
| 並列・バッチ評価 | 容易 | やや工夫が必要 | 自然に拡張しやすい(複数サンプル) |
| 連続空間への拡張 | 容易だが非効率 | 理論保証が複雑化 | GP-TSとして自然に拡張(本記事) |
トンプソンサンプリングの実務上の強みは、複数の候補を並行して評価するバッチ設定に自然に拡張できる点です。UCB系は「1点だけを決め打ちで選ぶ」ロジックのため、バッチ内の多様性を確保する工夫が別途必要になりますが、トンプソンサンプリングは事後分布から複数回独立にサンプリングするだけで、自然に多様な候補点の集合が得られます(それぞれ異なる乱数シードでサンプルすればよい)。この性質から、広告配信・レコメンドシステムなど大量のリクエストを並列に処理する実務のオンライン学習システムで広く採用されています(Chapelle & Li, 2011; Russo et al., 2018)。
まとめ
- 多腕バンディット問題における3手法(ε-greedy・UCB1・トンプソンサンプリング)を numpy でフルスクラッチ実装し、\(K=8\) 本のBernoulliバンディット・\(T=5000\) ステップ・200シード平均で累積リグレットを比較した結果、トンプソンサンプリング(81.25)がε-greedy(167.53)の約半分、UCB1(297.35)の約4分の1という最小の累積リグレットを達成した。
- トンプソンサンプリングの核心は、腕ごとの成功確率に対するBeta事後分布 \(\text{Beta}(1+s,1+f)\) を保持し、そこから1点サンプリングして最大の腕を選ぶという手続きにある。実測では試行回数0から500まで、事後分布の標準偏差が0.289から0.022まで単調に収縮する様子を確認した。
- ガウス過程と組み合わせたGP-Thompson Samplingは、連続空間の最適化にトンプソンサンプリングの考え方を拡張したものである。ベイズ最適化のテスト関数で1条件・30シードの追加検証を行った結果、本実験の低評価予算(15回)設定ではEIの方がGP-TSより平均regretが小さく(0.0146 対 0.0640)、頑健だった。ただしこれは低評価予算・低次元という条件下の結果であり、逐次的意思決定を大量に繰り返す設定ではトンプソンサンプリングが優位になりうることが知られている。
関連記事
- ベイズ最適化の基礎とPython実装 - 本記事のGP-TSは、この記事で解説したEI・UCB・PIと並ぶ第4の獲得関数として位置づけられます。
- ガウス過程回帰の理論とPython実装 - GP-TSの代理モデルであるガウス過程の理論的基盤を解説しています。
- ガウス過程回帰の実践:カーネル設計と応用 - Cholesky分解による数値安定な実装など、実務レベルのGP実装の詳細です。
- Hamiltonian Monte Carlo(HMC)の理論とPython実装 - 事後分布からのサンプリングという発想を、より高次元・高相関な分布に対して効率化する手法です。
- モンテカルロ最適化メソッド比較(CEM/SA/GA/MPPI/PSO) - サンプリングベースの最適化手法を横断的に比較するハブ記事です。
参考文献
- Thompson, W. R. (1933). On the likelihood that one unknown probability exceeds another in view of the evidence of two samples. Biometrika, 25(3/4), 285-294.
- Auer, P., Cesa-Bianchi, N., & Fischer, P. (2002). Finite-time analysis of the multiarmed bandit problem. Machine Learning, 47(2-3), 235-256.
- Chapelle, O., & Li, L. (2011). An empirical evaluation of Thompson sampling. NeurIPS 2011.
- Russo, D. J., Van Roy, B., Kazerouni, A., Osband, I., & Wen, Z. (2018). A tutorial on Thompson sampling. Foundations and Trends in Machine Learning, 11(1), 1-96.
- Chowdhury, S. R., & Gopalan, A. (2017). On kernelized multi-armed bandits. ICML 2017.
- Kandasamy, K., Krishnamurthy, A., Schneider, J., & Póczos, B. (2018). Parallelised Bayesian optimisation via Thompson sampling. AISTATS 2018.