Bayesiansk inferens

Klassisk statistik ger en punktskattning med ett konfidensintervall. Bayesiansk inferens ger en fördelning över parametervärden.

$\underbrace{p(\theta \mid D)}_{\text{posterior}} \propto \underbrace{p(D \mid \theta)}_{\text{likelihood}} \cdot \underbrace{p(\theta)}_{\text{prior}}$

Skillnaden märks mest när data är få. En A/B-test med 3 konverteringar av 10 besök ger punktskattningen 30 %. Det är inte fel, men det är inte heller särskilt informativt — den bayesianska posteriorn visar att allt mellan 10 % och 60 % är rimligt.

Den sista raden är den som gör bayesianska intervall lättare att kommunicera: de betyder faktiskt det folk redan tror att konfidensintervall betyder.

Konjugerade priorer ger posteriorn i sluten form — ingen sampling behövs.

Beta–Binomial är den mest användbara i praktiken. Tolkningen av priorn är konkret: Beta$(\alpha, \beta)$ motsvarar att ha sett $\alpha - 1$ lyckade och $\beta - 1$ misslyckade försök i förväg.

A/B-test bayesianskt löser flera praktiska problem på en gång:

1. Frågan blir «hur stor är sannolikheten att B är bättre än A?» — direkt tolkbart.
2. Man kan titta när som helst utan att inflatera felrisken, till skillnad från upprepade signifikanstest.
3. Man får en förväntad förlust av att välja fel, vilket är vad beslutet faktiskt beror på.

Punkt 2 är en verklig fördel, men den missförstås ofta: man slipper inte att data är brusiga, och att avbryta så snart posteriorn ser bra ut ger fortfarande dåliga beslut. Det man slipper är den specifika multiplicitetskorrektionen.

När data blir många spelar priorn nästan ingen roll — likelihooden dominerar. Bayesianska metoder är därför mest värda vid små dataset, hierarkiska strukturer och när förkunskap faktiskt finns.

Hierarkiska modeller är det starkaste argumentet: har du 50 skolor med olika många elever kan du låta skolorna dela en gemensam fördelning. Små skolor «lånar styrka» från de stora och dras mot gruppmedelvärdet — en partiell sammanslagning som varken behandlar alla lika eller var och en isolerad.

När sluten form inte finns används MCMC (NUTS i Stan eller PyMC) eller variationsinferens. Kontrollera alltid konvergensen: $\hat{R} < 1{,}01$ och tillräckligt effektivt stickprov.

import numpy as np
from math import lgamma, exp

# Beta–Binomial: posteriorn i sluten form
def posterior(alfa_prior, beta_prior, lyckade, forsok):
    return alfa_prior + lyckade, beta_prior + forsok - lyckade

def beta_kvantil(a, b, q, n=200_000, fro=0):
    x = np.random.default_rng(fro).beta(a, b, n)
    return float(np.quantile(x, q))

# 3 konverteringar av 10 — vad vet vi egentligen?
a, b = posterior(1, 1, 3, 10)
print(f"posterior Beta({a}, {b})")
print(f"  medel {a / (a + b):.3f}")
print(f"  95 %-intervall [{beta_kvantil(a, b, 0.025):.3f}, {beta_kvantil(a, b, 0.975):.3f}]")
# posterior Beta(4, 8)
#   medel 0.333
#   95 %-intervall [0.109, 0.612]     ← punktskattningen 0,30 döljer denna osäkerhet

# Priorns inflytande krymper med data
for n, k in ((10, 3), (100, 30), (1000, 300)):
    for namn, (ap, bp) in (("uniform", (1, 1)), ("stark 10 %", (10, 90))):
        a, b = posterior(ap, bp, k, n)
        print(f"  n={n:>4} {namn:<11} posteriormedel {a / (a + b):.4f}")
#  ↑ vid n=1000 är de två priorerna praktiskt taget omöjliga att skilja åt

# Bayesianskt A/B-test
def ab_test(k_a, n_a, k_b, n_b, prior=(1, 1), n_sim=200_000, fro=0):
    rng = np.random.default_rng(fro)
    pa = rng.beta(prior[0] + k_a, prior[1] + n_a - k_a, n_sim)
    pb = rng.beta(prior[0] + k_b, prior[1] + n_b - k_b, n_sim)
    b_battre = float((pb > pa).mean())
    return {
        "P(B > A)": round(b_battre, 4),
        "förväntad_lyft": round(float((pb - pa).mean()), 4),
        "95 %_intervall_lyft": [round(float(np.quantile(pb - pa, q)), 4) for q in (0.025, 0.975)],
        "förväntad_förlust_av_att_välja_B": round(float(np.maximum(pa - pb, 0).mean()), 5),
    }

print(ab_test(45, 500, 60, 500))
# {'P(B > A)': 0.9..., 'förväntad_lyft': 0.03, ...}

# Hierarkisk modell: små grupper lånar styrka från de stora
def partiell_sammanslagning(k, n, iterationer=200):
    """Empirisk Bayes: skatta gemensam prior ur data, dra grupperna mot den."""
    k, n = np.asarray(k, float), np.asarray(n, float)
    p_global = k.sum() / n.sum()
    styrka = 10.0
    for _ in range(iterationer):
        a, b = styrka * p_global, styrka * (1 - p_global)
        post = (a + k) / (a + b + n)
        p_global = float(post.mean())
    return {"rå": [round(float(x), 3) for x in k / n],
            "sammanslagen": [round(float(x), 3) for x in post],
            "gruppmedel": round(p_global, 3)}

print(partiell_sammanslagning(k=[1, 12, 48], n=[2, 40, 160]))
# {'rå': [0.5, 0.3, 0.3], 'sammanslagen': [0.348, 0.303, 0.301], 'gruppmedel': 0.317}
#  ↑ gruppen med 2 observationer dras kraftigt mot medelvärdet; den med 160 knappt alls