Python 4维球体上点的均匀分布

Posted 2023-02-16

【中文标题】Python 4维球体上点的均匀分布

【英文标题】：Python Uniform distribution of points on 4 dimensional sphere

【发布时间】：2013-03-30 14:38:00

【问题描述】：

我需要在 4 维球体上均匀分布点。我知道这不像选择 3 个角度和使用极坐标那么简单。

我使用 3 个维度

from random import random

u=random()
costheta = 2*u -1 #for distribution between -1 and 1
theta = acos(costheta)
phi = 2*pi*random

x=costheta
y=sin(theta)*cos(phi)
x=sin(theta)*sin(phi)

这给出了 x、y 和 z 的均匀分布。

如何获得 4 个维度的相似分布？

【问题讨论】：

如何在 N 球上随机生成均匀分布的点：en.wikipedia.org/wiki/…

等等，你想让这些点在一个球体上，但均匀分布在 x,y,z,（第 4 维）？这对我来说不算什么。我不认为均匀分布在球体上的点会映射到均匀分布在 4 空间中。

@SchighSchagh 所以你不能在 4 维中运行蒙特卡罗模拟？

@SameerPatel 这与蒙特卡洛或任何其他采样方法没有任何关系。这里有两个不同的空间，（一个是 R^4，另一个是 4 球体的表面），我们需要知道您希望对哪个空间进行均匀随机分布。

【参考方案1】：

A standard way，不过，或许not the fastest，是用Muller的方法在N-sphere上生成均匀分布的点：

import numpy as np
import matplotlib.pyplot as plt
import mpl_toolkits.mplot3d.axes3d as axes3d

N = 600
dim = 3

norm = np.random.normal
normal_deviates = norm(size=(dim, N))

radius = np.sqrt((normal_deviates**2).sum(axis=0))
points = normal_deviates/radius

fig, ax = plt.subplots(subplot_kw=dict(projection='3d'))
ax.scatter(*points)
ax.set_aspect('equal')
plt.show()

只需将 dim = 3 更改为 dim = 4 即可在 4 球体上生成点。

【讨论】：

采用均匀分布将产生一个球体，其采样点的密度朝向它将占据的立方体的角落更高。如何证明使用高斯分布也不会导致这个问题？如果我在脑海中想象它似乎是合理的，但它是真的吗，为什么？在指定大小时，我假设 numpy 的高斯分布对于每个坐标都是独立的。

那么除以零呢？

@GuillaumeChevalier：我承认我不明白细节，但也许你会找到this proof outline useful。关于被零除的问题：我相信除数为 0 的概率本身就是 0。在实际层面上，NumPy 通过发出 warning 并返回 NaN（和 Matplotlib跳过等于 NaN 的点）。

这是一个很好的证明！谢谢。

我们如何将points 写入txt 文件？我使用了np.savetxt('sphere.csv', points, delimiter=',')，但它将它们保存为三个大行而不是 3 列

【参考方案2】：

在4D空间中取一个坐标正态分布的点，计算其单位向量。这将在单位 4 球体上。

from random import random
import math
x=random.normalvariate(0,1)
y=random.normalvariate(0,1)
z=random.normalvariate(0,1)
w=random.normalvariate(0,1)
r=math.sqrt(x*x + y*y + z*z + w*w)
x/=r
y/=r
z/=r
w/=r
print (x,y,z,w)

【讨论】：

当然，这会在 4-sphere 上生成一个随机点，但是分布是否均匀？

@SchighSchagh：均匀分布？是的。

x,y,z,w 最初相对于 R^4 是均匀随机的，但随后它们经历了非线性变换，我仍然不清楚 OP 是否想要均匀随机相对于球体表面或相对于 R^4。编辑：你能指定关于你随机统一声明的内容，并证明它吗？

@SchighSchagh：糟糕，我忘了正常分发它们。见enwp.org/…。

除以零怎么样？

【参考方案3】：

我喜欢@unutbu 的回答，如果高斯采样真的创建了一个均匀间隔的球形分布（不像从立方体采样），但是为了避免在高斯分布上采样并且必须证明这一点，有一个简单的解决方案： strong>在球体上均匀分布采样（不是在立方体上）。

在均匀分布上生成点。 计算每个点的平方半径（避免平方根）。 丢弃积分： 丢弃平方半径大于 1 的点（因此，非平方半径大于 1）。 丢弃太接近零半径的点，以避免在下一步中与除法相关的数值不稳定。 对于保留的每个采样点，将采样点除以范数，以便将其重新归一化为单位半径。 清洗并重复因为丢弃的样品获得更多分数。

这显然适用于 n 维空间，因为半径在高维中始终是 L2 范数。

它的速度很快，可以避免在高斯分布上进行平方根和采样，但它不是矢量化算法。

【参考方案4】：

我找到了一个很好的 solution 用于从 N-dim 球体中采样。主要思想是：

如果 Y 是从不相关的多元正态分布中得出的，则 S = Y / ||Y|| 在单位 d 球上具有均匀分布。将 S 乘以 U^1/d，其中 U 在单位间隔 (0,1) 上具有均匀分布，从而在单位 d 维球中创建均匀分布。

这是执行此操作的 python 代码：

Y = np.random.multivariate_normal(mean=[0], cov=np.eye(1,1), size=(n_dims, n_samples))
Y = np.squeeze(Y, -1)
Y /= np.sqrt(np.sum(Y * sample_isotropic, axis=0))
U = np.random.uniform(low=0, high=1, size=(n_samples)) ** (1/n_dims)
Y *= distr * radius # in my case radius is one

这是我得到的球体：

【讨论】：

变量distr和sample_isotropic是什么？您能否修改代码示例以使其正常工作？