Mean-Square Riemann Integrals of a Jointly Gaussian Family are Jointly Gaussian

Let $(\Omega,\mathcal{F},P)$ be a probability space, let $a<b$ be real numbers, let $J$ be a set (possibly empty), and let $(X_j)_{j\in J}$ be a family of random variables on $(\Omega,\mathcal{F},P)$. Let $(H_t)_{t\in[a,b]}$ be a mean-square continuous family of square-integrable random variables, and suppose that the combined family consisting of all $X_j$ ($j\in J$) and all $H_t$ ($t\in[a,b]$), indexed by the disjoint union of $J$ and $[a,b]$, is jointly Gaussian. (When $J$ is empty this hypothesis is that $(H_t)_{t\in[a,b]}$ is jointly Gaussian; the index set remains nonempty since it contains a copy of $[a,b]$.)

Then, for every choice of versions of the mean-square Riemann integrals $\int_a^t H_u\,du$ ($t\in[a,b]$), which exist by Existence and Uniqueness of the Mean-Square Riemann Integral for Mean-Square Continuous Families, the enlarged family consisting of all $X_j$ ($j\in J$), all $H_t$ ($t\in[a,b]$), and all $\int_a^t H_u\,du$ ($t\in[a,b]$), indexed by the disjoint union of $J$ and two copies of $[a,b]$, is jointly Gaussian.