Let $\alpha>0$ be a positive parameter and consider the function

$$f(x) = \int_{\mathbb S^{n-1}} e^{-\alpha \left\lVert x-y \right\rVert^2} dS(y)$$ for $x \in \mathbb R^n.$ So, since this was asked, although we integrate over the unit sphere, the function "lives" on $\mathbb R^n.$

This function is clearly rotationally symmetric.

I would like to show that the global maxima are attained at one single radius $r$, only.

The rotational symmetry implies that we can consider it as a one-dimensional function by choosing $x=(x_1,0....,0)$, this way the exponent simplifies to $e^{-\alpha \left\lVert x-y \right\rVert^2}=e^{-\alpha (x_1-y_1)^2+1-y_1^2}.$

If anything is unclear about this question, then please let me know. I am happy to hear about any ideas how to approach this problem.

EDIT: Thanks to some interesting comments below, one can say that the global maximum is always attained at some radius $r \in [0,1]$ where for small $\alpha$ it seems to be attained close to zero and for large $\alpha$ it is attained closer to one.

The question remains however why is there only one radius at which the global maximum is attained? -In fact as George Lowther points out in the comments, for $\alpha \le n/2$ the unique maximum is attained at $r=0$ which leaves the case $\alpha >n/2$ when this does not hold true.

  • 4
    $\begingroup$ This function has no global minimum. It's positive everywhere and $\to 0$ at $\infty.$ Perhaps you meant "global maximum"? $\endgroup$ – zhw. Nov 1 '18 at 19:20
  • 2
    $\begingroup$ thank you, perhaps that is why nobody ever answered me :-) $\endgroup$ – Sascha Nov 1 '18 at 20:07
  • $\begingroup$ Is $\lVert x \rVert = r$ and is the radius of the $(n - 1)$-sphere equal to $1$? $\endgroup$ – Maxim Nov 3 '18 at 15:23
  • $\begingroup$ it is the unit sphere, so radius of the sphere is equal to $1$, but $x \in \mathbb R^n$. $\endgroup$ – Sascha Nov 3 '18 at 15:39
  • $\begingroup$ Then what is $r$? $\endgroup$ – Maxim Nov 3 '18 at 16:16

Let us write $\tilde{f}(x) = \frac{1}{|\mathbf{S}^{n-1}|}\int_{\mathbf{S}^{n-1}}e^{-\alpha|x-\omega|^2}dS(\omega)$. Since $\tilde{f}$ is radially symetric, we may assume $x = re$ where $r = |x|$ and $e = (1,0,\ldots,0)$. Then we may write

$$\tilde{f}(re) = \frac{1}{|\mathbf{S}^{n-1}|}\int_{\mathbf{S}^{n-1}}e^{-\alpha(r^2-2r\langle e,\omega\rangle+1)}dS(\omega) = \frac{e^{-\alpha(r^2+1)} }{|\mathbf{S}^{n-1}|}\int_{\mathbf{S}^{n-1}}e^{2\alpha r\langle e,\omega\rangle}dS(\omega).$$
Disregard the factor $e^{-\alpha}$ and define $f(r) = \frac{e^{-\alpha r^2} }{|\mathbf{S}^{n-1}|}\int_{\mathbf{S}^{n-1}}e^{2\alpha r\langle e,\omega\rangle}dS(\omega)$. We can observe that $X: \omega \mapsto \langle e,\omega\rangle\in [-1,1] $ is a random variable with probability density function $p_n(x) = c_n (1-x^2)^{\frac{n-3}{2}}$ on $[-1,1]. $($c_n$ is a positive constant such that $\int_{[-1,1]}p_n(x)dx =1$.)

Thus we have $f(r) = e^{-\alpha r^2}\int_{-1}^{1} e^{2\alpha r x}p_n(x)dx\;(*)$ and note that $M_X(t) = \int_{-1}^{1} e^{t x}p_n(x)dx$ is a moment-generating function(mgf) of $X$ , and $C_X(t) = \log M_X(t)$ is cumulant-generating function of $X$.

By taking logarithm on both sides of $(*)$, we maximize the following functional $$-\alpha r^2 + C_X(2\alpha r).$$

The first-order condition for this is:

$$r = C'_X(2\alpha r) = \frac{M'_X(2\alpha r)}{M_X(2\alpha r)}\quad \cdots (**).$$ Since $r=0$ is a trivial solution of $(**)$ we assume $r>0$.

What follows next is a probabilistic interpretation of this equation. Firstly, $$C'_X(t) = \frac{M'_X(t)}{M_X(t)} = E\left[X \frac{e^{tX}}{M_X(t)}\right] = \int x\cdot \frac{e^{tx}p(x)}{\int e^{ty}p(y)dy}dx$$ is an expected value of another random variable $X_t$ with the new density $\frac{e^{tx}p(x)}{\int e^{ty}p(y)dy}$. Thus $r = C'_X(2\alpha r) $ should be less than 1.

Secondly, $C''_X(t) = \frac{M''_X(t)}{M_X(t)}- \left(\frac{M'_X(t)}{M_X(t)}\right)^2$ in the similar way can be seen as variance of $X_t$ and thus is greater than $0$ (not zero since the distribution is not degenerate.) So $C'_X(t)$ is a strictly increasing function. Moreover, it should be that $C'_X(t) \uparrow 1$ as $t\to \infty$ since the distribution of $X_t$ converges to that of degenerate random variable $X\equiv 1$.

From the above argument, we see that there is unique $$\alpha= A(r) := \frac{(C'_X)^{-1}(r)}{2r}$$ for each $0<r<1$ satisfying $(**)$. The uniqueness of $r$'s corresponding to each $\alpha$ follows once if we show that $A(r)$ is strictly increasing on $(0,1)$. Since $C'_X$ is an increasing bijection from $(0,\infty)$ to $(0,1)$, we can substitute $r$ for $r = C'_X(u)$, and showing $1/A(r) = \frac{2C'_X(u)}{u}$ is decreasing in $u\in(0,\infty)$ is okay.

Our (almost) final claim is that $r \mapsto C'_X(r)/r$ is a decreasing function. Let us write by change of variable $u = (x+1)/2$, $$M_X(r) = e^{-r}\int_0^1 e^{2ru} \frac{u^{\frac{n-3}{2}}(1-u)^{\frac{n-3}{2}}}{B((n-1)/2,(n-1)/2)}du =: e^{-r}\Phi(2r).$$(Here, $B(x,y)$ is a beta function and $\Phi$ is a mgf of the beta distribution.)
We can calculate $\Phi(r)$ explicitly and omitting the procedure , we have

$$\Phi(r) =\sum_{k=0}^{\infty} \frac{a_k}{k!}r^k,\quad a_k := \prod_{j=0}^{k-1} \frac{\frac{n-1}{2}+j}{n-1+j}.$$ (https://en.wikipedia.org/wiki/Beta_distribution) From this, we have $$C'_X(r/2) = 2\left(\frac{\Phi'(r)}{\Phi(r)}-1/2\right). $$ Some more labor yields:

$$\Phi'(r) - \frac{1}{2}\Phi(r) = \sum_{k\geq 0} (a_{k+1}-a_k/2)\frac{r^k}{k!} = r \sum_{k=0}^{\infty}(a_{k+2}-a_{k+1}/2)\frac{r^k}{(k+1)!},$$ $$\frac{1}{4}\frac{C'_X(r/2)}{r/2} = \frac{\sum_{k=0}^{\infty}(a_{k+2}-a_{k+1}/2)\frac{r^k}{(k+1)!}}{\sum_{k=0}^{\infty} a_k\frac{r^k}{k!}}.$$

Using $\frac{a_{k+1}}{a_k} = \frac{1}{2}(1 + \frac{k}{n+k-1})$, we have $$\sum_{k=0}^{\infty}(a_{k+2}-a_{k+1}/2)\frac{r^k}{(k+1)!} = \sum_{k=0}^{\infty}\frac{a_{k+1}}{n+k}\frac{r^k}{k!}.$$ Finally, $$\frac{d}{dr}\left(\frac{1}{4}\frac{C'_X(r/2)}{r/2}\right) = \frac{\sum_{k=0}^{\infty}\frac{a_{k+2}}{n+k+1}\frac{r^k}{k!}\sum_{k=0}^{\infty} a_k\frac{r^k}{k!} -\sum_{k=0}^{\infty}\frac{a_{k+1}}{n+k}\frac{r^k}{k!}\sum_{k=0}^{\infty} a_{k+1}\frac{r^k}{k!} }{\left(\sum_{k=0}^{\infty} a_k\frac{r^k}{k!}\right)^2}.$$ Now the denominator is itself a power series whose $s$-th coefficient $c_s$is given by $$c_s = \sum_{j=0}^s \left(\frac{a_{j+2} a_{s-j}}{n+j+1}-\frac{a_{j+1} a_{s-j+1}}{n+j}\right)\frac{1}{j!(s-j)!}.$$ Using the recursion formula of $a_k$, $$\begin{multline} c_s= -\sum_{j=0}^s \frac{(s-2j)(n-1)+2(s-j)}{(n+j)(n+j+1)(n+2s-2j-1)}\frac{a_{j+1} a_{s-j+1}}{j!(s-j)!}\\ = -\frac{1}{2}\sum_{j=0}^s \left [\frac{(s-2j)(n-1)+2(s-j)}{(n+j)(n+j+1)(n+2s-2j-1)} +\frac{(2j-s)(n-1)+2j}{(n+2j-1)(n+s-j)(n+s-j+1)} \right]\frac{a_{j+1} a_{s-j+1}}{j!(s-j)!}. \end{multline}$$ And, $$\begin{multline}\frac{(s-2j)(n-1)+2(s-j)}{(n+j)(n+j+1)(n+2s-2j-1)} +\frac{(2j-s)(n-1)+2j}{(n+2j-1)(n+s-j)(n+s-j+1)} \\ =(n-1)(s-2j)\left[ \frac{1}{(n+j)(n+j+1)(n+2s-2j-1)} - \frac{1}{(n+2j-1)(n+s-j)(n+s-j+1)}\right] + \text{positive terms}\\ >\frac{(n-1)(s-2j)^2\{2j(s-j) +(n-1)(s-3) -4\}}{(n+j)(n+j+1)(n+2j-1)(n+s-j)(n+s-j+1)(n+2s-2j-1)} \geq 0, \end{multline}$$ if $s\geq 3, 1\leq j\leq s-1$.

In case $j=0$ or $j=s$, it is quite easy to check that $$\frac{s}{n(n+2s-1)} - \frac{s}{(n+s)(n+s+1)}\geq 0.$$ So, negativity of $c_s, s\geq 3$ is established. Some more labors yield: $$c_0 = 0, c_1 = c_2 = -\frac{1}{4n^2(n+2)}.$$

Since all $c_s \leq 0$, the denominator is less than 0, and the claim $\frac{d}{dr}\left(\frac{C'_X(r)}{r}\right)<0$ is proved.
Finally, observe that $\lim_{r\to 0^+}A(r) = \lim_{u\to 0^+}\frac{u}{2C'_X(u)} = \frac{1}{2C''_X(0)} = \frac{n}{2}.$ This is because $C''_X(0)$ is a variance of $X=\langle e,\omega\rangle=\omega_1$, and this value is equal to $$\frac{1}{|\mathbf{S}^{n-1}|}\int_{\mathbf{S}^{n-1}}\omega_1^2 dS(\omega)=\frac{1}{|\mathbf{S}^{n-1}|}\int_{\mathbf{S}^{n-1}}\frac{\omega_1^2 +\cdots \omega_n^2}{n} dS(\omega) = \frac{1}{n}.$$ Hence, $A((0,1)) = (\frac{n}{2}, \infty)$ and we conclude that:
(1) if $0<\alpha \leq \frac{n}{2}$, then the unique global maximizer is $r=0$.
(2) if $\alpha >\frac{n}{2}$, then the unique global maximizer is $r = A^{-1}(\alpha)>0$. Here, $0$ is not a maximizer even though it satisfies the FOC, since the second-order condition is not satisfied: $$(-\alpha r^2 + C_X(2\alpha r))''\Big|_{r=0}=2\alpha(\frac{2\alpha}{n}-1) >0. $$

  • $\begingroup$ If the dimension of $S^{n - 1}$ is $n - 1$, the exponent in $p_n(x)$ is off by $1/2$ ($p_2$ cannot be a constant). $\endgroup$ – Maxim Nov 6 '18 at 15:52
  • $\begingroup$ Thank you for your valuable comment. You're right. I should have been more careful. $\endgroup$ – Song Nov 7 '18 at 1:52

If $\mathcal S$ is the unit $(n - 1)$-sphere around the origin and $\lVert x \rVert = r$, then choosing hyperspherical coordinates s.t. $x_1 = r, \,y_1 = \cos \theta$ gives $$f(r) = \int_{\mathcal S} e^{-\alpha \lVert x - y \rVert^2} dS(y) = C_n \int_0^\pi e^{-\alpha (1 + r^2 - 2 r \cos \theta)} \sin^{n - 2} \theta \,d\theta = \\ C_n \sqrt \pi \,\Gamma {\left( \frac {n - 1} 2 \right)} (\alpha r)^{-n/2 + 1} e^{-\alpha (1 + r^2)} I_{n/2 - 1}(2 \alpha r), \\ f'(r) = 2 C_n \sqrt \pi \,\Gamma {\left( \frac {n - 1} 2 \right)} \alpha (\alpha r)^{-n/2 + 1} e^{-\alpha (1 + r^2)} (I_{n/2}(2 \alpha r) - r I_{n/2 - 1}(2 \alpha r)).$$ The goal is to prove that the ratio $\omega(r) = r I_{n/2 - 1}(2 \alpha r)/I_{n/2}(2 \alpha r)$ is monotonous in $r$. Since $$\omega'(r) = 2 \alpha r \left( 1 - \frac {I_{n/2 - 1}(2 \alpha r) I_{n/2 + 1}(2 \alpha r)} {I_{n/2}^2(2 \alpha r)} \right),$$ we need to prove $$I_\nu^2(z) > I_{\nu - 1}(z) I_{\nu + 1}(z), \quad\nu > 0, \,z > 0.$$ Multiplying the series expansions yields $$I_\nu^2(z) = \sum_{k \geq 0} \frac {\Gamma {\left( k + \nu + \frac 1 2 \right)}} {\sqrt \pi \,\Gamma(k + 1) \Gamma(k + \nu + 1) \Gamma(k + 2 \nu + 1)} z^{2 (k + \nu)}, \\ I_{\nu - 1}(z) I_{\nu + 1}(z) = \sum_{k \geq 0} \frac {\Gamma \!\left( k + \nu + \frac 1 2 \right)} {\sqrt \pi \,\Gamma(k + 1) \Gamma(k + \nu + 1) \Gamma(k + 2 \nu + 1)} \frac {k + \nu} {k + \nu + 1} z^{2 (k + \nu)}, \\ \frac {[z^{2 (k + \nu)}] I_\nu^2(z)} {[z^{2 (k + \nu)}] I_{\nu - 1}(z) I_{\nu + 1}(z)} = 1 + \frac 1 {k + \nu},$$ proving the inequality. Since $\omega(\infty) = \infty$, $\omega$ increases to infinity monotonously.

Finally, $\omega(0) = n/(2 \alpha)$. If $\omega(0) < 1$, there is a unique $r_0 > 0$ where $\omega(r_0) = 1$ and $f'(r_0) = 0$, therefore $f(r)$ increases on $[0, r_0]$ and decreases on $[r_0, \infty)$. If $\omega(0) \geq 1$, $f(r)$ decreases on $[0, \infty)$.

  • $\begingroup$ sorry, what precisely is this constant $C_n$? $\endgroup$ – Sascha Nov 5 '18 at 17:51
  • 1
    $\begingroup$ It's the integral of the area element over the remaining $n - 2$ hyperspherical angles, thus the surface area of a unit $(n - 2)$-sphere. $\endgroup$ – Maxim Nov 5 '18 at 18:12
  • $\begingroup$ thank you very much. $\endgroup$ – Sascha Nov 5 '18 at 18:13

This is certainly not a full solution.

First we would like to get an expression for $f(x)$. (Below I write $f(x)$ or $f_n(x)$ when I should probably should have written $f((x_1, 0,0,\ldots, 0)$.)

I am going to use the function $$ s(n,r) := \frac{2 \pi ^{\frac{n+1}{2}}}{\Gamma \left(\frac{n+1}{2}\right)} r^n, $$ the surface area of an n-dimensional sphere.

We can use $s$ to integrate $f$,

$$ f_n(x) = \int_{\theta=0}^\pi s (n-1 ,\sin(\theta)) \exp\left(-a \left((x-\cos (\theta))^2+\sin ^2(\theta)\right)\right)d\theta. $$

According to Mathematica, $\;f_2(x) = \frac{\pi e^{-a (x+1)^2} \left(e^{4 a x}-1\right)}{a x}$.

And $f_4(x) = \frac{\pi ^2 e^{-a (x+1)^2} \left(2 a x+e^{4 a x} (2 a x-1)+1\right)}{2 a^3 x^3}$.

For odd values of $n$, Mathematica gives Hypergeometric functions. $$f_3(x) = 2 \pi ^2 e^{-a \left(x^2+1\right)} \, _0\tilde{F}_1\left(;2;a^2 x^2\right)$$

In the Mathematica language, $f_3(x)$ is

"(2*Pi^2*Hypergeometric0F1Regularized[2, a^2*x^2])/E^(a*(1 + x^2))".

If you try to find the maximum of $f$, it seems that the maximum occurs at $x=0$ when $a<2$. One way to think about this is releasing a million drunkards randomly on the $S^{n-1}$ sphere and having them walk randomly. Initially, the most drunkards per unit volume occurs on the $S^{n-1}$ sphere, but as time goes on the radius of maximum drunkards moves toward the center.

Another way to look at it is as a convolution of the sphere with a Gaussian with $\sigma = \sqrt{1/(2 a)}$. When sigma is large, $f_n(x)$ is approximately Gaussian.

I was not able to find closed form solutions for $\mathrm{argmax}_x\, f_n(x)$. If we take the derivative when $n$ is even, set it equal to zero, and simplify, we get expressions of the form $$p_1(x,a) + \cosh(2 a x)p_2(x,a) +\sinh(2 a x)p_2(x,a) =0$$ where the $p_i$ are polynomials of degree of low degree with "nice" coefficients.

For example, for $n=10$, Mathematica says that the derivative of $f_{10}$ is zero if $$32 a^5 x^6 \sinh (2 a x)-32 a^5 x^5 \cosh (2 a x)-160 a^4 x^5 \cosh (2 a x)+240 a^4 x^4 \sinh (2 a x)+360 a^3 x^4 \sinh (2 a x)-840 a^3 x^3 \cosh (2 a x)-420 a^2 x^3 \cosh (2 a x)+1680 a^2 x^2 \sinh (2 a x)+210 a x^2 \sinh (2 a x)+945 \sinh (2 a x)-1890 a x \cosh (2 a x)=0.$$

Hope that helps.

  • $\begingroup$ yes, I feel that this does not tell us much about the location and number of global maxima so far. The hypergeometric function is still a somewhat opaque object. $\endgroup$ – Sascha Nov 3 '18 at 18:24

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.