8
$\begingroup$

Edit: I shall try to reformulate my question in order to make it -hopefully- more clear.

Let $X$ be a random variable that follows the $n$-dimensional Gaussian distribution. The probability density of $X$ is given by the following function:

$$ f_X(\mathbf{x;\mathbf{\mu}, \Sigma}) = \frac{1}{(2\pi)^{n/2}|\Sigma|^{1/2}} exp\{ -\frac{1}{2} (\mathbf{x}-\mathbf{\mu})^T \Sigma^{-1} (\mathbf{x}-\mathbf{\mu}) \}, $$

where $\mathbf{x},\mathbf{\mu}\in\Re^n$ and $\Sigma \in S_{++}^{n}$.

In addition, let $\mathcal{H}: \mathbf{x}^T\mathbf{w}=0$ be a hyperplane in the $n$-dimensional Euclidean space $\Re^n$, where $\mathbf{w}\in\Re^n$ and $b\in\Re$. The hyperplane $\mathcal{H}$ defines two half-spaces:

$$ \Omega_{+} = \{\mathbf{x} \in \Re^n | \mathbf{x}^T\mathbf{w}+b \geq 0 \}\\ \Omega_{-} = \{\mathbf{x} \in \Re^n | \mathbf{x}^T\mathbf{w}+b < 0 \} $$

What I would like to find out is what's happening when I try to integrate the above gaussian pdf over the region $\Omega_{+}$ (or $\Omega_{-}$), as it's well-known that -by definition- integrating over the whole $\Re^n$ gives 1. This is the spirit! And this is the reason of the existence if the (normalisation) term $\frac{1}{(2\pi)^{n/2}|\Sigma|^{1/2}}$.

Please correct me if I am wrong in something above.

Next, we could observe that if $b=0$, i.e., the hyperplane crosses the origin, then - due to the symmetry of the gaussian pdf - the integral should be equal to $1/2$. But how can I compute the value of the integral when $b\neq0$? Is there some clever trick to compute how the gaussian integral changes with respect to the bias term $b$?

Thanks in advance and apologies for being kinda "newbie" or not strict enough!... Every correction or helpful comment (either on my question or my notation) would be extremely appreciated!

$\endgroup$
8
  • 1
    $\begingroup$ Rotate it so that you get an integral over a half-space $\{ x_n \geqslant c\}$. Then you have $n-1$ full Gaussian integrals, and one $\int_c^\infty e^{-\xi^2/2}\,d\xi$. $\endgroup$ Nov 8, 2013 at 15:50
  • $\begingroup$ @DanielFischer Thank you very much for your instant answer. I will try it and I will return if I have other questions. It seems right, but I'm afraid that something is gonna go wrong! Thank you anyway! $\endgroup$ Nov 8, 2013 at 15:58
  • $\begingroup$ @DanielFischer, I didn't find the time to look at your approach yet, but because of the fact that I tried something similar, I remember that I stuck at the n-dimensional rotation matrix... Could you please clarify for me what exactly you mean by ``rotate it''? Rotate around what? What is the rotation matrix in n dimensions? Apologies for being newbie... $\endgroup$ Nov 10, 2013 at 17:47
  • $\begingroup$ I made some edits that I think that clarify my original question. This question could be concerned along with the following question math.stackexchange.com/questions/562332/… -- Thanks again! $\endgroup$ Nov 11, 2013 at 12:37
  • 1
    $\begingroup$ What is $b$? I imagine that instead of $0$ in the definitions of $\Omega_+$ and $\Omega_-$, you may have meant to have $b$. $\endgroup$
    – robjohn
    Nov 22, 2013 at 16:05

2 Answers 2

14
+50
$\begingroup$

To compute the integral

$$P = \frac{1}{(2\pi)^{n/2}\sqrt{\det \Sigma}} \int_{\Omega^+} \exp \left(-\frac12 (\mathbf{x}-\mu)^T \Sigma^{-1} (\mathbf{x}-\mu)\right)\,d\mathbf{x},$$

we will need several coordinate transforms. We start with a translation to get rid of the mean, $\mathbf{y} = \mathbf{x}-\mu$. Then $\Omega^+ = \{\mathbf{x} : \mathbf{x}^T\mathbf{w}+b \geqslant 0\}$ corresponds to $\Omega_1 = \{\mathbf{y} : \mathbf{y}^T\mathbf{w} + (\mu^T\mathbf{w} + b)\geqslant 0\}$, and we have

$$P = \frac{1}{(2\pi)^{n/2}\sqrt{\det \Sigma}} \int_{\Omega_1} \exp \left(-\frac12 \mathbf{y}^T \Sigma^{-1} \mathbf{y}\right)\,d\mathbf{y}.$$

Next, since $\Sigma$ is positive definite symmetric, there is an orthogonal matrix $U$ and a diagonal matrix $D$ with positive diagonal elements such that $\Sigma = U^TDU$. Then you have $\Sigma^{-1} = (U^TDU)^{-1} = U^{-1}D^{-1}(U^T)^{-1} = U^TD^{-1}U$ since $U^T = U^{-1}$. Then let $\mathbf{z} = U\mathbf{y}$, and $\mathbf{w}_1 = U\mathbf{w}$. With $\Omega_2 = \{\mathbf{z} : \mathbf{z}^T\mathbf{w}_1 + (\mu^T\mathbf{w}+b)\geqslant 0\}$, we obtain

$$P = \frac{1}{(2\pi)^{n/2}\sqrt{\det\Sigma}} \int_{\Omega_2} \exp \left(-\frac12 \mathbf{z}^T D^{-1} \mathbf{z} \right)\,d\mathbf{z},$$

since $\lvert \det U\rvert = 1$.

Now we do a bit of rescaling. Let $\mathbf{a} = \sqrt{D^{-1}}\mathbf{z}$, or $\mathbf{z} = \sqrt{D}\mathbf{a}$, $\mathbf{w}_2 = \sqrt{D}\mathbf{w}_1$, and $\Omega_3 = \{ \mathbf{a} : \mathbf{a}^T\mathbf{w}_2 + (\mu^T\mathbf{w}+b)\geqslant 0\}$. Since $\det \sqrt{D} = \sqrt{\det\Sigma}$, that yields

$$P = \frac{1}{(2\pi)^{n/2}} \int_{\Omega_3} \exp \left(-\frac12 \mathbf{a}^T\mathbf{a}\right)\,d\mathbf{a}.$$

Now find an orthogonal matrix $B$ such that $B\mathbf{w}_2 = \lVert\mathbf{w}_2\rVert\cdot e_n$, and let $\mathbf{b} = B\mathbf{a}$, $\Omega_4 = \{\mathbf{b} : \mathbf{b}^T(\lVert \mathbf{w}_2\rVert\cdot e_n) + (\mu^T\mathbf{w}+b)\geqslant 0\} = \{\mathbf{b} : \lVert\mathbf{w}_2\rVert b_n + (\mu^T\mathbf{w}+b)\geqslant 0\}$. That yields

$$P = \frac{1}{(2\pi)^{n/2}} \int_{\Omega_4} \exp\left(-\frac12 \mathbf{b}^T\mathbf{b}\right)\,d\mathbf{b}.$$

Now $\Omega_4 = \mathbf{R}^{n-1}\times [c,\infty)$, with $c = -\dfrac{(\mu^T\mathbf{w}+b)}{\lVert \mathbf{w}_2\rVert}$, and hence

$$P = \frac{1}{\sqrt{2\pi}} \int_c^\infty \exp\left(-\frac12 x^2\right)\,dx,$$

which can be evaluated in terms of the error function.

$\endgroup$
12
  • $\begingroup$ What can I say?! Thanks a lot! Great answer! Could you tell me if my approach below (using the circulant matrices) is good? You get what I wanted to do with the change of variables, don't you? Thanks again! $\endgroup$ Nov 25, 2013 at 21:40
  • $\begingroup$ I was afraid you'd ask. To be honest, I have no experience with circulant matrices, I can't tell you off-hand if the approach works. I have to read it very carefully to understand what's going on. Don't expect an answer to that soon. (Certainly not today, maybe tomorrow.) $\endgroup$ Nov 25, 2013 at 21:47
  • $\begingroup$ Whenever you can, of course! I'd just like to know if it's gonna work! If you have something, sometime, please let me know! Thanks again! $\endgroup$ Nov 25, 2013 at 21:51
  • 1
    $\begingroup$ We have $\det\sqrt{D}=\sqrt{\det\Sigma}$ since $\Sigma=U^TDU$, so $\det\Sigma=\det(U^T)\cdot\det D\cdot\det U=(\det U)^2\cdot\det D=\det D$, as $U$ is orthogonal, hence $\det U=\pm1$. Now $D=(\sqrt{D})^2$, whence $\det D=(\det\sqrt{D})^2$, and since everything is positive, $\det\sqrt{D}=\sqrt{\det D}=\sqrt{\det\Sigma}$. Then with the substitution $\mathbf{a}=(\sqrt{D})^{-1}\mathbf{z}$, the transformation formula in the short-cut version says $d\mathbf{z}=\lvert\det\sqrt{D}\rvert\cdot d\mathbf{a}$, and with the above, that cancels the $\sqrt{\det\Sigma}$. $\endgroup$ Nov 26, 2013 at 20:16
  • 1
    $\begingroup$ Yes, we have $\mathbf{w}_2 = \sqrt{D}U\mathbf{w}$, but unless $\Sigma$ is a multiple of the identity, you can't unfortunately express $\lVert\mathbf{w}_2\rVert$ that easily in terms of $\mathbf{w}$. However, $\lVert\mathbf{w}_2\rVert^2 = (\sqrt{D}U\mathbf{w})^T(\sqrt{D}U\mathbf{w}) = \mathbf{w}^TU^T\sqrt{D}\sqrt{D}U\mathbf{w} = \mathbf{w}^T\Sigma\mathbf{w}$, so you can write the norm as $\sqrt{\mathbf{w}^T\Sigma\mathbf{w}}$. $\endgroup$ Nov 26, 2013 at 23:04
0
$\begingroup$

There is something that could help... I would like some help on notation or whatever else needed, though!

Let $\mathbf{x}=R^T\mathbf{y}$. Then the hyperplane can be rewritten as follows:

$$ \mathcal{H}: \mathbf{w}^TR^T\mathbf{y}+b=0. $$

We define $R$ to be the circulant matrix generated by a vector $\mathbf{r}=(r_1,\dots,r_n)^T$. Then, I think that the following holds:

$$ \mathbf{w}^TR^T = (R\mathbf{w})^T = (W\mathbf{r})^T = \mathbf{r}^T W^T, $$

where $W$ is the circulant matrix generated by the vectors $\mathbf{w}$. Consequntly, if we demand

$$ \mathbf{r}^T W^T = \mathbf{v}^T, $$

then $\mathbf{r}^T = \mathbf{v}^T(W^T)^{-1}=\mathbf{v}^T(W^{-1})^T$. Thus,

$$ \mathbf{r} = W^{-1}\mathbf{v}. $$

If we set the vector $\mathbf{v}$ to be, for instance, $(1,0,\dots,0)^T$, then the change of variable we made will lead to a hyperplane as follows:

$$ \mathcal{H'}: y_1=c. $$

Am I right so far?

In order to find the matrix $R$, which constructs the change of variable ($\mathbf{x}=R^T\mathbf{y}$), we need to find the vector $\mathbf{r}$, which constructs the circulant matrix $R$. We can find $\mathbf{r}$ can be found from the following:

$$ \mathbf{r} = W^{-1}\mathbf{v}, $$

where $W$ is known and is generated as a circulant matrix from $\mathbf{w}$.

But $W$ is a circulant matrix, so it can be expressed as follows:

$$ W = Q\Lambda Q^{-1}, $$

where $Q$ is the matrix of the eigenvectors of $W$, while $\Lambda$ is the diagonal matrix of the corresponding eigenvalues of $W$. For a circulant matrix

$$ W= \begin{bmatrix} w_0 & w_{n-1} & \dots & w_{2} & w_{1} \\ w_{1} & w_0 & w_{n-1} & & w_{2} \\ \vdots & w_{1} & w_0 & \ddots & \vdots \\ w_{n-2} & & \ddots & \ddots & w_{n-1} \\ w_{n-1} & w_{n-2} & \dots & w_{1} & w_0 \\ \end{bmatrix}, $$

the $j$-th eigenvector is given by

$$ \mathbf{q}_j = (1,\omega_j, \omega_j^2, \dots, \omega_j^{n-1})^T, \: j=1,...,n-1 $$

where $\omega_j = exp\left(\frac{2\pi i j}{n}\right)$ are the $n$-th roots of unity, and $i$ the imaginary unit. Moreover, the corresponding $j$-th eigenvalue is given by $\lambda_j = w_0 + w_{n-1}\omega_j + \dots + c_1\omega_j^{n-1}$, $j=0,\dots,n-1$.

As a result we can find $\mathbf{r}$ as follows:

$$ \mathbf{r} = W^{-1}\mathbf{v} = Q\Lambda^{-1}Q^{-1}\mathbf{v}. $$

Thus, having the vectors $\mathbf{r}$ we can construct the circulant matrix $R$, and using the change of variable $\mathbf{x}=R^T\mathbf{y}$ we can go to a hyperplane of the form $y_n=c$, where the integration of simpler.

Am I right so far? Please give me some feedback! Thanks in advance!

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.