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

Материал из Википедии — свободной энциклопедии
Перейти к: навигация, поиск

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

Алгоритм Метрополиса-Гастингса позволяет семплировать любую функцию распределения. Он основан на создании цепи Маркова, то есть на каждом шаге алгоритма новое выбранное значение x^{t+1} зависит только от предыдущего x^t. Алгоритм использует вспомогательную функцию распределения Q(x'|x^t), зависящую от x^t, для которой генерировать выборку просто (например, нормальное распределение). На каждом шаге для этой функции генерируется случайное значение x'. Затем с вероятностью

u = \frac{P(x')Q(x^{t}|x')}{P(x^{t})Q(x'|x^{t})}

(или с вероятностью 1, если u > 1), выбранное значение принимается как новое: x^{t+1} = x', а иначе оставляется старое: x^{t+1} = x^{t}.

Например, если взять нормальную функцию распределения как вспомогательную функцию, то


Q( x'| x^t ) \sim N( x^t, \sigma^2 I).

Такая функция выдаёт новое значение в зависимости от значения на предыдущем шаге. Изначально алгоритм Метрополиса требовал, чтобы вспомогательная функция была симметрична: Q(x',x^t) = Q(x^t,x'), однако обобщение Гастингса снимает это ограничение.

Алгоритм[править | править вики-текст]

Пусть мы уже выбрали случайное значение x^t. Для выбора следующего значения сначала получим случайное значение x' для функции Q(x'|x^t). Затем найдем произведение a = a_1a_2, где

a_1 = \frac{P(x')}{P(x^t)}

является отношением вероятностей между промежуточным значением и предыдущим, а

a_2 = \frac{Q(x^t|x')}{Q(x'|x^t)}

это отношение между вероятностями пойти из x' в x^t или обратно. Если Q симметрична, то второй множитель равен 1. Случайное значение на новом шаге выбирается по правилу:


\begin{matrix}
\mbox{If } a \geq 1: &  \\
& x^{t+1} = x',
\end{matrix}

\begin{matrix}
\mbox{and if } a < 1: & \\
& x^{t+1} = \left\{
                   \begin{matrix}
                       x'\mbox{ with probability }a \\
                       x^t\mbox{ with probability }1-a.
                   \end{matrix}
            \right.
\end{matrix}

Алгоритм стартует из случайного значения x^0, и сначала прогоняется «вхолостую» некоторое количество шагов, чтобы «забыть» о начальном значении.

Лучше всего алгоритм работает тогда, когда форма вспомогательной функции близка к форме целевой функции P. Однако добиться этого априори зачастую невозможно. Для решения этой проблемы вспомогательную функцию настраивают в ходе подготовительной стадии работы алгоритма. Например, для нормального распределения настраивают его параметр \sigma^2 так, чтобы доля «принятых» случайных значений (то есть тех, для которых x^{t+1} = x') была близка к 60 %. Если \sigma^2 слишком мала, то значения будут получаться слишком близкими и доля принятых будет высока. Если \sigma^2 слишком велика, то с большой вероятностью новые значения будут выскакивать в зоны малой вероятности P, отчего доля принятых значений окажется низкой.