Bayesian Inference Foundations

In plain language

Bayesian inference is a method for updating beliefs in light of new evidence: you start with a "prior" — a probability distribution capturing what you already know about some unknown number (like a drug's event rate) — then multiply it by the probability of the data you just observed to get a "posterior," your revised belief. The key practical advantage is that the result is a genuine probability statement about the unknown number itself, for example "there is an 89% chance the true rate is above 0.25," which is what most people intend when they report a confidence interval but which frequentist statistics cannot technically say. In real-world evidence and health economics, Bayesian methods are especially useful when data are sparse (rare disease, small trials) and existing evidence can be formally incorporated, or when a cost-effectiveness model needs to propagate all sources of uncertainty coherently.

Core theorem: prior × likelihood → posterior

Bayesian inference rests on a single mechanistic rule. Before seeing data, the analyst encodes existing knowledge — from prior studies, expert opinion, theoretical constraints, or deliberate agnosticism — as a probability distribution p(θ) over the unknown parameter θ. This is the prior. The data D inform via the likelihood L(θ; D) = p(D | θ): the probability of observing exactly these data if θ were the true value. Bayes' theorem combines both:

posterior p(θ | D) ∝ prior p(θ) × likelihood p(D | θ)

The ∝ symbol ("proportional to") hides only a normalizing constant — a number that makes the posterior integrate to 1. That constant does not depend on θ and does not affect inference about θ. The posterior is the complete statement of post-data knowledge: it assigns a probability to every interval of plausible θ values.

In a finite conjugate update (e.g., the beta-binomial pair described below), the posterior is computed in closed form with no numerical integration. In all other cases — logistic or survival regressions with Bayesian priors, hierarchical multi-site models, models with correlated parameters — the posterior is approximated by Markov chain Monte Carlo (MCMC).

Credible interval vs confidence interval — the interpretation most people want

The single most practically important distinction between Bayesian and frequentist inference is in the interval that summarizes uncertainty. A 95% credible interval (CrI) is the interval [a, b] such that the posterior probability that θ ∈ [a, b] is exactly 0.95. The CrI is a direct probability statement about the parameter: "Given the prior and the data, I am 95% certain the true event rate is between 0.12 and 0.45."

This is precisely the interpretation most scientists, clinicians, and decision-makers apply when they see a confidence interval — and it is wrong for a frequentist confidence interval. A 95% CI is a statement about a procedure, not about a specific interval: if the study were repeated infinitely many times and a CI constructed each time, 95% of those intervals would contain the true fixed value of θ. For any single realized interval, θ either is or is not inside it — there is no 95% probability. The parameter θ is fixed and unknown; only the interval is random in the frequentist framework.

When sample sizes are large and priors are weak (non-informative), the numerical difference between a 95% CrI and a 95% CI is negligible. When samples are small, priors are informative, or parameters are on the boundary of their support, the two can diverge substantially. In rare-disease RWE settings with n < 50, the CrI is the only interval that correctly propagates all sources of uncertainty — prior, likelihood, and model structure — into a single coherent statement.

Conjugate priors and the beta-binomial example

A conjugate prior is a prior distribution family that, combined with a specific likelihood, yields a posterior in the same family — enabling exact closed-form updating. The canonical example for proportions and event rates is the beta-binomial pair:

Prior: Beta(α₀, β₀) [event rate p has prior mean α₀/(α₀+β₀)] Data: k events in n Bernoulli trials [likelihood: Binomial(n, p)] Posterior: Beta(α₀ + k, β₀ + n − k) [exact; no MCMC required]

The shape parameters α and β carry the interpretation of pseudo-counts: α₀ is the number of prior "events" and β₀ is the number of prior "non-events," giving the prior an effective sample size (ESS) of α₀ + β₀. The posterior simply adds observed counts to the pseudo-counts. The Jeffreys non-informative prior for a binomial proportion is Beta(0.5, 0.5), which contributes only 1 pseudo-observation. The link to the beta-distribution entry in this catalog: the beta-distribution entry explains how to fit Beta(α, β) from a published mean and SE (method of moments), which is exactly how informative priors are constructed from published literature.

Prior choice: informative, weakly informative, skeptical, and flat

The prior encodes what is known before the data are seen. Four classes arise in practice:

• Informative prior: encodes genuine external knowledge — a Beta(2, 8) prior on an event
• Weakly informative prior: rules out physically impossible or substantively absurd values
• Skeptical prior: deliberately weights against large treatment effects, raising the
• Flat / non-informative prior: places approximately equal mass across the parameter range.

Prior sensitivity analysis is mandatory reporting. Every Bayesian submission should include results under at least two or three prior specifications spanning from non-informative to the analyst's preferred informative prior. Robustness of the conclusion across this grid means the prior is not load-bearing; strong sensitivity means the data are insufficient to overcome prior beliefs — which is itself a substantive finding that must be disclosed. Regulatory agencies (FDA 2023 Bayesian guidance) and HTA bodies (NICE DSU Technical Support Documents) require prior sensitivity analyses for Bayesian submissions.

MCMC and convergence diagnostics

Most practical Bayesian models do not have conjugate closed forms. Markov chain Monte Carlo (MCMC) — Metropolis-Hastings, Gibbs sampling, and the modern No-U-Turn Sampler (NUTS) used in Stan, PyMC, and brms — approximates the posterior by constructing a Markov chain whose stationary distribution is the target posterior. After a burn-in (warmup) phase where the chain finds the high-probability region and tunes its step size, the retained draws are treated as samples from the posterior. Convergence diagnostics exist not to prove convergence (that cannot be proven) but to detect non-convergence: the R-hat (potential scale reduction factor) should be < 1.01 for every parameter; Effective Sample Size (ESS) should exceed 400 in the bulk and tail. Trace plots should look like well-mixed "hairy caterpillars" — no trends, no stuck chains, no chains drifting apart. Always report these diagnostics. A posterior from a non-converged chain is meaningless regardless of how sophisticated the model. Use Stan/brms (R) or PyMC (Python) — they implement NUTS and report diagnostics automatically.

Where Bayes earns its keep in RWE/HEOR

Bayesian inference is not universally superior to frequentist inference; it earns its additional modeling and computational cost in specific settings:

• Rare-disease historical borrowing: the beta-binomial conjugate update (or a MAP prior) lets
• Probabilistic sensitivity analysis in cost-effectiveness: PSA is already implicitly Bayesian
• Sequential safety monitoring: Bayesian adaptive designs accumulate posterior probability of
• Hierarchical shrinkage across sites: a Bayesian hierarchical model partially pools
• Bias priors in probabilistic bias analysis: unmeasured confounding can be modeled by

Frequentist-Bayes pragmatism (calibrated Bayes)

The tribal framing — "Bayesian vs frequentist" — is not useful in applied RWE/HEOR. Most production analyses use both paradigms at different layers. A Cox proportional hazards model is frequentist; adding Normal(0, σ²) priors on the regression coefficients and sampling via MCMC makes it Bayesian. Calibrated Bayes (Rubin 1984) is the pragmatic synthesis: use Bayesian machinery where it offers genuine advantages (informative priors, hierarchical structure, exact small-sample uncertainty, direct probability statements), and frequentist methods where they are adequate (large samples, no genuine prior, regulatory conventions requiring frequentist operating characteristics). The key sanity check: a Bayesian posterior interval should have approximately correct frequentist coverage when the model is well specified. When Bayesian and frequentist results agree, the prior is not carrying load. When they disagree, the prior is doing real work — and must be justified explicitly.

Pros, cons, and trade-offs

Pros of Bayesian inference: - Direct probability statements about parameters: "There is a 94% posterior probability the drug reduces 12-month mortality by more than 5 percentage points." - Natural, principled incorporation of prior information in small-sample and rare-disease settings — the prior ESS is a transparent, reportable quantity. - Coherent uncertainty propagation through complex hierarchical models without ad hoc variance adjustments. - No frequentist multiplicity penalty for pre-specified interim looks in adaptive designs. - The 95% credible interval has the interpretation practitioners routinely and incorrectly apply to confidence intervals.

Cons and limitations: - Prior specification adds a modeling choice that must be pre-specified, justified, and sensitivity-analyzed; a poorly chosen or post-hoc prior can bias inference in opaque ways. - Computational cost: MCMC can be slow for large RWE datasets or models with many parameters; convergence must be verified and can fail. - Frequentist operating characteristics (type-I error, power) must still be verified by simulation for regulatory submissions — the Bayesian framework does not guarantee error control unless the design is pre-specified. - Regulatory reviewers trained in frequentist methods may require additional justification; acceptance is increasing but uneven across jurisdictions.

When to use

• Small samples where prior information is genuine and documentable (rare disease, pediatric
• Hierarchical models with parameter sharing across sites, subgroups, or indications.
• Sequential adaptive designs where interim analyses are planned and a Bayesian stopping rule
• PSA frameworks where coherent posterior propagation of parameter uncertainty is preferred.
• When a direct probability statement about the parameter is the primary deliverable for a

When NOT to use — and when it is actively misleading

• Do not choose the prior to reach a desired result. Prior specification must be locked in
• Do not skip prior sensitivity analysis. Reporting only the result under one prior — even
• Do not interpret a high posterior probability as a causal claim. A 98% posterior
• Do not use flat priors on variance components in hierarchical models. A Uniform(0, ∞)
• Bayesian machinery does not fix a confounded design. Sophisticated MCMC sampling of a

Interpreting the output

From the worked example: prior Beta(2,8) (mean 0.2, equivalent to 10 pseudo-observations); observe 6 events in 20 trials; posterior Beta(8,22), mean 8/30 ≈ 0.267.

(1) Formal interpretation. The posterior Beta(8,22) is a probability distribution over the unknown event rate p. The posterior mean 8/30 ≈ 0.267 is a weighted compromise between the prior mean (2/10 = 0.2) and the data MLE (6/20 = 0.3), with weights proportional to the prior information (10 pseudo-observations) and the data (20 observations). The total information in the posterior is 8+22 = 30 pseudo-observations. A 95% credible interval from Beta(8,22) spans approximately [0.12, 0.45]: the posterior probability that the true event rate p lies in this interval is exactly 0.95 — a direct probability statement, not a repeated-sampling claim. This posterior is exact because of conjugacy; no MCMC is required for this model.

(2) Practical interpretation. Before seeing the new data, we believed the event rate was around 20%, with the strength of 10 observations. We then observed 30% in 20 patients. Our best updated estimate is approximately 26.7% — not 30%, because the prior pulled us back toward 20%; not 20%, because the new data pulled us upward. A decision-maker should note that the prior is doing real work: if its 10 pseudo-observations came from a comparable population, the shrinkage is appropriate. If the prior population differed (different disease stage, earlier era), the posterior is biased toward the prior, and a sensitivity analysis under Beta(0.5, 0.5) — contributing only 1 pseudo-observation — should accompany the primary result to show how much the conclusion depends on the informative prior.

Worked example

Scenario

A health outcomes analyst is studying an adverse event rate for a drug used in a rare pediatric disease. From published literature, the event rate is believed to be around 20%, based on evidence roughly equivalent to 10 patients. The analyst encodes this as a prior Beta(2,8) — 2 prior events and 8 prior non-events. A new registry study then observes 6 adverse events among 20 patients. The analyst wants to compute the updated (posterior) event rate using the exact beta-binomial conjugate formula, and confirm it lies between the prior belief (20%) and the data estimate (30%).

Dataset

Summary counts feeding the beta-binomial conjugate update. Each row represents one source of information: prior belief (encoded as pseudo-counts) and new observed data.

Steps

• Step 1 — Identify prior parameters. Prior Beta(2,8): alpha_0 = 2 prior events, beta_0 = 8 prior non-events. Prior mean = 2/(2+8) = 2/10 = 0.2.

• Step 2 — Read the new data. Observed k = 6 events in n = 20 trials. Data maximum-likelihood estimate = 6/20 = 0.3, which is higher than the prior mean of 0.2.

• Step 3 — Compute the posterior alpha parameter by adding observed events to the prior event count: alpha_post = 2+6 = 8.

• Step 4 — Compute the posterior beta parameter by adding observed non-events to the prior non-event count: n-k = 20-6 = 14, so beta_post = 8+14 = 22.

• Step 5 — Verify the posterior total pseudo-count (effective sample size): alpha_post + beta_post = 8+22 = 30, which equals the prior ESS (10) plus new data (20), confirming 10+20 = 30.

• Step 6 — Compute the posterior mean. Posterior mean = 8/30 ≈ 0.267. This lies between the prior mean (0.2) and the data MLE (0.3). The data (20 observations) outweigh the prior (10 pseudo-observations) by 2 to 1, pulling the posterior mean two-thirds of the way from 0.2 toward 0.3.

Result

Posterior Beta(8,22): prior mean 2/10 = 0.2, data MLE 6/20 = 0.3, posterior mean 8/30 ≈ 0.267. The posterior ESS = 8+22 = 30, composed of 10 prior pseudo-observations plus 20 real data observations. The 95% credible interval from Beta(8,22) is approximately [0.12, 0.45] — the posterior probability the true event rate lies in this range is exactly 0.95. Because the posterior mean (0.267) falls between the prior (0.2) and the data (0.3), it is correctly shrunk toward the prior in proportion to how informative the prior was relative to the new data.

Runnable example

import numpy as np
from scipy import stats

# ── Beta-binomial conjugate update (exact posterior, no MCMC) ─────────────────
alpha_0, beta_0 = 2, 8          # prior Beta(2,8): mean = 2/(2+8) = 0.2; ESS = 10
k, n            = 6, 20         # observe 6 events in 20 trials

alpha_post = alpha_0 + k        # 2 + 6 = 8
beta_post  = beta_0 + (n - k)  # 8 + (20 - 6) = 22

prior     = stats.beta(alpha_0, beta_0)
posterior = stats.beta(alpha_post, beta_post)

prior_mean = alpha_0 / (alpha_0 + beta_0)             # 2/10 = 0.200  (exact)
data_mle   = k / n                                     # 6/20 = 0.300  (exact)
post_mean  = alpha_post / (alpha_post + beta_post)     # 8/30 ≈ 0.267

print(f"Prior     Beta({alpha_0}, {beta_0}):  mean = {prior_mean:.4f}")
print(f"Data MLE: {k}/{n} = {data_mle:.4f}")
print(f"Posterior Beta({alpha_post}, {beta_post}): mean = {post_mean:.4f}")

# 95% credible interval — direct posterior probability statement
cri_lo, cri_hi = posterior.ppf(0.025), posterior.ppf(0.975)
print(f"95% credible interval: [{cri_lo:.3f}, {cri_hi:.3f}]")
print("Interpretation: P(event rate in [{:.3f}, {:.3f}] | data) = 0.95".format(cri_lo, cri_hi))

# Direct posterior probability query: P(p > 0.25 | data)
prob_q = 1 - posterior.cdf(0.25)
print(f"P(event rate > 0.25 | prior + data): {prob_q:.3f}")

# Prior sensitivity: Jeffreys Beta(0.5, 0.5) (non-informative)
a_jeff = 0.5 + k
b_jeff = 0.5 + (n - k)
post_jeff = stats.beta(a_jeff, b_jeff)
print(f"\nSensitivity — Jeffreys prior Beta(0.5,0.5):")
print(f"  Posterior Beta({a_jeff}, {b_jeff}): mean = {a_jeff/(a_jeff+b_jeff):.4f}")
print(f"  95% CrI: [{post_jeff.ppf(0.025):.3f}, {post_jeff.ppf(0.975):.3f}]")

# ── PyMC sketch: Bayesian logistic regression (non-conjugate) ─────────────────
# import pymc as pm
# import arviz as az
# with pm.Model():
#     # Weakly informative priors on log-odds scale
#     intercept = pm.Normal("intercept", mu=0, sigma=2)
#     beta_treat = pm.Normal("beta_treat", mu=0, sigma=1)   # skeptical prior
#     beta_age   = pm.Normal("beta_age",   mu=0, sigma=1)
#     p = pm.math.invlogit(intercept + beta_treat * treat + beta_age * age)
#     obs = pm.Bernoulli("obs", p=p, observed=y)
#     idata = pm.sample(2000, tune=2000, target_accept=0.9,
#                       chains=4, random_seed=1)
# # Convergence diagnostics — MANDATORY before interpreting posterior
# assert az.summary(idata)["r_hat"].max() < 1.01, "R-hat >= 1.01: chain not converged"
# print(az.summary(idata, var_names=["beta_treat"], hdi_prob=0.95))
# # Direct posterior probability: P(treatment increases rate)
# print("P(beta_treat > 0 | data):", float((idata.posterior["beta_treat"] > 0).mean()))

# ── Beta-binomial conjugate update (exact posterior, base R) ──────────────────
alpha_0 <- 2;  beta_0 <- 8    # prior Beta(2,8): mean = 2/(2+8) = 0.2; ESS = 10
k <- 6;  n <- 20              # observe 6 events in 20 trials

alpha_post <- alpha_0 + k              # 2 + 6 = 8
beta_post  <- beta_0 + (n - k)        # 8 + 14 = 22

prior_mean <- alpha_0 / (alpha_0 + beta_0)            # 2/10 = 0.200
data_mle   <- k / n                                    # 6/20 = 0.300
post_mean  <- alpha_post / (alpha_post + beta_post)   # 8/30 ≈ 0.267

cat(sprintf("Prior     Beta(%d,%d): mean = %.4f\n", alpha_0, beta_0, prior_mean))
cat(sprintf("Data MLE: %d/%d = %.4f\n", k, n, data_mle))
cat(sprintf("Posterior Beta(%d,%d): mean = %.4f\n", alpha_post, beta_post, post_mean))

# 95% credible interval — P(theta in CrI | data) = 0.95 exactly
cri <- qbeta(c(0.025, 0.975), alpha_post, beta_post)
cat(sprintf("95%% credible interval: [%.3f, %.3f]\n", cri[1], cri[2]))

# Direct posterior probability: P(p > 0.25 | data)
cat(sprintf("P(event rate > 0.25 | data): %.3f\n",
            1 - pbeta(0.25, alpha_post, beta_post)))

# Prior sensitivity: Jeffreys Beta(0.5, 0.5) — non-informative
a_jeff <- 0.5 + k;  b_jeff <- 0.5 + (n - k)
cat(sprintf("\nSensitivity check — Jeffreys prior Beta(0.5, 0.5):\n"))
cat(sprintf("  Posterior Beta(%.1f, %.1f): mean = %.4f\n",
            a_jeff, b_jeff, a_jeff / (a_jeff + b_jeff)))
cat(sprintf("  95%% CrI: [%.3f, %.3f]\n",
            qbeta(0.025, a_jeff, b_jeff), qbeta(0.975, a_jeff, b_jeff)))

# ── brms sketch: Bayesian logistic regression via Stan ────────────────────────
# library(brms)
#
# fit <- brm(
#   outcome ~ treatment + age + cci,
#   data   = rwe_df,
#   family = bernoulli(link = "logit"),
#   prior  = c(
#     prior(normal(0, 2), class = Intercept),
#     prior(normal(0, 1), class = b, coef = treatment),  # skeptical prior
#     prior(normal(0, 1), class = b)                      # weakly informative
#   ),
#   chains = 4, iter = 2000, warmup = 1000, seed = 1,
#   backend = "cmdstanr"   # faster; requires cmdstanr + CmdStan
# )
#
# # MANDATORY convergence check before interpreting results
# stopifnot(max(rhat(fit)) < 1.01)
# cat("All R-hat < 1.01: chains have converged.\n")
# summary(fit)  # posterior means, 95% CrI, Bulk ESS, Tail ESS for each parameter
#
# # Direct posterior probability: P(treatment coefficient > 0 | data)
# h <- hypothesis(fit, "treatment > 0")
# print(h)   # P+ = posterior probability treatment increases log-odds