I’m trying to extend the hierarchical regression model in this tutorial for a multi-dimensional dataset (dataset with multiple features).
I was trying to change the original model as shown below.
with pm.Model() as hierarchical_model: # Hyperpriors for group nodes mu_a = pm.Normal('mu_a', mu=0., sd=100**2) sigma_a = pm.HalfCauchy('sigma_a', 5) mu_b = pm.Normal('mu_b', mu=0., sd=100**2) sigma_b = pm.HalfCauchy('sigma_b', 5) # Intercept for each county, distributed around group mean mu_a # Above we just set mu and sd to a fixed value while here we # plug in a common group distribution for all a and b (which are # vectors of length n_counties). a = pm.Normal('a', mu=mu_a, sd=sigma_a, shape=n_counties) # Intercept for each county, distributed around group mean mu_a b = pm.Normal('b', mu=mu_b, sd=sigma_b, shape=(n_counties, 2)) # Model error eps = pm.HalfCauchy('eps', 5) radon_est = a[county_idx] + tt.dot(b[county_idx], data[["floor", "wave"]].values) # Data likelihood radon_like = pm.Normal('radon_like', mu=radon_est, sd=eps, observed=data.log_radon) with hierarchical_model: hierarchical_trace = pm.sample(draws=2000, n_init=1000)
However, I get an error when trying to perform the dot product between the random variable b and the predictors (X). Here, instead of single predictor “floor” I have selected the feature “wave” in order to test the model with multiple features.
Can someone please help me to figure out the correct way to implement this model