[proofplan]
We define the candidate marginal density by integrating the joint density over the $Y$-variable. Tonelli's theorem gives measurability of this candidate and identifies its total integral as the total mass of the joint density. Then, for an arbitrary Borel set $A\subset\mathbb R^m$, we compute $\mathbb P(X\in A)$ as the probability that $(X,Y)$ lies in $A\times\mathbb R^n$ and use Tonelli's theorem again to rewrite this as an integral of the candidate density over $A$. This proves that the candidate is a density of $X$, and the almost-everywhere uniqueness of Radon-Nikodym derivatives gives the final formulation for any other density $f_X$.
[/proofplan]
[step:Define the marginal candidate and verify its measurability]
Define the function
\begin{align*}
g:\mathbb R^m\to[0,\infty],\qquad g(x)=\int_{\mathbb R^n} f_{X,Y}(x,y)\,d\mathcal L^n(y).
\end{align*}
Since $f_{X,Y}$ is a density, it is $\mathcal B(\mathbb R^m)\otimes\mathcal B(\mathbb R^n)$-measurable and non-negative. By Tonelli's theorem for non-negative [measurable functions](/page/Measurable%20Functions) on the product [measure space](/page/Measure%20Space) $(\mathbb R^m,\mathcal B(\mathbb R^m),\mathcal L^m)\times(\mathbb R^n,\mathcal B(\mathbb R^n),\mathcal L^n)$, the section integral defining $g$ is $\mathcal B(\mathbb R^m)$-measurable.
[/step]
[step:Show that the marginal candidate has total mass one]
Using the standard identification $\mathcal L^{m+n}=\mathcal L^m\otimes\mathcal L^n$ on $\mathbb R^m\times\mathbb R^n$, Tonelli's theorem gives
\begin{align*}
\int_{\mathbb R^m} g(x)\,d\mathcal L^m(x)=\int_{\mathbb R^m}\int_{\mathbb R^n} f_{X,Y}(x,y)\,d\mathcal L^n(y)\,d\mathcal L^m(x).
\end{align*}
Again by Tonelli's theorem, the right-hand side equals the integral of the joint density over the product space:
\begin{align*}
\int_{\mathbb R^m} g(x)\,d\mathcal L^m(x)=\int_{\mathbb R^m\times\mathbb R^n} f_{X,Y}(x,y)\,d\mathcal L^{m+n}(x,y).
\end{align*}
Because $f_{X,Y}$ is a joint density of $(X,Y)$, this last integral is $1$. Hence
\begin{align*}
\int_{\mathbb R^m} g(x)\,d\mathcal L^m(x)=1.
\end{align*}
In particular, $g<\infty$ for $\mathcal L^m$-a.e. $x\in\mathbb R^m$, since a non-negative [measurable function](/page/Measurable%20Function) with finite integral cannot be infinite on a set of positive $\mathcal L^m$-measure.
[/step]
[step:Compute probabilities of events depending only on $X$]
Let $A\in\mathcal B(\mathbb R^m)$ be arbitrary. Define
\begin{align*}
H_A:\mathbb R^m\times\mathbb R^n\to[0,\infty),\qquad H_A(x,y)=\mathbb 1_A(x)f_{X,Y}(x,y).
\end{align*}
The function $H_A$ is non-negative and $\mathcal B(\mathbb R^m)\otimes\mathcal B(\mathbb R^n)$-measurable because $\mathbb 1_A:\mathbb R^m\to\{0,1\}$ is Borel-measurable and $f_{X,Y}$ is measurable. Since $(X,Y)$ has density $f_{X,Y}$,
\begin{align*}
\mathbb P(X\in A)=\mathbb P((X,Y)\in A\times\mathbb R^n)=\int_{A\times\mathbb R^n} f_{X,Y}(x,y)\,d\mathcal L^{m+n}(x,y).
\end{align*}
Using $\mathcal L^{m+n}=\mathcal L^m\otimes\mathcal L^n$ and applying Tonelli's theorem to $H_A$, we obtain
\begin{align*}
\int_{A\times\mathbb R^n} f_{X,Y}(x,y)\,d\mathcal L^{m+n}(x,y)=\int_{\mathbb R^m}\int_{\mathbb R^n}\mathbb 1_A(x)f_{X,Y}(x,y)\,d\mathcal L^n(y)\,d\mathcal L^m(x).
\end{align*}
Since $\mathbb 1_A(x)$ is independent of $y$, the inner integral is
\begin{align*}
\int_{\mathbb R^n}\mathbb 1_A(x)f_{X,Y}(x,y)\,d\mathcal L^n(y)=\mathbb 1_A(x)g(x).
\end{align*}
Therefore
\begin{align*}
\mathbb P(X\in A)=\int_{\mathbb R^m}\mathbb 1_A(x)g(x)\,d\mathcal L^m(x)=\int_A g(x)\,d\mathcal L^m(x).
\end{align*}
[guided]
We want to prove that $g$ is the density of $X$, so we must show that every Borel probability of $X$ is obtained by integrating $g$ over the corresponding Borel set. Fix an arbitrary Borel set $A\in\mathcal B(\mathbb R^m)$. The event $\{X\in A\}$ is the same as the event $\{(X,Y)\in A\times\mathbb R^n\}$, because no restriction is imposed on the $Y$-coordinate. Thus
\begin{align*}
\mathbb P(X\in A)=\mathbb P((X,Y)\in A\times\mathbb R^n).
\end{align*}
Since $(X,Y)$ has joint density $f_{X,Y}$ with respect to $\mathcal L^{m+n}$, the definition of density gives
\begin{align*}
\mathbb P((X,Y)\in A\times\mathbb R^n)=\int_{A\times\mathbb R^n} f_{X,Y}(x,y)\,d\mathcal L^{m+n}(x,y).
\end{align*}
To separate the $x$-variable from the $y$-variable, define
\begin{align*}
H_A:\mathbb R^m\times\mathbb R^n\to[0,\infty),\qquad H_A(x,y)=\mathbb 1_A(x)f_{X,Y}(x,y).
\end{align*}
This function is measurable because $\mathbb 1_A$ is Borel-measurable on $\mathbb R^m$ and $f_{X,Y}$ is product-measurable on $\mathbb R^m\times\mathbb R^n$. It is non-negative because both factors are non-negative. These are precisely the hypotheses needed to apply Tonelli's theorem for non-negative measurable functions. Using the identification $\mathcal L^{m+n}=\mathcal L^m\otimes\mathcal L^n$, Tonelli's theorem gives
\begin{align*}
\int_{A\times\mathbb R^n} f_{X,Y}(x,y)\,d\mathcal L^{m+n}(x,y)=\int_{\mathbb R^m}\int_{\mathbb R^n}\mathbb 1_A(x)f_{X,Y}(x,y)\,d\mathcal L^n(y)\,d\mathcal L^m(x).
\end{align*}
For each fixed $x\in\mathbb R^m$, the quantity $\mathbb 1_A(x)$ is constant as a function of $y$. Therefore it factors out of the inner integral:
\begin{align*}
\int_{\mathbb R^n}\mathbb 1_A(x)f_{X,Y}(x,y)\,d\mathcal L^n(y)=\mathbb 1_A(x)\int_{\mathbb R^n}f_{X,Y}(x,y)\,d\mathcal L^n(y).
\end{align*}
By the definition of $g$, this is
\begin{align*}
\int_{\mathbb R^n}\mathbb 1_A(x)f_{X,Y}(x,y)\,d\mathcal L^n(y)=\mathbb 1_A(x)g(x).
\end{align*}
Substituting this into the preceding display gives
\begin{align*}
\mathbb P(X\in A)=\int_{\mathbb R^m}\mathbb 1_A(x)g(x)\,d\mathcal L^m(x).
\end{align*}
Finally, the defining property of the indicator function converts this integral over all of $\mathbb R^m$ into an integral over $A$:
\begin{align*}
\mathbb P(X\in A)=\int_A g(x)\,d\mathcal L^m(x).
\end{align*}
This is exactly the density identity for the random vector $X$.
[/guided]
[/step]
[step:Conclude that $X$ is continuous and identify all versions of the density]
The preceding step proves that for every $A\in\mathcal B(\mathbb R^m)$,
\begin{align*}
\mathbb P(X\in A)=\int_A g(x)\,d\mathcal L^m(x).
\end{align*}
Thus the law $\mathbb P\circ X^{-1}$ is absolutely continuous with respect to $\mathcal L^m$, and $g$ is a density of $X$. Hence $X$ is a continuous random vector.
If $f_X:\mathbb R^m\to[0,\infty]$ is any other density of $X$ with respect to $\mathcal L^m$, then for every $A\in\mathcal B(\mathbb R^m)$,
\begin{align*}
\int_A f_X(x)\,d\mathcal L^m(x)=\mathbb P(X\in A)=\int_A g(x)\,d\mathcal L^m(x).
\end{align*}
By uniqueness of Radon-Nikodym derivatives up to almost-everywhere equality, $f_X=g$ $\mathcal L^m$-a.e. Therefore
\begin{align*}
f_X(x)=\int_{\mathbb R^n} f_{X,Y}(x,y)\,d\mathcal L^n(y)
\end{align*}
for $\mathcal L^m$-a.e. $x\in\mathbb R^m$. This is the marginal density formula.
[/step]
