Dimensions in Multivariate Normal Exa

Hi all,

I’m trying to fit a multivariate normal analogous to this example

https://docs.pymc.io/en/v3/pymc-examples/examples/case_studies/LKJ.html

Different from the example, I have several sets of observations for which I want to estimate independent parameters. I can’t seem to figure out how to add a dimension to the example to handle this (and I want to avoid an ugly loop thing). Normally I know how to work with the shapes but the observations being pairs somehow makes it different…

I’ve tried this (among other things)

t, trials  = pd.factorize(baseline_data['Trial Index'], sort=True)
with pm.Model(coords=coords) as model:
     chol, corr, stds = pm.LKJCholeskyCov('chol', n=2, eta=2.,
     sd_dist=pm.HalfCauchy.dist(500,shape=(len(trials), 2)),
                                   compute_corr=True)

     μ = pm.Normal('μ', 0., 2000., 
     testval=baseline_data['X'].mean(axis=0), shape=(len(trials), 2))
     obs = pm.MvNormal('obs', μ[t], chol=chol[t], observed=np.array([baseline_data_static['X'].values, 
     baseline_data_static['Y'].values]).T)
     trace = pm.sample(tune=1000, return_inferencedata=True, 
     init='advi+adapt_diag')

It results in "
IndexError: index 2 is out of bounds for axis 0 with size 2"

I have also tried different orders of the dimensions. And somewhere there is a version of the example with named dimensions but that didn’t help either…

Any suggestions are highly appreciated!

Many thanks
Jan

1 Like