Markov Chain Monte Carlo
Markov Chain Monte Carlo (MCMC) is a method for estimating probability distributions by generating samples from a Markov chain. The idea is to create a sequence of samples such that each sample is dependent only on the previous one, and each step of the sequence is chosen randomly.
Here’s a simple example: imagine you have a bag of 5 marbles — 3 red and 2 blue. To estimate the probability of picking a red marble, you can use MCMC.
Start by randomly picking a marble and record its color.
For each subsequent step, randomly choose a marble from the bag, with a preference for the color you just picked (since it’s the only information you have about the distribution).
Repeat step 2 many times, say 1000 times.
The ratio of the number of red marbles to the total number of marbles you picked is an estimate of the probability of picking a red marble from the bag.
MCMC is a widely used method in Bayesian statistics, where it is often used to estimate posterior distributions.
Another example is in optimization problems, MCMC can be used to find global minimum of a complex function where gradient descent or other optimization algorithms fail to work.
import pymc3 as pm
import numpy as np
# Define the data and model
data = np.random.normal(size=100)
with pm.Model() as model:
mu = pm.Normal("mu", mu=0, sigma=1)
sigma = pm.HalfNormal("sigma", sigma=1)
y = pm.Normal("y", mu=mu, sigma=sigma, observed=data)
# Sample from the posterior distribution
with model:
trace = pm.sample(draws=1000, chains=2, tune=500, random_seed=123)
In this example, we are trying to estimate the mean and standard deviation of a normal distribution using MCMC. The data variable is a sample of 100 values generated from a normal distribution. The model is defined using PyMC3, where we specify the prior distributions for mu and sigma. The y variable is observed data, and the model is fit to this data. Finally, we use the sample function to generate 1000 samples from the posterior distribution, using 2 chains and a tuning step of 500.