PyMC3 Normal с дисперсией на столбец - PullRequest
Новая
       35

PyMC3 Normal с дисперсией на столбец

1 голос
/ 25 марта 2019

Я пытаюсь определить переменную pymc3.Normal следующим образом: mu: 

import numpy as np
import pymc3 as pm

mx = np.array([[0.25 , 0.5  , 0.75 , 1.   ],    
               [0.25 , 0.333, 0.25 , 0.   ],
               [0.25 , 0.167, 0.   , 0.   ],
               [0.25 , 0.   , 0.   , 0.   ]])
epsilon = pm.Gamma('epsilon', alpha=10, beta=10)
p_ = pm.Normal('p_', mu=mx, shape = mx.shape, sd = epsilon)

Проблема состоит в том, что все случайные величины в p_ получают одинаковый std (epsilon).Я хотел бы, чтобы в первом ряду использовался epsilon1, во втором ряду epsilon2 и т. Д.

Как я могу это сделать?

1 Ответ

0 голосов
/ 01 апреля 2019

Для этого можно передать аргумент для параметра shape. Чтобы продемонстрировать это, давайте сделаем некоторые поддельные данные для передачи в соответствии с наблюдением, где мы используем фиксированные значения для эпсилона, которые мы можем сравнить с предполагаемыми.

Пример модели 

import numpy as np
import pymc3 as pm
import arviz as az

# priors
mu = np.array([[0.25 , 0.5  , 0.75 , 1.   ],    
               [0.25 , 0.333, 0.25 , 0.   ],
               [0.25 , 0.167, 0.   , 0.   ],
               [0.25 , 0.   , 0.   , 0.   ]])
alpha, beta = 10, 10

# fake data
np.random.seed(2019)

# row vector will use a different sd for each column
sd = np.random.gamma(alpha, 1.0/beta, size=(1,4))

# generate 100 fake observations of the (4,4) random variables
Y = np.random.normal(loc=mu, scale=sd, size=(100,4,4))

# true column sd's
print(sd)
# [[0.90055471 1.24522079 0.85846659 1.19588367]]

# mean sd's per column
print(np.mean(np.std(Y, 0), 0))
# [0.92028042 1.24437592 0.83383181 1.22717313]

# model
with pm.Model() as model:
    # use a (1,4) matrix to pool variance by columns
    epsilon = pm.Gamma('epsilon', alpha=10, beta=10, shape=(1, mu.shape[1]))

    p_ = pm.Normal('p_', mu=mu, sd=epsilon, shape=mu.shape, observed=Y)

    trace = pm.sample(random_seed=2019)

Это образцы хорошо, и дает следующее резюме

enter image description here

, который четко ограничивает истинные значения стандартных отклонений в HPD.

...