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

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

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

В конце концов, благодаря этому итеративному процессу в нашем саду остаются только те пути, которые логически согласуются как с нашими предыдущими знаниями, так и с наблюдаемыми данными. Эти оставшиеся ответвления представляют собой наиболее правдоподобные объяснения, которые соответствуют тому, что мы знаем. Тщательно сокращая и выбирая логически непротиворечивые пути, мы можем получить уточненный набор гипотез, которые лучше всего объясняют наблюдаемые данные.

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

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

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

Синий-белый-синий

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

Мы можем представить, как возможности расходятся, как ветви. Всего нужно рассмотреть 64 потенциальных пути (4^3=64). Только 3 согласуются как с нашим предположением о том, что находится в сумке, так и с данными, которые мы собрали.

Для каждой возможной гипотезы способы получения наблюдаемых данных:

Сочетание трех синих шариков и одного белого шарика обеспечивает наибольшую вероятность получения наблюдаемой последовательности. Эти первоначальные подсчеты шариков составляют нашу отправную точку, и впоследствии мы собираемся пересмотреть или обновить их. Назовем эти начальные подсчеты «предварительными».

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

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

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

ways = np.array([0, 3, 8, 9, 0])
ways / ways.sum()

# array([0.  , 0.15, 0.4 , 0.45, 0.  ])
  • p называется значением параметра. Это просто способ индексации возможных объяснений данных.
  • Относительная частота значения p, генерирующего данные, называется вероятностью. Он определяется путем рассмотрения всех возможных данных и удаления тех, которые не согласуются с наблюдаемыми данными.
  • Априорная правдоподобие любого конкретного p – это априорная вероятность.
  • Новая обновленная правдоподобие любого конкретного p — это апостериорная вероятность.

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

W L W W W L W L W

У нас есть 9 наблюдений (6 водных и 3 наземных), эта последовательность и есть наши данные.

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

  1. История данных: начните с построения повествования, объясняющего, как могли быть получены наблюдаемые данные. Этот процесс повествования помогает мотивировать и обеспечивать контекст для модели.
  2. Обновление: перейдите к обучению модели, введя доступные данные. Этот шаг включает в себя включение наблюдаемых данных в модель и обновление ее убеждений или параметров на основе этой новой информации.
  3. Оценка: Как и в случае любой статистической модели, важно оценивать и оценивать ее эффективность. Этот процесс оценки может выявить необходимость пересмотра или корректировки модели для повышения ее точности и полезности. Этот шаг включает в себя изучение соответствия модели данным, оценку ее прогностических способностей и рассмотрение любых ограничений или областей для улучшения.

Для нашей конкретной ситуации мы можем построить историю следующим образом:

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

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

Каждая пунктирная кривая представляет собой сплошную кривую предыдущего графика с переходом слева направо и сверху вниз. Всякий раз, когда наблюдается буква «W» (обозначающая воду), пик кривой правдоподобия смещается в сторону больших значений p, что указывает на более высокую вероятность большей доли воды. И наоборот, когда наблюдается буква «L» (обозначающая сушу), пик перемещается в противоположном направлении, что предполагает более высокую вероятность меньших значений p, что указывает на большую долю суши.

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

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

Компоненты модели

Переменные – это символы, которые могут иметь разные значения. Первая представляющая интерес переменная — p, представляющая долю воды на земном шаре. p не наблюдается (называется параметром), но его можно вывести из других переменных. Другими переменными являются наблюдаемые количества воды и суши.

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

  1. Каждый бросок не зависит от других,
  2. Вероятность получения воды остается постоянной для каждого броска.

Теория вероятностей предлагает окончательное решение, известное как биномиальное распределение.

Количество «воды» W и «земли» L распределяется биномиально, с вероятностью p «воды» при каждом броске.

W и L являются наблюдаемыми переменными.

stats.binom.pmf(6, n=9, p=0.5)

# 0.16406250000000003

Функция stats.binom.pmf() вычисляет функцию массы вероятности (PMF) биномиального распределения. PMF дает вероятность получения определенного количества успехов (в данном случае «6») в заданном количестве испытаний (указанном n) с фиксированной вероятностью успеха в каждом испытании (указанном p).

Итак, вызывая stats.binom.pmf(6, n=9, p=0.5), мы вычисляем вероятность получения ровно 6 успехов (или 6 наблюдений за водой) в 9 независимых попытках (вращениях), при этом вероятность успеха (получения воды) в каждой попытке равна 0,5 (50%).

probabilities = [stats.binom.pmf(6, n=9, p=(x+1)*0.1) for x in range(10)]
plt.plot(np.arange(0, 1., 0.1), probabilities, 'bo-')
plt.xlabel('Probability (p)')
plt.ylabel('Probability of 6 successes')
plt.title('Probability of 6 successes in 9 trials for different p values')
plt.show()

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

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

Теперь мы можем получить нашу модель.

W ∼ Биномиальный(N, p)

W в N испытаниях с вероятностью p.

p ∼ Равномерное(0, 1)

p имеет равномерную — плоскую — априорную оценку во всем возможном диапазоне.

Байесовская модель имеет возможность обновлять все предыдущие распределения до их логических следствий, что приводит к апостериорному распределению. Апостериорное распределение представляет собой распределение вероятностей параметров с учетом наблюдаемых данных. В данном конкретном случае мы обозначаем его как Pr(p|W, L), что означает вероятность каждого возможного значения p в зависимости от конкретных наблюдений W и L, которые мы видели.

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



В нем говорится, что вероятность любого конкретного значения p с учетом данных равна произведению относительной правдоподобности данных, зависящей от p, и априорной правдоподобности p, деленной на эту вещь Pr(W, L), которую я буду называть средней вероятностью данных.

Вероятность Pr(W, L) представляет собой среднюю вероятность наблюдаемых данных, которая вычисляется путем усреднения по предшествующему распределению. Его цель состоит в том, чтобы нормализовать или стандартизировать апостериорное распределение, гарантируя, что вероятности в пределах апостериорного распределения суммируются (или интегрируются) до 1.

E — это ожидание. Вычисление среднего по непрерывному распределению значений, например бесконечному количеству возможных значений p.

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

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

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

Аппроксимация сетки

Сетчатая аппроксимация — это простой метод кондиционирования в байесовском выводе. Вместо рассмотрения бесконечного числа значений параметров используется конечная сетка значений параметров для аппроксимации непрерывного апостериорного распределения.

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

Аппроксимация сетки нецелесообразна для сложного моделирования из-за проблем масштабируемости с увеличением числа параметров.

def posterior_grid_approx(grid_points=5, success=6, tosses=9):
    """
    """
    # define grid
    p_grid = np.linspace(0, 1, grid_points)

    # define prior
    prior = np.repeat(5, grid_points)  # uniform
    #prior = (p_grid >= 0.5).astype(int)  # truncated
    #prior = np.exp(- 5 * abs(p_grid - 0.5))  # double exp

    # compute likelihood at each point in the grid
    likelihood = stats.binom.pmf(success, tosses, p_grid)

    # compute product of likelihood and prior
    unstd_posterior = likelihood * prior

    # standardize the posterior, so it sums to 1
    posterior = unstd_posterior / unstd_posterior.sum()
    return p_grid, posterior
points = 20
w, n = 6, 9
p_grid, posterior = posterior_grid_approx(points, w, n)
plt.plot(p_grid, posterior, 'o-', label='success = {}\ntosses = {}'.format(w, n))
plt.xlabel('probability of water', fontsize=14)
plt.ylabel('posterior probability', fontsize=14)
plt.title('{} points'.format(points))
plt.legend(loc=0);

Квадратичная аппроксимация

Квадратичная аппроксимация является ценным подходом в байесовском выводе. Он основан на наблюдении, что при общих условиях область вокруг пика апостериорного распределения очень напоминает распределение Гаусса. Это приближение Гаусса упрощает апостериорное распределение, требуя только два параметра: среднее значение (центр) и дисперсию (разброс).

Термин «квадратичное приближение» возникает из-за того, что логарифм распределения Гаусса образует параболу, а парабола является квадратичной функцией. Следовательно, это приближение представляет собой логарифм апостериорного распределения с использованием параболической формы.

Используя приближение Гаусса, сложное апостериорное распределение может быть эффективно аппроксимировано более простым параболическим представлением.

data = np.repeat((0, 1), (3, 6))
with pm.Model() as normal_approximation:
    p = pm.Uniform('p', 0, 1)
    w = pm.Binomial('w', n=len(data), p=p, observed=data.sum())
    mean_q = pm.find_MAP()
    std_q = ((1/pm.find_hessian(mean_q, vars=[p]))**0.5)[0]
mean_q['p'], std_q

# (array(0.66666667), array([0.15713484]))
# shows the posterior mean value of p = 0.67
# the standard deviation of the posterior distribution 0.16
# Assuming the posterior is Gaussian, it is maximized at 0.67, 
# and its standard deviation is 0.16

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

norm = stats.norm(mean_q, std_q)
prob = .89
z = stats.norm.ppf([(1-prob)/2, (1+prob)/2])
pi = mean_q['p'] + std_q * z 
pi
# array([0.41553484, 0.91779849]) the lower and upper bounds of the confidence interval.

Используя stats.norm.ppf, мы вычисляем z-показатели, соответствующие нижнему и верхнему хвостам желаемого уровня вероятности. Часть [(1-prob)/2, (1+prob)/2] гарантирует, что мы вычислим z-показатели для двусторонней области.

Переменная pi представляет собой массив значений. Он вычисляет нижнюю и верхнюю границы доверительного интервала, умножая стандартное отклонение (std_q) на ранее рассчитанные z-значения и добавляя среднее значение (mean_q['p']).

# analytical calculation
w, n = 6, 9
x = np.linspace(0, 1, 100)
plt.plot(x, stats.beta.pdf(x , w+1, n-w+1),
         label='True posterior')

# quadratic approximation
plt.plot(x, stats.norm.pdf(x, mean_q['p'], std_q),
         label='Quadratic approximation')
plt.legend(loc=0, fontsize=13)

plt.title('n = {}'.format(n), fontsize=14)
plt.xlabel('Proportion water', fontsize=14)
plt.ylabel('Density', fontsize=14);

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

Цепь Маркова Монте-Карло

В случае сложных моделей, таких как многоуровневые модели с большим количеством параметров, сеточная аппроксимация и квадратичная аппроксимация могут оказаться неудовлетворительными. Сеточная аппроксимация становится непрактичной из-за вычислительных ограничений, а квадратичная аппроксимация может плохо соответствовать модели. Кроме того, в многоуровневых моделях часто отсутствует унифицированная функция для апостериорного распределения, что затрудняет поиск максимальной апостериорной оценки (MAP).

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

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

Источники

Все изображения и формулы, представленные в этом материале, взяты непосредственно из книги «Статистическое переосмысление, 2-е издание» Ричарда МакЭлрита, если не указано иное.

https://github.com/aloctavodia/Statistical-Rethinking-with-Python-and-PyMC3/blob/master/Chp_02.ipynb

Читать далее