使用`LKJCorr`优先级在PyMC3中修改的BPMF:使用`NUTS`的PositiveDefiniteError [英] Modified BPMF in PyMC3 using `LKJCorr` priors: PositiveDefiniteError using `NUTS`

问题描述

我以前实现了原始的贝叶斯概率矩阵分解(BPMF)pymc3中的a>模型. 请参阅我以前的问题以获取参考，数据源和问题设置.根据@twiecki对这个问题的回答，我使用相关关系矩阵的LKJCorr先验和标准差的统一先验实现了模型的变体.在原始模型中，协方差矩阵来自Wishart分布，但是由于pymc3的当前限制，无法从Wishart分布中正确采样.对一个松散相关问题的此答案为LKJCorr先验的选择提供了简洁的解释.新模型如下.

I previously implemented the original Bayesian Probabilistic Matrix Factorization (BPMF) model in pymc3. See my previous question for reference, data source, and problem setup. Per the answer to that question from @twiecki, I've implemented a variation of the model using LKJCorr priors for the correlation matrices and uniform priors for the standard deviations. In the original model, the covariance matrices are drawn from Wishart distributions, but due to current limitations of pymc3, the Wishart distribution cannot be sampled from properly. This answer to a loosely related question provides a succinct explanation for the choice of LKJCorr priors. The new model is below.

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


n, m = train.shape
dim = 10  # dimensionality
beta_0 = 1  # scaling factor for lambdas; unclear on its use
alpha = 2  # fixed precision for likelihood function
std = .05  # how much noise to use for model initialization

# We will use separate priors for sigma and correlation matrix.
# In order to convert the upper triangular correlation values to a
# complete correlation matrix, we need to construct an index matrix:
n_elem = dim * (dim - 1) / 2
tri_index = np.zeros([dim, dim], dtype=int)
tri_index[np.triu_indices(dim, k=1)] = np.arange(n_elem)
tri_index[np.triu_indices(dim, k=1)[::-1]] = np.arange(n_elem)

logging.info('building the BPMF model')
with pm.Model() as bpmf:
    # Specify user feature matrix
    sigma_u = pm.Uniform('sigma_u', shape=dim)
    corr_triangle_u = pm.LKJCorr(
        'corr_u', n=1, p=dim,
        testval=np.random.randn(n_elem) * std)

    corr_matrix_u = corr_triangle_u[tri_index]
    corr_matrix_u = t.fill_diagonal(corr_matrix_u, 1)
    cov_matrix_u = t.diag(sigma_u).dot(corr_matrix_u.dot(t.diag(sigma_u)))
    lambda_u = t.nlinalg.matrix_inverse(cov_matrix_u)

    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
    sigma_v = pm.Uniform('sigma_v', shape=dim)
    corr_triangle_v = pm.LKJCorr(
        'corr_v', n=1, p=dim,
        testval=np.random.randn(n_elem) * std)

    corr_matrix_v = corr_triangle_v[tri_index]
    corr_matrix_v = t.fill_diagonal(corr_matrix_v, 1)
    cov_matrix_v = t.diag(sigma_v).dot(corr_matrix_v.dot(t.diag(sigma_v)))
    lambda_v = t.nlinalg.matrix_inverse(cov_matrix_v)

    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,
        testval=np.random.randn(m, dim) * std)

    # Specify rating likelihood function
    R = pm.Normal(
        'R', mu=t.dot(U, V.T), tau=alpha * np.ones((n, m)),
        observed=train)

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

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

此重新实现的目标是生成一个可以使用NUTS采样器估算的模型.不幸的是，我在最后一行仍然遇到相同的错误:

The goal with this reimplementation was to produce a model that could be estimated using the NUTS sampler. Unfortunately, I'm still getting the same error at the last line:

PositiveDefiniteError: Scaling is not positive definite. Simple check failed. Diagonal contains negatives. Check indexes [   0    1    2    3    ...   1030 1031 1032 1033 1034   ]

我已经在此要点使复制错误变得简单.您所需要做的就是下载数据(要点中也有引用).

I've made all the code for PMF, BPMF, and this modified BPMF available in this gist to make it simple to replicate the error. All you need to do is download the data (also referenced in the gist).

推荐答案

您似乎正在将完整的精度矩阵传递给正态分布:

It looks like you are passing the complete precision matrix into the normal distribution:

mu_u = pm.Normal(
    'mu_u', mu=0, tau=beta_0 * lambda_u, shape=dim,
     testval=np.random.randn(dim) * std)

我假设您只想传递对角线值:

I assume you only want to pass the diagonal values:

mu_u = pm.Normal(
    'mu_u', mu=0, tau=beta_0 * t.diag(lambda_u), shape=dim,
     testval=np.random.randn(dim) * std)

对mu_u和mu_v的此更改是否为您解决?

Does this change to mu_u and mu_v fix it for you?

这篇关于使用`LKJCorr`优先级在PyMC3中修改的BPMF:使用`NUTS`的PositiveDefiniteError的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！