[proofplan]
We take the [Laplace transform](/page/Laplace%20Transform) of the state equation componentwise and use [integration by parts](/theorems/210) to convert the transform of $x'$ into $sX(s)$, with the boundary term vanishing because $x(0)=0$ and $\lim_{T\to\infty}e^{-sT}x(T)=0$. This gives the algebraic resolvent equation $(sI_n-A)X(s)=BU(s)$. Since $s \in \rho(A)$, the matrix $sI_n-A$ is invertible, so we solve for $X(s)$ and substitute into the transformed output equation.
[/proofplan]
[step:Transform the derivative term using the zero initial state and boundary condition]
For each component index $i \in \{1,\dots,n\}$, let $x_i: [0,\infty) \to \mathbb{C}$ denote the $i$-th component of $x$. Since $x$ is absolutely continuous on compact intervals, each $x_i$ is absolutely continuous on compact intervals and $x_i'$ exists for $\mathcal{L}^1$-a.e. $t \in [0,\infty)$.
Fix $T>0$. [Integration by parts](/theorems/2098) on $[0,T]$ gives
\begin{align*}
\int_0^{\!T}e^{-st}x_i'(t)\,d\mathcal{L}^1(t)=e^{-sT}x_i(T)-x_i(0)+s\int_0^{\!T}e^{-st}x_i(t)\,d\mathcal{L}^1(t).
\end{align*}
Because $x(0)=0$, we have $x_i(0)=0$. Because $\lim_{T\to\infty}e^{-sT}x(T)=0$ in $\mathbb{C}^n$, we have $\lim_{T\to\infty}e^{-sT}x_i(T)=0$. Passing to the limit as $T \to \infty$ therefore yields
\begin{align*}
\int_0^\infty e^{-st}x_i'(t)\,d\mathcal{L}^1(t)=s\int_0^\infty e^{-st}x_i(t)\,d\mathcal{L}^1(t).
\end{align*}
Since this holds for every component $i$, the vector identity is
\begin{align*}
\int_0^\infty e^{-st}x'(t)\,d\mathcal{L}^1(t)=sX(s).
\end{align*}
[guided]
The only analytic point in the proof is the transform of the derivative. We cannot simply replace the Laplace transform of $x'$ by $sX(s)$ without accounting for boundary terms, so we compute them explicitly.
For each component index $i \in \{1,\dots,n\}$, define $x_i: [0,\infty) \to \mathbb{C}$ to be the $i$-th coordinate function of $x$. Since $x$ is absolutely continuous on compact intervals, $x_i$ is absolutely continuous on compact intervals. Fix $T>0$. Applying integration by parts to the scalar functions $t \mapsto e^{-st}$ and $t \mapsto x_i(t)$ on the compact interval $[0,T]$ gives
\begin{align*}
\int_0^{\!T}e^{-st}x_i'(t)\,d\mathcal{L}^1(t)=e^{-sT}x_i(T)-x_i(0)+s\int_0^{\!T}e^{-st}x_i(t)\,d\mathcal{L}^1(t).
\end{align*}
Now we inspect the two boundary terms. The initial condition $x(0)=0$ implies $x_i(0)=0$. The assumed boundary decay
\begin{align*}
\lim_{T \to \infty} e^{-sT}x(T)=0
\end{align*}
in $\mathbb{C}^n$ implies, component by component, that
\begin{align*}
\lim_{T \to \infty} e^{-sT}x_i(T)=0.
\end{align*}
Finally, the definition of $X(s)$ gives
\begin{align*}
X(s)=\int_0^\infty e^{-st}x(t)\,d\mathcal{L}^1(t),
\end{align*}
so the $i$-th component of $X(s)$ is
\begin{align*}
\int_0^\infty e^{-st}x_i(t)\,d\mathcal{L}^1(t).
\end{align*}
Passing to the limit $T \to \infty$ in the integration-by-parts identity gives
\begin{align*}
\int_0^\infty e^{-st}x_i'(t)\,d\mathcal{L}^1(t)=s\int_0^\infty e^{-st}x_i(t)\,d\mathcal{L}^1(t).
\end{align*}
Because the same calculation holds for all $n$ components, it is exactly the vector identity
\begin{align*}
\int_0^\infty e^{-st}x'(t)\,d\mathcal{L}^1(t)=sX(s).
\end{align*}
[/guided]
[/step]
[step:Convert the state equation into a resolvent equation]
Since $x'(t)=Ax(t)+Bu(t)$ for $\mathcal{L}^1$-a.e. $t \in [0,\infty)$, multiplying by $e^{-st}$ and integrating with respect to $\mathcal{L}^1$ gives
\begin{align*}
\int_0^\infty e^{-st}x'(t)\,d\mathcal{L}^1(t)=\int_0^\infty e^{-st}Ax(t)\,d\mathcal{L}^1(t)+\int_0^\infty e^{-st}Bu(t)\,d\mathcal{L}^1(t).
\end{align*}
The matrices $A$ and $B$ are constant, so linearity of the vector-valued integral gives
\begin{align*}
sX(s)=AX(s)+BU(s).
\end{align*}
Rearranging the terms in $\mathbb{C}^n$ yields
\begin{align*}
(sI_n-A)X(s)=BU(s).
\end{align*}
[/step]
[step:Solve the resolvent equation for the state transform]
By the hypothesis $s \in \rho(A)$, the matrix $sI_n-A$ is invertible. Multiplying the identity
\begin{align*}
(sI_n-A)X(s)=BU(s)
\end{align*}
on the left by $(sI_n-A)^{-1}$ gives
\begin{align*}
X(s)=(sI_n-A)^{-1}BU(s).
\end{align*}
[/step]
[step:Transform the output equation and substitute the state transform]
The output equation is $y(t)=Cx(t)+Du(t)$ for $\mathcal{L}^1$-a.e. $t \in [0,\infty)$. Multiplying by $e^{-st}$ and integrating with respect to $\mathcal{L}^1$ gives
\begin{align*}
Y(s)=CX(s)+DU(s),
\end{align*}
again by linearity of the vector-valued integral and because $C$ and $D$ are constant matrices. Substituting
\begin{align*}
X(s)=(sI_n-A)^{-1}BU(s)
\end{align*}
therefore gives
\begin{align*}
Y(s)=C(sI_n-A)^{-1}BU(s)+DU(s).
\end{align*}
Factoring the common vector $U(s) \in \mathbb{C}^m$ on the right gives
\begin{align*}
Y(s)=\bigl(C(sI_n-A)^{-1}B+D\bigr)U(s).
\end{align*}
With
\begin{align*}
G(s)=C(sI_n-A)^{-1}B+D,
\end{align*}
this is precisely
\begin{align*}
Y(s)=G(s)U(s).
\end{align*}
[/step]
