I would like to choose arbitrarily two vectors $a$ and $b$ $\in \mathbb{S}^{n-1}$ and I would like to calculate the probability they are at least some distance $\delta$ apart. This probability should be $$\frac{|\mathbb{S}^{n-1}|-|\{ x \in \mathbb{S}^{n-1}: \|x-a\| \leq \delta \}|}{|\mathbb{S}^{n-1}|}$$ So (I think), my problem consists of finding $|\{ x \in \mathbb{S}^{n-1}: \|x-a\| \leq \delta \}|$ ($|\cdot|$ being the surface measure in $\mathbb{R}^n$). I suppose $a$ can be taken as the unit vector $(1, 0, \cdots, 0)$. By symmetry, the points $\{ x \in \mathbb{S}^{n-1}: \|x-a\| = \delta \}$ should be a "circle" with center somewhere on the line between $(0, \cdots, 0)$ and $a$. I need to find that center and the radius. I search for a point on the "circle" with only the first two components being non-zero: The center should have coordinates $(1-x, 0, \cdots, 0)$ and the point on the circle has coordinates $(1-x, y, 0, \cdots, 0)$. $x$ and $y$ need to verify $$x^2 + y^2 = \delta$$ and $$(1-x)^2 + y^2 = 1$$ This gives $x = \delta / 2$ and $y = (\sqrt{3}/4) \delta$. Going back to the surface we want to calculate: $$|\{ x \in \mathbb{S}^{n-1}: \|x-a\| \leq \delta \}| = |\{x \in \mathbb{S}^{n-1}\}\cap \{B((1-\delta / 2, 0, \cdots, 0), r = (\sqrt{3}/4) \delta)\}$$ This should just be the measure of the spherical cap of the sphere centered at $(1-\delta / 2, 0, \cdots, 0)$. To estimate this surface area, I wanted to use that surface area is the push-forward measure of a Gaussian measure on $\mathbb{R}^n$ under $x \rightarrow x/\|x\|$ (times the measure of $\mathbb{S}^{n-1}$ since the latter is a probability measure). So to estimate the above measure, we can also calculate the Gaussian measure of a cone (rays going from $0$ to points in $\{x \in \mathbb{S}^{n-1}: \|x-a\| \leq \delta\})$. Parametrizing this cone, I get $(t, \text{"disc of radius" }t*\frac{\sqrt{3}\delta}{4(1-\delta/2)}) = (t, disc)$ $$(2\pi)^{-n/2}\int_{(t, disc)}e^{-\|x\|^2/2} = (2\pi)^{-n/2}\int_0^{\infty}e^{-t^2/2}\int_{(disc)}e^{(x_2^2 + x_3^2 + \cdots + x_n^2)/2}dxdt$$ The last integral is over a radial function so $$ = (2\pi)^{-n/2}\int_0^{\infty}e^{-t^2/2}\omega_{n-2}\int_0^{t*\frac{\sqrt{3}\delta}{4(1-\delta/2)} }e^{-s^2/2}dsdt$$ We can take $\delta$ as very small, so $$= (2\pi)^{-n/2}\int_0^{\infty}e^{-t^2/2}\omega_{n-2}t*\frac{\sqrt{3}\delta}{4(1-\delta/2)}dt \leq (2\pi)^{-n/2}\omega_{n-2}\sqrt{3}\delta$$ EDIT: In the end, we have to multiply with $|\mathbb{S}^{n-1}|$ since the above is only the uniform probability on $\mathbb{S}^{n-1}$
Are these calculations correct? I am a bit confused about the term $(2\pi)^{-n/2}$, then $\delta$ could be quite large and we would still not occupy a large amount of the mass.
I don't think this has anything to do with the Gaussian integral.
Consider the case $n=3$ first. We want to compute the area of the spherical cap $C\subset S^2$ comprising all points $x\in S^2$ having distance $|x-e_3|\leq d$ from the point $e_3\in S^2$. To this end we transform $d$ into the corresponding spherical distance $\alpha$ related to $d$ by the formula $d=2\sin{\alpha\over2}$. The cap in question, call it $C_\alpha$, now consists of all points $x\in S^2$ satisfying $0\leq\theta\leq\alpha$, whereby $\theta$ is the complement of the geographical latitude on $S^2$, i.e., $\theta=0$ at the north pole $e_3$. Planes $z={\rm const.}$ partition $C_\alpha$ into "infinitesimal lampshades" of width ${\rm d}\theta$ and circumference $\omega_1\sin\theta$, whereby $\omega_1=2\pi$ is the one-dimensional measure of $S^1$. The area of such an infinitesimal lampshade is then given by $2\pi\sin\theta\>d\theta$, and the area of $C_\alpha$ computes to $$2\pi\int_0^\alpha\sin\theta\>d\theta=2\pi(1-\cos\alpha)=4\pi\,\sin^2{\alpha\over2}\ ,\tag{1}$$ hence is $\approx \pi d^2$ when $d\ll1$, as expected.
After this exercise we are ready for the general case. We want to compute the area of the spherical cap $C_\alpha\subset S^{n-1}$ comprising all points $x\in S^{n-1}$ having spherical distance $\leq \alpha$ from the point $e_n\in S^{n-1}$. Denote the $r$-dimensional area of the complete $r$-dimensional unit sphere by $\omega_r$ $(r\geq0)$. Then we obtain instead of $(1)$ the formula $${\rm area}_{n-1}(C_\alpha)=\omega_{n-2}\int_0^\alpha\sin^{n-2}\theta\>d\theta\ .$$ The probability you are interested in is then given by $$p(n-1,\alpha)={\omega_{n-2}\over\omega_{n-1}}\int_0^\alpha\sin^{n-2}\theta\>d\theta\ .$$ Depending on your needs you can compute the integral in elementary terms or in case $\alpha\ll1$ obtain good approximations.