PyMC3 Dirichlet распределение

Question

Я реализую модель линейной регрессии в pymc3, где неизвестный вектор весов ограничен функцией вероятностной массы, следовательно, смоделирован как распределение Дирихле, как в следующем коде:

with pm.Model() as model:
    #prior on precision of normal likelihood
    tau = pm.Gamma('tau', alpha=1, beta=1)

    phi = np.empty(ncountries, dtype=object)
    y = np.empty((nyears-1, ncountries), dtype=object)
    for icountry, country in enumerate(countries):
        #prior Dirichlet allocation for each country
        phi[icountry] = pm.Dirichlet('mix_{c}'.format(c=country),
                                     np.roll(mix, icountry),
                                     shape=ncountries)

        for iyear, year in enumerate(years[1:]):
            suffix = '_{y}-{c}'.format(y=year, c=country)
            previous_pop = Xs[iyear, :]
            #likelihood
            y[iyear, icountry] = pm.Normal('obs' + suffix,
                        mu=pm.Deterministic(
                            'mu' + suffix,
                            dot(phi[icountry], previous_pop)),
                        tau=tau,
                        observed=Ys[iyear, icountry])

После отбора пробы, запустив:

    start = pm.find_MAP()
    step = pm.Metropolis()
    nsteps = 1000
    trace = pm.sample(nsteps, step, start=start)

Я проанализировал трассировку переменных Дирихле и обнаружил, что их значения не добавляют к одной (ниже приведен пример):

array([[ 0.01029745,  0.00627394,  0.00996922, ...,  1.83955829,
     0.00962185,  0.01020659],
   [ 0.01029745,  0.00627394,  0.00996922, ...,  1.83955829,
     0.00962185,  0.01020659],
   [ 0.01029745,  0.00627394,  0.00996922, ...,  1.83955829,
     0.00962185,  0.01020659],
   ...,
   [ 0.02050308,  0.01685555,  0.01976797, ...,  1.92278065,
     0.03956622,  0.00473735],
   [ 0.01993214,  0.01632033,  0.01994876, ...,  1.92487858,
     0.04078728,  0.00481424],
   [ 0.01900882,  0.01528191,  0.02100671, ...,  1.92485693,
     0.0395159 ,  0.00524575]])

Я не знаком с переменными theano, и мне было трудно исследовать, как RV Дирихле выражается в pymc3... Я делаю что-то не так, или я должен просто нормализовать значения, возвращенные в трассировке, чтобы они суммировались в единицу?

Быстрое обновление похоже на функцию pm.find_MAP() использует своего рода оптимизацию градиентного спуска. При этом не учитывается ограничение, вытекающее из того факта, что вектор, представляющий ничью из распределения Дирихле, является функцией массы вероятности (ее значения должны быть положительными, а их сумма должна быть равна единице). Это ограничение, очевидно, также не применяется на этапе выборки алгоритма и вызывает проблемы сходимости, поскольку точность распределения правдоподобия падает до нуля.