在 python 中生成实心球体

Generating solid spheres in python

我想从 n 维实心球生成随机均匀样本。

我现在的方法是这样的

def sample_sphere(d, npoints):
    points = np.zeros((npoints, d))
    for i in range(npoints):
        r = np.random.rand() #random radius
        v = np.random.uniform(low= -1, high=1, size=d) #random direction
        v = v / np.linalg.norm(v)
        points[i] = r * v
    return points

但不出所料,这种方法在 0 附近给出了更高浓度的点,因为采样密度与体积不成正比:

如何统一采样?

球形分布有两种基本解释:

1。有界实心球体

在给定半径处得到一个点的概率由 shell 的体积给出,厚度为 dr 在该半径处:p(r) ~ r^D 直到一个常数。因此,径向 CDF 为 (r / R)^(D+1),其中 R 是外半径,D 是维数。球体必须有界,除非你想从所有 space.

select

给定 CDF 生成样本的标准方法是将其反转:

random_radius = R * np.pow(np.random.uniform(0, 1), 1 / D)

顺便说一句,选择正确统一方向的常用方法是使用高斯分布 (see here):

random_direction = np.random.normal(size=D)
random_direction /= np.linalg.norm(random_direction)

您可以而且应该对您的函数进行矢量化。不需要循环:

def sample_sphere(r, d, n):
    radii = r * np.pow(1 - np.random.uniform(0, 1, size=(n, 1)), 1 / d)
    directions = np.random.normal(size=(n, d))
    directions *= radii / np.linalg.norm(directions, axis=1, keepdims=True)
    return directions

注意 1 - np.random.uniform(...)。这是因为函数的范围是[0, 1),而我们想要(0, 1]。减法在不影响均匀性的情况下反转范围。

2。无限衰减分布

在这种情况下,您正在寻找类似具有球对称性的高斯分布。这种分布的径向 CDF 与 R^(D-1) 成比例。对于高斯分布,这在 2D 中称为 Rayleigh 分布,在 3D 中称为 Maxwell,在一般情况下称为 chi。这些可以分别通过 scipy.stats.rayleigh, scipy.stats.maxwell and scipy.stats.chi 对象直接获得。

在这种情况下,半径的意义不大,但可以解释为分布的标准偏差。

你函数的无限版本(也向量化),然后变成:

def sample_sphere(s, d, n):
    radii = scipy.stats.chi.rvs(df=d, scale=s, size=(n, 1))
    directions = np.random.normal(size=(n, d))
    directions *= radii / np.linalg.norm(directions, axis=1, keepdims=True)
    return directions