naver/roma

samping unit quaternions

ankurhanda opened this issue · 6 comments

The code to sample random unit quaternions in the repo

roma/roma/utils.py

Lines 57 to 58 in 34bd744

quat = torch.randn(size + (4,), dtype=dtype, device=device)
quat /= torch.norm(quat, dim=-1, keepdim=True)
isn't correct.

The correct way to sample is mentioned in Algorithm 2 in Kuffner at al. http://pajarito.materials.cmu.edu/documents/kuffner_james_2004_1.pdf

It might be worth adding that in. PyQuaternion implements it in the correct way: https://github.com/KieranWynn/pyquaternion/blob/master/pyquaternion/quaternion.py#L261

Hi,
Thank you for reviewing the code. Would you mind explaining me why it is incorrect, in your opinion?

The code samples 4D vectors from a multivariate normal distribution, which is isotropic.
Normalizing these vectors should lead to a uniform sampling of the unit quaternion hypersphere.

To complete my answer:
I am aware of the paper of J. Kuffner. Using the algorithm he suggested would allow to rely on simpler primitives like 'rand' instead of 'randn', but it would be at the cost of an increased complexity of the RoMa code.
I don't think that it would be worth it, but if you have strong arguments in favor of this algorithm, I would be interested to know them.

Thanks for your message. Do you know where I can read more about randn based unit quaternion being uniform on SO(3)? I guess one way to verify that would be to generate the quats via randn and then multiply a point [1, 0, 0] with the rotation and see if the new points mapped via this rotation are uniformly distributed on the sphere.

I understand that kuffner et al. might be slow due to sin/cos.

Also, I checked, eigen also implements it via the Kuffner et al. way
https://eigen.tuxfamily.org/dox/classEigen_1_1Quaternion.html#a7da87cda5567ff1e860782d638643868. I think you could perhaps add a flag to specify which way to generate the rotations then in roma?

I don't have explicit references in mind, but this is a classical technique to sample points uniformly on hyperspheres
(e.g. mentioned here: https://math.stackexchange.com/questions/1585975/how-to-generate-random-points-on-a-sphere).

The result is rather intuitive from a theoretical standpoint: a multivariate normal distribution is spherically symmetric ( $\exp(-x^2) \exp(-y^2) \exp(-z^2) \exp(-w^2) = \exp(-|q|^2)$, with $q=xi+yj+zj+w$), thus the distribution has to be spherically symmetric after normalization, and therefore to be uniform for the usual L2 metric.
This is just a theoretical argument however, and assessing numerically the quality of a random generator is quite challenging (for me at least).

One potential caveat of the current implementation is that one should in theory handle the case $|q|=0$. I don't think it is a big issue in practice, I never encountered a failure of the assertion

assert(torch.all( torch.abs(torch.norm(quat, dim=-1) - 1) < 1e-3 ))
.
Switching to a different generator as you suggest could remove this potential point of failure, I will try to have a look.

I changed the generation algorithm in v1.3.0.
Thank you for the feedback, and feel free to reopen the issue if needed.

Thank you!! I must say roma is a great library and I have been using it a lot for my work so I appreciate it.