# Reference Request: Integral of Gaussian over Unit Sphere

I am looking for a reference for integrals of the form $$$$\tag{1} \int_{S^{n-1}} \mathcal{N}_{\omega} ( \mu , \Sigma ) d \omega$$$$ where $$S^{n-1}$$ is the sphere in $$\mathbb{R}^n$$ and $$$$\mathcal{N}_x (\mu , \Sigma) = \frac{e^{- \frac{1}{2} ( x - \mu ) \cdot \Sigma^{-1} (x - \mu) }}{\sqrt{ ( 2 \pi )^n \det{\Sigma} }}$$$$ is the usual multivariate Gaussian distribution evaluated at $$x \in \mathbb{R}^n$$ with mean $$\mu \in \mathbb{R}^n$$ and covariance matrix $$\Sigma \in \mathbb{R}^{n \times n}$$. I am primarily concerned with the $$n = 2, 3$$ cases, but I find the case for general $$n$$ to be interesting as well.

There is a substantial simplification that can be made. By the spectral theorem and positive-definite nature of $$\Sigma$$ we may write $$\Sigma = V D^2 V^T$$ where $$V \in \mathrm{SO}_n (\mathbb{R})$$ and $$D = \mathrm{diag}(\sigma_1 , \cdots , \sigma_n)$$ is a diagonal matrix with positive, real values along the diagonal. Therefore, by spherical symmetry, I am essentially asking for a reference to the integral $$$$\tag{2} \int_{S^{n-1}} e^{ - \frac{1}{2} \| D^{-1} (\omega - \tilde{\mu}) \|^2 } d \omega ,$$$$ where $$\tilde{\mu} = V^T \mu$$.

Question: Does anyone have a reference or suggestion for how to calculate the closed-form expression for (2) and, hence, (1)?

Additional Thoughts: There are special cases that are simple enough to compute. For example, if $$\sigma_1 = \cdots = \sigma_n \equiv \sigma$$, then the integral in question becomes $$$$\frac{e^{- \frac{1 + |\mu|^2}{2 \sigma^2} }}{(2 \pi)^{n/2} \sigma} \int_{S^{n-1}} e^{\omega \cdot \mu / \sigma^2 } d \omega .$$$$ The previous integral may be evaluated using the formula $$$$\int_{S^{n-1}} e^{ \omega \cdot y } d \omega = (2 \pi)^{\frac{n}{2}} | y |^{1 - \frac{2}{n}} I_{\frac{n}{2} - 1} ( |y| ) ,$$$$ where $$I_{\nu}$$ is the hyperbolic Bessel function of order $$\nu$$. This formula may be found in, e.g., Loss and Lieb, Analysis, Section 7.11. Therefore, $$$$\tag{3} \int_{S^{n-1}} \mathcal{N}_{\omega} ( \mu , \sigma^2 ) d \omega = \frac{e^{- (1 + \mu^2)/(2 \sigma^2) }}{\sigma} \left( \frac{|\mu|}{\sigma^2} \right)^{1 - \frac{2}{n}} I_{\frac{n}{2} - 1} \left( \frac{|\mu|}{\sigma^2} \right) .$$$$

Another special case that I mention in the comments below is $$\mu = 0$$. As I point out there, the $$n = 3$$ case of (1) for $$\mu = 0$$ should be similar to the calculation of the normalization constant of the Kent distribution. For $$n = 2$$, (2) reduces to $$$$\tag{4} \int_0^{2 \pi} e^{- \frac{1}{2} \left( \frac{\cos^2{\theta}}{\sigma_1^2} + \frac{\sin^2{\theta}}{\sigma_2^2} \right) } d \theta = 2 \pi e^{- \frac{1}{4} \left( \frac{1}{\sigma_1^2} + \frac{1}{\sigma_2^2} \right) } I_0 \left( \frac{1}{4} \left( \frac{1}{\sigma_2^2} - \frac{1}{\sigma_1^2} \right) \right) .$$$$ The previous result follows from the cosine and sine reduction formula $$\cos^2{x} = (1 + \cos{(2x)})/2$$ and likewise for $$\sin^2{x}$$.

Considering that (3) and (4) have relatively simple expressions, I am hoping this is also the case for (2) with $$\tilde{\mu} \neq 0$$.

• What is $dx$? You're probably not interested in integrating w.r.t. the Lebesgue measure in $\mathbb{R}^n$.
– user140541
Apr 20, 2022 at 14:58
• No. I just fixed the notation to make it clear I am talking about integration over the unit sphere with respect to the usual spherical measure. Apr 20, 2022 at 15:10
• Also, I should note that another special case is $\mu = 0$. Then this integral is essentially the normalization constant for the Kent distribution en.wikipedia.org/wiki/Kent_distribution. However, this constant only seems to be known for $n = 2,3$ with the $n=3$ case written out explicitly in the Wikipedia article. Apr 20, 2022 at 15:15

For the $$n = 2$$ case, (2) in the OP can be evaluated as follows. Denote $$\tilde{\mu} = \mu$$. The integral in question may be reduced to $$$$\tag{1} e^{- \frac{1}{2} \left( \frac{\mu_1^2}{\sigma_1^2} + \frac{\mu_2^2}{\sigma_2^2} \right) } \int_0^{2 \pi} e^{ - \frac{1}{2} \left( \frac{\cos^2{\theta}}{\sigma_1^2} + \frac{\sin^2{\theta}}{\sigma_2^2} - 2 \frac{\mu_1}{\sigma_1^2} \cos{\theta} - 2 \frac{\mu_2}{\sigma_2^2} \sin{\theta} \right) } d \theta .$$$$ Using the double angle reduction formulas $$\cos^2{\theta} = ( 1 + \cos{(2 \theta)} ) / 2$$ and $$\sin^2{\theta} = ( 1 - \cos{(2 \theta)} ) / 2$$, we can simplify (1) to read $$$$\tag{2} e^{- \frac{1}{4} \left( \frac{2\mu_1^2 + 1}{\sigma_1^2} + \frac{2\mu_2^2 + 1}{\sigma_2^2} \right) } \int_0^{2 \pi} e^{ - \frac{1}{2} \left( \left( \frac{1}{2 \sigma_1^2} - \frac{1}{2 \sigma_2^2} \right) \cos{(2 \theta)} - 2 \frac{\mu_1}{\sigma_1^2} \cos{\theta} - 2 \frac{\mu_2}{\sigma_2^2} \sin{\theta} \right) } d \theta .$$$$ By the Jacobi-Anger expansion, we have $$$$e^{z \cos{2 \theta}} = I_0 (z) + 2 \sum_{k \geq 1} I_k (z) \cos{(2 k\theta)} ,$$$$ where $$z = - \frac{1}{4} \left( \frac{1}{\sigma_1^2} - \frac{1}{\sigma_2^2} \right)$$. Plugging this into (2), we have $$$$\tag{3} e^{- \frac{1}{4} \left( \frac{2\mu_1^2 + 1}{\sigma_1^2} + \frac{2\mu_2^2 + 1}{\sigma_2^2} \right) } \int_0^{2\pi} \left( I_0 (z) + 2 \sum_{k \geq 1} I_k (z) \cos{(2 k \theta)} \right) e^{\frac{\mu_1}{\sigma_1^2} \cos{\theta} + \frac{\mu_2}{\sigma_2^2} \sin{\theta} } d \theta ,$$$$ where, again, $$z = - \frac{1}{4} \left( \frac{1}{\sigma_1^2} - \frac{1}{\sigma_2^2} \right)$$. To proceed, note that $$$$y_1 \cos{\theta} + y_2 \sin{\theta} = \| y \| \cos{(\theta - \arctan{(y_2/y_1)})} ,$$$$ where $$y_1 = \mu_1 / \sigma_1^2$$ and $$y_2 = \mu_2 / \sigma_2^2$$. Plugging this into (3), changing variables, and using the cosine addition formula, we arrive at $$$$\int_{S^1} e^{- \frac{1}{2} \| D^{-1} ( \omega - \mu ) \|^2} d \omega = 2 \pi e^{- \frac{1}{4} \left( \frac{2\mu_1^2 + 1}{\sigma_1^2} + \frac{2\mu_2^2 + 1}{\sigma_2^2} \right) } \left( I_0 (z) I_0 (\|y\|) + 2 \sum_{k \geq 1} I_k (z) I_{2k} (\| y \|) \cos{\left(2 k \arctan{\left( \frac{y_2}{y_1} \right)}\right)} \right) , \tag{*}$$$$ where, again, $$z = - \frac{1}{4} \left( \frac{1}{\sigma_1^2} - \frac{1}{\sigma_2^2} \right)$$, $$y_1 = \mu_1 / \sigma_1^2$$ and $$y_2 = \mu_2 / \sigma_2^2$$.
We note how ($$\ast$$) agrees with some of the special cases outlined above. Indeed, for $$\sigma_1 = \sigma_2$$, $$z = 0$$ and all higher order Bessel functions in the sum in ($$\ast$$) are zero. When $$\mu_1 = \mu_2 = 0$$, the $$\arctan$$ isn't defined, so this special case isn't covered.