I am trying to establish that difference of two symmetrically doubly-truncated normals has finite variance(or that variance exists). Suppose $X,Y$ are two iid truncated normals with the following density function $$f(x,\sigma,a,b)=\frac{1}{\sqrt{2\pi}\sigma}\frac{e^{\frac{-x^2}{2}}}{\Phi(b)-\Phi(a)}$$ where a is the lower boundary of truncation, b is the upper boundary and $\Phi$ is the CDF of the standard normal and $\phi$ is the associated PDF. (See wikipedia https://en.wikipedia.org/wiki/Truncated_normal_distribution). Additionally,we have $a=-b$. Thus we have symmetrically truncated the distribution. Additionally, $\sigma=1$
My attempt: The MGF of the truncated normal is given by $e^{\frac{\sigma^2t^2}{2}}\left(\frac{\Phi(b-t)-\Phi(a-t)}{\Phi(b)-\Phi(a)} \right)$. Let $Z=X-Y$. Then the MGF of $Z$ is $$M_Z(t)=e^{\frac{2\sigma^2t^2}{2}}\left(\frac{\Phi(b-t)-\Phi(a-t)}{\Phi(b)-\Phi(a)} \right)^2$$ Then
\begin{equation*} \frac{d M_Z(t)}{dt}=e^{\sigma^2t^2}(2\sigma^2t)\left(\frac{\Phi(b-t)-\Phi(a-t)}{\Phi(b)-\Phi(a)} \right)^2+\left(e^{\sigma^2t^2}\left(\frac{\Phi(b-t)-\Phi(a-t)}{\Phi(b)-\Phi(a)} \right) \left(\frac{\phi(a)-\phi(b)}{\Phi(b)-\Phi(a)} \right) \right) \end{equation*}
And \begin{equation*} \begin{split} \frac{d^2 M_Z(t)}{dt^2}&=e^{\sigma^2t^2}(2\sigma^2t)^2\left(\frac{\Phi(b-t)-\Phi(a-t)}{\Phi(b)-\Phi(a)} \right)^2+e^{\sigma^2t^2}(2\sigma^2)\left(\frac{\Phi(b-t)-\Phi(a-t)}{\Phi(b)-\Phi(a)} \right)^2\\ &+2e^{\sigma^2t^2}(2\sigma^2t)\left(\frac{\phi(a)-\phi(b)}{\Phi(b)-\Phi(a)} \right)\left(\frac{\Phi(b-t)-\Phi(a-t)}{\Phi(b)-\Phi(a)} \right)\\ &+\frac{d}{dt}\left(e^{\sigma^2t^2}\left(\frac{\Phi(b-t)-\Phi(a-t)}{\Phi(b)-\Phi(a)} \right) \left(\frac{\phi(a)-\phi(b)}{\Phi(b)-\Phi(a)} \right) \right)+\text{higher order terms with t} \end{split} \end{equation*}
Evaluating $\frac{d2 M_Z(t)}{dt}$ and $\frac{d^2 M_Z(t)}{dt^2}$ at $t=0$ we get $\mathbb{E}(Z)=0$ and $\mathbb{E}(Z^2)=2\sigma^2$. (We used the fact that a=-b implies $\phi(a)=\phi(b)$). This implies that variance of $Z$ exists and is finite
Is the argument correct? I did 100 Monte Carlo rounds with 10000 draws in MATLAB and the the values are quite vary close. The below is the spread of the empirical variance with $a=-5$ and $b=5$
