Пример нескольких монет из одного монетного двора в PyMC

Пытаюсь изучить PyMC, перенеся некоторые модели из книги «Выполнение байесовского анализа данных» (Крушке). Один из основных примеров (из главы 9) состоит в том, чтобы предположить, что набор монет распределяется в соответствии с p~Bern(theta), где тета происходит из Beta распределения («чеканного двора») с фиксированными параметрами. Вот как я это закодировал (в PyMC2):

import pymc as pm
import numpy as np
from matplotlib import pyplot as plt
import seaborn as sbn
from pymc.Matplot import plot as mcplot

from pymc import Bernoulli, Beta, Gamma

flips = [[True, False, False, False],
         [False, False, False, True],
         [True, False, False, False],
         [False, False, False, False]]

mint = Beta('mint', alpha=2, beta=2)

coin0 = Bernoulli('coin0', p=mint, value=flips[0], observed=True)
coin1 = Bernoulli('coin1', p=mint, value=flips[1], observed=True)
coin2 = Bernoulli('coin2', p=mint, value=flips[2], observed=True)
coin3 = Bernoulli('coin3', p=mint, value=flips[3], observed=True)

mcmc = pm.MCMC([mint, coin0, coin1, coin2, coin3])
mcmc.sample(iter=10000, burn=1000)

mcmc.summary()
mint:

    Mean             SD               MC Error        95% HPD interval
    ------------------------------------------------------------------
    [[ 0.253]]       [[ 0.096]]       [[ 0.002]]       [ 0.074  0.439]


    Posterior quantiles:

    2.5             25              50              75             97.5
     |---------------|===============|===============|---------------|
    [[ 0.089]]       [[ 0.183]]      [[ 0.242]]     [[ 0.318]]    [[ 0.46]]

Кажется, это сработало, но мне интересно, как мне получить тета-значения для каждой монеты? Я предполагаю, что должны быть образцы, сгенерированные из апостериорного распределения для каждой монеты?


bayesian pymc
person Evan Zamir    schedule 09.03.2015    source источник

Ответы (1)


arrow_upward
3
arrow_downward

Не уверен, что вы подразумеваете под theta, так как в вашей модели нет тета. Вы имеете в виду вероятности монет (которые здесь представлены mint)? Вы указали одну вероятность для всех монет, а не 4 вероятности. Попробуйте изменить параметр mint на:

mint = Beta('mint', alpha=2, beta=2, size=4)

который будет определять векторный стохастик размера 4.

person Chris Fonnesbeck    schedule 11.03.2015
comment
Да, я имею в виду, что мне интересно найти тета для каждой отдельной монеты. Ваше предложение работает и имеет смысл. Единственная другая модификация, которая мне понадобилась, это установить p=mint[i] для каждой монеты. Спасибо! - person Evan Zamir; 11.03.2015