Given a normal random vector $$X\sim N(\mu,\Sigma)$$ for spd $\Sigma$, I'm interested in the covariance matrix $K=\mathrm{cov}(Y)$ of the variable $$Y = \mathrm{relu}(X)$$ where the relu is performed elementwise $Y_i = \mathrm{Max}(0,x_i)$, so $Y$ is distributed according to the rectified Gaussian distribution.

Given I know everything about $\Sigma$, how can I compute $K$?

The mean and variance of each $Y_i$ has been covered in other questions on this site, but the off-diagonal elements of $K$ seem pretty challenging to compute, and I haven't found anything on SO or elsewhere online about it.

I'm actually after the eigenvectors of $K$, so if anyone can relate the eigenvectors between $\Sigma$ and $K$ without directly computing $K$, that would be even more interesting.


Edit: Just to note there is a similar question asked here and thoroughly answered, but only in the scalar (or diagonal multivariate) case. For multi-dimensional $X$ with correlations, this seems much more challenging.

  • $\begingroup$ It may be possible, but efforts to deal with a truncated bivariate normal seem to suggest it may not be easy $\endgroup$
    – Henry
    Commented Mar 22, 2020 at 3:01
  • $\begingroup$ @MisterBlobfish I'm looking into a similar question and wondering whether you found a satisfactory answer in the end? What I'm looking for is to compute $\mathbb{E}(\text{relu}(X)\text{relu}(Y))$ when X, Y are dependent. $\endgroup$
    – smalldog
    Commented Feb 21, 2023 at 12:51
  • $\begingroup$ Related: math.stackexchange.com/questions/1963292/… $\endgroup$
    – ABIM
    Commented Feb 14 at 1:16

1 Answer 1


In the 2d case, let $$ \Sigma=\begin{pmatrix}\sigma_x^2 & r\sigma_x\sigma_y \\ r\sigma_x\sigma_y & \sigma_y^2\end{pmatrix}. $$ Then by Mathematica (or this paper, table 1), $\mathrm{E}[\max(X,0)\cdot\max(Y,0)] = \sigma_x\sigma_y f(r),$ where $$f(r) = \frac{2\sqrt{1-r^2}+r(\pi+2\arctan(\tfrac{r}{\sqrt{1-r^2}}))}{4\pi}.$$ We also have $E[\max(X,0)]=\frac{\sigma_x}{\sqrt{2\pi}}$.

Back to the general case, let $\Sigma \in \mathbb{R}^{n\times n}$ be the covariance matrix, and $\sigma$ be the diagonal of $\Sigma^{1/2}$. \begin{align}\mathrm{V}(\mathrm{relu}(X)) &=\mathrm{E}(\mathrm{relu}(X)\mathrm{relu}(X)^T) -\mathrm{E}(\mathrm{relu}(X))\mathrm{E}(\mathrm{relu}(X)^T) \\&=\sigma\sigma^T\circ f(\Sigma/(\sigma\sigma^T)) - \frac{\sigma\sigma^T}{2\pi}. \end{align}


We use the bounds $r/4+1/(2\pi)\le f(r)\le (1+r)/4$, we get \begin{align}\mathrm{V}(\mathrm{relu}(X)) &\le\sigma\sigma^T\circ (1+\Sigma/(\sigma\sigma^T))/4 - \frac{\sigma\sigma^T}{2\pi} \\&= \sigma\sigma^T(1/4-1/(2\pi)) + \Sigma/4. \end{align} and \begin{align}\mathrm{V}(\mathrm{relu}(X)) &\ge\sigma\sigma^T\circ (1/(2\pi)+\Sigma/(4\sigma\sigma^T)) - \frac{\sigma\sigma^T}{2\pi} \\&= \Sigma/4. \end{align}

Or in other words $$0 \le \mathrm{V}(\mathrm{relu}(X)) - \Sigma/4 \le 0.1 \sigma\sigma^T.$$

In my own experiments, $\Sigma/4$ seems like a pretty good approximation.

In all of the above we have assumed $\mu=0$. If this is not the case things get even more hairy.

I would assume, that if $(\mu_x, \mu_y)$ is in the positive quadrant, then most of the normal distribution will be within the non-zero region of the relu, so you can basically ignore the relu.

Meanwhile, if $(\mu_x, \mu_y)$ is in another quadrant, then only a very small part of the normal distribution will be in the non-zero region. Because of the sharp tails of the distribution, we can mostly focus on what happens on the border region of the $(\mu_x, \mu_y)$ quadrant and the positive one.


You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .