I’ve been having some issues reproducing the output from my model with a normal univariate prior distribution. In the model, I’m looking to introduce correlation between sate vector parameters in the prior, and so have been attempting to introduce MvNormal with a covariance or precision matrix.
I have noticed that using a prior MvNormal with an identity matrix as covariance (all zero off diagonals) and mu=np.ones(nx) is producing different output from a Normal with mu=1, sig=1, and shape=nx.
Statistically these should be equivalent. Using a NUTS sampler, both models are successfully converging on (slightly) different results (see image). I was wondering if there is a difference in the sampling or the RNG between univariate and multivariate distributions in pymc? Any help would be greatly appreciated.
