Hi, I totally agree with you on this point:
It is very common in Bayesian modeling to need to do further post-processing on the posterior samples in order to get meaningful and interpretable results, in many cases, these operations involve statistical or linear algebra.
Just want to add another popular case (in my opinion) is to do Decision Analysis with posterior distributions. I am no expert on this topic, but recently working on a similar problem in a real application.
The BDA3 book discussed Decision Analysis in chapter 9, which is referred in Stan’s document. I copy it here.
So basically, after the model has been built and we got the posterior inference. We may need to make some decisions from the inference result by doing some optimization on the posteriors with different options in a decision space.
Chapter 9 in BDA3 shows several examples how to do this, but there seems no available code examples yet (I searched for it before, but cannot find it).
One great example on Decision analysis using PyMC is in this blog post “Using Bayesian Decision Making to Optimize Supply Chains” by Thomas Wiecki & Ravin Kumar.
But I feel that it is kinds of difficult to work with posterior distributions to perform optimization. For example, in the above example (in cell 14), the authors only get some parts of posterior distribution (with [v:, 1]) for later optimization of the objective function.
supplier_yield_post_pred = pd.DataFrame({k: v[:, 1] for k, v in post_pred.items()})
I do not say this is wrong, but I wish that there is an easy way to work with posterior distributions to perform optimization in Decision Analysis.
So I think about a feature that could allow users to easily perform optimization with posterior distributions (xarray in this context). The optimization may use the scipy.optimize or use a optimization tool like Optuna.
It may difficult to do, but I think it could be a great feature to have 