Multivariable observables with error covariance


I am trying to construct an inference model where the N observables of multivariate normal with M dimension. Each data point comes with its own Error covariance, cov_err which we know, and an intrinsic covariance, cov_int which we want to learn.

cov = np.array([cov_int + cov_err[i] for i in range(N)])

yObs = pm.MvNormal('yObs', mu=mu, cov=cov, observed=obs)

but MvNormal does not support NxMxM type covariance. I was wondering if there is a fast solution for this problem assuming N ~ 100k, M ~ 10.

Thank you,

You are right MvNormal only support 2d cov matrix - you might find it easier to implement it as a new custom distribution.

Thank you for the quick response.

I tried a few different ways to implement it, they are just supper slow. I will try to speed them up, then. Is there any plan to add this functionality to MvNormal distribution?

You can also try MatrixNormal or KroneckerNormal, but I think neither is optimized to your usecase.