Thompson Sampling on Symmetric Alpha-Stable Bandits

1 Introduction

The multi-armed bandit (MAB) problem is a fundamental model in understanding the exploration-exploitation dilemma in sequential decision-making. The problem and several of its variants have been studied extensively over the years, and a number of algorithms have been proposed that optimally solve the bandit problem when the reward distributions are well-behaved, i.e. have a finite support, or are subexponential.

The most prominently studied class of algorithms are the Upper Confidence Bound (UCB) algorithms, that employ an "optimism in the face of uncertainty" heuristic [Auer et al.2002], which have been shown to be optimal (in terms of regret) in several cases [Capp´e et al.2013, Bubeck et al.2013]. Over the past few years, however, there has been a resurgence in interest in the Thompson Sampling (TS) algorithm [Thompson1933], that approaches the problem from a Bayesian perspective.

∗ Full version (with appendix) available at this link.

Rigorous empirical evidence in favor of TS demonstrated by [Chapelle and Li2011] sparked new interest in the theoretical analysis of the algorithm, and the seminal work of [Agrawal and Goyal2012,Agrawal and Goyal2013,Russo and Van Roy2014] demonstrated the optimality of TS when rewards are bounded in [0, 1] or are Gaussian. These results were extended in the work of [Korda et al.2013] to more general, exponential family reward distributions. The empirical studies, along with theoretical guarantees, have established TS as a powerful algorithm for the MAB problem.

However, when designing decision-making algorithms for complex systems, we see that interactions in such systems often lead to heavy-tailed and power law distributions, such as modeling stock prices [Bradley and Taqqu2003], preferential attachment in social networks [Mahanti et al.2013], and online behavior on websites [Kumar and Tomkins2010].

Specifically, we consider a family of extremely heavytailed reward distributions known as α-stable distributions. This family refers to a class of distributions parameterized by the exponent α, that include the Gaussian (α = 2), L´evy (α = 1/2) and Cauchy (α = 1) distributions, all of which are used extensively in economics [Frain2009], finance [Carr and Wu2003] and signal processing [Shao and Nikias1993].

The primary hurdle in creating machine learning algorithms that account for α-stable distributions, however, is their intractable probability density, which cannot be expressed analytically. This prevents even a direct evaluation of the likelihood under this distribution. Their heavy-tailed nature, additionally, often leads to standard algorithms (such as Thompson Sampling assuming Gaussian rewards), concentrating on incorrect arms.

In this paper, we create two algorithms for Thompson Sampling under symmetric α-stable rewards with finite means. Our contributions can be summarized as follows:

1. Using auxiliary variables, we construct a framework for posterior inference in symmetric α-stable bandits that leads to the first efficient algorithm for Thompson Sampling in this setting, which we call α-TS.
2. To the best of our knowledge, we provide the first finite-time polynomial bound on the Bayesian Regret of Thompson Sampling achieved by α-TS in this setting.
3. We improve on the regret by proposing a modified Thompson Sampling algorithm, called Robust α-TS,

that utilizes a truncated mean estimator, and obtain the first ˜O(N 1 1+ϵ ) Bayesian Regret in the α-stable setting. Our bound matches the optimal bound for α ∈ (1, 2) (within logarithmic factors).

1. Through a series of experiments, we demonstrate the proficiency of our two Thompson Sampling algorithms for α-stable rewards, which consistently outperform all existing benchmarks.

Our paper is organized as follows: we first give a technical overview of the MAB problem, the Thompson Sampling algorithm and α-stable distributions. Next, we provide the central algorithm α-TS and its analysis, followed by the same for the Robust α-TS algorithm. We then provide experimental results on multiple simulation benchmarks and finally, discuss the related work in this area prior to closing remarks.

2 Preliminaries

2.1 Thompson Sampling

The K-Armed Bandit Problem: In any instance of the Karmed bandit problem, there exists an agent with access to a set of K actions (or "arms"). The learning proceeds in rounds, indexed by t ∈ [1, T]. The total number of rounds, known as the time horizon T , is known in advance. The problem is iterative, wherein for each round t ∈ [T ]:

1. Agent picks armAt ∈ [K].
2. Agent observes reward rAt (t) from that arm.

For arm k ∈ [K], rewards come from a distribution Dk with mean µk = EDk [r]1 . The largest expected reward is denoted by µ∗ = maxk∈[K] µk, and the corresponding arm(s) is denoted as the optimal arm(s) k∗ . In our analysis, we will focus exclusively on the i.i.d. setting, that is, for each arm, rewards are independently and identically drawn from Dk, every time arm k is pulled.

To measure the performance of any (possibly randomized) algorithm we utilize a measure known as Regret R(T ), which, at any round T, is the difference of the cumulative mean reward of the algorithm against the expected reward of always playing an optimal arm.

Thompson Sampling (TS): Thompson Sampling [Thompson1933] proceeds by maintaining a posterior distribution over the parameters of the bandit arms. If we assume that for each arm k, the reward distribution Dk is parameterized by a (possibly vector) parameter θk that come from a set Θ with a prior probability distribution p(θk) over the parameters, the Thompson Sampling algorithm proceeds by selecting arms based on the posterior probability of the reward under the arms. For each round t ∈ [T], the agent:

1α-stable distributions with α ≤ 1 do not admit a finite first moment. To continue with existing measures of regret, we only consider rewards with α > 1.

1. Draws parameters ˆ θk(t) for each arm k ∈ [K] from the posterior distribution of parameters, given the previous rewards rk(t-1) = {r(1) k , r(2) k , ...} till round t-1 (note that the posterior distribution for each arm only depends on the rewards obtained using that arm). When t = 1, this is just the prior distribution over the parameters.

1. Given the drawn parameters ˆ θk(t) for each arm, chooses arm at with the largest mean reward over the posterior distribution.

1. Obtains reward rt after taking action At and updates the posterior distribution for arm At.

In the Bayesian case, the measure for performance we will utilize in this paper is the Bayes Regret (BR) [Russo and Van Roy2014], which is the expected regret over the priors. For any policy π, this is given by:

While the regret provides a stronger analysis, any bound on the Bayes Regret is essentially a bound on the expected regret, since if an algorithm admits a Bayes Regret of O(g(T )), then its Regret is also stochastically bounded by g(·) [Russo and Van Roy2014].

2.2 α-Stable Distributions

α-Stable distributions, introduced by L´evy [L´evy1925] are a class of probability distributions defined over R whose members are closed under linear transformations.

Definition 1 ( [Borak et al.2005]). Let X1 and X2 be two independent instances of the random variable X. X is stable if, for a1 > 0 and a2 > 0, a1X1 + a2X2 follows the same distribution as cX + d for some c > 0 and d ∈ R.

A random variable X ∼ Sα(β, µ, σ) follows an α-stable distribution described by the parameters α ∈ (0, 2] (characteristic exponent) and β ∈ [-1, 1] (skewness), which are responsible for the shape and concentration of the distribution, and parameters µ ∈ R (shift) and σ ∈ R+ (scale) which correspond to the location and scale respectively. While it is not possible to analytically express the density function for generic α-stable distributions, they are known to admit the characteristic function φ(x; α, β, σ, µ):

where Φα(x) is given by

For fixed values of α, β, σ and µ we can recover the density function from φ(·) via the inverse Fourier transform:

Algorithm 1 Chambers-Mallows-Stuck Generation

Most of the attention in the analysis of α-stable distributions has been focused on the stability parameter α, which is responsible for the tail "fatness". It can be shown that asymptotically, the tail behavior (x → ±∞) of X ∼ Sα(β, µ, σ) follows [Borak et al.2005]:

where α < 2 and Γ(·) denotes the Gamma function. The power-law relationship admitted by the density is responsible for the heaviness of said tails.

From Lemma 1 it follows that α-stable random variables only admit a finite mean for α > 1, and also admit a finite variance only when α = 2, which corresponds to the family of Gaussian distributions. To continue with existing measures of regret, we restrict our analysis hence to α-stable distributions only with α > 1. Note that for all our discussions, 1 < α < 2, hence, all distributions examined are heavytailed, with infinite variance. Additionally, we restrict ourselves to only symmetric (β = 0) distributions: asymmetric distributions do not allow a scaled mixture representation (which is the basis of our framework, see Section 3.1).

Sampling from α-Stable Densities

For general values of α, β, σ, µ, it is not possible to analytically express the density of α-stable distributions, and hence we resort to using auxiliary variables for sampling. The Chambers-Mallows-Stuck [Chambers et al.1976] algorithm is a straightforward method to generate samples from the density Sα(β, 1, 0) (for α ̸= 1) via a non-linear transformation of a uniform random variable V and an exponential random variable W, which can then be re-scaled to obtain samples from Sα(β, σ, µ) (Algorithm 1).

Products of α-Stable Densities

A central property of α-stable densities that will be utilized in the following section is the behavior of products of independent α-stable variables.

Lemma 2 ( [Borak et al.2005]). Let Z and Y > 0 be independent random variables such that Z ∼ Sγ(0, σ1, µ1) and Y ∼ Sδ(1, σ2, µ2). Then ZY 1/γ is stable with exponent γδ.

We now begin with the algorithm description and analysis.

3 α-Thompson Sampling

We consider the setting where, for an arm k, the corresponding reward distribution is given by Dk = Sα(0, σ, µk) where α ∈ (1, 2), σ ∈ R+ are known in advance, and µk is unknown 2 . We set a prior distribution over µk. We now derive the form of the prior distribution, and outline an algorithm for Bayesian inference.

3.1 Scale Mixtures of Normals

On setting γ = 2 (Gaussian) and β = α/2 < 1 in Lemma 2, the product distribution X = ZY 1/2 is stable with exponent α. This property is an instance of the general framework of scale mixtures of normals (SMiN) [Andrews and Mallows1974], which are described by the following:

This framework contains a large class of heavy-tailed distributions which include the exponential power law, Student'st and symmetric α-stable distributions [Godsill and Kuruoglu1999]. The precise form of the variable X depends on the mixing distribution pΛ(λ). For instance, when pΛ is the inverted Gamma distribution (the conjugate prior for a unknown variance Gaussian), the resulting pX follows a Student's-t distribution.

Bayesian inference directly from Sα(0, σ, µ) is difficult: the non-analytic density prevents a direct evaluation of the likelihood function, and the non-Gaussianity introduces difficulty in its implementation. However, the SMiN representation enables us to draw samples directly from Sα(0, σ, µ) using the auxiliary variable λ:

This sampling assists in inference since x is Gaussian conditioned on λ: given samples of λ, we can generate x from the induced Gaussian conditional distribution.

3.2 Posteriors for α-Stable Rewards

Let us examine approximate inference for a particular arm k ∈ [K]. At any time t ∈ [T], assume this arm has been pulled nk(t) times previously, and hence we have a vector of reward samples rk(t) = {r(1) k , ..., r(nk(t)) k } observed until time t. Additionally, assume we have nk(t) samples of an auxiliary variable λk(t) = {λ(1) k , ..., λ(nk(t)) i } where λk ∼ Sα/2(1, 1, 0).

Recall that rk ∼ Sα(0, σ, µk) for an unknown (but fixed) µk. From the SMiN representation, we know that rk is conditionally Gaussian given the auxiliary variable λk, that is p(rk|λk, µk) ∼ N(µk, λkσ2 ). We can then obtain the conditional likelihood as the following:

2 Typically, more general settings have been explored for TS in the Gaussian case, where the variance is also unknown. Our algorithm can also be extended to an unknown scale (σ) using an inverse Gamma prior, similar to [Godsill and Kuruoglu1999].

We can now assume the conjugate prior over µk, which is a normal distribution with mean µ0k and variance σ2 . We then obtain the posterior density for µk as (full derivation in the Appendix):

We know that ˆσ2 k (t) > 0 since λ(i) k > 0 ∀i as they are samples from a positive stable distribution (β = 1). Given rk(t) and µk, we also know that the individual elements of λk(t) are independent, which provides us with the following decomposition for the conditional density of λk(t),

Here, fα,β(·) is the density of a random variable following Sα(β, 1, 0). Our stepwise posterior sampling routine is hence as follows. At any time t, after arm k is pulled and we receive reward rnk(t) k , we set rk(t) = [rk(t - 1), rnk(t) k ]. Then for a fixed Q iterations, we repeat:

1. For i ∈ [1, nk(t)], draw λ(i) k ∼ p(λ(i) k |r(i) k , µk(t)).
2. Draw µk(t) ∼ p(µk|rk(t), λk(t)) .

Sampling from the conditional posterior of µk is straightforward since it is Gaussian. To sample from the complicated posterior of λ(i) k , we utilize rejection sampling.

3.3 Rejection Sampling for λ(i) k

Sampling directly from the posterior is intractable since it is not analytical. Therefore, to sample λ(i) k we follow the pipeline described in [Godsill and Kuruoglu1999]. We note that the likelihood of the mean-normalized reward v(i) k = r(i) k - µk(t) forms a valid rejection function since it is bounded:

Since v(i) k ∼N(0; λ(i) k σ2 ). Thus, we get the procedure:

1. Draw λ(i) k ∼ Sα/2(1, 1, 0) (using Algorithm 1).
2. Draw u ∼ U (0, (v(i) k √2π)-1 exp(-1/2)).
3. If u > p(v(i) k |0, λ(i) k σ2 ), reject λ(i) k and go to Step 1.

Combining all these steps, we can now outline our algorithm, α-Thompson Sampling (α-TS) as described in Algorithm 2.

It is critical to note that in Algorithm 2, we do not draw from the full vector of λk(t) at every iteration, but only from the last obtained reward. This is done to accelerate the inference process, and while it leads to a slower convergence of

Algorithm 2 α-Thompson Sampling

• 1: Input: Arms k ∈ [K], priors N(µ0k , σ2 ) for each arm.
• 2: Set Dk = 1, Nk = 0 for each arm k.
• 3: for For each iteration t ∈ [1, T] do
• 4: Draw ¯µk(t) ∼N ( µ0k+Nk Dk , σ2 Dk ) for each arm k.
• 5: Choose armAt = argmax k∈[K] ¯µk(t), and get reward rt.
• 6: for q ∈ [0, Q) do
• 7: Calculate v(t) At = rt - ¯µAt .
• 8: Draw λ(t) At following Section 3.3.

• 11: end for

• 13: end for

the posterior, we observe that it performs well in practice. Alternatively, one can re-sample λk(t) over a fixed window of the previous rewards, to prevent the runtime from increasing linearly while enjoying faster convergence.

3.4 Regret Analysis

In this section, we derive an upper bound on the finite-time Bayesian Regret (BR) incurred by the α-TS algorithm. We continue with the notation used in previous sections, and assume a K armed bandit with T maximum trials. Each arm k follows an α-stable reward Sα(0, σ, µk), and without loss of generality, let µ∗ = maxk∈[K] µi denote the arm with maximum mean reward.

Theorem 1 (Regret Bound). Let K > 1, α ∈ (1, 2), σ ∈ R+, µk:k∈[K] ∈ [-M,M]. For a K-armed stochastic bandit with rewards for each arm k drawn from Sα(0, σ, µk), we have, asymptotically, for ϵ chosen a priori such that ϵ → (α1)-,

Proof-sketch. We first utilize the characteristic function representation of the probability density to obtain the centered (1 + ϵ)th moment of the reward distribution for each arm. We use this moment to derive a concentration bound on the deviation of the empirical mean of α-stable densities. Next, we proceed with the decomposition of the Bayesian Regret in terms of upper-confidence bounds, as done in [Russo and Van Roy2014]. The complete proof is included in detail in the appendix.

We note the following: First, the only additional assumption we make on the reward distributions is that the true means are bounded, which is a standard assumption [Russo and Van Roy2014] and easily realizable in most application cases. Next, ϵ < α-1 must be chosen carefully to control the finite-time regret. As ϵ → α-1, we see that while the growth of T 2 1+ϵ decreases, the constants in the finite-time expression grow, and are not finite at ϵ = α - 1. This behavior arises

Algorithm 3 Robust α-Thompson Sampling

• 1: Input: Arms k ∈ [K], priors N(µ0k , σ2 ) for each arm.
• 2: Set Dk = 1, Nk = 0 for each arm k.
• 3: for For each iteration t ∈ [1, T] do
• 4: Draw ¯µk(t) ∼N ( µ0k+Nk Dk , σ2 Dk ) for each arm k.
• 5: Choose armAt = argmaxk ¯µk(t), and get reward rt.

• 7: Repeat steps 6-12 of Algorithm 2.
• 8: end for

from the non-compactness of the set of finite moments for α-stable distributions (see Appendix for detailed analysis).

Compared to the problem-independent regret bound of O(√KT log T ) for Thompson Sampling on the multi-armed Gaussian bandit problem demonstrated by [Agrawal and Goyal2013], our bound differs in two aspects: first, we admit a K 1 1+ϵ complexity on the number of arms, which contrasted with the Gaussian bandit is identical when ϵ → 1. Second, we have a superlinear dependence of order T 2 1+ϵ . The termof T 1 1+ϵ approaches the Gaussian case (√T) when ϵ → 1, leaving us with the additional sub-linear term (T 1 1+ϵ ) compared to the sub-logarithmic term (√log T ) from the Gaussian case.

In the next section, we address the issue of polynomial concentration by utilizing the more robust, truncated mean estimator instead of the empirical mean, and obtain a modified, robust version of α-TS.

3.5 Robust α-Thompson Sampling

Assume that for all arms k ∈ [K], µk ≤ M. Note that this assumption is equivalent to the boundedness assumption in the analysis of α-TS, and is a reasonable assumption to make in any practical scenario with some domain knowledge of the problem. Let δ ∈ (0, 1), ϵ ∈ (0, α - 1). Now, consider the estimator ˆ r∗ k (t) given by:

ˆ r∗ k (t) then describes a truncated mean estimator for an arm k (pulled nk(t) times), where a reward r(i) k at any trial i of the arm is discarded if it is larger than the bound. We choose such a form of the truncation since it allows us to obtain an exponential concentration for α-stable densities, which will eventually provide us tighter regret.

As can be seen, the corresponding Robust α-TS algorithm is identical to the basic α-TS algorithm, except for this step of rejecting a reward if it is greater than the bound. The algorithm is outlined in Algorithm 3.

Theorem 2 (Regret Bound). Let K > 1, α ∈ (1, 2), σ ∈ R+, µk:k∈[K] ∈ [ -M,M]. For a K-armed bandit with rewards for each arm k drawn from Sα(0, σ, µk), we have, for ϵ chosen a priori such that ϵ → (α - 1)- and truncated estimator given in Equation (11),

Proof-sketch. The truncated mean estimator can be used to obtain a tighter concentration bound on the empirical mean. This can be understood intuitively as by carefully rejecting values that are distant from the mean, our empirical mean is robust to outliers, and would require less samples to concentrate around the true mean. Using the concentration result, we follow an analysis identical to Theorem 1. The full proof is deferred to the Appendix for brevity.

We see that this bound is tight (upto logarithmic factors): it matches the problem-independent lower bound of O((KT) 1 1+ϵ ) [Bubeck et al.2013]. Algorithm 3's improvements increase as α decreases: the likelihood of obtaining confounding outliers increases as α → 1, and can perturb the posterior mean in the naive α-TS algorithm. In the next section, we discuss some experimental results that cement the value of α-TS for α-stable bandits.

4 Experiments

We conduct simulations with the following benchmarks - (i) an ε-greedy agent with linearly decreasing ε, (ii) Regular TS with Gaussian priors and a Gaussian assumption on the data (Gaussian-TS), (iii) Robust-UCB [Bubeck et al.2013] for heavy-tailed distributions using a truncated mean estimator, and (iv) α-TS and (v) Robust α-TS, both with Q(iterations for sampling) as 50.

Setting Priors for TS: In all the Thompson Sampling benchmarks, setting the priors are crucial to the algorithm's performance. In the competitive benchmarking, we randomly set the priors for each arm from the same range we use for setting the mean rewards.

Performance against Competitive Benchmarks

We run 100 MAB experiments each for all 5 benchmarks for α = 1.8 and α = 1.3, and K = 50 arms, and for each arm, the mean is drawn from [0, 2000] randomly for each experiment, and σ = 2500. Each experiment is run for T = 5000 iterations, and we report the regret averaged over time, i.e. R(t)/t at any time t.

In Figures 1A and 1B, we see the results of the regret averaged over all 100 experiments. We observe that for α = 1.3, there are more substantial variations in the regret (low α implies heavy outliers), yet both algorithms comfortably outperform the other baselines.

In the case of α = 1.8, the variations are not that substantial, the performance follows the same trend. It is important to note that regular Thompson Sampling (Gaussian-TS) performs competitively, although in our experiments, we observed that when K is large, the algorithm often concentrates on the incorrect arm, and subsequently earns a larger regret.

Intuitively, we can see that whenever arms have mean rewards close to each other (compared to the variance), ϵgreedy will often converge to the incorrect arm. Robust-UCB, however, makes very weak assumptions on the data distributions, and hence has a much larger exploration phase, leading to larger regret. Compared to α-TS, Robust α-TS is less affected by outlying rewards (as is more evident in α = 1.3 vs. α = 1.8) and hence converges faster.

Ablation Studies

We additionally run two types of ablation studies - first, we compare the performance of α-TS on the identical set up as before (same K and reward distributions), but with varying values of α ∈ (1, 2). We report the expected time-normalized regret averaged over 100 trials in Figure 1C, and observe that (i) as α increases, the asymptotic regret decreases faster, and (ii) as expected, for lower values of α there is a substantial amount of variation in the regret.

Secondly, we compare the effect of the sharpness of the priors. In the previous experiments, the priors are drawn randomly from the same range as the means, without any additional information. However, by introducing more information about the means through the priors, we can expect better performance. In this experiment, for each mean µk, we randomly draw the prior mean µ0k from [µk - δ, µk + δ], and observe the regret after T = 15K trials for δ from 50 to 1000. The results for this experiment are summarized in Figure 1D for K = 10 and σ = 25, and results are averaged over 25 trials each. We see that with uninformative priors, α-TS performs competitively, and gets better as the priors get sharper.

5 RelatedWork

A version of the UCB algorithm [Auer et al.2002] has been proposed in [Bubeck et al.2013] coupled with several robust mean estimators to obtain Robust-UCB algorithms with optimal problem-dependent (i.e. dependent on individual µks) regret when rewards are heavy-tailed. However, the optimal version of their algorithm has a few shortcomings that α-TS addresses: their algorithm requires the median-ofmeans estimator, which has an expensive space complexity of O(log T ) and time complexity of O(log log T ) per update. Secondly, there is no mechanism to incorporate prior information, which can be advantageous, as seen even with weak priors. Additionally, [Vakili et al.2013] introduce a deterministic exploration-exploitation algorithm, which achieves same order regret as Robust-UCB for heavy-tailed rewards.

There has been work in analysing Thompson Sampling for specific Pareto and Weibull heavy-tailed distributions in [Korda et al.2013], however, the Weibull and Pareto distributions typically have "lighter" tails owing to the existence of more higher order moments, and hence cannot typically be used to model very heavy tailed signals.

In related problems, [Yu et al.2018] provide a purely exploratory algorithm for best-arm identification under heavytailed rewards, using a finite (1 + ϵ)th moment assumption. Similarly, [Shao et al.2018,Medina and Yang2016] explore heavy-tailed rewards in the linear bandit setting.

6 Conclusion

In this paper, we first designed a framework for efficient posterior inference for the α-stable family, which has largely been ignored in the bandits literature owing to its intractable density function. We formulated the first polynomial problem-independent regret bounds for Thompson Sampling with α-stable densities, and subsequently improved the regret bound to achieve the optimal regret identical to the subGaussian case, providing an efficient framework for decisionmaking for these distributions.

Additionally, our intermediary concentration results provide a starting point for other machine learning problems that may be investigated in α-stable settings. There is ample evidence to support the existence of α-stability in various modeling problems across economics [Frain2009], finance [Bradley and Taqqu2003] and behavioral studies [Mahanti et al.2013]. With tools such as ours, we hope to usher scientific conclusions in problems that cannot make sub-Gaussian assumptions, and can lead to more robust empirical findings. Future work may include viewing more involved decision-making processes, such as MDPs, in the same light, leading to more robust algorithms.

Acknowledgements

We would like to thank Peter Krafft, Yan Leng, Dhaval Adjodah, Ziv Epstein, Michiel Bakker and Ishaan Grover for their comments on this manuscript, the anonymous reviewers for their suggestions, and MIT Quest for Intelligence for their generous support with computational resources.