How to compute the angle between two vectors expressed in the spherical coordinates?

Here is a recursive idea that might work. Let $s_1$, $s_2$,..., $s_{n-1}$ the sines of the angles of the first vector on $S^{n-1}$, $c_1$, ..., the cosines, and $s'_{\nu}$ and $c'_{\nu}$ those for the second one. With $p$ being the scalar product, we have $p=1$ for $n=1$, $p=c c'+s s'$ for $n=2$ and generically $$p(n)=c_1 c'_1 + s_1 s'_1 p(n-1)$$ where $p(n-1)$ is the scalar product of the vectors given by the angles $[\theta_2,...]$ and $[\theta'_2,...]$ on $S^{n-1}$(i.e. by removing the first angle of each sequence).

This, of course, still requires all sines and cosines, but my idea is to simplify $p(n-2)$ using the haversine formula and hope that this reduces the number of calls to trigoniometric functions.

Simple, alternative derivation of Ralf's formula, without the guesswork:

Angle between the vectors: $$\alpha = \arccos\left({{u \cdot v}\over{|u| |v|}}\right)$$ $$u \cdot v = \sum_i u_i v_i $$ in hyper-spherical coordinates (n+1 dimensions, hence n angles): $$u_i(n) = |u| \cos(\theta_i) \prod_{j=1}^{i-1}{\sin(\theta_j)} $$ except when i=n: $$u_n(n) = |u| \prod_{j=1}^{n}{\sin(\theta_j)} $$ and similar for v (I will use $\phi$ for the angles of v).

With unit vectors we combine the three formulae to: $$\cos(\alpha) = \sum_{i=1}^{n-1}{\cos(\theta_i)\cos(\phi_i) \prod_{j=1}^{i-1}{\sin(\theta_j)\sin(\phi_j)}} + \prod_{i=1}^{n}{\sin(\theta_j)\sin(\phi_j)} $$ which can be written as a recursive rule: $$\cos(\alpha)_n = \cos(\theta_{n-1})\cos(\phi_{n-1}) \cos(\alpha)_{n-1} + \prod_{j=1}^{n}{\sin(\theta_j)\sin(\phi_j)} $$

Which is similar to what Ralf's answer contains (cos/sin switched).