Блог

Я новичок в байесовском мире и PyMC3, и я борюсь с простой настройкой модели. В частности, как работать с настройкой, в которой «наблюдаемые» данные сами изменяются случайными переменными? В качестве примера, допустим, у меня есть набор 2d точек [Xi, Yi], которые образуют дугу вокруг окружности, центральную точку которой [Xc, Yc] я не знаю. Однако я ожидаю, что расстояния между точками и центром окружности, Ri, должны быть нормально распределены, примерно с известным радиусом, R. Поэтому я изначально думал, что могу назначить Xc и Yc равномерные значения (в некотором произвольно большом диапазоне), а затем пересчитать Ri в модели и назначить Ri как»наблюдаемые» данные для получения апостериорных оценок для Xc и Yc:

 import pymc3 as pm
import numpy as np

points = np.array([[2.95, 4.98], [3.28, 4.88], [3.84, 4.59], [4.47, 4.09], [2.1,5.1], [5.4, 1.8]])
Xi = points[:,0]
Yi = points[:,1]

#known [Xc,Yc] = [2.1, 1.8]
R = 3.3

with pm.Model() as Cir_model:
    
    Xc = pm.Uniform('Xc', lower=-20, upper=20)
    Yc = pm.Uniform('Yc', lower=-20, upper=20)
    
    Ri = pm.math.sqrt((Xi-Xc)**2   (Yi-Yc)**2)
    
    y = pm.Normal('y', mu=R, sd=1.0, observed=Ri)
    
    samples = pm.fit(random_seed=2020).sample(1000)
    
    pm.plot_posterior(samples, var_names=['Xc'])
    pm.plot_posterior(samples, var_names=['Yc']);
 

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

Эта модель на самом деле работает нормально, но есть несколько вещей, которые вы могли бы улучшить:

1. Использование переменной в качестве наблюдения не очень хорошо, поскольку вам следует подумать о том, что она делает с дистрибутивом, который вы устанавливаете. Это будет соответствовать распределению, но вам следует подумать о том, проводите ли вы двойной подсчет переменных в предыдущем и вероятностном. Хотя для этой игрушечной модели это не так уж и важно!
2. Вы используете pm.fit(...) , который использует вариационный вывод, но MCMC здесь подходит, поэтому замена всей этой строки на samples = pm.sample() works .
3. Предоставленные points вами данные почти точно расположены по кругу — эмпирическое стандартное отклонение составляет около 0,004, но стандартное отклонение, которое вы указываете в вероятности, составляет 1: около 250 раз больше истинного значения! Выборка из модели как есть позволяет центру точек находиться в двух разных местах:

 with pm.Model():
    Xc = pm.Uniform('Xc', lower=-20, upper=20)
    Yc = pm.Uniform('Yc', lower=-20, upper=20)
    Ri = pm.math.sqrt((Xi-Xc)**2   (Yi-Yc)**2)
    obs_sd = pm.HalfNormal('obs_sd', 1)
    y = pm.Normal('y', mu=R, sd=obs_sd, observed=Ri)
    
    samples = pm.sample()