Использование pymc3 для соответствия модели lomax - PullRequest
Новая
       5

Использование pymc3 для соответствия модели lomax

0 голосов
/ 23 мая 2018

У меня есть довольно простой пример, который не работает.Моя цель - построить модель Lomax, и, поскольку PyMC3 не имеет распределения Lomax, я использую тот факт, что экспонента, смешанная с гаммой, является Lomax (см. здесь ): 

import pymc3 as pm
from scipy.stats import lomax

# Generate artificial data with a shape and scale parameterization
data = lomax.rvs(c=2.5, scale=3, size=1000)

# if t ~ Exponential(lamda) and lamda ~ Gamma(shape, rate), then t ~ Lomax(shape, rate)
with pm.Model() as hierarchical:
    shape = pm.Uniform('shape', 0, 10)
    rate = pm.Uniform('rate', 0 , 10)
    lamda = pm.Gamma('lamda', alpha=shape, beta=rate)
    t = pm.Exponential('t', lam=lamda, observed=data)
    trace = pm.sample(1000, tune=1000)

Сводка: 

>>> pm.summary(trace)
           mean        sd  mc_error   hpd_2.5  hpd_97.5   n_eff      Rhat
shape  4.259874  2.069418  0.060947  0.560821  8.281654  1121.0  1.001785
rate   6.532874  2.399463  0.068837  2.126299  9.998271  1045.0  1.000764
lamda  0.513459  0.015924  0.000472  0.483754  0.545652  1096.0  0.999662

Я ожидаю, что оценки формы и скорости будут близки к 2,5 и 3 соответственно.Я пробовал различные неинформативные приоры для формы и скорости, включая pm.HalfFlat() и pm.Uniform(0, 100), но оба приводили к худшим припадкам.Есть идеи?

1 Ответ

0 голосов
/ 29 мая 2018

Понял: для получения lomax из смеси экспоненциально-гамма мне нужно указать lamda для каждого примера в наборе данных (lamda = pm.Gamma('lamda', alpha=shape, beta=rate, shape=len(data)).Это связано с тем, что модель предполагает, что каждый субъект в данных имеет свой собственный lamda_i, где lamda_i ~ Gamma(shape, rate) для каждого i.

...