[proofplan]
We split the proof according to whether the variance parameter is zero or positive. In the degenerate case $\sigma^2=0$, the definition of $\mathcal N(\mu,0)$ gives $X=\mu$ almost surely, so the expectation and variance follow directly. In the positive-variance case, we compute from the Gaussian density by translating and scaling to the standard normal density, then evaluate the first two standard normal moments using symmetry and [integration by parts](/theorems/210).
[/proofplan]
[step:Handle the degenerate Gaussian case]
Assume first that $\sigma^2=0$. By the degenerate-normal convention in the definition of $\mathcal N(\mu,0)$, the law of $X$ is the Dirac [probability measure](/page/Probability%20Measure) at $\mu$, equivalently $X=\mu$ $\mathbb P$-almost surely. Hence $X$ is integrable and
\begin{align*}
\mathbb E[X]=\int_\Omega X\,d\mathbb P=\int_\Omega \mu\,d\mathbb P=\mu\,\mathbb P(\Omega)=\mu.
\end{align*}
Since $\mathbb E[X]=\mu$, we also have $X-\mathbb E[X]=0$ $\mathbb P$-almost surely. Therefore
\begin{align*}
\operatorname{Var}(X)=\mathbb E[(X-\mathbb E[X])^2]=\int_\Omega 0\,d\mathbb P=0=\sigma^2.
\end{align*}
[/step]
[step:Reduce the positive-variance case to standard normal integrals]
Assume now that $\sigma^2>0$, and let $\sigma>0$ denote the positive square root of $\sigma^2$. By the density definition of $X\sim\mathcal N(\mu,\sigma^2)$, the law of $X$ has density
\begin{align*}
f_X:\mathbb R\to[0,\infty),\qquad x\mapsto \frac{1}{\sigma\sqrt{2\pi}}\exp\left(-\frac{(x-\mu)^2}{2\sigma^2}\right)
\end{align*}
with respect to one-dimensional [Lebesgue measure](/page/Lebesgue%20Measure) $\mathcal L^1$.
Define the standard normal density
\begin{align*}
\phi:\mathbb R\to[0,\infty),\qquad z\mapsto \frac{1}{\sqrt{2\pi}}\exp\left(-\frac{z^2}{2}\right).
\end{align*}
For every Borel-[measurable function](/page/Measurable%20Function) $g:\mathbb R\to\mathbb R$ for which the corresponding integral is absolutely integrable, the change of variables $x=\mu+\sigma z$ gives $d\mathcal L^1(x)=\sigma\,d\mathcal L^1(z)$ and transforms the domain $\mathbb R$ onto $\mathbb R$. Thus
\begin{align*}
\mathbb E[g(X)]=\int_{\mathbb R} g(x)f_X(x)\,d\mathcal L^1(x)=\int_{\mathbb R} g(\mu+\sigma z)\phi(z)\,d\mathcal L^1(z).
\end{align*}
[guided]
We want to compute moments of $X$, but the density of $X$ contains the location $\mu$ and the scale $\sigma$. The standard way to remove these parameters is to translate and rescale the integration variable. Since $\sigma^2>0$, there is a unique positive number $\sigma>0$ satisfying $\sigma^2$ equal to the variance parameter.
The density of $X$ is the map
\begin{align*}
f_X:\mathbb R\to[0,\infty),\qquad x\mapsto \frac{1}{\sigma\sqrt{2\pi}}\exp\left(-\frac{(x-\mu)^2}{2\sigma^2}\right).
\end{align*}
We also define the standard normal density as the map
\begin{align*}
\phi:\mathbb R\to[0,\infty),\qquad z\mapsto \frac{1}{\sqrt{2\pi}}\exp\left(-\frac{z^2}{2}\right).
\end{align*}
Now take a Borel-measurable function $g:\mathbb R\to\mathbb R$ such that the relevant integral is absolutely integrable. The substitution
\begin{align*}
x=\mu+\sigma z
\end{align*}
is a $C^1$ bijection from $\mathbb R$ to $\mathbb R$, and the one-dimensional change-of-variables formula gives
\begin{align*}
d\mathcal L^1(x)=\sigma\,d\mathcal L^1(z).
\end{align*}
Under this substitution, the exponent transforms as
\begin{align*}
\frac{(x-\mu)^2}{2\sigma^2}=\frac{z^2}{2}.
\end{align*}
Therefore
\begin{align*}
\int_{\mathbb R} g(x)f_X(x)\,d\mathcal L^1(x)=\int_{\mathbb R} g(\mu+\sigma z)\frac{1}{\sigma\sqrt{2\pi}}e^{-z^2/2}\sigma\,d\mathcal L^1(z).
\end{align*}
After cancelling the factor $\sigma$, this becomes
\begin{align*}
\mathbb E[g(X)]=\int_{\mathbb R} g(\mu+\sigma z)\phi(z)\,d\mathcal L^1(z).
\end{align*}
This reduction is the whole reason for introducing $\phi$: all dependence on $\mu$ and $\sigma$ has been moved into the simple affine expression $\mu+\sigma z$.
[/guided]
[/step]
[step:Evaluate the first two standard normal moments]
We record the two standard normal moment identities needed below:
\begin{align*}
\int_{\mathbb R} z\phi(z)\,d\mathcal L^1(z)=0.
\end{align*}
and
\begin{align*}
\int_{\mathbb R} z^2\phi(z)\,d\mathcal L^1(z)=1.
\end{align*}
For the first identity, the function $z\mapsto z\phi(z)$ is odd because $\phi(-z)=\phi(z)$. Its absolute integrability follows from
\begin{align*}
\int_{\mathbb R}|z|\phi(z)\,d\mathcal L^1(z)=2\int_0^\infty z\frac{1}{\sqrt{2\pi}}e^{-z^2/2}\,d\mathcal L^1(z)=\sqrt{\frac{2}{\pi}}<\infty,
\end{align*}
where the last equality uses the substitution $u=z^2/2$, so $d\mathcal L^1(u)=z\,d\mathcal L^1(z)$ on $(0,\infty)$. Hence the integral of the odd integrable function $z\phi(z)$ over $\mathbb R$ is zero.
For the second identity, note that $\phi$ is differentiable and $\phi'(z)=-z\phi(z)$ for every $z\in\mathbb R$. For $R>0$, [integration by parts](/theorems/2098) on $[-R,R]$ gives
\begin{align*}
\int_{-R}^{R} z^2\phi(z)\,d\mathcal L^1(z)=-\int_{-R}^{R}z\phi'(z)\,d\mathcal L^1(z).
\end{align*}
The integration-by-parts formula yields
\begin{align*}
-\int_{-R}^{R}z\phi'(z)\,d\mathcal L^1(z)=-R\phi(R)-R\phi(-R)+\int_{-R}^{R}\phi(z)\,d\mathcal L^1(z).
\end{align*}
Since $\phi(-R)=\phi(R)$ and $R\phi(R)\to0$ as $R\to\infty$, while $\int_{\mathbb R}\phi(z)\,d\mathcal L^1(z)=1$ by normalization of the standard normal density, the [monotone convergence theorem](/theorems/509) applied to the nonnegative functions $z^2\phi(z)\mathbb 1_{[-R,R]}(z)$ gives
\begin{align*}
\int_{\mathbb R} z^2\phi(z)\,d\mathcal L^1(z)=1.
\end{align*}
[/step]
[step:Compute the expectation and variance from the standard moments]
Apply the reduction formula with $g:\mathbb R\to\mathbb R$, $g(x)=x$. The absolute integrability follows from
\begin{align*}
\int_{\mathbb R}|\mu+\sigma z|\phi(z)\,d\mathcal L^1(z)\le |\mu|\int_{\mathbb R}\phi(z)\,d\mathcal L^1(z)+\sigma\int_{\mathbb R}|z|\phi(z)\,d\mathcal L^1(z)<\infty.
\end{align*}
Therefore
\begin{align*}
\mathbb E[X]=\int_{\mathbb R}(\mu+\sigma z)\phi(z)\,d\mathcal L^1(z)=\mu\int_{\mathbb R}\phi(z)\,d\mathcal L^1(z)+\sigma\int_{\mathbb R}z\phi(z)\,d\mathcal L^1(z)=\mu.
\end{align*}
Since $\mathbb E[X]=\mu$, the variance is
\begin{align*}
\operatorname{Var}(X)=\mathbb E[(X-\mu)^2].
\end{align*}
Apply the reduction formula with $g:\mathbb R\to\mathbb R$, $g(x)=(x-\mu)^2$. Its integrability follows from the second standard normal moment computed above. Thus
\begin{align*}
\operatorname{Var}(X)=\int_{\mathbb R}(\sigma z)^2\phi(z)\,d\mathcal L^1(z)=\sigma^2\int_{\mathbb R}z^2\phi(z)\,d\mathcal L^1(z)=\sigma^2.
\end{align*}
This proves both asserted identities in the positive-variance case, and together with the degenerate case completes the proof.
[/step]
