Using a Multivariate Normal as a prior

Questions
1

I was wondering how to go about using a the multivariate normal as a prior distribution. Suppose I have a simple linear regression problem where I want to sample from a line’s gradient m and y-intercept c, where prior on these two parameters is a correlated normal distribution with a known mean and covariance matrix. If the parameters were uncorrelated, and I had observed data with a Gaussian likelihood and known standard deviation, I would just do, e.g.

meanm = 0.
meanc = 5.
sigmam = 1.
sigmac = 3.

data = ...  # some data defined at points x
x = ...       # the points at which x is defined

with pm.Model() as model:
    m = pm.Normal('m', mu=meanm, sd=sigmam)
    c = pm.Normal('c', mu=meanc, sd=sigmac)

    obs = pm.Normal('obs', mu=m*x + c, sd=1., obs=data)
    trace = pm.sample()

But, how would I do something similar with the MvNormal distribution? Can I somehow unpack the m and c values from it to define my model to pass to the likelihood?

2

I think I’ve answered my own question with some simple experimentation. I can do, e.g.:

with pm.Model() as model:
    var = pm.MvNormal('var', mu=[0., 5.], cov=np.array(([1., 0.5], [0.5, 1.0])), shape=2)
    m = var[0]
    c = var[1]
    model = pm.Normal('model', mu=m*x + c, sd=1., observed=data)
    trace = pm.sample()

and trace['var'] contains an array sampled over both variables.

2 Likes