I am trying to create an observed variable in PyMC3 which has a multivariate normal distribution. I want the covariance matrix whose elements are other random variables. As an example, consider the following code:
import pymc3 as pm
with pm.Model() as model:
a = pm.Normal('a', mu=0, sigma=10)
b = pm.Normal('a', mu=0, sigma=10)
c = pm.Normal('a', mu=0, sigma=10)
# I recognize that the matrix is not positive definite
# with this parameterization. This is just a toy example.
# The main point is that I want the elements of the matrix
# to be random variables.
cov = [[a, b]
[b, c]]
data = pm.MvNormal('data', mu=[0, 0],
cov=cov, observed=obs)
This does not work. Neither does using np.array(cov)
. I imagine the solution is to use Theano tensors somehow. I am unable to figure out how to use them.
I’d appreciate any help with this. Thanks.
You can use
import theano.tensor as tt
tt.concatenate([a, b, b, c]).reshape((2, 2))
Or alternatively tt.stack
.
If you want priors on covariance matrices, you probably want to look into pm.LKJCholeskyCov
instead.
1 Like
Using tt.stack
worked. Thank you. Now, I am getting the following future warning:
/usr/local/lib/python3.6/dist-packages/theano/tensor/basic.py:6611: FutureWarning: Using a non-tuple sequence for multidimensional indexing is deprecated; use `arr[tuple(seq)]` instead of `arr[seq]`. In the future this will be interpreted as an array index, `arr[np.array(seq)]`, which will result either in an error or a different result.
result[diagonal_slice] = x
Do you know why this warning is being thrown and how I can remedy this?
Hi,
Yeah it’s a classical and harmless warning. Usually we just get rid of it with:
import warnings
warnings.simplefilter(action="ignore", category=FutureWarning)
4 Likes
Can someone provide an example of what this would look like using tt.stack? I tried tt.concatenate(…) like this answer suggests but I’m getting a “TypeError: Join cannot handle arguments of dimension 0”. Since the OP replied using “tt.stack worked”, I’m wondering if @awareeye ran into the same problem?