Sign of product of standard normal random variables The question is mainly in the title:

Given two standard normal random variables with correlation $\rho$, what is the distribution of sign of their product?

I understand that when $\rho=0$, we have two iid standard normal random variables and therefore, they take positive and negative values independently with probability $\frac12$.
But I don't know what to do if $\rho\ne0$. We can take $\rho>0$ without loss of generality (because if $X$ and $Y$ are std normal with correlation $\rho$, then $-X$ and $Y$ are std normal with correlation $-\rho$). But I could not proceed further.
 A: Consider $X,Y\sim N(0,1)$ with correlation $\rho$. Then (Nadarajah & Pogány, 2016; Gaunt, 2018) their product is variance-gamma distributed:
$$ XY \sim \text{VG}(1,\rho,\sqrt{1-\rho^2},0). $$
Its PDF is
$$ f_{XY}(x) = \frac{1}{\pi\sqrt{1-\rho^2}}\exp\left(\frac{\rho x}{1-\rho^2}\right) K_0\left(\frac{|x|}{1-\rho^2}\right), $$
where $K_0$ is the modified Bessel function of the second kind of order $0$.
Thus
$$ P(\text{sign}(XY) = -1) = \int_{-\infty}^0 f_{XY}(x)\,dx. $$
You may be able to evaluate this using your favorite computer algebra system. Unfortunately, this exceeds the standard computation time for WolframAlpha. COOLSerdash notes that the integral evaluates to a nice round $\frac{\arccos\rho}{\pi}$.
Alternatively, you could map the above parameterization to the one employed in the VarianceGamma package for R and use the functions in there, if all you are interested in is numerical results.
A: The correlation alone is not enough to be able to derive the probability distribution for the sign.
See the example below for a case where $\rho = 2/\pi$ and the sign of $XY$ is always positive.

A: As long as $(X,Y)$ is standard bivariate normal with correlation $\rho$, the probability that $XY$ is positive or negative can be found using the well-known result for the positive quadrant probability $$P(X>0,Y>0)=\frac14+\frac1{2\pi}\sin^{-1}\rho \tag{1}$$
(This is likely discussed here before but I cannot quite find the question.)
You have
\begin{align}
P(XY>0)&=P(X>0,Y>0)+P(X<0,Y<0)
\\&=P(X>0,Y>0)+P(-X>0,-Y>0)
\end{align}
Because $(-X,-Y)$ has the same distribution as $(X,Y)$, this probability is just $$P(XY>0)=2P(X>0,Y>0)$$
Similarly,
\begin{align}
P(XY<0)&=P(X>0,Y<0)+P(X<0,Y>0)
\\&=P(X>0,-Y>0)+P(-X>0,Y>0)
\end{align}
Again, $(X,-Y)$ and $(-X,Y)$ have the same distribution, so
$$P(XY<0)=2P(X>0,-Y>0)$$
And since $(X,-Y)$ is bivariate normal with correlation $-\rho$, we have from $(1)$ that
$$P(X>0,-Y>0)=\frac14-\frac1{2\pi}\sin^{-1}\rho$$
