Derivative of joint probability wrt mean For two jointly normally distributed random variables $X$ and $Y$, with $X\sim N(0,1)$ and $Y\sim N(\mu,1)$, I'd like to establish how $\Pr(X>a, Y<b)$ changes with $\mu$. That is, establish the sign of $$\frac{d }{d \mu} \int_{-\infty}^b \int_a^\infty  f(x,y) dx \, dy$$
where $f(x,y)$ notes the joint normal PDF of $X$ and $Y$.
Numerical simulations tell me that this is negative independently of $b,a$ and $\rho_{xy}$ $-$which I find graphically intuitive$-$, but I can't show it mathematically.
 A: Suppose $$(X,Y) \sim \operatorname{BiNormal}((0,0), (1,\sigma^2), \rho),$$ where $\sigma > 0$ is the standard deviation of $Y$, and $-1 \le \rho \le 1$ is the correlation.  Then $$f_{X,Y}(x,y) = \frac{\exp\left(-\frac{1}{1-\rho^2}\left(\frac{x^2}{2}-\frac{\rho x y}{\sigma}+\frac{y^2}{2\sigma^2}\right)\right)}{2\pi\sqrt{1-\rho^2}\sigma}.$$  A scaling transformation of $Y$, namely $W = Y/\sigma$, results in $(X,W)$ being bivariate normal with variances equal to $1$ and the same correlation $\rho$, so we can consider instead the probability $$\Pr[(X > a) \cap (Y < b)] = \Pr[(X > a) \cap (W < b/\sigma)] = \int_{w=-\infty}^{b/\sigma} \int_{x=a}^\infty f_{X,W}(x,w) \, dx \, dw.$$  It is clear that the integrand is independent of $\sigma$, and the probability depends only on the scale parameter $\sigma$ through the upper limit.  Explicitly, let $$\begin{align*} g(w) = \int_{x=a}^\infty f_{X,W}(x,w) \, dx &= \int_{x=a}^\infty \frac{e^{-(x-\rho w)^2/(2(1-\rho^2))}e^{-w^2/2}}{2\pi \sqrt{1-\rho^2}} \, dx \\ &= \frac{e^{-w^2/2}}{\sqrt{2\pi}} \left(1 - \Phi\left( \frac{a - \rho w}{\sqrt{1-\rho^2}}\right) \right) \end{align*}$$ where $\Phi$ is the CDF of the (univariate) standard normal distribution.  Consequently, $$\frac{d}{d\sigma}\int_{w=-\infty}^{b/\sigma} g(w) \, dw = -g(b/\sigma) \frac{b}{\sigma^2}.$$  Because $g(w) > 0$ for all $w \in \mathbb R$ and $\sigma^2 > 0$, it follows that the sign of the derivative of the probability in question is solely determined by the sign of $b$; namely, $$\frac{d}{d\sigma}\left[\Pr[(X > a) \cap (Y < b)]\right] > 0$$ if and only if $b < 0$.

Update.  Since we are instead looking at $$(X,Y) \sim \operatorname{BiNormal}((0,\mu), (1,1), \rho),$$ the same technique applies: we perform a location transformation of the form $W = Y - \mu$ to get $$\Pr[(X > a) \cap (Y < b)] = \Pr[(X > a) \cap (W < b-\mu)] = \int_{w=-\infty}^{b-\mu} \int_{x=a}^\infty f_{X,W}(x,w) \, dx \, dw.$$  And again, $g(w)$ is the same function as before.  The only difference (no pun intended) is that $$\frac{d}{d\mu} \int_{w=-\infty}^{b-\mu} g(w) \, dw = -g(b - \mu).$$  And since we established the positivity of $g$, it follows that the derivative $$\frac{d}{d\mu}\left[\Pr[(X > a) \cap (Y < b)]\right] < 0$$ for all $a, b, \rho$.
