Я изучаю пакет PyMC3, где меня интересовала реализация пакета в сценарии, в котором у меня есть несколько разных сигналов, и каждый сигнал имеет разные амплитуды.
Однако я застрял в том, какой тип приоров мне нужно использовать, чтобы реализовать в нем PyMC3 и реализовать распределение вероятности. Пример сценария показан на следующем изображении:
Я пытался реализовать его здесь, но каждый раз получаю сообщение об ошибке:
pymc3.exceptions.SamplingError: Bad initial energy
Мой код
## Signal 1:
with pm.Model() as model:
# Parameters:
# Prior Distributions:
# BoundedNormal = pm.Bound(pm.Exponential, lower=0.0, upper=np.inf)
# c = BoundedNormal('c', lam=10)
# c = pm.Uniform('c', lower=0, upper=300)
alpha = pm.Normal('alpha', mu = 0, sd = 10)
beta = pm.Normal('beta', mu = 0, sd = 1)
sigma = pm.HalfNormal('sigma', sd = 1)
mu = pm.Normal('mu', mu=0, sigma=1)
sd = pm.HalfNormal('sd', sigma=1)
# Observed data is from a Multinomial distribution:
# Likelihood distributions:
# bradford = pm.DensityDist('observed_data', logp=bradford_logp, observed=dict(value=S1, loc=mu, scale=sd, c=c))
# observed_data = pm.Beta('observed_data', mu=mu, sd=sd, observed=S1)
observed_data = pm.Beta('observed_data', alpha=alpha, beta=beta, mu=mu, sd=sd, observed=S1)
with model:
# obtain starting values via MAP
startvals = pm.find_MAP(model=model)
# instantiate sampler
# step = pm.Metropolis()
step = pm.HamiltonianMC()
# step = pm.NUTS()
# draw 5000 posterior samples
trace = pm.sample(start=startvals, draws=1000, step=step, tune=500, chains=4, cores=1, discard_tuned_samples=True)
# Obtaining Posterior Predictive Sampling:
post_pred = pm.sample_posterior_predictive(trace, samples=500)
print(post_pred['observed_data'].shape)
plt.title('Trace Plot of Signal 1')
pm.traceplot(trace, var_names=['mu', 'sd'], divergences=None, combined=True)
plt.show(block=False)
plt.pause(5) # Pauses the program for 5 seconds
plt.close('all')
pm.plot_posterior(trace, var_names=['mu', 'sd'])
plt.title('Posterior Plot of Signal 1')
plt.show(block=False)
plt.pause(5) # Pauses the program for 5 seconds
plt.close('all')
Дополнительные вопросы
Я также изучал идею реализации проверки пригодности и фильтра Калмана при использовании других дистрибутивов, кроме гауссовского, так что, если вы есть время, я был бы признателен, если бы вы могли взглянуть на них ?. Оба вопроса можно найти здесь:
Ссылка на тесты соответствия: Тест соответствия
Ссылка на фильтр Калмана: Фильтр Калмана
Редактировать 1
Скажем, у меня есть около 5 сигналов, и я хотел бы реализовать байесовский интерфейс, чтобы увидеть разницу в PDF-файлах сигналов. Как мне подойти к этой проблеме? Нужно ли мне создавать несколько моделей и получать их апостериорное распределение? Как на этом изображении:
Если мне нужно получить апостериорное распределение, могу ли я использовать следующий код?
# Obtaining Posterior Predictive Sampling:
post_pred = pm.sample_posterior_predictive(trace, samples=500)
Edit 2
Если у меня несколько сигналов, могу ли я реализовать это таким образом, чтобы увидеть изменения в alpha и beta во всех сигналах?
observed_data_S1 = pm.Beta('observed_data_S1', alpha=alpha[0], beta=beta[0], observed=S1[0])
observed_data_S2 = pm.Beta('observed_data_S2', alpha=alpha[1], beta=beta[1], observed=S2[0])
observed_data_S3 = pm.Beta('observed_data_S3', alpha=alpha[2], beta=beta[2], observed=S3[0])
observed_data_S4 = pm.Beta('observed_data_S4', alpha=alpha[3], beta=beta[3], observed=S4[0])
observed_data_S5 = pm.Beta('observed_data_S5', alpha=alpha[4], beta=beta[4], observed=S5[0])
observed_data_S6 = pm.Beta('observed_data_S6', alpha=alpha[5], beta=beta[5], observed=S6[0])
Редактировать 3:
как построить несколько трасс на одном графике? потому что я смотрю на несколько сигналов и подумываю объединить все альфа и бета в один график.