Для более сложных распределений вероятностей, вам могут потребоваться более продвинутые методы для генерации выборок, чем методы, описанные в Common Pseudorandom Number Generation Methods. Такие распределения возникают, например, в байесовском анализе данных и в больших комбинаторных задачах симуляции цепи Маркова Monte Carlo (MCMC). Альтернативой является построение марковской цепи со стационарным распределением, равным целевому дискретизирующему распределению, с использованием состояний цепи для генерации случайных чисел после начального периода горения, в котором распределение состояний сходится к цели.
Алгоритм Metropolis-Hastings черпает выборки из распределения, которое известно только до константы. Случайные числа генерируются из распределения с функцией плотности вероятностей, которая равна или пропорциональна функции предложения.
Чтобы сгенерировать случайные числа:
Примите начальное значение x (t).
Нарисуйте выборку, y (t), из q распределения предложений (y | x (t)).
Примите y (t) как следующую x выборки (t + 1) с r вероятностей (x (t), y (t)) и сохраните x (t) как следующую x выборки (t + 1) с вероятностью 1 - r (x (t), y (t)), где:
Увеличьте t → t + 1 и повторите шаги 2 и 3, пока не получите желаемое количество выборок.
Сгенерируйте случайные числа с помощью метода Metropolis-Hastings с mhsample функция. Чтобы эффективно производить выборки с помощью алгоритма Metropolis-Hastings, очень важно выбрать хорошее распределение предложений. Если трудно найти эффективное распределение предложений, используйте выборку среза (slicesample) или гамильтониан Монте-Карло (hmcSampler) вместо этого.
В образцах, когда трудно найти эффективное распределение предложений Metropolis-Hastings, алгоритм дискретизации срезов не требует явной спецификации. Алгоритм выборки среза черпает выборки из области под функцией плотности с помощью последовательности вертикальных и горизонтальных шагов. Во-первых, он выбирает высоту случайным образом от 0 до функции плотности f (x). Затем он выбирает новое значение x случайным образом путем выборки из горизонтального «среза» плотности выше выбранной высоты. Аналогичный алгоритм дискретизации среза используется для многомерного распределения.
Если задана функция f (x), пропорциональная функции плотности, то сделайте следующее, чтобы сгенерировать случайные числа:
Примите начальное значение x (t) в области f (x).
Нарисуйте действительное значение y равномерно из (0, f (x (t))), тем самым определив горизонтальный «срез» как S = {x: y < f (x)}.
Найдите интервал I = (L, R) вокруг x (t), который содержит все или большую часть S «среза».
Нарисуйте новую x точки (t + 1) в этом интервале.
Увеличьте t → t + 1 и повторите шаги 2-4, пока не получите желаемое количество выборок.
Выборка среза может сгенерировать случайные числа из распределения с произвольной формой функции плотности при условии, что доступна эффективная числовая процедура для нахождения интервала I = (L, R), который является «срезом» плотности.
Сгенерируйте случайные числа с помощью метода выборки среза с slicesample функция.
Выборка Metropolis-Hastings и срез может производить цепи MCMC, которые медленно смешиваются и долго сходятся к стационарному распределению, особенно в среднемерных и высоко-размерных задачах. Используйте основанную на градиенте гамильтоновую выборку Монте-Карло (HMC) для ускорения выборки в этих ситуациях.
Чтобы использовать отбор проб HMC, необходимо задать журнал f(x) (до аддитивной константы) и его градиент. Можно использовать численный градиент, но это приводит к более медленной выборке. Все переменные выборки должны быть без ограничений, что означает, что логарифмические f(x) и его градиент четко определены для всех реальных x. Чтобы выборка ограниченные переменные, преобразуйте эти переменные в без ограничений таковых перед использованием HMC-сэмплера.
Алгоритм квантования HMC вводит случайный «вектор импульса» z и определяет плотность z соединений и x «вектор положения» следующим P(x,z) = f(x)g(z). Цель состоит в том, чтобы выбрать из этого распределения соединений и затем игнорировать значения z - маргинальное распределение x имеет необходимую f(x) плотности.
Алгоритм HMC присваивает Гауссову плотность с ковариацией матричной M («большой матрицы»), чтобы z:
Затем он задает «энергетическую функцию» как
с U(x) = - логарифмируйте f(x) «потенциальную энергию» и K (z) = zTM-1z/2 «кинетическая энергия». Плотность соединений определяется P(x,z) ∝ exp {-E(x,z) }.
Чтобы сгенерировать случайные выборки, алгоритм HMC:
Принимает начальное значение x вектора положения.
Генерирует выборку вектора импульса: z ∼ g(z).
Развивает (x, z) состояния на некоторое количество фиктивного времени, τ к новому состоянию, (x’,z’) используя «уравнения движения»:
Если бы уравнения движения могли быть решены точно, энергия (и, следовательно, плотность) оставалась бы постоянной: E(x,z) = E(x’,z’). На практике уравнения движений должны быть решены численно (обычно с использованием так называемой чехарды интегрирования) и энергия не сохранена.
Принимает x’ как следующую выборку с вероятностью p acc = min (1, exp {E(x,z) – E(x’,z’)}), и сохраняет x как следующую выборку с вероятностью 1 - p acc.
Повторяет шаги 2-4 до тех пор, пока не сгенерирует желаемое количество выборок.
Чтобы использовать выборку HMC, создайте дискретизатор, используя hmcSampler функция. После создания сэмплера можно вычислить оценки точек MAP (MAP-a-posteriori), настроить сэмплер, отобрать выборки и проверить диагностику сходимости. Пример этого рабочего процесса см. в Байесовской линейной регрессии с использованием гамильтонового Монте-Карло.