Markov Chain Monte Carlo (MCMC): Steekproeven uit complexe kansverdelingen

Sommige kansverdelingen zijn gewoon te complex om er rechtstreeks mee te werken.

Als je echte data modelleert, loopt de wiskunde vaak vast voordat je ergens komt. Het gebeurt vaak dat integralen op papier beheersbaar lijken, maar onhandelbaar worden zodra je een paar latente variabelen toevoegt. Dit komt vooral voor bij Bayesiaanse inferentie, waar de posteriorverdeling je prior beliefs combineert met geobserveerde data — en het resultaat iets is dat je niet met een simpele formule kunt samenvatten.

Het basisidee van Markov Chain Monte Carlo is dat MCMC in plaats van rechtstreeks de wiskunde te doen, de verdeling via simulatie verkent. Het trekt steekproeven die de vorm ervan weerspiegelen zonder die ooit volledig te hoeven uitrekenen.

In dit artikel behandel ik de kernconcepten achter MCMC, loop ik door de meest gebruikte algoritmen en laat ik je zien hoe je het in Python inzet.

Heb je een opfriscursus Python-wiskunde nodig? Lees onze blogposts Demystifying Mathematical Concepts for Deep Learning om te zien hoe wiskunde wordt toegepast in Numpy.

Wat is Markov Chain Monte Carlo?

Markov Chain Monte Carlo (MCMC) is een familie van algoritmen die steekproeven genereren uit kansverdelingen — zelfs wanneer die verdelingen te complex zijn om er direct mee te werken.

De naam bestaat uit twee delen. De Markov-keten bepaalt hoe het algoritme door mogelijke toestanden beweegt. Elke stap hangt alleen af van waar je nu bent, niet van de volledige geschiedenis van hoe je daar bent gekomen. Het Monte Carlo-deel betekent dat je willekeurige steekproeven gebruikt om grootheden van interesse te schatten.

Samen bouwt MCMC een keten van willekeurige steekproeven die na verloop van tijd de vorm van je doelverdeling weerspiegelen. Het is in de eerste plaats een steekproeftechniek. Je lost de wiskunde niet exact op, je benadert die via simulatie.

Waarom Markov Chain Monte Carlo nodig is

Het probleem met verdelingen van real-world data is dat ze lang niet zo netjes zijn als in de leerboeken.

In Bayesiaanse inferentie probeer je vaak een posteriorverdeling te berekenen — de geüpdatete kans op je modelparameters na het zien van data. De formule lijkt op papier eenvoudig: vermenigvuldig de prior met de likelihood en deel dan door de marginale likelihood. Die laatste term vereist integreren over elke mogelijke parameterwaarde. In hoge dimensies is die integraal praktisch niet te berekenen.

Het wordt alleen maar erger naarmate je model groeit. Als je meer parameters of meer latente variabelen toevoegt, loopt de exacte berekening dood. Je komt dit tegen in allerlei veelvoorkomende scenario’s:

• Bayesiaanse hiërarchische modellen: Waar parameters afhankelijk zijn van andere parameters, waardoor geneste verdelingen ontstaan die niet netjes te integreren zijn
• Hoge-dimensionale posteriors: Waar het aantal parameters rastergebaseerde benaderingen computationeel onhaalbaar maakt
• Niet-conjugate priors: Waar de prior en de likelihood niet samenkomen tot een herkenbare verdeling

MCMC is in deze scenario’s een goede uitweg. In plaats van de verdeling te berekenen, trekt het steekproeven eruit. Die steekproeven bevatten alles wat je nodig hebt, zonder ooit de integraal op te lossen.

De twee ideeën achter MCMC

MCMC combineert twee ideeën die op zichzelf eenvoudig zijn, maar krachtig in combinatie. Ik loop ze met je door.

Markov-ketens

Een Markov-keten is een reeks toestanden waarbij elke stap alleen afhangt van waar je nu bent.

Waar je eerder was doet er niet toe. Alleen je huidige toestand bepaalt waar je vervolgens naartoe gaat. Deze "geheugenloze" eigenschap — formeel de Markov-eigenschap genoemd — maakt de wiskunde beheersbaar en het algoritme praktisch.

De keten beweegt stap voor stap door toestanden en met de juiste inrichting komt die uiteindelijk in een stationaire verdeling terecht — een stabiel patroon waarbij de kans om in een bepaalde toestand te zijn niet meer verandert. Met die stabiele verdeling werkt MCMC precies.

Montecarlo-methoden

Montecarlo-methoden gebruiken willekeurige steekproeven om grootheden te schatten die lastig direct te berekenen zijn.

Het idee is om genoeg willekeurige steekproeven uit een verdeling te trekken en vervolgens het gemiddelde, de variantie of een andere eigenschap te schatten door simpelweg naar de steekproeven te kijken. Hoe meer steekproeven je trekt, hoe dichter je schattingen bij de echte waarden komen.

Op zichzelf vereisen Montecarlo-methoden dat je direct uit de verdeling kunt steekproeven — precies het probleem dat we proberen op te lossen. Markov-ketens pakken dat deel aan.

Hoe Markov Chain Monte Carlo werkt

MCMC is een lus met bij elke stap een eenvoudige beslissing.

1. Begin met een initiële toestand: Kies een startpunt in je parameter-ruimte, zelfs als het een ruwe gok is
2. Stel een nieuwe toestand voor: Gebruik een voorstelmechanisme om te suggereren waar je vervolgens heen gaat
3. Beslis of je die accepteert: Vergelijk de voorgestelde toestand met de huidige en accepteer of verwerp op basis van hoe waarschijnlijk elk is onder de doelverdeling
4. Herhaal: Voer deze lus duizenden keren uit en verzamel onderweg steekproeven
5. Gebruik de steekproeven: Schat gemiddelden, credible intervals of elke andere grootheid van interesse op basis van de verzamelde steekproeven

De acceptatie-/verwerpstap is waar de “magie” gebeurt.

Door betere toestanden vaker te accepteren dan slechtere, beweegt de keten naar gebieden met hoge waarschijnlijkheid — zonder ooit de volledige verdeling te hoeven berekenen.

Vroege steekproeven zijn afhankelijk van je startpunt, dus die gooi je weg. Na genoeg iteraties “vergeet” de keten het startpunt en weerspiegelen de overblijvende steekproeven de echte vorm van je doelverdeling.

De doelverdeling en stationaire verdeling

MCMC is gebouwd rond het doel om steekproeven te genereren uit een doelverdeling waar je niet direct uit kunt steekproeven.

De doelverdeling is waar je iets over wilt leren — meestal een posteriorverdeling bij Bayesiaanse inferentie. Je kent de vorm tot op een normaliseringsconstante, maar je kunt die constante niet direct berekenen. MCMC heeft die niet nodig.

Elk MCMC-algoritme is zo ontworpen dat zijn Markov-keten een stationaire verdeling heeft die overeenkomt met de doelverdeling. Een stationaire verdeling is de verdeling waarin de keten na genoeg stappen uitkomt.

Laat de keten doorlopen en hij begint steekproeven te produceren die eruitzien als trekkingen uit je doelverdeling. De integraal wordt omzeild.

Populaire Markov Chain Monte Carlo-algoritmen

Er zijn een handvol MCMC-algoritmen die je in de praktijk ziet. Ze volgen allemaal dezelfde kernlus, maar verschillen in hoe ze nieuwe toestanden voorstellen en hoe ze informatie over de doelverdeling gebruiken.

Metropolis-algoritme

Het Metropolis-algoritme is het eenvoudigste MCMC-algoritme en degene waarmee het allemaal begon.

Bij elke stap stelt het een nieuwe toestand voor door willekeurige ruis aan de huidige toe te voegen. Als de voorgestelde toestand een hogere waarschijnlijkheid heeft onder de doelverdeling, wordt die altijd geaccepteerd. Is die lager, dan wordt die geaccepteerd met een kans evenredig aan de verhouding van de twee waarschijnlijkheden — anders blijft de keten staan.

Dit acceptatie-/verwerpmechanisme betekent dat de keten meer tijd doorbrengt in gebieden met hoge waarschijnlijkheid zonder ooit de volledige verdeling te berekenen.

Het Metropolis-algoritme gebruikt een symmetrische voorstelverdeling, wat betekent dat een stap in elke richting even waarschijnlijk is. Het gaat vaak haperen naarmate modellen groter worden.

Metropolis-Hastings-algoritme

Het Metropolis-Hastings (MH)-algoritme generaliseert Metropolis door asymmetrische voorstelverdelingen toe te staan.

MH past de acceptatiekans aan om te corrigeren voor het feit dat sommige voorstellen waarschijnlijker zijn dan andere. Je kunt het voorstel afstemmen op de vorm van je doel, wat leidt tot betere verkenning en snellere convergentie.

De meeste moderne MCMC-methoden zijn uitbreidingen van MH of gebouwd op dezelfde principes. Als je Metropolis-Hastings begrijpt, snap je de basis van het vakgebied.

Gibbs sampling

Gibbs sampling werkt één variabele tegelijk bij in plaats van in één keer een nieuwe toestand voor alle parameters voor te stellen.

Bij elke stap trekt het voor elke variabele een steekproef uit de conditionele verdeling — de verdeling van die variabele gegeven de huidige waarden van alle andere. Als je alle variabelen hebt doorlopen, heb je één volledige iteratie voltooid.

Dit vermijdt volledig de acceptatie-/verwerpstap, omdat elke conditionele trekking altijd wordt geaccepteerd. Het is handig wanneer de volledige gezamenlijke verdeling lastig te bemonsteren is maar de conditionelen wel tractable zijn, wat vaak voorkomt in Bayesiaanse hiërarchische modellen.

Hamiltonian Monte Carlo

Hamiltonian Monte Carlo (HMC) is het eerste algoritme dat moderne Bayesiaanse inferentie op schaal praktisch maakte.

In plaats van willekeurig nieuwe toestanden voor te stellen, gebruikt HMC gradiëntinformatie van de doelverdeling om toestanden voor te stellen die ver van de huidige positie liggen maar toch waarschijnlijk geaccepteerd worden. Het beweegt veel beter door de parameter-ruimte dan random-walk-methoden. Er zijn minder verworpen voorstellen en een betere verkenning van hoge-dimensionale verdelingen.

Random-walk-methoden zoals Metropolis schalen niet mee als het aantal parameters toeneemt. HMC heeft dat probleem in veel mindere mate.

HMC is de motor achter Stan, een van de meest gebruikte probabilistische programmeerplatformen. De No-U-Turn Sampler (NUTS), een adaptieve uitbreiding van HMC die in PyMC wordt gebruikt, haalt de noodzaak weg om stapgrootte en aantal stappen handmatig af te stemmen.

Markov Chain Monte Carlo in de Bayesiaanse statistiek

Als er één gebied is waar MCMC de grootste impact heeft gehad, is het Bayesiaanse inferentie.

Bayesiaanse statistiek draait om de posteriorverdeling, de geüpdatete kans op je modelparameters na het zien van data. Die berekenen betekent de prior met de likelihood vermenigvuldigen en normaliseren. Die normalisatiestap vereist een integraal die zelden tractable is.

MCMC haalt deze stap volledig weg. Je evalueert gewoon de ongenormaliseerde posterior op een willekeurig punt en laat de keten de rest doen.

Hier is een eenvoudig voorbeeld. Stel dat je de bias van een munt schat. Je begint met de priorovertuiging dat de munt waarschijnlijk eerlijk is en observeert vervolgens een reeks worpen. Voor een simpel muntmodel heeft de posterior een gesloten vorm. Als je een hiërarchische structuur toevoegt — dus de bias over honderd munten tegelijk schatten — wordt het onmogelijk om te berekenen.

Met MCMC laat je de keten lopen, verzamel je steekproeven uit de posterior en gebruik je die om te berekenen wat je nodig hebt.

Burn-in, convergentie en mixing begrijpen

Deze drie concepten verwarren data scientists die nieuw zijn met MCMC. Als je ze verkeerd aanpakt, krijg je wel resultaten, maar weet je niet waarom ze onbetrouwbaar zijn.

Burn-in

Wanneer een Markov-keten start, heeft hij geen idee waar de gebieden met hoge waarschijnlijkheid van je doelverdeling liggen.

Die vroege steekproeven worden beïnvloed door je startpunt, niet door de doelverdeling. Burn-in is de praktijk om ze weg te gooien. Je laat de keten eerst een aantal iteraties lopen, gooit die steekproeven weg en houdt alleen wat overblijft nadat de keten de tijd heeft gehad om een goed startgebied te vinden.

Er is geen universele regel voor hoe lang de burn-in moet zijn. Het hangt af van je model, je startpunt en hoe goed de keten mixt. In de praktijk diagnosticeer je het visueel met trace plots in plaats van vooraf een vast aantal te kiezen.

Convergentie

Convergentie betekent dat de keten niet langer wordt beïnvloed door het startpunt en nu steekproeven trekt die de doelverdeling weerspiegelen.

Een keten die niet is geconvergeerd, produceert bevooroordeelde steekproeven. Het gemiddelde dat je eruit berekent, komt niet overeen met het echte posteriorgemiddelde. In plaats daarvan weerspiegelt het waar de keten toevallig vastzat.

Convergentie beoordeel je achteraf met diagnostiek. Meerdere ketens vanaf verschillende startpunten laten lopen en controleren of ze overeenkomen is een van de betrouwbaarste manieren om mislukte convergentie te ontdekken.

Mixing

Een keten die convergeert maar slecht mixt, blijft problematisch.

Mixing beschrijft hoe goed de keten de doelverdeling verkent. Een goed mixende keten beweegt vrij rond, bezoekt zowel hoge als lage waarschijnlijkheidsregio’s en produceert steekproeven die grofweg onafhankelijk van elkaar zijn. Een slecht mixende keten blijft vele iteraties in één gebied voordat hij zich verplaatst en produceert sterk gecorreleerde steekproeven die de volledige verdeling niet representeren.

Slechte mixing zie je vaak terug als een trace plot die lijkt op een langzame, kronkelende rivier in plaats van een ruisende horizontale band. Zie je dat, dan heeft je sampler tuning nodig — een betere voorstelverdeling of een geheel ander algoritme.

Vergelijkingsgrafiek van mixing

Hoe MCMC-resultaten te evalueren

Ik laat je nu vier manieren zien om MCMC te evalueren en leg uit wanneer je welke gebruikt.

Trace plots

Een trace plot toont de bemonsterde waarde van een parameter bij elke iteratie. Het is het eerste waar je naar moet kijken na het draaien van MCMC.

Een gezonde trace plot ziet eruit als witte ruis rond een stabiel gemiddelde. Je zou geen trends, lange vlakke stukken of langzame drift moeten zien. Als je ziet dat de keten ronddwaalt of lang in één gebied blijft hangen, is dat een mixingprobleem en zijn je steekproeven onbetrouwbaar.

Trace plots gevisualiseerd

Autocorrelatie

MCMC-steekproeven zijn nooit volledig onafhankelijk. Elke steekproef wordt beïnvloed door de vorige. Autocorrelatie meet hoe sterk steekproeven over iteraties heen gecorreleerd zijn.

Hoge autocorrelatie betekent dat je steekproeven minder informatie bevatten dan hun aantal suggereert. Tweeduizend gecorreleerde steekproeven kunnen je dezelfde informatie geven als tweehonderd onafhankelijke. De meeste MCMC-bibliotheken bevatten autocorrelatieplots zodat je kunt zien hoe snel de correlatie afneemt naarmate steekproeven verder uit elkaar liggen.

Autocorrelatieplots gevisualiseerd

Effectieve steekproefgrootte

Effectieve steekproefgrootte (ESS) vertaalt die autocorrelatie naar een praktisch getal: aan hoeveel onafhankelijke steekproeven je keten gelijkstaat.

Als je 5.000 steekproeven hebt getrokken maar de ESS is 200, werk je met de statistische power van 200 onafhankelijke trekkingen. Een lage ESS betekent dat je je keten langer moet laten lopen, de sampler moet tunen, of beide. De meeste practitioners mikken op een ESS van minstens een paar honderd per parameter voordat ze hun schattingen vertrouwen.

Convergentiediagnostiek

Als je meerdere ketens draait, kun je formeel testen of ze naar dezelfde verdeling zijn geconvergeerd. De Gelman-Rubin-diagnostiek, gerapporteerd als R-hat, vergelijkt variantie binnen elke keten met variantie tussen ketens.

Een R-hat dicht bij 1,0 betekent dat de ketens overeenkomen — een goed teken. Waarden boven 1,01 of 1,05 (afhankelijk van de drempel die je bibliotheek gebruikt) duiden erop dat de ketens niet zijn geconvergeerd en je meer iteraties nodig hebt. De meeste moderne bibliotheken, zoals PyMC, berekenen R-hat automatisch en geven een waarschuwing als die te hoog is.

Markov Chain Monte Carlo in Python

Python heeft een paar bibliotheken voor MCMC, elk met een eigen filosofie.

• PyMC: De meest toegankelijke optie voor de meeste data scientists. Je definieert je model met Python-syntaxis en PyMC regelt de sampling onder de motorkap met NUTS (een adaptieve versie van HMC). Het is dé keuze voor Bayesiaans modelleren in Python.
• CmdStanPy / PyStan: Python-interfaces naar Stan, een toegewijde probabilistische programmeertaal. Stan compileert je model naar C++ voordat het wordt uitgevoerd, wat het snel maakt en geschikt voor complexe hiërarchische modellen. Het nadeel is een steilere leercurve.
• NumPy: Geen dedicated MCMC-bibliotheek, maar je kunt algoritmen zoals Metropolis-Hastings vanaf nul implementeren. Handig om te begrijpen wat er onder de motorkap gebeurt voordat je naar tools op hoger niveau gaat.

Voor de meeste praktische werkzaamheden begin je bij PyMC. Dat is wat ik ga gebruiken, dus als je meedoet, installeer de bibliotheek dan eerst:

pip install pymc

Om het simpel te houden blijf ik bij een eenvoudig voorbeeld: de bias van een munt schatten op basis van een reeks worpen.

Stap 1: Definieer het model

import pymc as pm
import numpy as np

# 1 = heads, 0 = tails
observed_flips = np.array([1, 0, 1, 1, 0, 1, 1, 1, 0, 1])

with pm.Model() as coin_model:
    # Prior: we believe the coin is probably fair
    bias = pm.Beta("bias", alpha=2, beta=2)

    # Likelihood: observed flips given the bias
    flips = pm.Bernoulli("flips", p=bias, observed=observed_flips)

De pm.Beta-prior codeert een zwakke overtuiging dat de munt eerlijk is. De pm.Bernoulli-likelihood koppelt het model aan de geobserveerde data.

Stap 2: Start de sampler

with coin_model:
    trace = pm.sample(2000, tune=1000, return_inferencedata=True)

Uitvoer van het draaien van de sampler

tune stelt het aantal burn-in-stappen in — die steekproeven worden weggegooid. sample trekt na tuning 2000 posteriorsteekproeven per keten.

Stap 3: Inspecteer de posterior

import arviz as az

az.plot_trace(trace, var_names=["bias"])
az.summary(trace, var_names=["bias"])

Model-trace plot en samenvattingsresultaten

az.summary() geeft je het posteriorgemiddelde, de standaardafwijking en R-hat voor elke parameter. Als R-hat dicht bij 1,0 ligt, zijn de ketens geconvergeerd. az.plot_trace() toont de trace en posteriorverdeling naast elkaar voor elke parameter.

Voor deze dataset — 7 kop uit 10 worpen — is het posteriorgemiddelde 0,642 met een standaardafwijking van 0,124. Dit weerspiegelt het bewijs in de data terwijl het dicht bij de fair-coin-prior blijft. R-hat is 1,00 en ESS ligt ruim boven 2000, dus de ketens zijn geconvergeerd en de steekproeven zijn betrouwbaar.

Veelgemaakte fouten met MCMC

MCMC is makkelijk te draaien, maar ook makkelijk verkeerd te gebruiken. Dit zijn de fouten die het vaakst voorkomen.

• Te vroeg aannemen dat er is geconvergeerd: Dat de sampler zonder fouten draaide, betekent niet dat die is geconvergeerd. Controleer altijd R-hat en trace plots voordat je je resultaten vertrouwt. Een keten kan stabiel lijken en toch vastzitten in één regio van de parameter-ruimte.
• Burn-in negeren: Vroege steekproeven zijn bevooroordeeld richting je startpunt, niet je doelverdeling. Gooi ze weg. De meeste bibliotheken doen dit automatisch via de tune-parameter, maar controleer of je die steekproeven niet per ongeluk in je analyse opneemt.
• Slechte mixing: Een keten die te langzaam beweegt, produceert sterk gecorreleerde steekproeven die minder informatie bevatten dan hun aantal doet vermoeden. Als je trace plot eruitziet als een langzame drift in plaats van een ruisende horizontale band, moet je je voorstelverdeling tunen of een beter algoritme gebruiken. Overschakelen van Metropolis naar NUTS lost dit vaak direct op.
• Te weinig steekproeven: Een lage effectieve steekproefgrootte (ESS) betekent dat je schattingen onbetrouwbaar zijn, zelfs als de keten technisch is geconvergeerd. Als je ESS een klein deel is van je totale trekkingen, laat de keten dan langer lopen of verbeter de mixing voordat je conclusies trekt.

MCMC vs andere samplingmethoden

MCMC is niet de enige manier om een verdeling te benaderen. Zo verhoudt het zich tot de alternatieven.

• Rejection sampling trekt voorstellen uit een eenvoudigere verdeling en accepteert ze op basis van hoe goed ze bij het doel passen. Het is exact wanneer het werkt, maar acceptatiepercentages storten in bij hoge dimensionaliteit.
• Importance sampling weegt steekproeven uit een voorstelverdeling om verwachtingen onder de doelverdeling te benaderen. Het is snel en werkt goed in lage dimensies, maar de variantie van de gewichten explodeert naarmate de dimensionaliteit toeneemt. In de praktijk dragen een paar steekproeven bijna al het gewicht, waardoor schattingen onbetrouwbaar worden.
• Variational inference (VI) past een eenvoudigere verdeling aan om de doelverdeling te benaderen door een divergatiemaat te minimaliseren. VI is veel sneller dan MCMC en schaalt naar grote datasets, maar introduceert benaderingsfouten. De gepaste verdeling kan posteriorstructuur missen die MCMC wel zou oppikken.

Hier is de korte versie:

MCMC vergeleken met alternatieven

MCMC is de juiste keuze wanneer nauwkeurigheid belangrijker is dan snelheid. Als je moet schalen naar grote datasets of inferentie in real time nodig hebt, kan variational inference de nauwkeurigheidsafruil waard zijn.

Conclusie

MCMC is zo’n tool die er van buiten intimiderend uitziet, maar veel begrijpelijker wordt zodra je snapt wat het daadwerkelijk doet — een keten van steekproeven opbouwen die geleidelijk de vorm weerspiegelt van een verdeling die je niet direct kunt uitrekenen.

Het is ook veel makkelijker te begrijpen als je het opsplitst: Markov-ketens en Montecarlo-methoden.

De rol ervan in de Bayesiaanse statistiek is nauwelijks te overschatten. Posteriorverdelingen die anders onbereikbaar zouden zijn, worden oplosbaar zodra je een betrouwbare sampler hebt. Daarom staat MCMC centraal in probabilistische programmeerbibliotheken als PyMC en Stan.

Maar voordat je in de implementatie duikt, moet je de intuïtie goed hebben. Begrijp waarom de keten moet inbranden, wat mixing echt betekent en hoe je een trace plot leest. De code zelf is het makkelijke deel, want Python-bibliotheken verbergen alle abstracties achter eenvoudige functieaanroepen.

Als je vaardig wilt worden in machine learning, schrijf je dan in voor onze Machine Learning Scientist in Python-track. 85 uur aan materiaal helpt je om in 2026 job-ready te zijn.