How to characterize rotations in $\mathbb{R}^n$?

Question

I am studying the performance of an optimizer algorithm to find the $$ \textrm{argmin}_{x\in \mathbb{R}^n} f(x) \text{ where } f : \mathbb{R}^n \rightarrow \mathbb{R} $$

I would like to test how the performance changes when I rotate the coordinate system. I need to do this in a systematic manner so that each rotation has "same probability".

In the $n=2$ case this is relatively easy, since any 2-dimensional rotational matrix can be written as:

\begin{pmatrix} cos(\phi) & sin(-\phi) \\ sin(\phi) & cos(\phi) \end{pmatrix}

The rotation is characterized by a single number $\phi$ so I can sample the $\phi$ uniformly (equidistantly) on $[-\pi; \pi]$.

In the $n=3$ case things get more difficult. There are 3 axes by which I can rotate, but the composition of these rotations is not commutative. So I could get biased results when I sample these 3 rotations uniformly.

I was thinking about using the Euler's rotation theorem:

in 3D space, any two Cartesian coordinate systems with a common origin are related by a rotation about some fixed axis

This means that the rotations in 3d are characterized by an axis and an angle.

So, I could sample all possible axes of some fixed length and all possible rotations around such axis. But this approach is very nontrivial to generalize to arbitrary dimensions.

How can I sample n-dimensional rotations uniformly?

Answer 1

The space of all rotations in $\Bbb R^n$ is the Lie group $SO(n)$, and it may help to read the wikipedia article on them. They are usually represented as matrices, and can also be viewed as an ordered set of basis vectors that are mutually orthogonal. In view of this, it is easy to generate such a set: roughly speaking, you apply Gram-Schmidt to a random set of vectors.

To be more precise: generate a random unit vector $v_1\in\Bbb R^n$. The easiest way to do this is to take advantage of the rotation invariance of the normal distribution in $n$ dimensions: if $x_1,\dots,x_n\sim{\cal N}(0,1)$, then $u_1=(x_1,\dots,x_n)$ is distributed in a rotationally invariant way, and $v_1=\frac {u_1}{|u_1|}$ is a uniformly chosen unit vector.

Next, we choose a random vector orthogonal to $v_1$. You can do this by generating another random normal vector $u_2$, then subtract the part in the direction of $v_1$, and normalize to get $v_2$.

Continue in this way, orthogonalizing random vectors with respect to all previous vectors, to get $v_1,\dots,v_{n-1}$. However, $v_n$ requires a bit of special handling. There are only two possible unit vectors that are orthogonal to all previous vectors, since there is only one dimension remaining. If we want a rotation, which is to say, not a reflection, then our choice should be $v_n=v_1\times v_2\times\dots\times v_{n-1}$, which sets the orientation correctly and gives us an element of $SO(n)$. If reflections don't matter, then the same method as before will give you an element of $O(n)$, the orthogonal group.

As for computational complexity, this requires $n^2$ samples from a normal distribution, as well as roughly $n^2/2$ orthogonalizations of dimension $n$ vectors, so it is complexity $O(n^3)$.

Answer 2

You could do it inductively as follows:

1) Pick a unit vector $u$ uniformly at random in $\mathbb{R}^n$. This vector will be the image of $e_1$ under the rotation.

2) Use Gram-Schmidt on the basis $u,e_2,...,e_n$ to get a new orthonormal basis $u,w_1,...,w_{n-1}$. Change the sign on $w_1$ if necessary so that you have an oriented basis, ie, so that the determinant of the matrix $[u, w_1, ..., w_{n-1}]$ is positive. You don't actually have to calculate this determinant- because of the way Gram-Schmidt works the determinant equals that of $[u, e_2, ..., e_n]$, which is the sign of the first component of $u$. (This process will fail if the first component of $u$ equals $0$, but that is a probability zero event). Let $U$ denote the matrix $[u, w_1, w_2, ..., w_{n-1}]$.

3) Pick a rotation $A$ of $\mathbb{R}^{n-1}$ at random. Then let $\bar{A}$ be the $n \times n$ matrix whose top left entry is $1$, whose bottom right square equals $A$, and the rest of whose entries are zero. Schematically:

$$\bar{A} = \left(\begin{matrix} 1 & 0 \\ 0 & A \end{matrix}\right)$$

Then $\bar{A}U$ is a random rotation of $\mathbb{R}^n$.

By the way, I don't know much about computational complexity so if this is way too computationally inefficient, I apologize.