# Tensorflow probability - Bayesian hierarchical model **URL:** \ **Category:** General Discussion\ **Tags:** models, help\_request\ **Created:** [August 9, 2021, 5:07pm UTC](https://discuss.ai.google.dev/t/tensorflow-probability-bayesian-hierarchical-model/28194 "2021-08-09T17:07:40Z")\ **Posts on this page:** 6\ **Page:** 1
**Author:** ![Majid\_Mohammadi](https://d1dlmcr85iqnpo.cloudfront.net/user_avatar/discuss.ai.google.dev/majid_mohammadi/32/8625_2.png) [@Majid\_Mohammadi](https://discuss.ai.google.dev/u/Majid_Mohammadi)\ **Post date:** [August 9, 2021, 5:07pm UTC](https://discuss.ai.google.dev/t/tensorflow-probability-bayesian-hierarchical-model/28194/1 "2021-08-09T17:07:40Z")
I have a simple hierarchical model, currently developed in PyMC3 as follows: ``` a_b = np.array([[3, 4, 6, 1, 5, 2, 9, 7], [1, 2, 8, 4, 5, 3, 9, 6], [2, 2, 3, 1, 5, 5, 9, 8], [2, 1, 8, 2, 9, 3, 8, 8], [2, 4, 9, 1, 4, 3, 5, 5], [1, 2, 9, 1, 3, 5, 5, 4]]) a_w = np.array([ [7, 6, 4, 9, 5, 8, 1, 3], [9, 8, 2, 5, 4, 5, 1, 3], [8, 8, 5, 9, 5, 5, 1, 2], [8, 9, 2, 8, 1, 8, 2, 2], [8, 6, 1, 9, 6, 7, 4, 4], [9, 8, 1, 9, 7, 5, 5, 6]]) aw_tc = np.sum(a_w, axis=1) ab_tc = np.sum(a_b, axis=1) dm_no, c_no = a_w.shape with pm.Model() as BayesianBWM: w_star = pm.Dirichlet("w_star", a=0.01*np.ones(c_no)) gamma_agg = pm.Gamma("beta_agg",alpha=0.01, beta=0.01) for i in range(dm_no): globals()[f"w{i}"] = pm.Dirichlet('w'+str(i), a=gamma_agg*w_star) globals()[f"aw{i}"] = pm.Multinomial('a_w'+str(i),n=aw_tc[i], p=globals()[f"w{i}"], observed=a_w[i,:]) globals()[f"ab{i}"] = pm.Multinomial('a_b'+str(i),n=ab_tc[i], p= 1 / globals()[f"w{i}"], observed=a_b[i,:]) ``` As you see, we have two observed matrices, a\_w and a\_b, and for each row of these matrices, we learn a w, and finally based on these w’s, we learn w\_star. I have difficulties implementing this in Tensorflow JointDistributionSequential. I want to develop it purely in TF probability and not using other high level libraries such as Edward2. Any help is highly appreciated! ---
**Author:** ![Bhack](https://d1dlmcr85iqnpo.cloudfront.net/user_avatar/discuss.ai.google.dev/bhack/32/6388_2.png) [@Bhack](https://discuss.ai.google.dev/u/Bhack)\ **Post date:** [August 11, 2021, 3:38pm UTC](https://discuss.ai.google.dev/t/tensorflow-probability-bayesian-hierarchical-model/28194/3 "2021-08-11T15:38:57Z")
/cc @Christopher_Suter ---
**Author:** ![Christopher\_Suter](https://d1dlmcr85iqnpo.cloudfront.net/user_avatar/discuss.ai.google.dev/christopher_suter/32/7551_2.png) [@Christopher\_Suter](https://discuss.ai.google.dev/u/Christopher_Suter)\ **Post date:** [August 11, 2021, 4:29pm UTC](https://discuss.ai.google.dev/t/tensorflow-probability-bayesian-hierarchical-model/28194/4 "2021-08-11T16:29:53Z")
The code below is the spiritual TFP equivalent to your PyMC model. I used JointDistributionCoroutine instead of Sequential, since it’s a bit easier to write for loops therein. There’s a bug in the turnkey sampling routine at present; I’ve raised the issue with the team and will circle back when we have a solution. Meanwhile, this is the rough idea! ```python a_b = np.array([ [3, 4, 6, 1, 5, 2, 9, 7], [1, 2, 8, 4, 5, 3, 9, 6], [2, 2, 3, 1, 5, 5, 9, 8], [2, 1, 8, 2, 9, 3, 8, 8], [2, 4, 9, 1, 4, 3, 5, 5], [1, 2, 9, 1, 3, 5, 5, 4]]).astype(np.float32) a_w = np.array([ [7, 6, 4, 9, 5, 8, 1, 3], [9, 8, 2, 5, 4, 5, 1, 3], [8, 8, 5, 9, 5, 5, 1, 2], [8, 9, 2, 8, 1, 8, 2, 2], [8, 6, 1, 9, 6, 7, 4, 4], [9, 8, 1, 9, 7, 5, 5, 6]]).astype(np.float32) aw_tc = np.sum(a_w, axis=1) ab_tc = np.sum(a_b, axis=1) dm_no, c_no = a_w.shape @tfd.JointDistributionCoroutine def model(): w_star = yield tfd.Dirichlet(concentration=1. * np.ones(c_no, dtype=np.float32), name='w_star') beta_agg = yield tfd.Gamma(concentration=2., rate=2., name='beta_agg') for i in range(dm_no): w = yield tfd.Dirichlet(concentration=beta_agg * w_star, name=f'w{i}') aw = yield tfd.Multinomial(total_count=aw_tc[i], probs=w, name=f'aw{i}') ab = yield tfd.Multinomial(total_count=ab_tc[i], probs = 1 / w, name=f'ab{i}') pin_values = dict(([ (f'aw{i}', a_w[i, :]) for i in range(dm_no) ] + [ (f'ab{i}', a_b[i, :]) for i in range(dm_no) ])) pinned = model.experimental_pin(**pin_values) results = tfp.experimental.mcmc.windowed_adaptive_hmc( n_draws=100, joint_dist=pinned, num_leapfrog_steps=10, n_chains=100, num_adaptation_steps=300) ``` ---
**Author:** ![Christopher\_Suter](https://d1dlmcr85iqnpo.cloudfront.net/user_avatar/discuss.ai.google.dev/christopher_suter/32/7551_2.png) [@Christopher\_Suter](https://discuss.ai.google.dev/u/Christopher_Suter)\ **Post date:** [August 11, 2021, 8:40pm UTC](https://discuss.ai.google.dev/t/tensorflow-probability-bayesian-hierarchical-model/28194/6 "2021-08-11T20:40:49Z")
The fix is in, and should be in tfp-nightly tomorrow: [Fix bug in windowed\_mcmc related to constrained distributions having … · tensorflow/probability@956d09f · GitHub](https://github.com/tensorflow/probability/commit/956d09f0a9c61abc2cb9c375339f935aa60b6275) After that, the code I snipped above should work. Please let us know if you have any more issues! ---
**Author:** ![Majid\_Mohammadi](https://d1dlmcr85iqnpo.cloudfront.net/user_avatar/discuss.ai.google.dev/majid_mohammadi/32/8625_2.png) [@Majid\_Mohammadi](https://discuss.ai.google.dev/u/Majid_Mohammadi)\ **Post date:** [August 12, 2021, 4:59am UTC](https://discuss.ai.google.dev/t/tensorflow-probability-bayesian-hierarchical-model/28194/7 "2021-08-12T04:59:00Z")
Thanks a lot for taking the time, Christopher! I first updated tfp-nightly, executed the code, and got the following error! I do not know if I must stay longer so that the fix in tfp-nightly will be available! Here is the error (only the error and not the traceback): ValueError: Input `event_shape` does not match `event_shape_in` ([7] vs [8]). ---
**Author:** ![Christopher\_Suter](https://d1dlmcr85iqnpo.cloudfront.net/user_avatar/discuss.ai.google.dev/christopher_suter/32/7551_2.png) [@Christopher\_Suter](https://discuss.ai.google.dev/u/Christopher_Suter)\ **Post date:** [August 12, 2021, 8:29pm UTC](https://discuss.ai.google.dev/t/tensorflow-probability-bayesian-hierarchical-model/28194/9 "2021-08-12T20:29:58Z")
HI @Majid_Mohammadi – the error you see is what the commit I mentioned should fix. I think you just tried it a little too soon. The new tfp-nightly package goes out around 08:00 UTC. I think you will also need a couple of amendments to the code I wrote: ```auto w = yield tfd.Dirichlet(concentration=beta_agg[..., tf.newaxis] * w_star, name=f'w{i}') ``` will make sure the shapes align correctly. And you’ll also want to wrap the hmc call with tf.function: ```auto results = tf.function(lambda: tfp.experimental.mcmc.windowed_adaptive_hmc( n_draws=100, joint_dist=pinned, num_leapfrog_steps=10, n_chains=100, num_adaptation_steps=300))() ``` This will ensure the sampler is fully compiled and should reduce run time.