Bayesian correlation coefficient

Questions v5
1

Hi, I was attempting to infer Pearson-like correlation coefficients using this PYMC3-based notebook as basis, but I cannot figure out how to initiate the precision matrix for PYMC5. The core code is this:

import theano.tensor as T

def precision(sigma, rho):
    C = T.alloc(rho, 2, 2)
    C = T.fill_diagonal(C, 1.)
    S = T.diag(sigma)
    return T.nlinalg.matrix_inverse(S.dot(C).dot(S))


def analyze(data):
    with pm.Model() as model:
        # priors might be adapted here to be less flat
        mu = pm.Normal('mu', mu=0., tau=0.000001, shape=2, testval=np.mean(data, axis=1))
        sigma = pm.Uniform('sigma', lower=1e-6, upper=1000., shape=2, testval=np.std(data, axis=1))
        rho = pm.Uniform('r', lower=-1., upper=1., testval=0)

        prec = pm.Deterministic('prec', precision(sigma, rho))
        mult_n = pm.MvNormal('mult_n', mu=mu, tau=prec, observed=data.T)

    return model

I will try to figure out how in the meantime and might answer myself, I think it is the type of functionality someone else might find this useful.

2

The advice is usually to use a LKJCholeskyCov as a prior, because every other form of covariance prior tends to sample poorly.

There are a couple of examples in the gallery that use it.

This one is still in PyMC v4 but should be similar with the latest release: LKJ Cholesky Covariance Priors for Multivariate Normal Models — PyMC example gallery

1 Like
3

Thank so much for referring me to that link, I was unaware of this. Thanks @ricardoV94 ! Marking as solved as well

4

Hello, I have a related question and hope it is ok to ask it here. Is there a good way to get a posterior distribution for a Pearson’s correlation coefficient, if I don’t want to assume normality for the variables?