From 1D gaussian to 2D gaussian I read this:

The Gaussian kernel for dimensions higher than one, say N, can be described as a regular product of N one-dimensional kernels. Example: 
g2D(x,y,$\sigma_1^2 + \sigma_2^2$) = g1D(x,$\sigma_1^2$)g2D(y,$\sigma_2^2$)

saying that the product of two 1 dimensional gaussian functions with variances $\sigma_1^2$ and $\sigma_2^2$ is equal to a two dimensional gaussian function with the sum of the two variances as its new variance.
I tried to deduce this by using:
g1D(x,$\sigma1^2$)g2D(y,$\sigma2^2$) = $\frac{1}{\sqrt{2\pi}\sigma_1}e^{\frac{-x^2}{2\sigma_1^2}}\frac{1}{\sqrt{2\pi}\sigma_2}e^{\frac{-y^2}{2\sigma_2^2}}$ = $\frac{1}{2\pi\sigma_1\sigma_2}e^{-(\frac{x^2}{2\sigma_1^2}+\frac{y^2}{2\sigma_2^2})}$
but I fail to obtain 
$\frac{1}{2\pi(\sigma_1^2 + \sigma_2^2)}e^{\frac{-(x^2+y^2)}{2\sigma_1^2 + 2\sigma_2^2}}$
which is equal to g2D(x,y,$\sigma_1^2 + \sigma_2^2$).
Someone know how to get there?
 A: Your first expression,
$G(x,\sigma_1^2) G(y,\sigma_2^2) = \frac{1}{\sqrt{2\pi}\sigma_1}e^{\frac{-x^2}{2\sigma_1^2}}\frac{1}{\sqrt{2\pi}\sigma_2}e^{\frac{-y^2}{2\sigma_2^2}}$ = $\frac{1}{2\pi\sigma_1\sigma_2}e^{-(\frac{x^2}{2\sigma_1^2}+\frac{y^2}{2\sigma_2^2})}$
is a correct anisotropic 2D Gaussian. The expression you are trying to obtain,
$\frac{1}{2\pi(\sigma_1^2 + \sigma_2^2)}e^{\frac{-(x^2+y^2)}{2\sigma_1^2 + 2\sigma_2^2}}$
is not (it’s an isotropic Gaussian). There is no way you can go from the one to the other. However, one can incorrectly match the two expressions by creating a new coordinate system that is a scaled version of the original one, in such a way that the Gaussian is isotropic in the new coordinate system. The other answer does this. But that is not useful.
It seems that OP is confused about the combination of the variances. The convolution of two 1-dimensional Gaussian functions with variances $\sigma_1^2$ and $\sigma_2^2$ is equal to a 1-dimensional Gaussian function with variance $\sigma_1^2 + \sigma_2^2$. There is no such expectation for the multiplication of Gaussians (in fact, when multiplying them, assuming the same orientation and the same mean, the resulting variance is smaller, not larger. 
A: Hint: You see in the exponential that we have
$$ \exp \left ( - \frac{x^2 + y^2}{2( \sigma_1^2 + \sigma_2^2 ) }  \right )=\exp \left (  \frac{-x^2 }{2( \sigma_1^2 + \sigma_2^2 ) } \right )\exp \left (  \frac{- y^2}{2( \sigma_1^2 + \sigma_2^2 ) }  \right )  $$
The normalization hanging around should be
$$\frac{1}{2 \pi ( \sigma_1^2 + \sigma_2^2 ) } $$
as you've seen...but under the change of variables
$$ z_1 = \frac{ \sigma_1x}{\sqrt{\sigma_1^2 + \sigma_2^2}}  \quad \& \quad z_2 = \frac{ \sigma_2y}{\sqrt{\sigma_1^2 + \sigma_2^2}}  $$
What do you obtain? Also note that
$$ 2 \pi = \sqrt{2\pi} \sqrt{2 \pi} $$
Edit: We see that
$$\frac{1}{2 \pi ( \sigma_1^2 + \sigma_2^2 ) }  \int_{\mathbb{R}^2} \exp \left ( - \frac{x^2 + y^2}{2( \sigma_1^2 + \sigma_2^2 ) }  \right ) dxdy = $$
$$= \left[ \frac{1}{\sqrt{2 \pi ( \sigma_1^2 + \sigma_2^2 )} }  \int_{\mathbb{R}} \exp \left ( - \frac{x^2 }{2( \sigma_1^2 + \sigma_2^2 ) }  \right ) dx \right] \left [ \frac{1}{\sqrt{2 \pi ( \sigma_1^2 + \sigma_2^2 )} }   \int_{\mathbb{R}} \exp \left ( - \frac{y^2 }{2( \sigma_1^2 + \sigma_2^2 ) }  \right )dy  \right]$$
$$ =\left[ \frac{1}{\sqrt{2 \pi }\sigma_1 }  \int_{\mathbb{R}} \exp \left ( - \frac{z_1^2 }{2 \sigma_1^2   }  \right ) dz_1 \right] \left [ \frac{1}{\sqrt{2 \pi }\sigma_2 }   \int_{\mathbb{R}} \exp \left ( - \frac{z_2^2 }{2\sigma_2^2  }  \right )dz_2  \right] $$
