Hi,
There are a couple of small modifications to your code that should fix the problem, however, I’d recommend taking advantage of ArviZ and xarray as it is shown in this notebook. If you have any trouble with the notebook please ask here below, I hope the model will be similar enough and the explanations clear enough.
for i in range(groups):
alpha = trace_rsp_rt['α'][:, i]
beta = trace_rsp_rt['β'][:, i]
mu = alpha + beta * rt
# there may be broadcasting issues requiring to use rt[None, :]
# xarray would handle broadcasting automatically ass seen in the notebook
ax_rsp_rt[i].plot(rt, mu.mean(), c='k', label= f'rsp = {alpha:.2f} + {beta:.2f} * rt')
az.plot_hpd(rt, mu, credible_interval=0.98, color='k', ax=ax_rsp_rt[i])
ax_rsp_rt[i].legend()
# combining pyplot and object based commands can yield unexpected results
ax.set_xlim(1.2, 1.8)
ax.set_ylim(0.6, 1)
On a somewhat tangential topic, note that μ
in trace_rsp_rt
is the distribution of the mean, not the posterior predictive distribution which is what should be used for predicting and for posterior predictive checks (comparing the predictions of the model with observed data).
EDIT: The notebook is not yet available in https://docs.pymc.io/ but will be updated soon.