Gaussian integrals over a half-space

Question

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!

Answer 1

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.

Answer 2

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!