Showing that $P (|X_1| < \Gamma , |X_2| < \Gamma)$ is increasing in $|\rho|$ Assume $(X_1,X_2 )^T$ is mean $0$ bivariate normal distributed with covariance matrix $\Sigma = \left (\begin{matrix} 1 & \rho \\ \rho & 1 \end{matrix} \right)$ and let $\Gamma > 0$ a positive constant. Then i would like to show that
$P (|X_1| < \Gamma , |X_2| < \Gamma)$ is increasing in $|\rho|$.
Any tips? I already tried to simply use the integral representation of the probability, but could not show it.
 A: This follows from Mehler's Formula, in the form 
$$f(x_1,x_2)=\phi(x_1)\phi(x_2)\sum_{k=0}^\infty\frac {\rho^n}{n!} He_n(x_1)He_n(x_2),$$ 
where $f$ is the joint density of $(X_1,X_2)$,  $\phi$ is the marginal density of the $X_i$, and $He_n(x)$ is a "probabilist's" Hermite polynomial.  Integrating over $[-\Gamma,\Gamma]\times[-\Gamma,\Gamma]$ yields an expression of form $$
I(\Gamma,\rho) := P(|X_1|<\Gamma, |X_2|<\Gamma) = 
\sum_{k=0}^\infty\frac {\rho^n}{n!}\left(\int_{-\Gamma}^\Gamma \phi(x)He_n(x)\,dx\right)^2,\tag 1$$
all of whose terms are non-negative.  In fact $He_n$ is an odd polynomial if $n$ is odd, so the odd terms in (1) vanish, so we see $I(\Gamma,\rho)$ is given by a power series in $\rho^2$ with non-negative coefficients.
A: Note that the pdf $f(x_1,x_2)$ satisfies $\frac{\partial f}{\partial \rho} > 0$.  From there, it suffices to observe that
$$
\frac{d P}{d \rho} = \frac{d}{d\rho} \iint_R f(x_1,x_2) dx_1 dx_2= 
\ \iint_R \frac{\partial f}{\partial \rho}\,(x_1,x_2) \,dx_1 dx_2 < 0,
$$
where $R$ denotes the rectangle $[-\Gamma,\Gamma]\times [-\Gamma,\Gamma]$.  After computing $\frac{\partial f}{\partial \rho}$, you should find that one of the resulting integrals can be solved as an iterated integral; solve the inner-integral with the substitution $u_2 = x_2^2$.
A: Let $f_{\rho}$ be the density of $(X_1,X_2)$.
An interesting identity I stumbled upon in this connection is 
$$\frac{\partial}{\partial\rho}f_{\rho}(x,y)=\frac{\partial^2}{\partial x\partial y}f_{\rho}(x,y)$$
This is part of a general result due to Plackett (1954).
Then differentiating under the integral sign,
\begin{align}
\frac{d}{d\rho}P(|X_1|<c,|X_2|<c)&=\int_{-c}^c\int_{-c}^c \frac{\partial}{\partial\rho}f_{\rho}(x,y)\,dx\,dy
\\&=\int_{-c}^c\int_{-c}^c \frac{\partial^2}{\partial x\partial y}f_{\rho}(x,y)\,dx\,dy
\\&=\int_{-c}^c \frac{\partial}{\partial x}\left(\int_{-c}^c \frac{\partial}{\partial y} f_{\rho}(x,y)\,dy\right)dx
\\&=\int_{-c}^c \frac{\partial}{\partial x} f_{\rho}(x,c)\,dx-\int_{-c}^c \frac{\partial}{\partial x} f_{\rho}(x,-c)\,dx
\\&=f_{\rho}(c,c)-f_{\rho}(-c,c)-f_{\rho}(c,-c)+f_{\rho}(-c,-c)
\\&=2\left(f_{\rho}(c,c)-f_{\rho}(c,-c)\right)
\end{align}
That is, for some $c>0$ we have 
$$\frac{d}{d\rho}P(|X_1|<c,|X_2|<c)=\frac1{\pi\sqrt{1-\rho^2}}\left[e^{-c^2/(1+\rho)}-e^{-c^2/(1-\rho)}\right]$$
Note that $\rho > 0\implies -\frac1{1+\rho}>-\frac1{1-\rho}$, so that right hand side of the above equation is positive for $\rho\in (0,1)$. Similarly it is negative for $\rho\in (-1,0)$. Hence the probability is increasing in $|\rho|$.
