Anyone who has studied pseudo random numbers should know that the following will not generate points uniformly on the sphere. Here is a demonstration of the fact:

$ ipython -wthread

In [1]: from numpy import *

In [2]: from enthought.mayavi import mlab

In [3]: p = random.random(10000)*2*pi

In [4]: t = random.random(10000)*pi

In [5]: x, y, z = sin(t)*cos(p), sin(t)*sin(p), cos(t)

In [6]: g = mlab.points3d(x, y, z, z, mode='point')

In [7]: g.actor.property.point_size = 2

The point size is made larger so you can see the points more clearly. Here is the picture it produces:

Notice the clustering at the poles. OTOH, the following works nicely,

In [3]: p = random.random(10000)*2*pi

In [10]: mlab.clf()

In [12]: t = arccos(random.random(10000)*2.0 - 1.0)

In [13]: x, y, z = sin(t)*cos(p), sin(t)*sin(p), cos(t)

In [14]: g1 = mlab.points3d(x, y, z, z, mode='point')

In [15]: g1.actor.property.point_size = 2

This is the picture this one produces:

As you can see the transformation has worked and generates what appears to be a uniform distribution on the surface of a sphere. Visual testing isn't good enough uniformness but that isn't the point here.

## 2 comments:

Hi,

another way is to sample Gaussian vectors of dimension d, with mean 0_d, and covariance matrix I_d, and to normalise them:

x = randn(n)

u = x/sqrt(sum(x))

Code is more compact, and looks faster, but I have tried to compare both approaches.

Nicolas CHOPIN, ENSAE

oops, sorry, code is not correct:

x = randn(d,n)

u = x/sqrt(sum(x**2,0))

if d==2:

plot(u[0,:],u[1,:])

Post a Comment