I am trying to implement a simple Dirichlet model with partial pooling to estimate changes to preferences between two samples of count data (think vote intentions between two polls as an example). I can estimate the preferences for candidates A,B and C in two polls using the simple models describe below. I try to implement a hierarchical model and I was expecting that it would show the means for a candidate in each poll would shrink towards the mean for the candidate in both polls… However in my case, there is no change to the means. I am assuming I have incorrectly setup the model.
I was following similar ideas as to what is described here: https://stats.stackexchange.com/questions/44144/multinomial-dirichlet-model-with-hyperprior-distribution-on-the-concentration-pa
A reproducible example is:
#data, row = "poll", candidate = column
counts=np.array([[317,120,563],[308,85,607]])
shape=counts.shape
totals=np.sum(counts,axis=1)
#Non-pooled model
with pm.Model() as simple_model:
theta = pm.Dirichlet("theta", a=np.ones(shape),shape=shape)
results = pm.Multinomial("results", n=totals, p=theta, observed=counts)
with simple_model:
simple_model = pm.sample(draws=1000)
#check means:
#pm.summary(simple_model)
#Hierarchical model
with pm.Model() as hierarchical_model:
beta = pm.Dirichlet("beta", a=np.ones(shape[1]),shape=shape[1])
lambd = pm.Exponential('lambd',np.ones(shape[1]),shape=shape[1])
theta = pm.Dirichlet("theta", beta*lambd,shape=shape)
results = pm.Multinomial("results", n=totals, p=theta, observed=counts)
with hierarchical_model:
hierarchical_model_trace = pm.sample(draws=10000)
#check means:
#pm.summary(hierarchical_model_trace)