Thank you Junpeng!
That really helped! I think I’m almost there.
doh forgot to include in my description that g is observed (g=np.array([19.2,20.8,10.0,13.0]) - incorporating this prior info is the whole motivation whoops.
Even though it makes sense for H to be stochastic I’ve kept H as determined for now (cos i’m trying to replicate the example):
with pm.Model() as Maher_model:
theta = pm.MvNormal('theta', mu=mu0, cov=V0*500, shape = 6)
#H0 = pm.Normal('H0', mu = H, sd=5, shape=(4, 6))
Ymu = tt.dot(H, theta)
funcY = pm.MvNormal('funcY', mu=Ymu, cov=sigma, observed=g)
samples3 = pm.fit(random_seed=RANDOM_SEED).sample(10000)
It seems replicate the estimates for theta answers OK:
mu = (7.4, 4.4, 4.9, 6.2, 6.5, 7.0) with fairly low SDs
but there is something I can’t understand. I would have expected that when I scaled up the covariance V0 that then the posterior means would change to not consider mu0 as much and also they would have a large standard deviation??? For some reason it doesn’t seem to be influencing very much at all. Am I missing something?