PyMC3 передает стохастическую ковариационную матрицу в pm.MvNormal () - PullRequest
Новая
       108

PyMC3 передает стохастическую ковариационную матрицу в pm.MvNormal ()

0 голосов
/ 09 ноября 2018

Я пытался приспособить простую двухмерную гауссову модель к данным наблюдений, используя PyMC3. 

import numpy as np
import pymc3 as pm

n = 10000;
np.random.seed(0)
X = np.random.multivariate_normal([0,0], [[1,0],[0,1]], n);

with pm.Model() as model:
    # PRIORS
    mu = [pm.Uniform('mux', lower=-1, upper=1), 
          pm.Uniform('muy', lower=-1, upper=1)]
    cov = np.array([[pm.Uniform('a11', lower=0.1, upper=2), 0],
                    [0, pm.Uniform('a22', lower=0.1, upper=2)]])

    # LIKELIHOOD
    likelihood = pm.MvNormal('likelihood', mu=mu, cov=cov, observed=X)

with model:
    trace = pm.sample(draws=1000, chains=2, tune=1000)

Хотя я могу сделать это в 1D, передав sd в pm.Normal, у меня есть некоторые проблемы с передачей ковариационной матрицы в pm.MvNormal.

Где я иду не так?

1 Ответ

0 голосов
/ 14 ноября 2018

Объекты распространения PyMC3 не являются простыми числовыми объектами или массивами . Вместо этого они являются узлами в графике анонимных вычислений и часто требуют операций из pymc3.math или theano.tensor, чтобы манипулировать ими. Более того, помещать объекты PyMC3 в массивы numy не нужно, поскольку они уже многомерны.

Оригинальная модель

Придерживаясь намерения вашего кода, рабочая версия будет выглядеть примерно так: 

import numpy as np
import pymc3 as pm
import theano.tensor as tt

N = 10000
np.random.seed(0)
X = np.random.multivariate_normal(np.zeros(2), np.eye(2), size=N)

with pm.Model() as model:
    # use `shape` argument to define tensor dimensions
    mu = pm.Uniform('mu', lower=-1, upper=1, shape=2)

    # diagonal values on covariance matrix
    a = pm.Uniform('a', lower=0.1, upper=2, shape=2)

    # convert vector to a 2x2 matrix with `a` on the diagonal
    cov = tt.diag(a)

    likelihood = pm.MvNormal('likelihood', mu=mu, cov=cov, observed=X)

Альтернативная модель

Я полагаю, что приведенный вами пример - всего лишь игрушка для решения проблемы. Но на всякий случай я упомяну, что основное преимущество использования многомерной нормали (моделирование ковариации между параметрами) теряется при ограничении ковариационной матрицы диагональю. Кроме того, теория априоров для ковариационных матриц хорошо разработана, поэтому стоит задуматься о существующих решениях. В частности, есть пример PyMC3, использующий априор LKJ для ковариационных матриц .

Вот простое применение этого примера в этом контексте: 

with pm.Model() as model_lkj:
    # use `shape` argument to define tensor dimensions
    mu = pm.Uniform('mu', lower=-1, upper=1, shape=2)

    # LKJ prior for covariance matrix (see example)
    packed_L = pm.LKJCholeskyCov('packed_L', n=2,
                                 eta=2., sd_dist=pm.HalfCauchy.dist(2.5))
    # convert to (2,2)
    L = pm.expand_packed_triangular(2, packed_L)

    likelihood = pm.MvNormal('likelihood', mu=mu, chol=L, observed=X)
...