Метод Монте‑Карло на марковских цепях (MCMC): выборка из сложных распределений вероятностей

Некоторые распределения вероятностей слишком сложны, чтобы работать с ними напрямую.

При моделировании реальных данных математика часто рушится ещё до того, как вы получите что-то полезное. Нередко интегралы выглядят приемлемо на бумаге, но становятся неразрешимыми, как только вы добавляете пару скрытых переменных. Это особенно характерно для байесовского вывода, где апостериорное распределение объединяет априорные представления с наблюдаемыми данными — и результат уже нельзя выразить простой формулой.

Базовая идея метода Монте‑Карло на марковских цепях в том, что вместо прямой работы с формулами MCMC исследует распределение через моделирование. Он генерирует выборки, отражающие его форму, не вычисляя её полностью.

В этой статье я разберу ключевые концепции MCMC, пройду по наиболее распространённым алгоритмам и покажу, как применить их в Python.

Нужен повтор математики для Python? Прочитайте наши публикации Demystifying Mathematical Concepts for Deep Learning, чтобы увидеть применение математики в Numpy.

Что такое метод Монте‑Карло на марковских цепях?

Markov Chain Monte Carlo (MCMC) — это семейство алгоритмов, которые генерируют выборки из распределений вероятностей, даже когда эти распределения слишком сложны для прямой работы с ними.

Название распадается на две части. Марковская цепь управляет тем, как алгоритм перемещается по возможным состояниям. На каждом шаге всё зависит только от текущего положения, а не от всей истории попадания в него. Часть Монте‑Карло означает использование случайной выборки для оценки интересующих величин.

Вместе MCMC строит цепь случайных выборок, которые со временем отражают форму целевого распределения. Прежде всего это техника выборки. Вы не решаете математику точно, а приближаете её с помощью моделирования.

Зачем нужен метод Монте‑Карло на марковских цепях

Проблема реальных распределений данных в том, что они далеки от учебниковых.

В байесовском выводе вы часто пытаетесь вычислить апостериорное распределение — обновлённые вероятности параметров модели после наблюдения данных. Формула на бумаге проста: перемножьте априорное и правдоподобие, затем разделите на маргинальное правдоподобие. Последний член требует интегрирования по всем возможным значениям параметров. В высоких размерностях этот интеграл практически невозможно вычислить.

С ростом модели становится только хуже. При добавлении большего числа параметров или скрытых переменных точное вычисление заходит в тупик. С этим вы столкнётесь во множестве типичных ситуаций:

• Байесовские иерархические модели: когда параметры зависят от других параметров, образуя вложенные распределения, которые нельзя аккуратно проинтегрировать
• Высокомерные апостериоры: когда число параметров делает приближения на сетке вычислительно неприемлемыми
• Неконъюгированные априоры: когда априор и правдоподобие не объединяются в узнаваемое распределение

В этих случаях MCMC — удачный обходной путь. Вместо вычисления распределения он берёт из него выборки. И в этих выборках есть всё необходимое — без решения интеграла.

Две идеи, лежащие в основе MCMC

MCMC объединяет две идеи, по отдельности достаточно простые, но мощные в связке. Давайте пройдёмся по ним.

Марковские цепи

Марковская цепь — это последовательность состояний, где каждый шаг зависит только от текущего положения.

Где вы были раньше, не важно. Только текущее состояние определяет, куда вы пойдёте дальше. Это свойство «без памяти» — формально называемое свойством Маркова — делает математику управляемой, а алгоритм — практичным.

Цепь перемещается по состояниям шаг за шагом и при правильной настройке в итоге выходит на стационарное распределение — устойчивый режим, в котором вероятность нахождения в любом состоянии перестаёт меняться. Именно с таким устойчивым распределением и работает MCMC.

Методы Монте‑Карло

Методы Монте‑Карло используют случайную выборку для оценки величин, которые сложно вычислить напрямую.

Идея в том, чтобы взять достаточно случайных выборок из распределения и оценить его среднее, дисперсию или другие свойства, просто анализируя выборки. Чем больше вы выбираете, тем ближе оценки к истинным значениям.

Сами по себе методы Монте‑Карло требуют уметь напрямую выбирать из распределения — а это как раз наша проблема. С этим справляются марковские цепи.

Как работает метод Монте‑Карло на марковских цепях

MCMC — это цикл с простым решением на каждом шаге.

1. Начните с начального состояния: выберите отправную точку в пространстве параметров, даже если это грубая догадка
2. Предложите новое состояние: используйте механизм предложения, чтобы указать, куда двигаться дальше
3. Решите, принять ли его: сравните предлагаемое состояние с текущим и примите или отклоните его в зависимости от того, насколько вероятны оба под целевым распределением
4. Повторяйте: прогоните этот цикл тысячи раз, по пути собирая выборки
5. Используйте выборки: оцените средние, правдоподобные интервалы и любые другие интересующие величины по собранным данным

Именно на шаге принятия/отклонения происходит «магия».

Принимая лучшие состояния чаще, чем худшие, цепь тяготеет к областям высокой вероятности — не вычисляя распределение полностью.

Ранние выборки зависят от точки старта, поэтому их отбрасывают. После достаточного числа итераций цепь «забывает» начальную точку, и оставшиеся выборки отражают истинную форму целевого распределения.

Целевое распределение и стационарное распределение

MCMC нацелен на генерацию выборок из целевого распределения, из которого нельзя выбрать напрямую.

Целевое распределение — это то, что вы хотите изучить, обычно апостериор в байесовском выводе. Вы знаете его форму с точностью до нормирующей константы, но не можете напрямую вычислить эту константу. MCMC она не нужна.

Каждый алгоритм MCMC устроен так, чтобы стационарное распределение его марковской цепи совпадало с целевым. Стационарное распределение — это то, к которому цепь приходит после достаточного числа шагов.

Продолжайте работу цепи, и она начнёт производить выборки, выглядящие как прямые выборки из целевого распределения. Интеграл обходится стороной.

Популярные алгоритмы MCMC

На практике используется несколько алгоритмов MCMC. Все они следуют одному и тому же базовому циклу, но различаются тем, как предлагают новые состояния и как используют информацию о целевом распределении.

Алгоритм Метрополиса

Алгоритм Метрополиса — самый простой алгоритм MCMC и с него всё началось.

На каждом шаге он предлагает новое состояние, добавляя к текущему случайный шум. Если предложенное состояние имеет более высокую вероятность под целевым распределением, оно всегда принимается. Если ниже — принимается с вероятностью, пропорциональной отношению двух вероятностей; иначе цепь остаётся на месте.

Механизм принятия/отклонения заставляет цепь проводить больше времени в областях высокой вероятности, не вычисляя распределение целиком.

Алгоритм Метрополиса использует симметричное распределение предложения, то есть одинаково вероятно сделать шаг в любом направлении. По мере роста моделей он склонен «ломаться».

Алгоритм Метрополиса — Гастингса

Алгоритм Метрополиса — Гастингса (MH) обобщает Метрополиса, позволяя использовать несимметричные распределения предложения.

MH корректирует вероятность принятия с учётом того, что одни предложения вероятнее других. Вы можете подстроить предложение под форму целевого распределения, что ведёт к лучшему исследованию пространства и более быстрой сходимости.

Большинство современных методов MCMC — это расширения MH или построены на тех же принципах. Поэтому, поняв Метрополиса — Гастингса, вы поймёте фундамент области.

Сэмплирование Гиббса

Сэмплирование Гиббса обновляет по одной переменной за раз вместо предложения нового состояния сразу для всех параметров.

На каждом шаге оно выбирает каждую переменную из её условного распределения — распределения этой переменной при фиксированных текущих значениях остальных. Когда вы прошли по всем переменным, это завершает одну полную итерацию.

Шаг принятия/отклонения полностью отсутствует, так как каждое условное вытягивание всегда принимается. Это удобно, когда из полного совместного распределения выбирать трудно, но условные распределения tractable — что часто бывает в байесовских иерархических моделях.

Гамильтоново Монте‑Карло

Hamiltonian Monte Carlo (HMC) — первый алгоритм, сделавший современный байесовский вывод практически масштабируемым.

Вместо случайного предложения новых состояний HMC использует градиентную информацию целевого распределения, чтобы предлагать состояния, находящиеся далеко от текущего положения, но всё ещё с высокой вероятностью принятия. Он куда лучше перемещается по пространству параметров, чем методы случайного блуждания. Отказов меньше, исследование высокоразмерных распределений — лучше.

Методы случайного блуждания вроде Метрополиса плохо масштабируются с ростом числа параметров. У HMC этой проблемы в той же степени нет.

HMC — двигатель Stan, одной из самых популярных платформ вероятностного программирования. Алгоритм No‑U‑Turn Sampler (NUTS), адаптивное расширение HMC, используемое в PyMC, снимает необходимость вручную настраивать размер шага и число шагов.

MCMC в байесовской статистике

Если есть область, где MCMC оказал наибольшее влияние, то это байесовский вывод.

Байесовская статистика сосредоточена на апостериорном распределении — обновлённой вероятности параметров модели после наблюдения данных. Его вычисление означает умножить априор на правдоподобие и нормировать. Этап нормировки требует интеграла, который редко бывает разрешимым.

MCMC полностью снимает этот шаг. Вы просто вычисляете ненормированный апостериор в произвольной точке и отдаёте остальное цепи.

Вот простой пример. Допустим, вы оцениваете смещение монеты. Вы начинаете с априорного убеждения, что монета, скорее всего, честная, затем наблюдаете последовательность бросков. Для простой модели монеты апостериор имеет замкнутый вид. Если добавить иерархию — то есть оценивать смещение сразу для сотни монет, — точное вычисление становится невозможным.

С MCMC вы запускаете цепь, собираете выборки из апостериора и по ним считаете нужные вам величины.

Что такое разогрев, сходимость и смешивание

Эти три понятия часто путают новички в MCMC. Ошибка здесь приведёт к результатам, ненадёжность которых вы не поймёте.

Разогрев (burn‑in)

Когда марковская цепь стартует, она не знает, где находятся области высокой вероятности вашего целевого распределения.

Ранние выборки зависят от стартовой точки, а не от целевого распределения. Разогрев — это практика их отбрасывания. Сначала вы запускаете цепь на определённое число итераций, выбрасываете эти выборки и сохраняете только то, что идёт после того, как цепь нашла хорошую область.

Универсального правила длительности разогрева нет. Всё зависит от вашей модели, стартовой точки и того, насколько хорошо цепь смешивается. На практике его диагностируют визуально по трейсплотам, а не выбирают фиксированное число заранее.

Сходимость

Сходимость означает, что на выборки перестала влиять стартовая точка, и теперь они отражают целевое распределение.

Цепь, не сошедшаяся, даёт смещённые выборки. Среднее, посчитанное по ним, не совпадёт с истинным апостериорным средним, а отразит область, в которой цепь «застряла».

Сходимость оценивают постфактум с помощью диагностик. Один из самых надёжных способов — запускать несколько цепей из разных стартовых точек и проверять, согласуются ли они.

Смешивание (mixing)

Цепь, которая сошлась, но плохо смешивается, всё ещё проблемна.

Смешивание описывает, насколько хорошо цепь исследует целевое распределение. Хорошо смешивающаяся цепь свободно перемещается, посещая области как высокой, так и низкой вероятности, и даёт выборки, примерно независимые друг от друга. Плохо смешивающаяся цепь надолго застревает в одной области и производит сильно коррелированные выборки, не отражающие всего распределения.

Плохое смешивание часто видно по трейсплоту, который выглядит как медленно текущая река, а не как шумная горизонтальная лента. Увидели такое — настраивайте сэмплер: улучшайте распределение предложения или вовсе меняйте алгоритм.

График сравнения смешивания

Как оценивать результаты MCMC

Далее я покажу четыре способа оценки MCMC и объясню, когда каждый использовать.

Трейсплоты (trace plots)

Трейсплот показывает значение параметра на каждой итерации. Это первое, на что стоит посмотреть после запуска MCMC.

Здоровый трейсплот выглядит как белый шум вокруг устойчивого среднего. Не должно быть трендов, длинных плато или медленного дрейфа. Если цепь блуждает или надолго застревает в одной области, это проблема смешивания, и выборки ненадёжны.

Визуализация трейсплотов

Автокорреляция

Выборки MCMC никогда не полностью независимы. Каждая зависит от предыдущей. Автокорреляция измеряет, насколько сильно выборки связаны между итерациями.

Высокая автокорреляция означает, что ваши выборки несут меньше информации, чем их количество. Две тысячи коррелированных выборок могут дать тот же объём информации, что двести независимых. В большинстве библиотек MCMC есть графики автокорреляции, на которых видно, как быстро корреляция спадает по мере увеличения лага.

Визуализация графиков автокорреляции

Эффективный размер выборки

Эффективный размер выборки (ESS) переводит автокорреляцию в практическую величину: сколько независимых выборок эквивалентна ваша цепь.

Если вы получили 5 000 выборок, но ESS равен 200, то фактически у вас статистическая мощность 200 независимых вытягиваний. Низкий ESS означает, что цепь нужно гонять дольше, настраивать сэмплер или и то и другое. На практике стремятся к ESS хотя бы в несколько сотен на параметр, прежде чем доверять оценкам.

Диагностики сходимости

При запуске нескольких цепей можно формально проверить, сошлись ли они к одному и тому же распределению. Диагностика Гельмана — Рубина, отчётная метрика — R‑hat, сравнивает дисперсию внутри цепей с дисперсией между цепями.

R‑hat, близкий к 1.0, означает согласованность цепей — хороший знак. Значения выше 1.01 или 1.05 (в зависимости от порога в библиотеке) указывают, что цепи не сошлись, и нужно больше итераций. Современные библиотеки вроде PyMC автоматически считают R‑hat и выводят предупреждения при завышенных значениях.

MCMC в Python

В Python есть несколько библиотек для MCMC, каждая со своей философией.

• PyMC: самый доступный вариант для большинства дата‑сайентистов. Вы определяете модель на Python, а PyMC берёт сэмплирование на себя с помощью NUTS (адаптивной версии HMC). Это основной инструмент для байесовского моделирования в Python.
• CmdStanPy / PyStan: интерфейсы Python к Stan, специализированному языку вероятностного программирования. Stan компилирует модель в C++ перед запуском, что делает его быстрым и подходящим для сложных иерархических моделей. Обратная сторона — более крутая кривая обучения.
• NumPy: отдельной библиотеки для MCMC нет, но можно реализовать алгоритмы вроде Метрополиса — Гастингса с нуля. Полезно, чтобы понять, что происходит «под капотом», прежде чем перейти к более высокоуровневым инструментам.

Для большинства практических задач стоит начинать с PyMC. Я буду использовать его, так что, если вы следуете за мной, сначала установите библиотеку:

pip install pymc

Чтобы упростить, я возьму лёгкий пример — оценку смещения монеты по последовательности бросков.

Шаг 1. Определите модель

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)

Априор pm.Beta кодирует слабое убеждение, что монета честная. Правдоподобие pm.Bernoulli связывает модель с наблюдаемыми данными.

Шаг 2. Запустите сэмплер

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

Вывод работы сэмплера

tune задаёт число шагов разогрева — эти выборки отбрасываются. sample после настройки извлекает по 2000 апостериорных выборок на цепь.

Шаг 3. Изучите апостериор

import arviz as az

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

Трейсплот модели и сводка результатов

az.summary() выдаёт апостериорное среднее, стандартное отклонение и R‑hat для каждого параметра. Если R‑hat близок к 1.0, цепи сошлись. az.plot_trace() строит трейсплот и апостериорное распределение рядом для каждого параметра.

Для этого набора — 7 орлов из 10 бросков — апостериорное среднее равно 0.642 со стандартным отклонением 0.124. Это отражает данные и остаётся близким к априору честной монеты. R‑hat равен 1.00, а ESS значительно выше 2000, так что цепи сошлись и выборки надёжны.

Типичные ошибки при работе с MCMC

MCMC легко запустить, но также легко использовать неверно. Вот самые частые ошибки.

• Слишком раннее предположение о сходимости: то, что сэмплер отработал без ошибок, не означает сходимость. Всегда проверяйте R‑hat и трейсплоты, прежде чем доверять результатам. Цепь может выглядеть стабильной и при этом быть «застрявшей» в одном регионе пространства параметров.
• Игнорирование разогрева: ранние выборки смещены к стартовой точке, а не к целевому распределению. Отбрасывайте их. Большинство библиотек делают это автоматически через параметр tune, но проверьте, что вы случайно не включили эти выборки в анализ.
• Плохое смешивание: цепь, которая движется слишком медленно, даёт сильно коррелированные выборки, несущие меньше информации, чем их количество. Если трейсплот выглядит как медленный дрейф, а не шумная горизонтальная лента, настройте распределение предложения или используйте более удачный алгоритм. Переход с Метрополиса на NUTS часто решает проблему.
• Слишком мало выборок: низкий эффективный размер выборки (ESS) делает оценки ненадёжными, даже если цепь формально сошлась. Если ESS — малая доля от общего числа вытягиваний, гоняйте цепь дольше или улучшайте смешивание, прежде чем делать выводы.

MCMC и другие методы выборки

MCMC — не единственный способ приблизить распределение. Вот как он сравнивается с альтернативами.

• Отбор с отклонением (rejection sampling) предлагает варианты из более простого распределения и принимает их в зависимости от того, насколько хорошо они соответствуют целевому. Когда работает — даёт точный результат, но в высоких размерностях доля принятий рушится.
• Важностное сэмплирование (importance sampling) взвешивает выборки из распределения предложения, чтобы аппроксимировать математические ожидания под целевым. Быстро и хорошо работает в низких размерностях, но дисперсия весов взрывается с ростом размерности. На практике несколько выборок несут почти весь вес, из‑за чего оценки ненадёжны.
• Вариационный вывод (VI) подбирает более простое распределение для аппроксимации целевого, минимизируя меру расхождения. VI намного быстрее MCMC и масштабируется на большие наборы данных, но вносит погрешность аппроксимации. Подогнанное распределение может упустить структуру апостериора, которую MCMC уловит.

Коротко:

MCMC против альтернатив

MCMC — правильный выбор, когда точность важнее скорости. Если нужно масштабироваться на большие данные или выполнять вывод в реальном времени, вариационный вывод может оправдать компромисс по точности.

Заключение

MCMC — из тех инструментов, что кажутся пугающими снаружи, но становятся понятными, как только вы осознаёте, что он на самом деле делает: строит цепь выборок, постепенно отражающую форму распределения, которое напрямую не вычислить.

Его гораздо легче понять, если разложить на части: марковские цепи и методы Монте‑Карло.

Его роль в байесовской статистике трудно переоценить. Апостериоры, до которых не дотянуться иначе, становятся решаемыми, как только у вас есть надёжный сэмплер. Поэтому MCMC — в ядре библиотек вероятностного программирования вроде PyMC и Stan.

Но прежде чем браться за реализацию, важно правильно понять интуицию. Разберитесь, зачем цепи нужен разогрев, что такое смешивание и как читать трейсплот. Сам код — простая часть, так как библиотеки Python скрывают всю абстракцию за простыми вызовами функций.

Если вы хотите стать профессионалом в машинном обучении, запишитесь на наш Machine Learning Scientist in Python трек. 85 часов материалов помогут вам быть готовыми к работе в 2026 году.