Correlation between two nodes of a single layer MLP for joint-Gaussian input

Let's say you have a jointly gaussian vector random variable $\mathbf{x}$, with mean $\mathbf{M}$ and covariance $\mathbf{S}$. I now transform each scalar element of $\mathbf{x}$ , say $x_j$, with a sigmoid:

$$y_j = 1/(1+\exp(-x_j))$$

I am interested in the expectation between two variables $y_j$ and $y_j'$ of the resulting distribution, that is

$$E\{y_j y_j'\}$$

where $E\{\}$ is the expectation operator. Note:

• The whole PDF of $y$ can be computed just by applying a change of variable. Unfortunately the integrals leading to the expectations are intractable (to the best of my knowledge).

• I'm looking for closed-form approximations, no Markov Chain-Monte Carlo, no variational stuff. That is, approximations to the variable change, to the expectation integrals, to the sigmoid, to the resulting PDF, that allow computing $E\{y_j y_j'\}$.

• Some dead ends: Taylor, often used in papers on the topic, is inaccurate by the mile. Gradshteyn and Ryzhik does not seem to contain the integrals.

• What do you mean by "expectation between two variables"? Mar 2 '11 at 15:20
• @user3509, you can get an easy lower bound by noting that $x^{-1}$ is convex and if $X$ is marginally Gaussian, then $e^{-X}$ is lognormal and so $\mathbb{E}(Y_i Y_j) > (1 + \mathbb{E} e^{-X_i} + \mathbb{E} e^{-X_j} + \mathbb{E} e^{-X_i-X_j})^{-1}$. The expectations on the right hand side are easy to compute. Mar 2 '11 at 15:53
• Stationary phase (using a linear approximation to find the critical point) looks like it will work, but the resulting expressions--although they will be algebraic combinations of $\mathbf{M}$ and $\mathbf{S}$--will be extremely messy. How far are you prepared to go with this and what's the purpose of a "closed form" solution (when perhaps a numerical algorithm will in practice be more efficient and accurate)?
– whuber
Mar 2 '11 at 16:40
• @bayer, I mean E{y_1,y_2}=\int\int s(x_1)s(x_2) p(x_1,x_2) dx_1 dx_2 where y=s(x) is the sigmoid function, p(x,y) is the joint-PDF (gaussian). The integral is over RxR. Mar 2 '11 at 17:05
• @user3509 TeX is enabled in comments as well as questions and replies: surround it with $\$$, as in "\$$E\{y_1,y_2\}=\int\int s(x_1)s(x_2) p(x_1,x_2) dx_1 dx_2$\$$" gives E\{y_1,y_2\}=\int\int s(x_1)s(x_2) p(x_1,x_2) dx_1 dx_2 – whuber Mar 2 '11 at 17:12 1 Answer The question really concerns pairs of normal variates. Let's call them x_1 and x_2 with means \mu_i, standard deviations \sigma_i, and correlation \rho. Whence their joint pdf is$$\frac{1}{2 \pi \sqrt{1 - \rho^2} \sigma_1 \sigma_2} e^{-\frac{1}{1-\rho^2} \left(\frac{(x_1 - \mu_1)^2}{2 \sigma_1^2} + \frac{(x_2 - \mu_2)^2}{2 \sigma_2^2} - \frac{\rho (x_1 - \mu_1)(x_2 - \mu_2)}{\sigma_1 \sigma_2}\right)} dx_1 dx_2\text{.}$$Let f(x_1,x_2) be the product of this with the y_i (as functions of the x_i). The first component of the gradient of \log(f) is$$\frac{\partial \log(f)}{\partial x_1} = \frac{1}{1 + e^{x_1}} + \frac{\rho(\mu_2 - x_2) \sigma_1 + (x_1 - \mu_1)\sigma_2}{(\rho^2-1)\sigma_1^2 \sigma_2},$$with a similar expression for the second component (via the symmetry achieved by exchanging the subscripts 1 and 2). There will be a unique global maximum, which we can detect by setting the gradient to zero. This pair of nonlinear equations has no closed form solution. It is rapidly found by a few Newton-Raphson iterations. Alternatively, we can linearize these equations. Indeed, through second order, the first component equals$$\frac{1}{2} + x_1\left(\frac{-1}{4} + \frac{1}{(\rho^2-1)\sigma_1^2}\right) + \frac{-\rho x_2 \sigma_1 + \rho \mu_2 \sigma_1 - \mu_1 \sigma_2}{(\rho^2 -1)\sigma_1^2 \sigma_2}.$$This gives a pair of linear equations in (x_1, x_2), which therefore do have a closed form solution, say \hat{x}_i(\mu_1, \mu_2, \sigma_1, \sigma_2, \rho), which obviously are rational polynomials. The Jacobian at this critical point has 1,1 coefficient$$\frac{e^\hat{x_1}\left(2 - (\rho^2-1)\sigma_1^2 + 2\cosh(\hat{x_1})\right)}{(1+e^\hat{x_1})^2(\rho^2-1)\sigma_1^2},$$1,2 and 2,1 coefficients$$\frac{\rho}{\sigma_1 \sigma_2(1 - \rho^2)},$$and 2,2 coefficient obtained from the 1,1 coefficient by symmetry. Because this is a critical point (at least approximately), we can substitute$$e^\hat{x_1} = \frac{(\rho^2-1)\sigma_1^2 \sigma_2}{(\mu_2 - \hat{x_2})\rho \sigma_1 + (\hat{x_1} - \mu_1)\sigma_2} - 1

and use that also to compute $\cosh(\hat{x_1}) = \frac{e^\hat{x_1} - e^{-\hat{x_1}}}{2}$, with a similar manipulation for $e^\hat{x_2}$ and $\cosh(\hat{x_2})$. This enables evaluation of the Hessian (the determinant of the Jacobian) as a rational function of the parameters.

The rest is routine: the Hessian tells us how to approximate the integral as a binormal integral (a saddlepoint approximation). The answer equals $\frac{1}{2\pi}$ times a rational function of the five parameters: that's your closed form (for what it's worth!).

• I am a bit concerned over the linearization but it is definitely worth trying. Do you know of some available literature on stationary phase aproximation of integrals?. Mar 2 '11 at 19:20
• I was concerned about the linearization, too, but (a) the series actually is linear through second order and (b) extensive plotting of the pdf for means between -3 and 3, sds from 0.1 to 3, and correlations from -1 to 1, indicate it should generally work well (the linearized critical point appears quite close to the global maximum). I learned this technique from an article in The American Statistician perhaps some 10-15 years ago, but I can't locate it at this moment.
– whuber
Mar 2 '11 at 19:27
• @user3509 I found the TAS reference and linked to it in the text: it goes by the name "saddlepoint approximation." I apologize for confusing it with stationary phase (they are essentially the same thing: one with a real variable, the other with an imaginary variable).
– whuber
Mar 2 '11 at 19:42
• Yes, that was exactly what I could not understand I was looking for the complex exponential somewhere. Thank you very much. Mar 2 '11 at 20:13
• +1 for the nice answer and the link to the article. Did not know that it is feasible to use such approximations. I would have resorted to numerical integration without even trying to get the closed form solution. Mar 3 '11 at 11:29