# Using a Multivariate Normal as a prior

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?

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
c = var
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