Compound distributions allow us to compose probability distributions together. In other words, it allows us to have parameters in our model that are themselves random variables. This facilitates including a hierarchical structure (e.g. like in a BHM) in our models.
A few weeks ago I learned that SymPy supports compound distributions.
Here is an example. Suppose
\[X \sim \text{Normal}(2, 4)\]
\[L \sim \text{Exponential}(100)\]
\[N \sim \text{Normal}(X, L)\]
In SymPy we can represent this model as follows:
"""Represent a compound distribution in SymPy.
A compound distribution has random parameters: the mean and the scale of ``N``
are themselves random variables.
"""
from __future__ import annotations
from sympy.stats import Exponential, Normal
from sympy.stats.compound_rv import CompoundDistribution
from sympy.stats.crv_types import NormalDistribution
def compound_distribution() -> CompoundDistribution:
"""Build the compound distribution used in the post.
Returns the ``CompoundDistribution`` for
.. math::
X \\sim \\mathrm{Normal}(2, 4), \\quad
L \\sim \\mathrm{Exponential}(100), \\quad
N \\sim \\mathrm{Normal}(X, L).
"""
X = Normal("X", 2, 4)
L = Exponential("L", 100)
N = NormalDistribution(X, L)
return CompoundDistribution(N)
def main() -> None:
"""Console-script entry point: build the compound distribution."""
print(compound_distribution())
if __name__ == "__main__":
main()I’m looking forward to seeing more advancements in SymPy for statistics. Mathematical statistics can get particularly difficult when performed by hand. Having specialized classes with handy methods for distributions or random variables could be a huge benefit to bringing in reproducible mathematical derivations of moments (e.g. expectations) and distributions.