Байесовская вероятностная матричная факторизация (BPMF) с PyMC3: PositiveDefiniteError с использованием `NUTS`

Question

Я реализовал алгоритм Bayesian Probabilistic Matrix Factorization, используя pymc3 в Python. Я также реализовал его предшественник, вероятностную матричную факторизацию (PMF). See my previous question для ссылки на данные, используемые здесь.

У меня возникли проблемы с отображением образцов MCMC с использованием пробоотборника NUTS. Я инициализирую параметры модели, используя MAP из PMF, и гиперпараметры, использующие гауссовские случайные дроби, посыпанные вокруг 0. Однако при настройке объекта шага для сэмплера я получаю PositiveDefiniteError. Я проверил, что оценка MAP от PMF является разумной, поэтому я ожидаю, что она имеет какое-то отношение к тому, как инициализируются гиперпараметры. Вот модель PMF:

import pymc3 as pm 
import numpy as np 
import pandas as pd 
import theano 
import scipy as sp 

data = pd.read_csv('jester-dense-subset-100x20.csv')  
n, m = data.shape 
test_size = m/10 
train_size = m - test_size 

train = data.copy() 
train.ix[:,train_size:] = np.nan # remove test set data 
train[train.isnull()] = train.mean().mean() # mean value imputation 
train = train.values 

test = data.copy() 
test.ix[:,:train_size] = np.nan # remove train set data 
test = test.values  

# Low precision reflects uncertainty; prevents overfitting 
alpha_u = alpha_v = 1/np.var(train) 
alpha = np.ones((n,m)) * 2 # fixed precision for likelihood function 
dim = 10 # dimensionality 

# Specify the model. 
with pm.Model() as pmf: 
    pmf_U = pm.MvNormal('U', mu=0, tau=alpha_u * np.eye(dim), 
         shape=(n, dim), testval=np.random.randn(n, dim)*.01) 
    pmf_V = pm.MvNormal('V', mu=0, tau=alpha_v * np.eye(dim), 
         shape=(m, dim), testval=np.random.randn(m, dim)*.01) 
    pmf_R = pm.Normal('R', mu=theano.tensor.dot(pmf_U, pmf_V.T), 
         tau=alpha, observed=train) 

    # Find mode of posterior using optimization 
    start = pm.find_MAP(fmin=sp.optimize.fmin_powell) 

А вот BPMF:

n, m = data.shape 
dim = 10 # dimensionality 
beta_0 = 1 # scaling factor for lambdas; unclear on its use 
alpha = np.ones((n,m)) * 2 # fixed precision for likelihood function 

logging.info('building the BPMF model') 
std = .05 # how much noise to use for model initialization 
with pm.Model() as bpmf: 
    # Specify user feature matrix 
    lambda_u = pm.Wishart(
     'lambda_u', n=dim, V=np.eye(dim), shape=(dim, dim), 
     testval=np.random.randn(dim, dim) * std) 
    mu_u = pm.Normal(
     'mu_u', mu=0, tau=beta_0 * lambda_u, shape=dim, 
     testval=np.random.randn(dim) * std) 
    U = pm.MvNormal(
     'U', mu=mu_u, tau=lambda_u, shape=(n, dim), 
     testval=np.random.randn(n, dim) * std) 

    # Specify item feature matrix 
    lambda_v = pm.Wishart(
     'lambda_v', n=dim, V=np.eye(dim), shape=(dim, dim), 
     testval=np.random.randn(dim, dim) * std) 
    mu_v = pm.Normal(
     'mu_v', mu=0, tau=beta_0 * lambda_v, shape=dim, 
     testval=np.random.randn(dim) * std) 
    V = pm.MvNormal(
     'V', mu=mu_v, tau=lambda_v, shape=(m, dim), 
     testval=np.random.randn(m, dim) * std) 

    # Specify rating likelihood function 
    R = pm.Normal(
     'R', mu=theano.tensor.dot(U, V.T), tau=alpha, 
     observed=train) 

# `start` is the start dictionary obtained from running find_MAP for PMF. 
for key in bpmf.test_point: 
    if key not in start: 
     start[key] = bpmf.test_point[key] 

with bpmf: 
    step = pm.NUTS(scaling=start) 

В последней строке, я получаю следующее сообщение об ошибке:

PositiveDefiniteError: Scaling is not positive definite. Simple check failed. Diagonal contains negatives. Check indexes [ 0 2 ... 2206 2207 ] 

Как я понимаю, я могу 't используйте find_MAP с моделями, у которых есть гиперпорции, такие как BPMF. Вот почему я пытаюсь инициализировать значения MAP из PMF, которые используют точечные оценки параметров на U и V, а не параметризованные гиперприорности.

Answer 1

К сожалению, распределение Уишарта нефункционально. Я недавно добавил предупреждения здесь: https://github.com/pymc-devs/pymc3/commit/642f63973ec9f807fb6e55a0fc4b31bdfa1f261e

Смотрите здесь для более дискуссий по этому хитрому распределению: https://github.com/pymc-devs/pymc3/issues/538

Вы могли бы подтвердить, что это источник, фиксируя ковариационную матрицу. Если это так, я бы попытался использовать предыдущий дистрибутив JKL: https://github.com/pymc-devs/pymc3/blob/master/pymc3/examples/LKJ_correlation.py