В сегодняшнем упражнении по кодированию мы узнаем, как сопоставить параметры модели с данными (с планками погрешностей) и получить наиболее вероятное распределение параметров моделирования, которое лучше всего объясняет данные, называемое апостериорным распределением. Мы сделаем это в рамках байесовской модели, которая является очень мощным подходом, поскольку позволяет нам учитывать априорные знания и неопределенности, а также обновлять наши представления о параметрах модели по мере того, как мы наблюдаем больше данных. Мы выберем апостериорное распределение параметров модели, используя простой и общий метод Монте-Карло с цепями Маркова (MCMC), известный как алгоритм Метрополиса-Гастингса.

Вы можете найти сопроводительный Код Python на GitHub.

Прежде чем мы начнем, ниже показано, как выглядит запуск нашего алгоритма подбора модели:

Байесовская структура

Мы будем использовать байесовскую структуру для подбора модели и оценки параметров. Предположим, что есть модель с набором параметров θ и задан набор данных D. Мы пытаемся найти апостериорное распределение P(θ|D) для параметров модели с использованием теоремы Байеса:

где

  • P(θ|D) – это апостериорнаявероятность того, что параметры имеют значение θ, при условии, что данные D верны.
  • P(D|θ) — вероятность данных D, при условии, что θ правда.
  • P(θ) – это априорнаявероятность того, что набор параметров θ true (до просмотра каких-либо данных).
  • P(D) — это маргинализация, то есть вероятность того, что D истинна.

Это все, что нам нужно для байесовского вывода. Простая, но мощная идея.

Цепь Маркова Метрополис-Гастингс Монте-Карло

Алгоритм Metropolis-Hastings MCMC будет случайным образом выбирать апостериорное распределение. Цепь Маркова — это стохастический процесс, в котором каждое состояние в цепи зависит только от предыдущего состояния. Алгоритм будет:

  1. Нарисуйте случайное значение для θиз предыдущего распределения, θprev
  2. Для i = 1 … N, где N — длина цепи Маркова:
  • Предложите новое значение θprop для добавления в цепочку, добавив случайное возмущение к θprev с помощью распространение предложения
  • Оцените апостериорные вероятности Pprop и Pprev параметров θprop и θ пред.
  • Нарисуйте случайное число U из равномерного распределения от 0 до 1.
  • Если U ‹min(1, Pprop/Pprev), добавьте θprop к цепи Маркова и установите θprev=θprop, иначе добавьте (еще одну копию) θprev в цепочку.

3. Отрезать начало цепочки (область «пригара»), где параметры модели далеки от хорошего соответствия.

С точки зрения кода это выглядит так:

Априорное распределение может, например, быть равномерным распределением с минимальным и максимальным значением, если о параметрах известно немногое.

Распределение предложений (в функции «propose()») — это случайное отклонение, которое добавляется к каждому параметру как способ обхода пространства параметров. Например, шаг возмущения может быть случайным значением, взятым из гауссовой кривой. В этом случае важно установить стандартное отклонение Гаусса. Имеет смысл, чтобы возмущение не было слишком большим или слишком маленьким, поэтому в этом примере мы установили его равным 2% от диапазона параметров. Точный размер шага не имеет большого значения, но он может повлиять на скорость принятия/отклонения предложения в алгоритме MCMC. Если параметры ограничены, нам также необходимо рассмотреть случай, когда предлагаемый набор параметров θprop выходит за границы; в таком случае мы можем отразить значение через границу, чтобы вернуться в пределы.

Наконец, в основе алгоритма лежит вычисление апостериорного значения. И наиболее важной частью этого шага является расчет вероятности, которая зависит от подбираемой модели и параметров (подробнее об этом в следующем разделе). Поскольку все, что имеет значение в алгоритме MCMC, — это отношение апостериорных значений Pprop/Pprev, нам не нужно вычислять маргинализацию P(D), поскольку она аннулируется в соотношение. Кроме того, если априорные значения являются однородными/плоскими, мы также можем игнорировать их вычисление, потому что они также отменяются.

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

Пример: подгонка параметров орбиты экзопланеты к фиктивному набору данных о радиальной скорости

Мы разработаем конкретный пример использования MCMC Метрополис-Гастингс для подбора (фиктивного) набора данных измерений лучевой скорости экзопланеты для восстановления орбитальных параметров модели.

Телескоп измеряет «колебание» (лучевую скорость) звезды из-за экзопланеты, вращающейся вокруг нее. Эта лучевая скорость v(t), измеренная в несколько моментов времени, t, зависит от параметров орбиты:

где

Данные измерений vᵢ в моменты времени tᵢс (гауссовскими) ошибкамиσᵢ,а также с неизвестным звездным дрожанием s (еще один параметр, который мы также используем в нашем анализе),вероятность заданного набора параметров θ=[V, s, K, ϖ, e, P,χ]:

Запуск кода позволяет визуализировать алгоритм MCMC Метрополиса-Хастингса в режиме реального времени, поскольку он используется для подбора (фиктивного) набора данных измерений радиальной скорости экзопланеты. Он отображает модели, которые выбирают апостериорное распределение параметров, а также гистограммы окончательного апостериорного распределения вместе со значениями «наземной истины» (вертикальные красные линии). Мы видим, что алгоритм MCMC достаточно хорошо восстанавливает реальные значения, поскольку все красные линии находятся вблизи пиковых значений гистограммы.

Алгоритм MCMC довольно прост и эффективен и может применяться ко всем видам задач подбора модели. Ниже приведен пример значений апостериорных параметров, полученных в результате подгонки набора космологических данных в Pellejero-Ibañez et al. (2009), демонстрируя возможные корреляции между некоторыми подгоночными параметрами.

Ниже показана расширенная идея MCMC Metropolis-Hastings, называемая Вложенная выборка, которая используется для сопоставления набора данных о радиальных скоростях двух планет. Вложенная выборка имеет то преимущество, что она имеет несколько частиц, одновременно пересекающих пространство параметров (лучшее покрытие), а также оценивает показатель доказательства модели, который можно использовать для выбора между двумя конкурирующими моделями.

Загрузите код Python на GitHub для нашего алгоритма Metropolis-Hastings MCMC, чтобы визуализировать подбор параметров в режиме реального времени и, возможно, изменить его для применения к вашим собственным задачам моделирования. Наслаждаться!

Повышение уровня кодирования

Спасибо, что являетесь частью нашего сообщества! Перед тем, как ты уйдешь:

  • 👏 Хлопайте за историю и подписывайтесь на автора 👉
  • 📰 Смотрите больше контента в публикации Level Up Coding
  • 💰 Бесплатный курс собеседования по программированию ⇒ Просмотреть курс
  • 🔔 Подписывайтесь на нас: Twitter | ЛинкедИн | "Новостная рассылка"

🚀👉 Присоединяйтесь к коллективу талантов Level Up и найдите прекрасную работу