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