We square an integral, but why change a variable? 
If we square an integral, we also change the integration variable in one of the integrals. But why is this actually correct?

For example, say I have the following:
Solve $\int_{-\infty}^\infty e^{-x^2} dx$. Let $I=\int_{-\infty}^\infty e^{-x^2} dx$, so
\begin{align}
I^2 &=\bigg(\int_{-\infty}^\infty e^{-x^2} dx\bigg)^2\\
&= \bigg(\int_{-\infty}^\infty e^{-x^2} dx\bigg) \times 
\bigg(
\underbrace{
\int_{-\infty}^{\infty} e^{-y^2}dy
}_{\text{Why}?}
\bigg) 
\\
&=\int_{-\infty}^\infty\bigg(\int_{-\infty}^\infty e^{-(x^2+y^2)} dx\bigg)dy 
\end{align}
But why is the following wrong:
\begin{align}
I^2 &=\bigg(\int_{-\infty}^\infty e^{-x^2} dx\bigg)^2\\
&= \bigg(\int_{-\infty}^\infty e^{-x^2} dx\bigg) \times \bigg(\int_{-\infty}^{\infty} e^{-x^2}dx\bigg) \\
&=\int_{-\infty}^\infty\bigg(\int_{-\infty}^\infty e^{-x^2-x^2} dx\bigg)dx \\
&=\int_{-\infty}^\infty\bigg(\int_{-\infty}^\infty e^{-2x^2} dx\bigg)dx \qquad ?
\end{align}

The short answer is that in your second to last identity you are neglecting cross terms in your multiplication.
A simple example may be used to demonstrate the error, and we will use a strict summation instead of an integration for clarity.  Consider the sum
\begin{equation}
  \sum_{x=1}^3 x = 1 + 2 +3 = 6.
\end{equation}
Now, squaring the summation yields
\begin{align}
  \left(\sum_{x=1}^3x\right)^2 &= (1+2+3)^2 = (1+2+3)\times(1+2+3)
\end{align}
In order to properly calculate this quantity (long-hand) requires one to completely distribute including cross terms:
\begin{align}
  &\ {\color{white}+}\ 1\times 1 + 1\times 2 + 1\times 3 \notag \\
  (1+2+3)\times(1+2+3) =& +2\times 1 + 2\times 2 + 2\times 3 = 36. \notag \\
  & + 3\times 1 + 3\times 2 + 3\times 3
\end{align}
However, if the cross terms were neglected one would obtain
\begin{equation}
  (1+2+3)\times(1+2+3) \ne (1^2+2^2+3^2) = 1\times 1 + 2\times 2 + 3\times 3 = 14.
\end{equation}
Written another way, we conclude that
\begin{equation}
  \left(\sum_{x=1}^3x\right)^2 = \left(\sum_{x=1}^3x\right)\left(\sum_{x=1}^3x\right) \ne \sum_{x=1}^3x^2.
\end{equation}
However, as noted in the comments, the $x$ in the summation (or the in the definite integral) is just a dummy variable, and may be replaced with another symbol, such as $y$:
\begin{equation}
  \left(\sum_{x=1}^3x\right)^2 = \left(\sum_{x=1}^3x\right)\left(\sum_{y=1}^3y\right).
\end{equation}
The reason for changing the dummy variable is that we may now say that the summation over $x$ and $y$ are independent, and so we may rearrange the summations:
\begin{equation}
  \left(\sum_{x=1}^3x\right)^2 = \left(\sum_{x=1}^3x\right)\left(\sum_{y=1}^3y\right) = \sum_{x=1}^3\sum_{y=1}^3x\times y.
\end{equation}
Again, integration is just summation, and so the same fact holds true for integrals.  Furthermore, it doesn't matter what the function inside the summation/integral is.  Thus, we may write
\begin{equation}
  \left(\int_a^b f(x)\ dx\right)^2 = \left(\int_a^b f(x)\ dx\right)\left(\int_a^b f(y)\ dy\right) = \int_a^b\int_a^b f(x)\times f(y)\ dx\ dy.
\end{equation}
Using $f(x)=e^{-x^2}$ results in your example problem.
