Thompson Sampling With Poisson Outcomes

Why Model Event Rates?

Some clinical outcomes are counts accumulated over a period of observation: infections, hospital visits, or other events. Participants may contribute different amounts of observation time, so comparing raw counts can be misleading. A participant observed for twice as long has more opportunity to experience an event, even if the underlying event rate is the same. A Poisson model represents counts relative to exposure, while a Gamma prior gives a convenient way to update uncertainty about each treatment’s event rate.

In Thompson sampling, allocation is based on posterior uncertainty, not just on which treatment has the most favorable estimated rate. For an outcome where more events are desirable, the relevant quantity is the posterior probability that an arm has the highest rate. For an adverse event, where fewer events are desirable, it is the probability that an arm has the lowest rate. These probabilities become the randomization probabilities for the next participant, subject to the trial’s specified allocation rule.

Gamma-Poisson Model and Posterior Update

For treatment arm $k$, let $Y_k$ be the observed event count and $E_k$ its known exposure, such as person-years. Let $\lambda_k$ denote the arm’s event rate per unit of exposure. The model assumes that, conditional on the rate and exposure, the count has a Poisson distribution with mean $E_k\lambda_k$. The prior below uses the shape-rate parameterization of the Gamma distribution.

Definition: The Gamma distribution with shape $a>0$ and rate $b>0$, written $\operatorname{Gamma}(a,b)$, has density $$ f(\lambda)=\frac{b^a}{\Gamma(a)}\lambda^{a-1}e^{-b\lambda},\qquad \lambda>0. $$ Its mean is $a/b$ and its variance is $a/b^2$. If $Y_k\mid\lambda_k,E_k\sim\operatorname{Poisson}(E_k\lambda_k)$ and $\lambda_k\sim\operatorname{Gamma}(a_k,b_k)$, then the posterior is $\lambda_k\mid Y_k,E_k\sim\operatorname{Gamma}(a_k+Y_k,b_k+E_k)$. Here $b_k$ and $E_k$ are expressed in the same exposure units.

The rate parameter $b_k$ is not a treatment event rate: it has units of exposure, so that $a_k/b_k$ has units of events per exposure. This distinction matters when interpreting the prior and when entering exposure into the update. A count is added to the shape; exposure is added to the Gamma rate parameter.

Why the Update Has This Form

For one arm, the Poisson likelihood as a function of $\lambda$ is proportional to $(E\lambda)^Y e^{-E\lambda}$. Terms that do not depend on $\lambda$ can be omitted when forming the posterior density. Multiplying by the Gamma prior gives the kernel below.

\begin{aligned} L(\lambda\mid Y,E)\,\pi(\lambda) &\propto \lambda^Y e^{-E\lambda}\lambda^{a-1}e^{-b\lambda}\\ &= \lambda^{a+Y-1}e^{-(b+E)\lambda}. \end{aligned}

This is the kernel of a Gamma density with shape $a+Y$ and rate $b+E$, proving the stated update. For several arms, use the corresponding count and exposure for each arm. If the prior distributions are independent across arms and the arm-specific counts are conditionally independent given their rates, the posterior distributions remain independent across arms. That independence is what permits separate posterior draws for each treatment in the basic allocation rule.

The assumptions are substantive. Exposure must be measured on a comparable and appropriate scale; the rate is assumed constant over the modeled observation period; and, conditional on the rate and exposure, counts are assumed to follow the Poisson model. The basic model also treats a participant’s event count as conditionally independent of other participants’ counts. If these assumptions do not reasonably describe the data, the posterior probabilities can be precise calculations under a poor model rather than reliable summaries of treatment uncertainty.

Worked Example: Updating a Rate After Observed Events

Suppose, hypothetically, that one arm has prior $\lambda\sim\operatorname{Gamma}(2,2)$, with exposure measured in person-years. In a completed group, there are 5 events over 10 person-years. Use the Poisson-Gamma model to update the rate and find the posterior mean.

The data are $Y=5$ and $E=10$ person-years. The prior shape and rate are $a=2$ and $b=2$ person-years. The observations are treated as one count over the stated total exposure, under the model’s constant-rate and Poisson assumptions.

The update is $$ \lambda\mid Y,E\sim\operatorname{Gamma}(a+Y,b+E) =\operatorname{Gamma}(2+5,2+10) =\operatorname{Gamma}(7,12). $$

The posterior mean rate is $$ \frac{7}{12}=0.5833 $$ events per person-year. This is a posterior mean, not the observed count divided by exposure, which is $5/10=0.5000$ events per person-year. The prior contributes information through both its shape and rate, so the posterior mean need not equal the crude observed rate.

From Posterior Rates to Allocation Probabilities

After updating each arm, draw one rate from each arm’s posterior distribution. For a higher-is-better outcome, assign the next participant to the arm whose draw is largest. For a lower-is-better outcome, assign the participant to the arm whose draw is smallest. Repeating this sampling rule defines a probability of assignment for each arm.

Definition: For $K$ treatment arms with posterior rates $\lambda_1,\ldots,\lambda_K$, the higher-is-better Thompson probability for arm $k$ is $$ \pi_k=P(\lambda_k>\lambda_j\text{ for every }j\ne k\mid\text{trial data}). $$ For a lower-is-better outcome, it is $$ \pi_k=P(\lambda_k<\lambda_j\text{ for every }j\ne k\mid\text{trial data}). $$ These probabilities sum to 1 when posterior ties have probability 0, as they do for independent continuous Gamma posteriors.

The draws represent uncertainty about the unknown rates. They are not simulated participant outcomes, and $\pi_k$ is not the probability that a participant assigned to arm $k$ will experience an event. Nor is it simply an ordering of posterior means. An arm can have a lower posterior mean but still have a meaningful chance of being best because its posterior is more uncertain.

Worked Example: Allocating Toward the Higher Rate

Consider a hypothetical trial with two arms and a higher-is-better event rate. Each arm has an independent $\operatorname{Gamma}(1,1)$ prior. Arm A has 0 events over 1 person-year; arm B has 0 events over 2 person-years. Assume the counts satisfy the Poisson model and both exposures are measured in person-years.

The posteriors are $$ \lambda_A\mid\text{data}\sim\operatorname{Gamma}(1,2),\qquad \lambda_B\mid\text{data}\sim\operatorname{Gamma}(1,3). $$ A Gamma distribution with shape 1 is exponential, with the displayed rate parameter. Thus the posterior densities are $2e^{-2x}$ and $3e^{-3x}$ for $x>0$.

The probability that A has the higher rate is $$ \begin{aligned} P(\lambda_A>\lambda_B\mid\text{data}) &=\int_0^\infty P(\lambda_B<x\mid\text{data})\,2e^{-2x}\,dx\\ &=\int_0^\infty (1-e^{-3x})\,2e^{-2x}\,dx\\ &=1-\frac{2}{5}=0.6000. \end{aligned} $$ The probability that B has the higher rate is $1-0.6000=0.4000$. Therefore the next participant is assigned to A with probability 0.6000 and to B with probability 0.4000. This result follows from the posterior distributions, not from comparing the two posterior means alone.

Worked Example: Allocating Toward the Lower Rate

Suppose, hypothetically, that a lower event rate is preferable because the outcome is an adverse event. Use independent $\operatorname{Gamma}(1,1)$ priors. Arm A has 0 events over 3 person-years and arm B has 0 events over 1 person-year. Assume the same constant-rate Poisson model and exposure units for both arms.

The posterior distributions are $$ \lambda_A\mid\text{data}\sim\operatorname{Gamma}(1,4),\qquad \lambda_B\mid\text{data}\sim\operatorname{Gamma}(1,2). $$ To find the probability that A has the lower rate, integrate the probability that $\lambda_A=x$ while B’s rate exceeds $x$: $$ \begin{aligned} P(\lambda_A<\lambda_B\mid\text{data}) &=\int_0^\infty 4e^{-4x}P(\lambda_B>x\mid\text{data})\,dx\\ &=\int_0^\infty 4e^{-4x}e^{-2x}\,dx\\ &=\frac{4}{6}=\frac{2}{3}=0.6667. \end{aligned} $$ Thus, under this lower-is-better rule, A receives the next participant with probability 0.6667, and B with probability $1-0.6667=0.3333$. The direction of the comparison must match the clinical meaning of the outcome.

Calculation, Diagnostics, and Alternatives

For two independent exponential posteriors with shape 1 and rate parameters $b_A$ and $b_B$, a useful exact calculation follows by integration. Since $P(\lambda_B<x)=1-e^{-b_Bx}$,

P(\lambda_A>\lambda_B) =\int_0^\infty (1-e^{-b_Bx})b_Ae^{-b_Ax}\,dx =1-\frac{b_A}{b_A+b_B} =\frac{b_B}{b_A+b_B}.

More generally, the best-arm probability is an integral over all the other arms’ posterior distributions, and it usually has no simple closed form. It can be estimated by repeatedly drawing one rate from each arm’s posterior and recording which arm is highest or lowest. The estimated probability is the fraction of draws on which that arm wins. Monte Carlo error means this fraction is an approximation; it is not an additional source of clinical outcome data.

Several checks help reveal a mismatch between the model and the trial:

Check exposure accounting. Confirm that counts and exposures use compatible units and that the exposure for each arm is the total exposure contributing to its count. Adding a count where exposure belongs, or omitting exposure, changes the posterior rate.
Assess the Poisson assumption. If counts vary much more than the Poisson model predicts, or events cluster within participants or sites, a simple Poisson model may understate uncertainty. A posterior predictive check can compare observed count patterns with counts generated under the fitted model.
Review the constant-rate assumption. If risk changes substantially over time or participants have systematically different rates, one common arm-specific rate may not represent the outcome process.
Inspect prior influence. With little exposure or few events, the prior can materially affect the posterior and thus the allocation probabilities. Check whether plausible prior choices change the allocation conclusions.
Verify the direction of benefit. Using the highest-rate probability for an adverse-event outcome reverses the intended preference, even when every calculation is otherwise correct.

This approach differs from a binary-outcome Beta-Binomial model: the Poisson-Gamma update uses counts and exposure to learn an event rate, rather than treating each observation as a single success or failure. It also differs from a plug-in greedy rule, which selects the arm with the most favorable point estimate and does not randomize according to uncertainty. A normal-outcome model instead describes a continuous measurement and its sampling variation. Each model answers a different data-generating question; the outcome scale and sampling process should guide the choice.

Gamma-Poisson conjugacy makes the update transparent, but it does not make the assumptions automatic. In particular, the model does not by itself account for dependence, changing rates, or heterogeneity among participants. If these features are important, the simple independent arm-specific posteriors may produce misleading best-rate probabilities. Allocation probabilities are model-based decision quantities, not guarantees about future outcomes.

Common Mistake

A frequent error is to add the number of events to the Gamma rate parameter and the exposure to the shape. That reverses the update. The likelihood contributes $\lambda^Y$, so the event count increases the power of $\lambda$ and therefore the shape. It contributes $e^{-E\lambda}$, so exposure increases the exponential rate parameter. The correct update is always $\operatorname{Gamma}(a+Y,b+E)$ under this parameterization.

Key takeaway: With a Gamma shape-rate prior and a Poisson count over known exposure, events update the shape and exposure updates the rate parameter. Thompson sampling then allocates according to the posterior probability of the clinically preferred extreme rate: highest when higher is better, lowest when lower is better.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions.

A prior is $\operatorname{Gamma}(3,4)$. The arm observes 2 events over 5 person-years. What is the posterior distribution?
In the first worked allocation example, why is the probability for B the complement of the probability for A?
For an adverse event, should Thompson sampling favor the arm with the highest or lowest posterior rate? State what probability determines its allocation.
Explain why a Poisson count model needs exposure information when comparing arms with different observation times.

Solutions. (1) Add the count to the shape and exposure to the rate: $\operatorname{Gamma}(3+2,4+5)=\operatorname{Gamma}(5,9)$. (2) The two independent continuous posterior rates tie with probability 0, so exactly one of the two arms has the higher rate; their probabilities sum to 1. (3) It should favor the lowest rate, with allocation probability $P(\lambda_k<\lambda_j\text{ for every }j\ne k\mid\text{trial data})$. (4) Under the model, the count mean is exposure multiplied by rate. Different exposures therefore imply different expected counts even at the same rate; accounting for exposure is necessary to compare rates rather than raw counts.