Марковские цепи Монте-Карло

Марковские цепи Монте-Карло (МЦМК, англ. Markov chain Monte Carlo, MCMC) — это класс алгоритмов для семплирования (получения случайных выборок) из вероятностных распределений, прямое моделирование которых затруднено. Метод основан на построении марковской цепи, стационарное распределение которой совпадает с целевым распределением. После достижения стационарности (сходимости) состояния цепи используются как выборки из целевого распределения. МЦМК широко применяются в байесовской статистике, машинном обучении, физике, биоинформатике и других областях, где требуется численное интегрирование или оценка сложных многомерных распределений.

История

Истоки метода МЦМК лежат в работе над Манхэттенским проектом в 1940-х годах. В 1947 году Джон фон Нейман, Станислав Улам и Николас Метрополис предложили метод Монте-Карло для моделирования нейтронных цепных реакций. Однако ключевой шаг был сделан в 1953 году, когда Метрополис, Ариэль Розенблут, Маршалл Розенблут и Августа Теллер (в соавторстве с Эдвардом Теллером) опубликовали статью, в которой описали алгоритм, позже названный алгоритмом Метрополиса. Этот алгоритм позволял генерировать выборки из распределения Больцмана для систем многих частиц, используя марковскую цепь.

В 1970 году Уилфред Кейт Гастингс обобщил алгоритм Метрополиса, что привело к созданию алгоритма Метрополиса — Гастингса, который стал основой для большинства современных МЦМК-методов. В 1984 году Стюарт Геман и Дональд Геман применили МЦМК к задачам восстановления изображений, а в 1990-х годах метод получил широкое распространение в байесовской статистике благодаря работам Адриана Смита, Алана Гельмана и других. Развитие вычислительной техники в конце XX века сделало МЦМК практичным инструментом для анализа сложных моделей.

Основные принципы

Марковская цепь

Марковская цепь — это последовательность случайных величин \(X_1, X_2, \dots\), в которой распределение следующего состояния зависит только от текущего состояния (свойство марковости). В МЦМК цепь строится так, чтобы её стационарное распределение \(\pi(x)\) совпадало с целевым распределением \(p(x)\), из которого требуется семплировать. Для этого переходное ядро \(T(x \to x')\) должно удовлетворять условию детального баланса: \[ p(x) T(x \to x') = p(x') T(x' \to x). \] Это гарантирует, что \(p(x)\) является стационарным распределением цепи.

Сходимость

После начального этапа («разогрева», англ. burn-in) цепь достигает стационарного состояния. Состояния, полученные до сходимости, обычно отбрасываются. Для оценки сходимости используются графики трассировки (trace plots), автокорреляционные функции и статистические тесты (например, критерий Гельмана — Рубина). Сходимость может быть медленной для многомерных или сильно коррелированных распределений.

Основные алгоритмы

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

Наиболее общий и широко используемый алгоритм МЦМК. На каждом шаге:

1. Предлагается новое состояние \(x'\) из предложного распределения \(q(x' | x)\).
2. Вычисляется вероятность принятия:

\[ \alpha = \min\left(1, \frac{p(x') q(x | x')}{p(x) q(x' | x)}\right). \]

1. С вероятностью \(\alpha\) цепь переходит в \(x'\), иначе остаётся в \(x\).

Если \(q\) симметрично (\(q(x' | x) = q(x | x')\)), формула упрощается до \(\alpha = \min(1, p(x')/p(x))\) (алгоритм Метрополиса).

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

Частный случай алгоритма Метрополиса — Гастингса, где каждое новое состояние получается последовательным семплированием из условных распределений каждой переменной при фиксированных остальных. Для многомерного распределения \(p(x_1, x_2, \dots, x_d)\) на каждом шаге:

1. Выбирается переменная \(x_i\) (циклически или случайно).
2. Новое значение \(x_i'\) семплируется из \(p(x_i | x_{\setminus i})\), где \(x_{\setminus i}\) — все остальные переменные.

Семплирование Гиббса не требует настройки предложного распределения и эффективно для моделей с простыми условными распределениями (например, в иерархических байесовских моделях).

Гамильтонов МЦМК (HMC)

Метод, использующий градиенты логарифма целевого распределения для эффективного исследования пространства состояний. Вводится вспомогательная переменная «импульс», и цепь строится с помощью уравнений гамильтоновой динамики. HMC особенно эффективен для непрерывных распределений большой размерности, так как уменьшает случайное блуждание и автокорреляцию. Требует настройки параметров (шаг интегрирования, число шагов).

Срезовое семплирование (Slice sampling)

Алгоритм, который адаптивно выбирает шаги, не требуя предложного распределения. Идея: семплировать равномерно из области под графиком целевого распределения. Для этого:

1. Выбирается высота \(y \sim U(0, p(x))\).
2. Семплируется \(x'\) равномерно из множества \(\{x: p(x) \geq y\}\).

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

Применение

Байесовская статистика

МЦМК является основным инструментом для апостериорного вывода в байесовских моделях. Позволяет вычислять маргинальные распределения, моменты, доверительные интервалы и предсказательные распределения для моделей с произвольными априорными распределениями и функциями правдоподобия. Примеры: байесовская линейная регрессия, иерархические модели, модели со скрытыми переменными.

Машинное обучение

Используется в обучении вероятностных графических моделей (например, латентное размещение Дирихле, LDA), в байесовских нейронных сетях, для оценки сложных интегралов в вариационных методах. В задачах с большими данными применяются приближённые методы, такие как стохастический градиентный МЦМК.

Физика и химия

МЦМК применяется для моделирования систем многих частиц (метод Монте-Карло в статистической физике), расчёта термодинамических свойств, моделирования фазовых переходов, в квантовой химии (метод Монте-Карло для квантовых систем).

Биоинформатика

Используется в филогенетике (оценка эволюционных деревьев), анализе последовательностей ДНК/белков, моделировании пространственной структуры белков.

Проблемы и ограничения

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

Состояния марковской цепи коррелированы, что снижает эффективную численность выборки (ESS). Для уменьшения автокорреляции применяют прореживание (thinning) — сохранение каждого k-го состояния. Однако прореживание не всегда улучшает оценку, так как выбрасывает информацию.

Сходимость

Гарантировать сходимость за конечное время невозможно. Для многомодальных распределений (с несколькими изолированными пиками) цепь может застревать в одном моде, не исследуя другие. Для решения этой проблемы используются методы, такие как параллельный тюнинг (parallel tempering) или реплики с обменом.

Выбор предложного распределения

Для алгоритма Метрополиса — Гастингса эффективность сильно зависит от выбора предложного распределения \(q\). Слишком узкое распределение приводит к медленному исследованию пространства, слишком широкое — к высокой частоте отказов. На практике часто используют адаптивные методы, которые корректируют \(q\) в процессе работы.

Сравнение с другими методами

МЦМК отличается от методов прямого семплирования (например, обратного преобразования) тем, что не требует знания аналитической формы распределения. По сравнению с методами важности (importance sampling) МЦМК не требует подбора предложного распределения, близкого к целевому, но может быть медленнее для низкоразмерных задач. Для многомерных распределений МЦМК часто превосходит методы Монте-Карло с независимыми выборками, но уступает вариационным методам по скорости (хотя обычно точнее).

Интересные факты

• Алгоритм Метрополиса был назван в честь Николаса Метрополиса, хотя его роль в разработке была скорее административной.
• Термин «Марковские цепи Монте-Карло» впервые появился в 1980-х годах, хотя сам метод использовался ранее.
• В 2020 году метод МЦМК был применён для анализа распространения COVID-19 в эпидемиологических моделях.
• Для задач с очень большими размерностями (тысячи и миллионы переменных) МЦМК часто заменяют на вариационные методы, такие как автоэнкодеры.

Источники

• Метрополис Н., Розенблут А., Розенблут М., Теллер А., Теллер Э. «Уравнения состояния для жидкостей и твердых тел с помощью быстрых вычислительных машин». Журнал химической физики, 1953.
• Гастингс У. К. «Методы Монте-Карло с выборкой с помощью марковских цепей». Биометрика, 1970.
• Геман С., Геман Д. «Стохастическая релаксация, распределения Гиббса и байесовская реставрация изображений». Труды IEEE по анализу изображений и машинному интеллекту, 1984.
• Гельман А., Карлин Дж. Б., Стерн Х. С., Дансон Д. Б., Вехтари А., Рубин Д. Б. «Байесовский анализ данных». 3-е издание, CRC Press, 2013.
• Брукс С., Гельман А., Джонс Г., Мэн С.-Л. (ред.). «Справочник по марковским цепям Монте-Карло». CRC Press, 2011.