Probability Theory on Product Spaces
Too this point, our discussion of conditional probability theory has been a bit hampered by its abstraction. Fortunately, both the theory and application of conditional probability theory become much more straightforward when applied to the product spaces that dominate practice. In this chapter, we will work through the application of conditional probability theory to product spaces step by step, with a specific emphasis the probability density function representations and practical computation.
1 Product Probability Spaces
Before considering conditional probability theory, let’s begin by discussing the probabilistic objects that we can derive from the structure of an arbitrary product space.
1.1 Product \sigma-Algebras
The construction of product \sigma-algebras from component \sigma-algebras neatly parallels the construction of topologies.
Component \sigma-algebras, \mathcal{X}_{i}, defines measurable component subsets. When each component space is equipped with a component \sigma-algebra, we can construct rectangle subsets from these measurable pieces, \mathsf{x} = \times_{i \in 1:I} \mathsf{x}_{i} with \mathsf{x}_{i} \in \mathcal{X}_{i}, Because these particular rectangle subsets are compatible with all of the component-wise \sigma-algebras at the same time, they are natural candidates for measurable subsets over the full product space.
The only problem with these candidates is they do not form a consistent \sigma-algebra. Specifically, their unions, intersections, and complements define non-rectangular subsets that aren’t neatly related to the component \sigma-algebras. The solution is to just include all of these successors, generating a \sigma-algebra from the rectangular subsets.
Although the resulting \sigma-algebra is well behaved when we have only a finite number of component spaces, it becomes a bit problematic when we consider infinitely many component spaces. Fortunately we can avoid these issues in the same way that we did when building up product topologies: by starting from cylinder subsets instead of more general rectangle subsets.
A product \sigma-algebra, \mathcal{X}^{1:I}, is generated by taking cylinder subsets, \mathsf{x} = \times_{i \in 1:I} \mathsf{x}_{i} with \mathsf{x}_{i} \in \mathcal{X}_{i} and \mathsf{x}_{i} = X_{i} for all but a finite number of i \in 1:I, and then adding their unions, intersections, and complements. Unsurprisingly, this also known as a cylinder \sigma-algebra.
One useful feature of cylinder \sigma-algebras is that they are always compatible with all of the projection functions that we can define over a product space. More formally, given the component indices \mathsf{i} \, \vec{\subset} \, 1:I, the projection function \varpi_{ \mathsf{i} } is always ( \mathcal{X}^{1:I}, \mathcal{X}^{ \mathsf{i} } )-measurable. Consequently, we never rarely have to worry about measurability on product spaces so long as the component spaces are well behaved.
1.2 Independent Product Measures and Probability Distributions
Product \sigma-algebras are not the only probabilistic structure that we can derive from component probabilistic structure. We can also derive entire measures and probability distributions.
Consider each component space being equipped with not only a component \sigma-algebra \mathcal{X}_{i}, but also a component measure \mu_{i} : \mathcal{X}_{i} \rightarrow [0, \infty]. Each of these component measures defines allocations to component measurable subsets. In turn, these component allocations define a notion of allocation to measurable cylinder subsets.
Specifically, given the measurable cylinder subset, \mathsf{x} = \times_{i \in 1:I} \mathsf{x}_{i} with \mathsf{x}_{i} \in \mathcal{X}_{i}, we can define the product allocations \begin{align*} \mu^{1:I} ( \mathsf{x} ) &= \mu^{1:I} \left( \times_{i \in 1:I} \mathsf{x}_{i} \right) \\ &\equiv \prod_{i = 1}^{I} \mu_{i} ( \mathsf{i} ). \end{align*}
To define the allocations for arbitrary measurable product subsets, we can apply Carathéodory’s extension theorem as we did when constructing the Lebesgue measure in Chapter 4, Section 5.2. When each component measure is \sigma-finite, these derived allocations define a unique product measure \mu^{1:I} : \mathcal{X}^{1:I} \rightarrow [0, \infty]. Because each component measure acts independently of the others, this derived measure is often referred to as an independent product measure.
If each component measure is not only \sigma-finite but also a probability distribution, then the derived measure will satisfy all of the Kolmogorov axioms. In other words, component probability distributions define an independent product probability distribution over the corresponding product space. Given how cumbersome this name is, however, shorthands like independent product distribution or even just independent distribution are common.
One nice feature of this construction is that it does not lose any information about the component spaces. In particular, we can always recover any of the component measures by pushing an independent product measure forward along the corresponding projection function, \mu_{i} = ( \varpi_{i} )_{*} \mu.
We have already encountered an independent product measure in Chapter 4, Section 5.3. Multivariate real spaces \mathbb{R}^{I} are built up from I component real lines. Fixing the component metric on each of these real lines defines \sigma-finite component Lebesgue measures, and the independent product of these component Lebesgue measures formally defines the multivariate Lebesgue measure.
1.3 The Fubini-Tonelli Theorem
One of most useful properties of independent product measures is that their integrals can be computed entirely from component-wise integrals.
Before going into more detail, however, let’s briefly discuss notation. The integral of a sufficiently well-behaved function g : X^{1:I} \rightarrow \mathbb{R}. with respect to an independent product measure is most compactly denoted \mathbb{I}_{ \mu^{1:I} } [ g ] = \int \mu^{1:I}( \mathrm{d} x_{1:I} ) \, g( x_{1:I} ) Often, however, it can be more useful to list the components out explicitly using the integral notation \mathbb{I}_{ \mu^{1:I} } [ g ] = \int \mu^{1:I}( \mathrm{d} x_{1}, \ldots, \mathrm{d} x_{I} ) \, g( x_{1}, \ldots, x_{I} ). This is especially useful when we denote the component spaces not with indices, but rather different variables names, such as \mathbb{I}_{ \mu^{1:3} } [ g ] = \int \mu^{1:3}( \mathrm{d} x, \mathrm{d} y, \mathrm{d} z ) \, g( x, y, z ).
With notation settled, let’s now consider two measurable component spaces X_{1} and X_{2} that are equipped with the \sigma-finite measures \mu_{1} and \mu_{2}, respectively. The Fubini-Tonelli theorem (Folland 1999) shows that the integral of an integrable, measurable function g : X_{1} \times X_{2} \rightarrow \mathbb{R} with respect to the independent product measure, \mathbb{I}_{ \mu^{1,2} }[ g ] = \int \mu^{1,2}( \mathrm{d} x_{1}, \mathrm{d} x_{2} ) \, g( x_{1}, x_{2} ), can be computed by nesting or iterating component integrals in either order, \begin{align*} \mathbb{I}_{ \mu^{1,2} }[ g ] &= \int \mu_{1} ( \mathrm{d} x_{1} ) \left[ \int \mu_{2} ( \mathrm{d} x_{2} ) \, g( x_{1}, x_{2} ) \right] \\ &= \int \mu_{2} ( \mathrm{d} x_{2} ) \left[ \int \mu_{1} ( \mathrm{d} x_{1} ) \, g( x_{1}, x_{2} ) \right]. \end{align*}
By repeatedly applying this result, we can decompose integrals with respect to independent product measures on arbitrary product spaces into any sequence of component-wise integrals. More formally, for any permutation of the component indices, ( s(1), \ldots, s(I) ) we can compute the independent product integral with the nested component integrals \mathbb{I}_{ \mu^{1:I} }[ g ] = \int \mu_{s(1)} ( \mathrm{d} x_{s(1)} ) \, \ldots \int \mu_{s(I)} ( \mathrm{d} x_{s(I)} ) \, g( x_{1}, \ldots, x_{I} ).
To make this a bit more concrete, consider the product space \mathbb{Z}^{I} built up from I integer spaces, each of which is naturally equipped with a counting measure \chi_{i}. The independent product of these component counting measures, \chi^{1:I} defines a product measure over \mathbb{Z}^{I}.
Because each \chi_{i} is \sigma-finite, integration with respect to \chi^{1:I} can be implemented by nested component integrals. Moreover, in this case the component integrals can be evaluated directly by summation.
Consequently, \chi^{1:I} integrals can be evaluated by any sequence of nested summations, such as \begin{align*} \int \chi^{1:3}( \mathrm{d} x_{1}, \mathrm{d} x_{2}, \mathrm{d} x_{3} ) \, &g( x_{1}, x_{2}, x_{3} ) \\ &= \sum_{ x_{1} } \sum_{ x_{2} } \sum_{ x_{3} } g( x_{1}, x_{2}, x_{3} ) \\ &= \sum_{ x_{1} } \left[ \sum_{ x_{2} } \left[ \sum_{ x_{3} } g( x_{1}, x_{2}, x_{3} ) \right] \right]. \end{align*}
Another common example is a multivariate Lebesgue measure over \mathbb{R}^{I}. Because the component Lebesgue measures \lambda_{i} are all \sigma-finite, we can implement multivariate Lebesgue integrals by nesting one-dimensional Lebesgue integrals. When the integrand is sufficiently well-behaved, we can then evaluate these component Lebesgue integrals as Riemann integrals.
For instance, the integral of f : \mathbb{R}^{3} \rightarrow \mathbb{R} can be calculated with any sequence of nested Riemann integrals, such as \begin{align*} \int \lambda^{1:3}( \mathrm{d} &x_{1}, \mathrm{d} x_{2}, \mathrm{d} x_{3} ) \, g( x_{1}, x_{2}, x_{3} ) \\ &= \int \lambda_{2} ( \mathrm{d} x_{2} ) \left[ \int \lambda_{1} ( \mathrm{d} x_{1} ) \left[ \int \lambda_{3} ( \mathrm{d} x_{3} ) \, g( x_{1}, x_{2}, x_{3} ) \right] \right] \\ &= \int \mathrm{d} x_{2} \, \int \mathrm{d} x_{1} \, \int \mathrm{d} x_{3} \, g( x_{1}, x_{2}, x_{3} ) \end{align*} or \begin{align*} \int \lambda^{1:3}( \mathrm{d} &x_{1}, \mathrm{d} x_{2}, \mathrm{d} x_{3} ) \, g( x_{1}, x_{2}, x_{3} ) \\ &= \int \lambda_{3} ( \mathrm{d} x_{3} ) \left[ \int \lambda_{2} ( \mathrm{d} x_{2} ) \left[ \int \lambda_{1} ( \mathrm{d} x_{1} ) \, g( x_{1}, x_{2}, x_{3} ) \right] \right] \\ &= \int \mathrm{d} x_{3} \, \int \mathrm{d} x_{2} \, \int \mathrm{d} x_{1} \, g( x_{1}, x_{2}, x_{3} ). \end{align*}
The Fubini-Tonelli theorem does not require that the component spaces are homogeneous. Consider, for example, the product of a rigid real line and an integer space, X = \mathbb{R} \times \mathbb{Z}. The two component spaces are naturally equipped with a Lebesgue measure and counting measure, respectively. The independent product of these two component measures defines a kind of mixed measure over X.
Integrals of sufficiently well-behaved functions with respect to this mixed measure can implemented by iterating summation and Riemann integration, \begin{align*} \mathbb{I}_{ \mu^{1,2} }[ g ] &= \int \mu_{1} ( \mathrm{d} x_{1} ) \left[ \int \mu_{2} ( \mathrm{d} x_{2} ) \, g( x_{1}, x_{2} ) \right] \\ &= \int \lambda ( \mathrm{d} x_{1} ) \left[ \int \chi ( \mathrm{d} x_{2} ) \, g( x_{1}, x_{2} ) \right] \\ &= \int \mathrm{d} x_{1} \, \left[ \sum_{ x_{2} } g( x_{1}, x_{2} ) \right] \end{align*} or, equivalently, \begin{align*} \mathbb{I}_{ \mu^{1,2} }[ g ] &= \int \mu_{2} ( \mathrm{d} x_{2} ) \left[ \int \mu_{1} ( \mathrm{d} x_{1} ) \, g( x_{1}, x_{2} ) \right] \\ &= \int \chi ( \mathrm{d} x_{1} ) \left[ \int \lambda ( \mathrm{d} x_{2} ) \, g( x_{1}, x_{2} ) \right] \\ &= \sum_{ x_{1} } \left[ \int \mathrm{d} x_{2} \, g( x_{1}, x_{2} ) \right]. \end{align*}
1.4 General Product Measures and Probability Distributions
Once we have defined a cylinder \sigma-algebra, we are not limited to only independent product measures. We can always define more general product measures axiomatically. Specifically, any countably-addictive function \nu : \mathcal{X}^{1:I} \rightarrow [0, \infty] defines a measure over the product space X^{1:I}, while while any countably-additive function \pi : \mathcal{X}^{1:I} \rightarrow [0, 1]. defines a probability distribution.
Once we introduce general product measures, however, we have to be weary of the potential for terminological confusion. Here, I have used “product measure” to refer to any measure over a product space. In particular, “product” refers to the structure of the space and not the structure of the measure. If one were to assume the latter, however, then “product measure” would more naturally refer to the specific independent product measures that we derived in the previous section.
One way to avoid this potential confusion is to be a bit more descriptive, referring to arbitrary product measures as joint product measures, or simply joint measures. The adjective “joint” here implies that the component spaces cannot necessarily be treated separately; it also nicely complements the use of “independent” for those exceptional measures that can be decomposed into separate component contributions.
Equivalently, I will refer to arbitrary product probability distributions as joint product probability distributions, joint probability distributions, or even just joint distributions for short.
Let’s pause briefly to review. Given a finite collection of measurable spaces, \{ ( X_{1}, \mathcal{X}_{1} ), ..., ( X_{I}, \mathcal{X}_{I} ) \}, we can immediately construct a product space X^{1:I} and a cylinder \sigma-algebra, \mathcal{X}^{1:I}. In other words, measurable component spaces define a measurable product space.
Once we’ve constructed the measurable product space ( X^{1:I}, \mathcal{X}^{1:I} ), we can then define arbitrary joint measures and probability distributions. At the same time, a finite collection of measure spaces \{ ( X_{1}, \mathcal{X}_{1}, \mu_{i} ), ..., ( X_{I}, \mathcal{X}_{I}, \mu_{i} ) \}, defines not only a product space and cylinder \sigma-algebra, but also an independent product measure. Component measure spaces also define a product measure space.
In addition to deriving independent product measures from explicit component measures, we can also define them axiomatically as special joint measures. From this more abstract perspective, an independent product measure is any measure that satisfies \mu( \mathsf{x} ) = \sum_{i \in 1:I} (\varpi_{i})_{*} \mu( \mathsf{x}_{i} ) for all cylinder subsets \mathsf{x} = \times_{i \in 1:I} \mathsf{x}_{i} with measurable components, \mathsf{x}_{i} \in \mathcal{X}_{i}.
1.5 Component Means
If a component space can be embedded into a rigid real line, then the corresponding projection function defines a measurable function that can be integrated. More formally, given a projection function, \varpi_{i} : X^{1:I} \rightarrow X_{i}, and an embedding function, \iota : X_{i} \rightarrow \mathbb{R}, we can construct the measurable, real-valued function \iota \circ \varpi_{i} : X^{1:I} \rightarrow \mathbb{R}. A joint probability distribution, \pi : \mathcal{X}^{1:I} \rightarrow [0, 1], then defines a corresponding expectation value \mathbb{E}_{ \pi } \! \left[ \iota \circ \varpi_{i} \right].
Using the transformation properties of expectation values that we discussed in Chapter 7, we can also write \mathbb{E}_{ \pi } \! \left[ \iota \circ \varpi_{i} \right] as a component-wise expectation value, \mathbb{E}_{ \pi } \! \left[ \iota \circ \varpi_{i} \right] = \mathbb{E}_{ \pi } \! \left[ ( \varpi_{i} )^{*} \iota \right] = \mathbb{E}_{ (\varpi_{i})_{*} \pi } \! \left[ \iota \right]. From this perspective, the expectation value is just the mean of the pushforward distribution (\varpi_{i})_{*} \pi!
When taking the embedding for granted, \mathbb{E}_{ \pi } \! \left[ \varpi_{i} \right] \equiv \mathbb{E}_{ \pi } \! \left[ \iota \circ \varpi_{i} \right], we often refer to the expectation of a projection function as a component mean. From this component mean, we can then build up more general component moments and cumulants accordingly.
The convention for denoting a component mean varies wildly across different communities. For example, it’s not uncommon to encounter a notation like \mathbb{E}_{\pi}[ x_{i} ], where the projection function is implied by the component variable.
1.6 Expanding Product \sigma-Algebras
Let’s end this section on a technical note that can be safely skipped by those not interested in more rigorous details.
Cylinder \sigma-algebras inherit many of the nice features that are shared by all of the component \sigma-algebras. For example, if the component \sigma-algebras are all Hausdorff, then the cylinder \sigma-algebra will also be Hausdorff.
That said, not all useful properties are preserved by this construction. In particular, even if all of the component \sigma-algebras are complete, the product \sigma-algebra might not be. To facilitate the construction of complete product measures, we have to expand the initial cylinder \sigma-algebra to include cylinder subsets where only some of the components are measurable.
This raises the question of how we should define the action of an independent product measure on a cylinder subset of mixed measurability. Typically we just set the allocations to zero, so that the new subsets are all null subsets. In this case, the new subsets have have no impact on practical applications. The difference between a nominal cylinder \sigma-algebra and an expanded one can then be safely ignored.
The distinction, however, is important in more theoretical calculations that leverage properties like completeness. Consequently, readers moving on to more technical treatments of probability theory should prepare themselves for this subtlety.
2 Conditional Probability Theory on Product Spaces
At this point, we have discussed the explicit construction of independent product measures, but only the implicit existence of more general, joint product measures. How do we construct explicit joint product measures in practice? With our good friend, conditional probability theory.
Applying conditional probability theory to the projection functions on a product space allows us to build up joint measures and joint distributions from more manageable pieces.
2.1 Decomposing Joint Probability Distributions
Recall that any subsequence of component indices, \mathsf{i} \, \vec{\subset} \, 1:I, defines a projection function that maps the full product space into a thinned product space, \varpi_{\mathsf{i}} : X^{1:I} \rightarrow X^{\mathsf{i}}. If each component space X_{i} is equipped with a component \sigma-algebra, then we can construct a cylinder \sigma-algebra over both X^{1:I} and X^{\mathsf{i}}.
If the component \sigma-algebras are all Hausdorff, then these cylinder \sigma-algebra will also be Hausdorff. This allows us to also define a subspace \sigma-algebra over the cross sections \varpi_{\mathsf{i}}^{-1}( x_{\mathsf{i}} ) = X^{\mathsf{i}^{c}} \times \{ x_{\mathsf{i}} \}. In fact, the subspace \sigma-algebra over any cross section is isomorphic to the complementary cylinder \sigma-algebra \mathcal{X}^{\mathsf{i}^{c}}. Each cross section X^{\mathsf{i}^{c}} \times \{ x_{\mathsf{i}} \} can still be treated as a copy of \mathcal{X}^{\mathsf{i}^{c}} after we have introduced measurable structure!
Because the projection function \varpi_{\mathsf{i}} is measurable with respect to all of these \sigma-algebras, it disintegrates any joint probability distribution (Figure 1) into a marginal probability distribution over X^{\mathsf{i}}, ( \varpi_{ \mathsf{i} } )_{*} \pi ( \mathsf{x}_{\mathsf{i}} ), and a conditional probability kernel, \pi^{ \varpi_{\mathsf{i}} }( \mathsf{x} \mid x_{ \mathsf{i} } ). The conditional probability kernel consists of conditional probability distributions that concentrate on each cross section. If we interpret the conditional probability distributions as not just concentrating but actually being defined on these cross sections, then we can also write the kernel in terms of cross section variables, \pi^{ \varpi_{\mathsf{i}} }( (\mathsf{x}_{\mathsf{i}^{c}} )_{x_{\mathsf{i}}} \mid x_{ \mathsf{i} } ).
In general, the behavior of these conditional probability distributions will be heterogeneous. Although each cross section X^{\mathsf{i}^{c}} \times \{ x_{\mathsf{i}} \} can be identified with \mathcal{X}^{\mathsf{i}^{c}}, the probabilistic behavior will can with x_{\mathsf{i}}. Once we consider conditional probability theory, we have to be careful to distinguish between each cross section.
The formal notation that we have used to this point takes nothing for granted. It is also, however, sufficiently dense that it can be difficult to parse. If we’re conditioning with only projection functions, however, then some of this notation is redundant and amenable to streamlining.
Specifically, the variable to the right of the conditioning bar completely determines both the relevant projection function and the corresponding cross section. In this case, we can safely overload our notation a bit, writing ( \varpi_{ \mathsf{i} } )_{*} \pi ( \mathsf{x}_{\mathsf{i}} ) = \pi( \mathsf{x}_{\mathsf{i}} ) and \pi^{\varpi_{\mathsf{i}}}( (\mathsf{x}_{\mathsf{i}^{c}} )_{x_{\mathsf{i}}} \mid x_{ \mathsf{i} } ) = \pi( \mathsf{x}_{\mathsf{i}^{c}} \mid x_{ \mathsf{i} } ). Under the assumption that we are conditioning with only projection functions, the arguments fully differentiate each probabilistic object without any risk of ambiguity. While these particular disintegrations dominate practical applications, they are not completely universal. Consequently, we do have to be on the lookout for the occasional exception.
When working with only a few component spaces, this streamlined notation becomes even easier to parse in integral notation. For instance, given four component spaces ( X_{1}, X_{2}, X_{3}, X_{4} ), the conditional expectations defined by the disintegration of a joint measure with respect to the projection function \varpi_{(2, 4)} : X^{(1, 2, 3, 4)} \rightarrow X^{(2, 4)} can be written as \int \pi( \mathrm{d} x_{1}, \mathrm{d} x_{3} \mid x_{2}, x_{4} ) \, g( x_{1}, x_{3} ). Similarly, the law of total expectation can be written as \begin{align*} \int \pi( \mathrm{d} x_{1}, &\mathrm{d} x_{2}, \mathrm{d} x_{3}, \mathrm{d} x_{4} ) \, g( x_{1}, x_{2}, x_{3}, x_{4} ) \\ &= \int \pi( \mathrm{d} x_{2}, \mathrm{d} x_{4} ) \, \int \pi( \mathrm{d} x_{1}, \mathrm{d} x_{3} \mid x_{2}, x_{4} ) \, g( x_{1}, x_{2}, x_{3}, x_{4} ), \end{align*} or even just \pi( \mathrm{d} x_{1}, \mathrm{d} x_{2}, \mathrm{d} x_{3}, \mathrm{d} x_{4} ) = \pi( \mathrm{d} x_{2}, \mathrm{d} x_{4} ) \, \pi( \mathrm{d} x_{1}, \mathrm{d} x_{3} \mid x_{2}, x_{4} ).
2.2 Constructing Joint Probability Distributions
Conditional probability theory can also be used to build up joint distributions rather than tearing them down.
Given the component indices \mathsf{i} \, \vec{\subset} \, 1:I, we can first engineer an appropriate marginal distribution \pi( \mathsf{x}_{ \mathsf{i} } ). over the thinned product space X^{ \mathsf{i} }. Next, we can construct a conditional probability kernel, \pi( \mathsf{x}_{ \mathsf{i}^{c} } \mid x_{ \mathsf{i} } ) over the cross sections X^{ \mathsf{i}^{c} } \times \{ x_{ \mathsf{i} } \}. The conditional probability kernel “lifts” the marginal distribution into a joint distribution over the full product space.
Because both X^{ \mathsf{i} } and X^{ \mathsf{i}^{c} } \times \{ x_{ \mathsf{i} } \} contain fewer degrees of freedom than X^{1:I}, engineering probabilistic objects over these spaces will often be more manageable than trying to work with the full product space all at once. Moreover, when X^{ \mathsf{i} } is still too complex we can always decompose it further. This will be the subject of the next chapter.
3 Product Probability Density Functions
In order to realize all of this probability theory in a practical application, we need to be able to implement the probabilistic calculations. To do that, we’ll need to consider probability density functions. Fortunately, the component structure of a product space provides all of the ingredients we need to construct useful joint and conditional probability density functions.
3.1 Product Reference Measures
Theoretically, we can always define a reference measure directly over a given product space. As with most product structures, however, it will typically be more useful to build up a product reference measure from component pieces. If each component space X_{i} is equipped with a \sigma-finite reference measure \nu_{i} that admits practical integration, then their independent product \nu^{1:I} defines a natural reference measure over the corresponding product space.
Provided that each of the component reference measures is \sigma-finite, the Fubini-Tonelli theorem will allow us to evaluate integrals with respect to \nu^{1:I} by iterating the integration operations defined by each \nu_{i}. As we discussed in Section Blah, for instance, integrals with respect to a product counting measure can be evaluated by nesting discrete sums, while integrals with respect to a product of Lebesgue measures can often be evaluated by nesting one-dimensional Riemann integrals.
Another powerful feature of independent product reference measures is their disintegrations with respect to projection functions are particularly well-behaved. For example, the component indices, \mathsf{i} \, \vec{\subset} \, 1:I, allows us to construct not only the reference measure \nu^{1:I} over X^{1:I}, but also the reference measure \nu^{ \mathsf{i} } over X^{ \mathsf{i} }. Together, the thinned reference measure \nu^{ \mathsf{i} } and the projection function \varpi_{\mathsf{i}} : X^{1:I} \rightarrow X^{\mathsf{i}}. disintegrate the full reference measure \nu^{1:I} into a collection of conditional measures \nu^{\varpi_{\mathsf{i}}}( \mathsf{x}_{\mathsf{i}^{c}} \mid x_{ \mathsf{i} } ) defined over each of the cross sections \varpi_{\mathsf{i}}^{-1}( x_{\mathsf{i}} ) = X^{\mathsf{i}^{c}} \times \{ x_{\mathsf{i}} \}.
These conditional measures behave like independent product measures over the complementary components, \nu^{\varpi_{\mathsf{i}}}( \mathsf{x}_{\mathsf{i}^{c}} \mid \tilde{x}_{ \mathsf{i} } ) \cong \nu^{ \mathsf{i}^{c} } ( \mathsf{x}_{\mathsf{i}^{c}} ) for all \tilde{x}_{ \mathsf{i} } \in X^{ \mathsf{i} }. This means that we can also use the Fubini-Tonelli theorem to evaluate conditional integrals as nested component integrals.
If we can implement integration on each component space, then independent product reference measures allow us to integrate on all of the full product space, any thinned product space, and any cross section.
3.2 Joint, Marginal, And Conditional Density Functions
Once we have an independent product reference measure, the construction of probability density functions, and consequently the evaluation of expectation values, is immediate.
For any probability distribution \pi defined over the full product space, we can construct the joint probability density function, or simply joint density function, using the independent product reference measure, p( x_{1}, \ldots, x_{I} ) = \frac{ \mathrm{d} \pi }{ \mathrm{d} \nu^{1:I} } ( x_{1}, \ldots, x_{I} ). Given this joint density function, we can then compute joint expectation values by nesting component integrals, for example \begin{align*} \mathbb{E}_{ \pi }[ g ] &= \int \nu_{1}( \mathrm{d} x_{1} ) \bigg[ \ldots \\ & \quad\quad \int \nu_{i}( \mathrm{d} x_{i} ) \bigg[ \ldots \\ & \quad\quad\quad\quad \int \nu_{I}( \mathrm{d} x_{I} ) \, p( x_{1}, \ldots, x_{I} ) \, g( x_{1}, \ldots, x_{I} ) \\ & \hspace{18mm} \ldots \bigg] \\ & \hspace{15mm} \ldots \bigg]. \end{align*} Provided that the joint expectation value is well-defined, any ordering of the component integrals will yield the same answer.
Similarly, for any \mathsf{i} \, \vec{\subset} \, 1:I we can construct a marginal probability density function, or simply marginal density function, p( x_{ \mathsf{i} } ) = \frac{ \mathrm{d} ( \varpi_{\mathsf{i}} )_{*} \pi }{ \mathrm{d} \nu^{ \mathsf{i} } } ( x_{ \mathsf{i} } ), On the cross sections \varpi_{\mathsf{i}}^{-1}( x_{\mathsf{i}} ) we can define conditional probability density functions, or simply conditional density functions, p( x_{\mathsf{i}^{c}} \mid x_{ \mathsf{i} } ) = \frac{ \mathrm{d} \pi^{\varpi_{\mathsf{i}}} }{ \mathrm{d} \nu^{\varpi_{\mathsf{i}}} } ( x_{\mathsf{i}^{c}} \mid x_{ \mathsf{i} } ). The computation of marginal and conditional expectation values once again follows from the Fubini-Tonelli theorem and the component integration operations.
This is all a bit easier to follow when we consider explicit components. To that end, let’s combine four rigid real lines into the product space X^{1:I} = \mathbb{R}^{4}. Because each component space is equipped with a specific metric, we can construct a \sigma-finite Lebesgue reference measure \lambda_{i} on each of them. The independent product of these component reference measures defines the joint reference measure \lambda^{1:4}.
This allows us to construct joint density functions, p( x_{1}, x_{2}, x_{3}, x_{4} ) = \frac{ \mathrm{d} \pi }{ \mathrm{d} \nu^{1:I} } ( x_{1}, x_{2}, x_{3}, x_{4} ), and evaluate joint expectation values with nested Riemann integrals. For instance we might take \mathbb{E}_{ \pi }[ g ] = \int \mathrm{d} x_{1} \, \int \mathrm{d} x_{2} \, \int \mathrm{d} x_{3} \, \int \mathrm{d} x_{4} \, p( x_{1}, x_{2}, x_{3}, x_{4} ) \, g( x_{1}, x_{2}, x_{3}, x_{4} ), but any ordering of the component integrals is equally valid.
Given the component indices \mathsf{i} = ( 2, 4 ), with \mathsf{i}^{c} = ( 1, 3 ), we can disintegrate this joint density function into a marginal density function, p( x_{2}, x_{4} ), and conditional density functions, p( x_{1}, x_{3} \mid x_{2}, x_{4} ). Marginal and conditional expectations are, once again, given by nesting the appropriate component integrals.
3.3 Relating Probability Density Functions
In practice, we rarely if ever construct joint, marginal, and conditional density functions separately. Instead, it’s almost always more convenient derive them from each other.
3.3.1 Decomposition
As we saw in Chapter 9, for example, a marginal density function can be derived from the a joint density function by integrating over the cross sections, p( x_{ \mathsf{i} } ) \overset{ \nu^{\mathsf{i}} }{=} \int \nu^{\varpi_{\mathsf{i}}} ( x_{\mathsf{i}^{c}} ) p( x_{\mathsf{i}^{c}}, x_{ \mathsf{i} } ). This is no longer an abstract result, however, as the conditional integrals can now be evaluated with the Fubini-Tonelli theorem.
At the same time, the product rule that we derived in Chapter 9, Section 5.2 relates the joint, marginal, and conditional probability density functions together, p( x_{\mathsf{i}^{c}}, x_{ \mathsf{i} } ) \overset{ \nu^{1:I} }{=} p( x_{\mathsf{i}^{c}} \mid x_{ \mathsf{i} } ) \, p( x_{ \mathsf{i} } ). This allows us to derive any one of these objects from the other two.
Given a joint density function, we can derive a marginal density function by integrating over the cross sections and then construct conditional density functions from point-wise division, p( x_{\mathsf{i}^{c}} \mid x_{ \mathsf{i} } ) \overset{ \nu^{1:I} }{=} \frac{ p( x_{\mathsf{i}^{c}}, x_{ \mathsf{i} } ) }{ p( x_{ \mathsf{i} } ). }
To demonstrate, let’s consider the product space X^{1:2} = \mathbb{R}^{2} equipped with a multivariate Lebesgue reference measure. This allows us to define a joint distribution through the joint density function (Figure 2 (a)) p( x_{1}, x_{2} ) = \frac{ 1 }{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \left[ -\frac{1}{2} Q( x_{1}, x_{2} ) \right] where Q ( x_{1}, x_{2} ) = \frac{ \sigma_{2}^{2} \, x_{1}^{2} - 2 \, \rho \, \sigma_{1} \, \sigma_{2} \, x_{1} \, x_{2} + \sigma_{1}^{2} \, x_{2}^{2} }{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ). }
The marginal density function over X_{1} is given by integrating the joint density function over the level sets of \varpi_{1}, which are just lines where x_{1} is fixed (Figure 2 (b)), p ( x_{1} ) = \int \mathrm{d} x_{2} \, p( x_{1}, x_{2} ). After some calculus, which I’ve sequestered to Appendix A.1, this becomes a normal density function (Figure 2 (c)) p ( x_{1} ) = \text{normal} ( x_{1} \mid 0, \sigma_{1} ).
Similarly, the conditional density functions become normal density functions of their own (Figure 2 (d)). \begin{align*} p ( x_{2} \mid x_{1} ) &= \frac{ p ( x_{1}, x_{2} ) }{ p ( x_{1} ) } \\ \text{normal} \left( x_{2} \biggm| \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1}, \sigma_{2} \sqrt{ 1 - \rho^{2} } \right). \end{align*}
Because of the symmetry of the joint density function, the same calculations apply when deriving the marginal density over X_{2} and the corresponding conditional density functions. The result there just swaps all of the component indices (Figure 3) \begin{align*} p ( x_{2} ) &= \text{normal} ( x_{2} \mid 0, \sigma_{2} ) \\ p ( x_{1} \mid x_{2} ) &= \text{normal} \left( x_{1} \biggm| \rho \frac{ \sigma_{1} }{ \sigma_{2} } x_{2}, \sigma_{1} \sqrt{ 1 - \rho^{2} } \right). \end{align*}
3.3.2 Composition
These relationships can also be used to construct joint density functions from more manageable pieces. Any marginal density function p( x_{ \mathsf{i} } ) over the thinned product space X^{ \mathsf{i} } and collection of conditional density functions p( x_{\mathsf{i}^{c}} \mid x_{ \mathsf{i} } ) over the cross sections X^{\mathsf{i}^{c}} \times \{ x_{ \mathsf{i} } \} together define a joint density function over the full product space by point-wise multiplication, p( x_{\mathsf{i}^{c}}, x_{ \mathsf{i} } ) \overset{ \nu^{1:I} }{=} p( x_{\mathsf{i}^{c}} \mid x_{ \mathsf{i} } ) \, p( x_{ \mathsf{i} } ).
Consider, for example, the mixed product space X^{1:2} = \mathbb{R} \times \mathbb{Z}. Taking p( x_{1} ) \propto x_{1}^{2} \, (1 - x_{1})^{2} and p( x_{2} \mid x_{1} ) \propto x_{1}^{ x_{2} } \, ( 1 - x_{1} )^{10 - x_{2} } gives the joint density function (Figure 4) \begin{align*} p( x_{1}, x_{2} ) &= p( x_{2} \mid x_{1} ) \, p( x_{1} ) \\ &\propto x_{1}^{ 2 + x_{2} } \, ( 1 - x_{1} )^{ 12 - x_{2} }. \end{align*}
3.4 The Robustness of Logarithmic Density Functions
On a computer, evaluating a joint probability density function defined by a product of density functions can be surprisingly awkward. In particular, the multiplication of density function evaluations, each of which can result in a number close to zero, is prone to numerical underflow.
One convenient way to avoid these numerical issues, is to not work with probability density functions but rather work with their composition with the natural logarithm function, p( x_{1}, \ldots, x_{I} ) \mapsto \log \circ \, p( x_{1}, \ldots, x_{I} ). The logarithm of the joint density function is then given by adding together the logarithm of the conditional density functions and marginal density function, \log \circ \, p( x_{1:I} ) \overset{ \nu^{1:I} }{=} \log \circ \, p( x_{ \mathsf{i}^{c} } \mid x_{ \mathsf{i} } ) + \log \circ \, p( \mathsf{i} ). Because addition is a much more numerically stable operation than multiplication, \log \circ \, p( x_{1:I} ) is much easier to evaluate accurately on a computer.
Once we’ve computed this logarithmic density function, we can always recover the joint density function whenever needed by applying the exponential function, p( x_{1}, \ldots, x_{I} ) = \exp \left( \log \circ \, p( x_{1}, \ldots, x_{I} ) \right).
3.5 Conditioning and Partial Evaluation
In theory, the product rule allows us to immediately condition a joint probability density function on any of the cross sections defined by the component indices \mathsf{i}, p( x_{\mathsf{i}^{c}} \mid x_{ \mathsf{i} } ) \overset{ \nu^{1:I} }{=} \frac{ p( x_{\mathsf{i}^{c}}, x_{ \mathsf{i} } ) }{ p( x_{ \mathsf{i} } ). } In practice, however, the application of this equation can be limited by our ability to calculate the marginal density function in the denominator, p( x_{ \mathsf{i} } ) \overset{ \nu^{\mathsf{i}} }{=} \int \nu^{\varpi_{\mathsf{i}}} ( x_{\mathsf{i}^{c}} ) p( x_{\mathsf{i}^{c}}, x_{ \mathsf{i} } ), for each x_{ \mathsf{i} } \in X^{ \mathsf{i} }.
That said, when we bind x_{ \mathsf{i} } to a particular value \tilde{x}_{ \mathsf{i} }, and fix out attention to a single cross section, the denominator reduces to a single number, p( \tilde{x}_{ \mathsf{i} } ) \in \mathbb{R}^{+}. In this case, the particular conditional density function is equal to the partial evaluation of the joint density function up to proportionality, p( x_{\mathsf{i}^{c}} \mid \tilde{x}_{ \mathsf{i} } ) \overset{ \nu^{1:I} }{=} \frac{ p( x_{\mathsf{i}^{c}}, \tilde{x}_{ \mathsf{i} } ) }{ p( \tilde{x}_{ \mathsf{i} } ). } \propto p( x_{\mathsf{i}^{c}}, \tilde{x}_{ \mathsf{i} } ).
Consequently, we don’t need to be able to evaluate the marginal density function in order to derive unnormalized conditional density functions. This ends up being a surprisingly powerful shortcut in many practical applications.
Consider, for instance, a joint density function over the binary product space X_{1} \times X_{2} defined by p( x_{1}, x_{2} ) \overset{ \nu^{1:2} }{=} p( x_{2} \mid x_{1} ) \, p( x_{1} ). Because we specify p( x_{2} \mid x_{1} ) in the definition of the joint density function, any application that needs these conditional density functions is immediately satisfied.
The same, however, cannot be said for any application that needs an opposing conditional density function, p( x_{1} \mid \tilde{x}_{2} ). This requires being able to calculate \begin{align*} p( x_{1} \mid \tilde{x}_{2} ) &= \frac{ p( x_{1}, \tilde{x}_{2} ) }{ p( \tilde{x}_{2} ) } \\ &= \frac{ p( x_{1}, \tilde{x}_{2} ) }{ \int \nu_{1}( \mathrm{d} x_{1}) p( x_{1}, \tilde{x}_{2} ). } \end{align*} If we cannot evaluate the normalizing integral, then we cannot construct p( x_{1} \mid \tilde{x}_{2} ). That said, if the normalization is irrelevant, then we can skip the normalizing integral entirely, p( x_{1} \mid \tilde{x}_{2} ) \overset{ \nu^{1:I} }{\propto} p( x_{1}, \tilde{x}_{2} ).
One circumstance where the marginal density function can’t be neglected is when we want to compare two conditional distributions to each other. The Radon-Nikodym derivative between two conditional distributions is given by the ratio of their conditional density functions. Using the product rule, we can write this ratio as \begin{align*} \frac{ p( x_{1} \mid \tilde{x}_{2} ) }{ p( x_{1} \mid \tilde{x}'_{2} ) } &\overset{ \nu^{1:2} }{=} \frac{ p( x_{1}, \tilde{x}_{2} ) }{ p( \tilde{x}_{2} ) } \frac{ p( \tilde{x}'_{2} ) }{ p( x_{1}, \tilde{x}'_{2} ) } \\ &\overset{ \nu^{1:2} }{=} \frac{ p( x_{1}, \tilde{x}_{2} ) }{ p( x_{1}, \tilde{x}'_{2} ) } \frac{ p( \tilde{x}'_{2} ) }{ p( \tilde{x}_{2} ). } \end{align*} Here the ratio of the marginal density function evaluations is critical to an accurate comparison.
4 Bayes’ Theorem
Given a subsequence of component indices, we can decompose a joint distribution two different ways. The conditional probability kernel and marginal distribution from each of these decompositions are intimately related.
4.1 General Form
The component indices \mathsf{i} \, \vec{\subset} \, 1:I define the conditional probability kernel \pi^{\varpi_{\mathsf{i}}}( \mathsf{x}_{\mathsf{i}^{c}} \mid x_{ \mathsf{i} } ) and the marginal probability distribution ( \varpi_{\mathsf{i}} )_{*} \pi.
At the same time, the complementary indices \mathsf{i}^{c} \, \vec{\subset} \, 1:I define the conditional probability kernel \pi^{\varpi_{ \mathsf{i}^{c} }}( \mathsf{x}_{\mathsf{i}} \mid x_{ \mathsf{i}^{c} } ) and the marginal probability distribution ( \varpi_{ \mathsf{i}^{c} } )_{*} \pi.
Each conditional probability distribution in the conditional probability kernel of the first decomposition is defined on a cross section, X^{\mathsf{i}^{c}} \times \{ x_{\mathsf{i}} \}. Each of these cross sections can be identified with the thinned product space X^{\mathsf{i}^{c}}. This thinned product space, however, is exactly the space over which the marginal probability distribution in the second decomposition is defined. Consequently, we can define Radon-Nikodym derivatives between these objects, \frac{ \mathrm{d} \pi^{ \varpi_{\mathsf{i}} } }{ \mathrm{d} ( \varpi_{ \mathsf{i}^{c} } )_{*} \pi } ( x_{ \mathsf{i}^{c} } \mid x_{ \mathsf{i} } ).
Similarly, we can construct Radon-Nikodym derivatives between each conditional probability distribution in the second decomposition and the marginal probability distribution in the first decomposition, \frac{ \mathrm{d} \pi^{ \varpi_{ \mathsf{i}^{c} } } }{ \mathrm{d} ( \varpi_{ \mathsf{i} } )_{*} \pi } ( x_{ \mathsf{i} } \mid x_{ \mathsf{i}^{c} } ).
In order for the two decompositions to be consistent, and define the same joint probability distribution, these Radon-Nikodym derivatives have to be almost always equal, \frac{ \mathrm{d} \pi^{ \varpi_{\mathsf{i}} } }{ \mathrm{d} ( \varpi_{ \mathsf{i}^{c} } )_{*} \pi } ( x_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}} ) \overset{ \pi }{=} \frac{ \mathrm{d} \pi^{ \varpi_{ \mathsf{i}^{c} } } }{ \mathrm{d} ( \varpi_{ \mathsf{i} } )_{*} \pi } ( x_{ \mathsf{i} } \mid x_{ \mathsf{i}^{c} } ).
Although not immediately obvious, the identification between these Radon-Nikodym derivatives results in a particularly interpretable consequence. Consider the conditional expectations \int \pi^{ \varpi_{ \mathsf{i}^{c} } } ( \mathrm{d} x_{ \mathsf{i} } \mid \tilde{x}_{ \mathsf{i}^{c} } ) \, g ( x_{ \mathsf{i} } ). Using Radon-Nikodym derivatives, we can compute these conditional expectations as a complementary marginal expectation, \int \pi^{ \varpi_{ \mathsf{i}^{c} } } ( \mathrm{d} x_{ \mathsf{i} } \mid \tilde{x}_{ \mathsf{i}^{c} } ) \, g ( x_{ \mathsf{i} } ) = \int ( \varpi_{ \mathsf{i} } )_{*} \pi ( \mathrm{d} x_{ \mathsf{i} } ) \, \frac{ \mathrm{d} \pi^{ \varpi_{ \mathsf{i}^{c} } } }{ \mathrm{d} ( \varpi_{ \mathsf{i} } )_{*} \pi } ( x_{ \mathsf{i} } \mid \tilde{x}_{ \mathsf{i}^{c} } ) \, g ( x_{ \mathsf{i} } ). Using our consistency identity, we can almost always swap the Radon-Nikodym derivatives to give \int \pi^{ \varpi_{ \mathsf{i}^{c} } } ( \mathrm{d} x_{ \mathsf{i} } \mid \tilde{x}_{ \mathsf{i}^{c} } ) \, g ( x_{ \mathsf{i} } ) \overset{ \pi }{=} \int ( \varpi_{ \mathsf{i} } )_{*} \pi ( \mathrm{d} x_{ \mathsf{i} } ) \, \frac{ \mathrm{d} \pi^{ \varpi_{\mathsf{i}} } }{ \mathrm{d} ( \varpi_{ \mathsf{i}^{c} } )_{*} \pi } ( \tilde{x}_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}} ) \, g ( x_{ \mathsf{i} } ).
In other words, binding x_{ \mathsf{i}^{c} } to the particular value \tilde{x}_{ \mathsf{i}^{c} } “updates” the marginal expectation \int ( \varpi_{ \mathsf{i} } )_{*} \pi ( \mathrm{d} x_{ \mathsf{i} } ) \, g ( x_{ \mathsf{i} } ) into the conditional expectation \int ( \varpi_{ \mathsf{i} } )_{*} \pi ( \mathrm{d} x_{ \mathsf{i} } ) \, \frac{ \mathrm{d} \pi^{ \varpi_{\mathsf{i}} } }{ \mathrm{d} ( \varpi_{ \mathsf{i}^{c} } )_{*} \pi } ( \tilde{x}_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}} ) \, g ( x_{ \mathsf{i} } ) with the insertion of \frac{ \mathrm{d} \pi^{ \varpi_{\mathsf{i}} } }{ \mathrm{d} ( \varpi_{ \mathsf{i}^{c} } )_{*} \pi } ( \tilde{x}_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}}. ) This conditional expectation is also known as a posterior expectation.
Equivalently, we can write this integral relationship out in short hand, \pi^{ \varpi_{ \mathsf{i}^{c} } } ( \mathrm{d} x_{ \mathsf{i} } \mid \tilde{x}_{ \mathsf{i}^{c} } ) \overset{ \pi }{=} \frac{ \mathrm{d} \pi^{ \varpi_{\mathsf{i}} } }{ \mathrm{d} ( \varpi_{ \mathsf{i}^{c} } )_{*} \pi } ( \tilde{x}_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}} ) \, ( \varpi_{ \mathsf{i} } )_{*} \pi ( \mathrm{d} x_{ \mathsf{i} } ).
This relationship between marginal and conditional expectation values is known as Bayes’ Rule or Bayes’ Theorem. The Reverend Thomas Bayes first derived a special case of this relationship which was then posthumously published in 1763 (Bayes 1763).
4.2 Probability Density Function Form
This general form of Bayes’ Theorem is, admittedly, a bit abstract. Fortunately, it becomes much more straightforward in terms of probability density functions.
Assuming consistent reference measures, the component indices \mathsf{i} induce the probability density function decomposition p( x_{\mathsf{i}}, x_{ \mathsf{i}^{c} } ) \overset{ \nu^{1:I} }{=} p( x_{\mathsf{i}} \mid x_{ \mathsf{i}^{c} } ) \, p( x_{ \mathsf{i}^{c} } ), while the complementary indices \mathsf{i}^{c} give p( x_{\mathsf{i}}, x_{ \mathsf{i}^{c} } ) \overset{ \nu^{1:I} }{=} p( x_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}}) \, p( x_{ \mathsf{i} } ).
Bayes’ Theorem relates all of these density functions together in a single equation for a posterior density function, \begin{align*} \int \pi^{ \varpi_{ \mathsf{i}^{c} } } ( \mathrm{d} x_{ \mathsf{i} } \mid \tilde{x}_{ \mathsf{i}^{c} } ) \, g ( x_{ \mathsf{i} } ) &\overset{ \pi }{=} \int ( \varpi_{ \mathsf{i} } )_{*} \pi ( \mathrm{d} x_{ \mathsf{i} } ) \, \frac{ \mathrm{d} \pi^{ \varpi_{\mathsf{i}} } }{ \mathrm{d} ( \varpi_{ \mathsf{i}^{c} } )_{*} \pi } ( \tilde{x}_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}} ) \, g ( x_{ \mathsf{i} } ) \\ \int \nu^{ \mathsf{i} } ( \mathrm{d} x_{ \mathsf{i} } ) \, p ( x_{ \mathsf{i} } \mid \tilde{x}_{ \mathsf{i}^{c} } ) \, g ( x_{ \mathsf{i} } ) &\overset{ \pi }{=} \int \nu^{ \mathsf{i} } ( \mathrm{d} x_{ \mathsf{i} } ) \, p ( x_{ \mathsf{i} } ) \, \frac{ \mathrm{d} \pi^{ \varpi_{\mathsf{i}} } }{ \mathrm{d} ( \varpi_{ \mathsf{i}^{c} } )_{*} \pi } ( \tilde{x}_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}} ) \, g ( x_{ \mathsf{i} } ) \\ \int \nu^{ \mathsf{i} } ( \mathrm{d} x_{ \mathsf{i} } ) \, p ( x_{ \mathsf{i} } \mid \tilde{x}_{ \mathsf{i}^{c} } ) \, g ( x_{ \mathsf{i} } ) &\overset{ \pi }{=} \int \nu^{ \mathsf{i} } ( \mathrm{d} x_{ \mathsf{i} } ) \, p ( x_{ \mathsf{i} } ) \, \frac{ p ( \tilde{x}_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}} ) }{ p ( \tilde{x}_{ \mathsf{i}^{c} } ) }\, g ( x_{ \mathsf{i} } ), \end{align*} or p ( x_{ \mathsf{i} } \mid \tilde{x}_{ \mathsf{i}^{c} } ) \overset{ \nu^{1:I} }{=} p ( x_{ \mathsf{i} } ) \, \frac{ p ( \tilde{x}_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}} ) }{ p ( \tilde{x}_{ \mathsf{i}^{c} } ) }.
Conveniently, we can immediately derive this form of Bayes’ Theorem directly from the consistency of the two conditional decompositions, \begin{align*} p( x_{\mathsf{i}}, x_{ \mathsf{i}^{c} } ) &\overset{ \nu^{1:I} }{=} p( x_{\mathsf{i}}, x_{ \mathsf{i}^{c} } ) \\ p( x_{\mathsf{i}} \mid x_{ \mathsf{i}^{c} } ) \, p( x_{ \mathsf{i}^{c} } ) &\overset{ \nu^{1:I} }{=} p( x_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}}) \, p( x_{ \mathsf{i} } ) \\ p( x_{\mathsf{i}} \mid x_{ \mathsf{i}^{c} } ) &\overset{ \nu^{1:I} }{=} \frac{ p( x_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}}) }{ p( x_{ \mathsf{i}^{c} } ) } \, p( x_{ \mathsf{i} } ). \end{align*}
Bayes’ Theorem becomes somewhat trivial to implement when we need only an unnormalized posterior density function. As we saw in Section 3.5, an unnormalized conditional density function is given by just partially evaluating the joint density function, \begin{align*} p( x_{\mathsf{i}}, \tilde{x}_{ \mathsf{i}^{c} } ) &\overset{ \nu^{1:I} }{=} \frac{ p( \tilde{x}_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}}) }{ p( \tilde{x}_{ \mathsf{i}^{c} } ) } \, p( x_{ \mathsf{i} } ) \\ &\overset{ \nu^{1:I} }{ \propto } p( \tilde{x}_{ \mathsf{i}^{c} } \mid x_{\mathsf{i}}) \, p( x_{ \mathsf{i} } ) \\ &\overset{ \nu^{1:I} }{ \propto } p( \tilde{x}_{ \mathsf{i}^{c} } , x_{\mathsf{i}} ). \end{align*} Indeed, this is how Bayes’ Theorem is most often implemented in practice.
5 Convolution
All of this product space probability theory can also be applied to probabilistic systems that don’t necessarily have any composite structure of their own.
Consider, for example, the probability space (X, \mathcal{X}, \pi). Even if the ambient space X does not have any composite structure, we can always introduce product structure by combining it with itself. Specifically, we can construct the power space X^{2} = X \times X' which is naturally equipped with the projection functions \begin{alignat*}{6} \varpi :\; & X^{2} & &\rightarrow& \; & X & \\ & (x, x') & &\mapsto& & x &, \end{alignat*} and \begin{alignat*}{6} \varpi' :\; & X^{2} & &\rightarrow& \; & X & \\ & (x, x') & &\mapsto& & x' &. \end{alignat*}
Once we have defined X^{2} and its projection functions, we can introduce a conditional probability kernel over the level sets of \varpi, \pi( \mathrm{d} x' \mid x ). Any choice of conditional probability kernel lifts the initial probability distribution \pi to a joint probability distribution over the power space, \pi( \mathrm{d} x' \mid x ) \, \pi( \mathrm{d} x ) = \pi( \mathrm{d} x, \mathrm{d} x').
By construction, pushing this joint distribution forward along \varpi returns the original probability distribution, (\varpi)_{*} \pi ( \mathrm{d} x ) = \pi( \mathrm{d} x ). Pushing it forward along the second projection function \varpi', however, defines an entirely new probability distribution over X, (\varpi')_{*} \pi ( \mathrm{d} x' ) = \int \pi( \mathrm{d} x' \mid x ) \, \pi( \mathrm{d} x ) \ne \pi( \mathrm{d} x' ).
This process of lifting on one projection function and then pushing forward along the other, transforming an initial probability distribution into a new one, is known as convolution. In this context, the conditional probability kernel \pi( \mathrm{d} x' \mid x ) is referred to as a convolution kernel.
Consider, for example, the integers X = \mathbb{I} equipped with a probability distribution defined by the Poisson density function p( x ) = \text{Poisson}( x \mid \lambda ) = \frac{ \lambda^{x} e^{-\lambda} }{ x! }. and the somewhat awkward convolution kernel defined by the conditional density functions p( x' \mid x ) = \frac{ \eta^{x' - x} e^{-\eta} }{ (x' - x)! } I[ x \le x' ].
The convolution of this initial probability distribution with the convolution kernel is a new probability distribution. After some algebra, which is worked out in Appendix A.2, the new probability distribution is defined by another Poisson density function, \begin{align*} p ( x' ) &= \sum_{x = 0}^{ \infty } p( x' \mid x ) \, p (x) \\ &= \text{Poisson}( x' \mid \lambda + \eta ). \end{align*} In this case, the convolution translates the intensity of the Poisson parameter by \eta.
For a second demonstration, consider a rigid real line X = \mathbb{R} equipped with a probability distribution defined by the normal density function \begin{align*} p( x ) &= \text{normal} ( x \mid \mu, \sigma ) \\ &= \frac{1}{ \sqrt{ 2 \, \pi \, \sigma^{2} } } \exp \left[ -\frac{1}{2} \frac{ ( x - \mu )^{2} }{ \sigma^{2} } \right]. \end{align*}
Given a convolution kernel defined by the conditional density functions \begin{align*} p( x' \mid x ) &= \text{normal} ( x' \mid \eta - x, \tau ) \\ &= \frac{1}{ \sqrt{ 2 \, \pi \, \tau^{2} } } \exp \left[ -\frac{1}{2} \frac{ ( x' - \eta + x )^{2} }{ \tau^{2} } \right], \end{align*} we can construct the convolution \begin{align*} p( x' ) &= \int \mathrm{d} x \, p( x, x' ) \\ &= \text{normal} \left( x' \mid \mu + \eta, \sqrt{ \sigma^{2} + \tau^{2} } \right). \end{align*} The full derivation can be found in Appendix A.3.
In this case, the convolution operation takes in initial normal density function, translates the location parameter by \eta, and then inflates the square of the scale.
When we take the limit \tau \rightarrow 0, the convolution becomes a pure translation of the location parameter, p( x' ) = \text{normal} \left( x' \mid \mu + \eta, \sqrt{ \sigma^{2} + \tau^{2} } \right). Indeed, this is exactly the pushforward of a normal density function along the translation operation that we discussed in Chapter 7, Section 4.3.2.1.
The convolution kernel in this limit acts like a singular delta function that shifts an input x to an output x + \eta. This limiting convolution is also known as a linear convolution.
6 Conclusion
All of the abstraction of conditional probability theory becomes concrete when we apply it to product spaces. Fortunately, all of the applications that we will consider in this book will be on product spaces. The material in this chapter will be the foundation on which we will build those applications.
That said, we have not yet completed our discussion of decomposing joint probability distributions. Once we have decomposed an initial joint distribution into a conditional probability kernel and marginal distribution, we can iterate and decompose the marginal distribution into simpler pieces. These recursive decompositions will be the subject of the next chapter.
Appendix: “Explicit” Calculations
Per tradition, I have have sequestered the more involved calculations in this chapter to this appendix where they cannot harm anyone. At least not without consent.
A.1 Density Function Disintegration
Consider the binary product space X^{1:2} = \mathbb{R}^{2} equipped with a two-dimensional Lebesgue measure and the joint distribution defined by the joint density function p( x_{1}, x_{2} ) = \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \left[ -\frac{1}{2} Q( x_{1}, x_{2} ) \right] where Q ( x_{1}, x_{2} ) = \frac{ \sigma_{2}^{2} \, x_{1}^{2} - 2 \, \rho \, \sigma_{1} \, \sigma_{2} \, x_{1} \, x_{2} + \sigma_{1}^{2} \, x_{2}^{2} }{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ). }
The marginal density function over X_{1} is given by integrating over x_{2}, \begin{align*} p ( x_{1} ) &= \int \mathrm{d} x_{2} \, p( x_{1}, x_{2} ) \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \int \mathrm{d} x_{2} \, \exp \left[ -\frac{1}{2} Q( x_{1}, x_{2} ) \right]. \end{align*}
In order to evaluate this integral, we need to isolate all of the x_{2} terms in Q ( x_{1}, x_{2} ). This requires completing the square for x_{2} in the numerator, \begin{align*} \sigma_{2}^{2} \, x_{1}^{2} - 2 \, \rho \, \sigma_{1} \, \sigma_{2} \, x_{1} \, x_{2} + \sigma_{1}^{2} \, x_{2}^{2} &= \sigma_{1}^{2} \, \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} + \sigma_{2}^{2} \, x_{1}^{2} - \rho \, \sigma_{2}^{2} \, x_{1}^{2} \\ &= \sigma_{1}^{2} \, \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} + \sigma_{2}^{2} \, x_{1}^{2} \, (1 - \rho^{2} ). \end{align*}
Then \begin{align*} Q ( x_{1}, x_{2} ) &= \frac{ \sigma_{2}^{2} \, x_{1}^{2} - 2 \, \rho \, \sigma_{1} \, \sigma_{2} \, x_{1} \, x_{2} + \sigma_{1}^{2} \, x_{2}^{2} }{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } \\ &= \frac{ \sigma_{1}^{2} \, \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} + \sigma_{2}^{2} \, x_{1}^{2} \, (1 - \rho^{2} ) }{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } \\ &= \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } + \frac{ x_{1}^{2} }{ \sigma_{1}^{2} }. \end{align*}
This allows us to write the joint density function as \begin{align*} p( x_{1}, x_{2} ) &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \left[ -\frac{1}{2} Q( x_{1}, x_{2} ) \right] \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \left[ -\frac{1}{2} \left( \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } + \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \right) \right] \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \Bigg[ -\frac{1}{2} \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } \Bigg] \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg]. \end{align*}
With this factorization, the marginal density function reduces to a normal integral which we can evaluate in closed form, \begin{align*} p ( x_{1} ) &= \int \mathrm{d} x_{2} \, p( x_{1}, x_{2} ) \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \int \mathrm{d} x_{2} \, \exp \Bigg[ -\frac{1}{2} \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } \Bigg] \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg] \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg] \int \mathrm{d} x_{2} \, \exp \Bigg[ -\frac{1}{2} \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } \Bigg] \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg] \sqrt{ 2 \, \pi \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } \\ &= \frac{1}{ \sqrt{2 \, \pi \, \sigma_{1}^{2} } } \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg]. \end{align*}
After all of the that, the marginal density function is just a normal density function, \begin{align*} p ( x_{1} ) &= \frac{1}{ \sqrt{2 \, \pi \, \sigma_{1}^{2} } } \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg] \\ &= \text{normal} ( x_{1} \mid 0, \sigma_{1} ). \end{align*}
We can use this same factorization of the joint density function to simplify the calculation of the conditional density functions, \begin{align*} p ( x_{2} \mid x_{1} ) &= \frac{ p ( x_{1}, x_{2} ) }{ p ( x_{1} ) } \\ &= \frac{ \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \Bigg[ -\frac{1}{2} \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } \Bigg] \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg] }{ \frac{1}{ \sqrt{2 \, \pi \, \sigma_{1}^{2} } } \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg] } \\ &= \frac{1}{ \sqrt{ 2 \, \pi \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \Bigg[ -\frac{1}{2} \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } \Bigg]. \end{align*}
Perhaps surprisingly, each conditional density function is also a normal density function, \begin{align*} p ( x_{2} \mid x_{1} ) &= \frac{1}{ \sqrt{ 2 \, \pi \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \Bigg[ -\frac{1}{2} \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } \Bigg] \\ &= \text{normal} \left( x_{2} \biggm| \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1}, \sigma_{2} \sqrt{ 1 - \rho^{2} } \right). \end{align*}
Because of the symmetry of the joint density function, the same calculations apply when deriving the other marginal density function and corresponding conditional density functions. The result just swaps all of the indices, \begin{align*} p ( x_{2} ) &= \text{normal} ( x_{2} \mid 0, \sigma_{2} ) \\ p ( x_{1} \mid x_{2} ) &= \text{normal} \left( x_{1} \biggm| \rho \frac{ \sigma_{1} }{ \sigma_{2} } x_{2}, \sigma_{1} \sqrt{ 1 - \rho^{2} } \right). \end{align*}
Finally, this is an example where some pattern matching can point to the marginal and conditional density functions directly, without having to evaluate a single integral. Examining p( x_{1}, x_{2} ) = \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \Bigg[ -\frac{1}{2} \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } \Bigg] \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg], we might notice that both of the exponentiated quadratics define normal density functions, \begin{align*} p( x_{1}, x_{2} ) &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma_{1}^{2} \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \Bigg[ -\frac{1}{2} \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } \Bigg] \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg] \\ &= \frac{1}{ \sqrt{ 2 \, \pi \, \sigma_{2}^{2} \, (1 - \rho^{2} ) } } \exp \Bigg[ -\frac{1}{2} \frac{ \left( x_{2} - \rho \frac{ \sigma_{2} }{ \sigma_{1} } x_{1} \right)^{2} }{ \sigma_{2}^{2} \, (1 - \rho^{2} ) } \Bigg] \frac{1}{ \sqrt{ 2 \, \pi \, \sigma_{1}^{2} } } \exp \Bigg[ -\frac{1}{2} \frac{ x_{1}^{2} }{ \sigma_{1}^{2} } \Bigg] \\ &= \text{normal} \left( x_{1} \biggm| \rho \frac{ \sigma_{1} }{ \sigma_{2} } x_{2}, \sigma_{1} \sqrt{ 1 - \rho^{2} } \right) \, \text{normal} ( x_{1} \mid 0, \sigma_{1} ). \end{align*}
In other words, by isolating x_{2} we have incidentally manipulated the joint density function into the desired conditional decomposition p( x_{1}, x_{2} ) = p( x_{2} \mid x_{1} ) \, p( x_{1} ).
This approach isn’t always useful; in particular, it requires being able to identify normalized density functions by eye. That said, it can save time when working with similar calculations.
A.2 Discrete Convolution
We begin with the space X = \mathbb{I} equipped with a counting reference measure and a probability distribution defined by the Poisson density function p( x ) = \text{Poisson}( x \mid \lambda ) = \frac{ \lambda^{x} e^{-\lambda} }{ x! }.
Then we’ll consider the convolution kernel defined by the conditional density functions p( x' \mid x ) = \frac{ \eta^{x' - x} e^{-\eta} }{ (x' - x)! } I[ x \le x' ]. Before deriving the convolution, however, let’s first verify that this somewhat awkward looking object is, in fact, a conditional probability kernel. To do that, we’ll need to verify that each conditional distribution is properly normalized for all values of \eta > 0.
The normalizations are given by explicit summations, \begin{align*} \sum_{x' = 0}^{ \infty } p( x' \mid x ) &= \sum_{x' = 0}^{ \infty } \frac{ \eta^{x' - x} e^{-\eta} }{ (x' - x)! } I[ x \le x' ] \\ &= \sum_{x' = x}^{ \infty } \frac{ \eta^{x' - x} e^{-\eta} }{ (x' - x)! }. \end{align*} Changing the index of the summation to k = x' - x gives \begin{align*} \sum_{x' = 0}^{ \infty } p( x' \mid x ) &= \sum_{k = 0}^{ \infty } \frac{ \eta^{k} e^{-\eta} }{ k! } \\ &= e^{-\eta} \sum_{k = 0}^{ \infty } \frac{ \eta^{k} }{ k! } \\ &= e^{-\eta} e^{\eta} \\ &= 1, \end{align*} as required.
With the validity of the convolution kernel verified, we can proceed to the convolution itself, \begin{align*} p ( x' ) &= \sum_{x = 0}^{ \infty } p( x' \mid x ) \, p (x) \\ &= \sum_{x = 0}^{ \infty } \frac{ \eta^{x' - x} e^{-\eta} }{ (x' - x)! } I[ x \le x' ] \, \frac{ \lambda^{x} e^{-\lambda} }{ x! } \\ &= e^{-(\lambda + \eta)} \sum_{x = 0}^{ \infty } I[ x \le x' ] \, \frac{ 1 }{ x! \, (x' - x)! } \eta^{x' - x} \, \lambda^{x} \\ &= \frac{ e^{-(\lambda + \eta)} }{ x'! } \sum_{x = 0}^{ x' } \frac{ x'! }{ x! \, (x' - x)! } \eta^{x' - x} \, \lambda^{x}. \end{align*}
Conveniently, the summation is immediately given by the binomial theorem, ( a + b )^{K} = \sum_{k = 0}^{K} a^{k} \, b^{K - k}. Consequently, \begin{align*} p ( x' ) &= \frac{ e^{-(\lambda + \eta)} }{ x'! } \sum_{x = 0}^{ x' } \frac{ x'! }{ x! \, (x' - x)! } \eta^{x' - x} \, \lambda^{x} \\ &= \frac{ e^{-(\lambda + \eta)} }{ x'! } (\lambda + \eta)^{x'}. \end{align*} This, however, is just another Poisson density function, p ( x' ) = \frac{ e^{-(\lambda + \eta)} }{ x'! } (\lambda + \eta)^{x'} = \text{Poisson}( x' \mid \lambda + \eta ).
A.3 Continuous Convolution
Consider a rigid real line X = \mathbb{R} equipped with a Lebesgue reference measure and a probability distribution defined by the normal density function \begin{align*} p( x ) &= \text{normal} ( x \mid \mu, \sigma ) \\ &= \frac{1}{ \sqrt{ 2 \, \pi \, \sigma^{2} } } \exp \left[ -\frac{1}{2} \frac{ ( x - \mu )^{2} }{ \sigma^{2} } \right]. \end{align*}
The convolution kernel defined by the conditional density functions \begin{align*} p( x' \mid x ) &= \text{normal} ( x' \mid \eta - x, \tau ) \\ &= \frac{1}{ \sqrt{ 2 \, \pi \, \tau^{2} } } \exp \left[ -\frac{1}{2} \frac{ ( x' - \eta + x )^{2} }{ \tau^{2} } \right]. \end{align*} lifts the initial probability distribution into a joint distribution over X^{2} defined by the joint density function \begin{align*} p( x, x' ) &= p( x' \mid x ) \, p( x ) \\ &= \text{normal} ( x \mid \mu, \sigma ) \, \text{normal} ( x' \mid \eta - x, \tau ) \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma^{2} \, \tau^{2} } } \exp \left[ -\frac{1}{2} Q(x, x') \right], \end{align*} where Q(x, x') = \frac{ ( x - \mu )^{2} }{ \sigma^{2} } + \frac{ ( x' - \eta + x )^{2} }{ \tau^{2} }.
The convolution of the initial probability distribution with the convolution kernel is given by projecting this joint density function to the auxiliary space, \begin{align*} p( x' ) &= \int \mathrm{d} x \, p( x, x' ) \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma^{2} \, \tau^{2} } } \int \mathrm{d} x \, \exp \left[ -\frac{1}{2} Q(x, x') \right]. \end{align*}
To evaluate this integral, we’ll need to isolate x with some – surprise! – completion of squares, \begin{align*} Q(x, x') &= \frac{ ( x - \mu )^{2} }{ \sigma^{2} } + \frac{ ( x' - \eta + x )^{2} }{ \tau^{2} } \\ &= \frac{ x^{2} - 2 \, \mu \, x + \mu^{2} }{ \sigma^{2} } + \frac{ x^{2} - 2 \, (x' - \eta) \, x + (x' - \eta)^{2} }{ \tau^{2} } \\ &= \frac{1}{ \sigma^{2} \, \tau^{2} } \bigg[ ( \sigma^{2} + \tau^{2} ) \, x^{2} - 2 ( \tau^{2} \, \mu + \sigma^{2} \, (x' - \eta) ) \, x + \tau^{2} \, \mu^{2} + \sigma^{2} ( x' - \eta)^{2} \bigg] \\ &= \frac{1}{\sigma^{2} \, \tau^{2} } \bigg[ \quad ( \sigma^{2} + \tau^{2} ) \, \left( x - \frac{ \tau^{2} \, \mu + \sigma^{2} \, (x' - \eta) }{ \sigma^{2} + \tau^{2} } \right)^{2} \\ &\hspace{16mm} + \tau^{2} \, \mu^{2} + \sigma^{2} ( x' - \eta)^{2} - \frac{ \left( \tau^{2} \, \mu + \sigma^{2} \, (x' - \eta) \right)^{2} }{ \sigma^{2} + \tau^{2} } \bigg] \\ &=\quad \frac{1}{ \sigma^{2} \, \tau^{2} } \bigg[ ( \sigma^{2} + \tau^{2} ) \, \left( x - \frac{ \tau^{2} \, \mu + \sigma^{2} \, (x' - \eta) }{ \sigma^{2} + \tau^{2} } \right)^{2} \bigg] \\ &\quad+ \frac{1}{ \sigma^{2} \, \tau^{2} } \bigg[ \quad \frac{ \tau^{4} \, \mu^{2} + \sigma^{2} \, \tau^{2} \, \mu^{2} + \sigma^{4} ( x' - \eta)^{2} + \sigma^{2} \, \tau^{2} \, ( x' - \eta)^{2} }{ \sigma^{2} + \tau^{2} } \\ &\hspace{20mm} - \frac{ \tau^{4} \, \mu^{2} + 2 \tau^{2} \, \sigma^{2} \, \mu \, (x' - \eta) + \sigma^{4} (x' - \eta)^{2} }{ \sigma^{2} + \tau^{2} } \hspace{15mm} \bigg] \\ &=\quad \frac{1}{ \sigma^{2} \, \tau^{2} } \bigg[ ( \sigma^{2} + \tau^{2} ) \, \left( x - \frac{ \tau^{2} \, \mu + \sigma^{2} \, (x' - \eta) }{ \sigma^{2} + \tau^{2} } \right)^{2} \bigg] \\ &\quad+ \frac{1}{ \sigma^{2} \, \tau^{2} } \bigg[ \frac{ \sigma^{2} \, \tau^{2} \, \mu^{2} - 2 \, \tau^{2} \, \sigma^{2} \, \mu \, (x' - \eta) + \sigma^{2} \, \tau^{2} \, ( x' - \eta)^{2} }{ \sigma^{2} + \tau^{2} } \bigg] \\ &=\quad \frac{\sigma^{2} + \tau^{2} }{ \sigma^{2} \, \tau^{2} } \left( x - \frac{ \tau^{2} \, \mu + \sigma^{2} \, (x' - \eta) }{ \sigma^{2} + \tau^{2} } \right)^{2} \\ &\quad+ \frac{ \mu^{2} - 2 \, \mu \, (x' - \eta) + ( x' - \eta)^{2} }{ \sigma^{2} + \tau^{2} } \\ &=\quad \frac{\sigma^{2} + \tau^{2} }{ \sigma^{2} \, \tau^{2} } \left( x - \frac{ \tau^{2} \, \mu + \sigma^{2} \, (x' - \eta) }{ \sigma^{2} + \tau^{2} } \right)^{2} \\ &\quad+ \frac{ ( x ' - (\mu + \eta) )^{2} }{ \sigma^{2} + \tau^{2} }. \end{align*}
This allows us to reduce the convolution operation to a normal integral, \begin{align*} p( x' ) &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma^{2} \, \tau^{2} } } \int \mathrm{d} x \, \exp \left[ -\frac{1}{2} Q(x, x') \right] \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma^{2} \, \tau^{2} } } \int \mathrm{d} x \, \; \exp \left[ -\frac{1}{2} \frac{\sigma^{2} + \tau^{2} }{ \sigma^{2} \, \tau^{2} } \left( x - \frac{ \tau^{2} \, \mu + \sigma^{2} \, (x' - \eta) }{ \sigma^{2} + \tau^{2} } \right)^{2} \right] \\ &\hspace{29.75mm} \cdot \exp \left[ -\frac{1}{2} \frac{ ( x ' - (\mu + \eta) )^{2} }{ \sigma^{2} + \tau^{2} } \right] \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma^{2} \, \tau^{2} } } \, \exp \left[ -\frac{1}{2} \frac{ ( x ' - (\mu + \eta) )^{2} }{ \sigma^{2} + \tau^{2} } \right] \\ &\hspace{13mm} \int \mathrm{d} x \, \exp \left[ -\frac{1}{2} \frac{\sigma^{2} + \tau^{2} }{ \sigma^{2} \, \tau^{2} } \left( x - \frac{ \tau^{2} \, \mu + \sigma^{2} \, (x' - \eta) }{ \sigma^{2} + \tau^{2} } \right)^{2} \right] \\ &= \frac{1}{ 2 \, \pi \, \sqrt{ \sigma^{2} \, \tau^{2} } } \, \exp \left[ -\frac{1}{2} \frac{ ( x ' - (\mu + \eta) )^{2} }{ \sigma^{2} + \tau^{2} } \right] \sqrt{ 2 \, \pi \, \frac{ \sigma^{2} \, \tau^{2} }{\sigma^{2} + \tau^{2} } } \\ &= \frac{1}{ \sqrt{ 2 \, \pi \, ( \sigma^{2} + \tau^{2} ) } } \, \exp \left[ -\frac{1}{2} \frac{ ( x ' - (\mu + \eta) )^{2} }{ \sigma^{2} + \tau^{2} } \right] \\ &= \text{normal} \left( x' \mid \mu + \eta, \sqrt{ \sigma^{2} + \tau^{2} } \right). \end{align*}
Acknowledgements
I thank jd for helpful comments.
A very special thanks to everyone supporting me on Patreon: Alessandro Varacca, Alex D, Alexander Noll, Amit, Andrea Serafino, Andrew Mascioli, Andrew Rouillard, Andrés Castro Araújo, Ara Winter, Ari Holtzman, Austin Rochford, Aviv Keshet, Avraham Adler, Ben Matthews, Ben Swallow, Benoit Essiambre, boot, Brendan Galdo, Bryan Chang, Cameron Smith, Canaan Breiss, Cat Shark, Cathy Oliveri, Charles Naylor, Chase Dwelle, Chris Jones, Christina Van Heer, Christopher Mehrvarzi, Colin Carroll, Colin McAuliffe, Damien Mannion, dan mackinlay, Dan W Joyce, Dan Waxman, Dan Weitzenfeld, Daniel Hammarström, Danny Van Nest, David Burdelski, Dr. Jobo, Dr. Omri Har Shemesh, Dylan Maher, Dylan Spielman, Ebriand, Ed Cashin, Eric LaMotte, Erik Banek, Eugene O’Friel, Felipe González, Felipe Vaca, Fergus Chadwick, Francesco Corona, Geoff Rollins, Glenn Williams, Granville Matheson, Guilherme Marthe, Hamed Bastan-Hagh, haubur, Hector Munoz, Horace Guy, hs, Hugo Botha, Håkan Johansson, Ian Costley, idontgetoutmuch, Ignacio Vera, Ilaria Prosdocimi, iris pisscauldron, Isaac Vock, jacob pine, Jair Andrade, James C, James Hodgson, James Wade, Janek Berger, Jarrett Byrnes, Jason Pekos, Jason Wong, jd, Jeff Burnett, Jeff Dotson, Jeff Helzner, Jeffrey Erlich, Jerry Lin , Jesper Fischer Ehmsen, Jessica Graves, Joe Sloan, John Flournoy, Jonathan H. Morgan, Jonathon Vallejo, Josh Knecht, JU, Julian Lee, Justin Bois, Karim Naguib, Karim Osman, Karsten Skogsholm, Konstantin Shakhbazov, Kristian Gårdhus Wichmann, Kádár András, Lars Barquist, lizzie , LOU ODETTE, Mads Christian Hansen, Marek Kwiatkowski, Mark Donoghoe, Markus P., Daniel Edward Marthaler, Matthieu LEROY, Mattia Arsendi, Matěj, Maurits van der Meer, Max, Michael Colaresi, Michael DeJesus, Michael DeWitt, Michael Dillon, Michael Lerner, Mick Cooney, MisterMentat , Márton Vaitkus, N Sanders, Nathaniel Burbank, Nicholas Cowie, Nick S, Octavio Medina, Ole Rogeberg, Olivier Ma, Paolo, Pat LS, Patrick Kelley, Patrick Boehnke, Pau Pereira Batlle, Pieter van den Berg , ptr, quasar, Ramiro Barrantes Reynolds, Raúl Peralta Lozada, Rex, Riccardo Fusaroli, Richard Nerland, Rob Davies, Robert Frost, Robert Goldman, Robert kohn, Robin Taylor, Ryan Gan, Ryan Grossman, Ryan Kelly, S Hong, Sean Wilson, Sergiy Protsiv, Seth Axen, shira, Simon Duane, Simon Lilburn, Simon Steiger, Simone, Spencer, sssz, Stefan Lorenz, Stephen Lienhard, Steve Forrest, Steve Harris, Steven Forrest, Stew Watts, Stone Chen, Stoyan Georgiev, Susan Holmes, Svilup, Tate Tunstall, Tatsuo Okubo, Teresa Ortiz, Theodore Dasher, Thomas Siegert, Thomas Vladeck, Tobychev , Tomas Capretto, Tony Wuersch, Virgile Andreani, Virginia Fisher, Vitalie Spinu, Vladimir M, VO2 Maximus Decius, Will Farr, Will Lowe, Will Wen, William A Vauter, yolhaj, yureq, and Zach A.
References
License
The text and figures in this chapter are copyrighted by Michael Betancourt and licensed under the CC BY-NC 4.0 license.