I am working on a hierarchical multinomial model with partial pooling (over some common hyper priors). Basically something similar to a multinomial with random effect on panel data.

My dataset has 4 predictors (x1,x2,x3,x4) and the y has 3 possible states (0,1,2).

Furthermore the dataset has N total individuals (N=10930) for which I am trying to forecast the different betas, and in a simple logit model I am able to map the right betas to the right person in the “mu” variable using an array that maps every row in my dataset to their “owner”.

Now my approach is something like this:

```
model = pm.Model()
with model:
mu_common_prior = pm.Normal('common_mu', mu=0., sd=100**2, shape=(4, 3))
sd_common_prior = pm.HalfCauchy('common_sd', 5, shape=(4, 3))
beta = pm.Normal('beta', mu=mu_prior, sd=sd_prior, shape=(N, 4, 3))
mu = pm.math.dot(beta[individual_indexes], df[['x1', 'x2', 'x3', 'x4']].values)
p = T.nnet.softmax(mu)
y = pm.Categorical('y', p=p, observed=df['y'].values)
res = pm.fit(20000, callbacks=[CheckParametersConvergence()])
trace = res.sample(1000)
pm.traceplot(trace)
```

but unfortunately the pm.dot call fails miserably complaining about the shapes of my data.

`shapes (43720,3) and (10930,4) not aligned: 3 (dim 1) != 10930 (dim 0)`

any thoughts on what is wrong with my model? thanks