python之在多维环中均匀采样，无排斥

这个question中的算法告诉我们如何有效地从多维球中采样。有没有一种方法可以从多维环中同样有效地采样，即有 r1<r<r2

我希望那个缩放函数的修改不是太复杂 r*(gammainc(s2/2,n/2).^(1/n))./sqrt(s2)是可能的。 (平庸免责声明:甚至还没有想出原始缩放函数的代数/几何)。

原始matlab代码复制粘贴:

function X = randsphere(m,n,r) 
 
% This function returns an m by n array, X, in which  
% each of the m rows has the n Cartesian coordinates  
% of a random point uniformly-distributed over the  
% interior of an n-dimensional hypersphere with  
% radius r and center at the origin.  The function  
% 'randn' is initially used to generate m sets of n  
% random variables with independent multivariate  
% normal distribution, with mean 0 and variance 1. 
% Then the incomplete gamma function, 'gammainc',  
% is used to map these points radially to fit in the  
% hypersphere of finite radius r with a uniform % spatial distribution. 
% Roger Stafford - 12/23/05 
 
X = randn(m,n); 
s2 = sum(X.^2,2); 
X = X.*repmat(r*(gammainc(s2/2,n/2).^(1/n))./sqrt(s2),1,n); 

等效的 python 代码和来自 Daniel's answer 的演示:

import numpy as np 
from scipy.special import gammainc 
from matplotlib import pyplot as plt 
def sample(center,radius,n_per_sphere): 
    r = radius 
    ndim = center.size 
    x = np.random.normal(size=(n_per_sphere, ndim)) 
    ssq = np.sum(x**2,axis=1) 
    fr = r*gammainc(ndim/2,ssq/2)**(1/ndim)/np.sqrt(ssq) 
    frtiled = np.tile(fr.reshape(n_per_sphere,1),(1,ndim)) 
    p = center + np.multiply(x,frtiled) 
    return p 
 
fig1 = plt.figure(1) 
ax1 = fig1.gca() 
center = np.array([0,0]) 
radius = 1 
p = sample(center,radius,10000) 
ax1.scatter(p[:,0],p[:,1],s=0.5) 
ax1.add_artist(plt.Circle(center,radius,fill=False,color='0.5')) 
ax1.set_xlim(-1.5,1.5) 
ax1.set_ylim(-1.5,1.5) 
ax1.set_aspect('equal') 

请您参考如下方法:

最后一个方法here ⁽¹⁾适用于任意维度球体:

在球体上随机选取一个点:
- 生成 N 个高斯随机变量 x1,x2..xN
- 获取 x[i] 的范数

 L = Sqrt(x1*x1 + x2*x2 + .. + xn*xn) 
 ux1 = x1 / L 
 ux2 = x2 / L 
 ... 

那么向量 ux[i] 在表面 S^N-1

在环中提供均匀分布:
- 在范围内生成统一的随机数

R_NPow = RandomUniform(R_Inner^N, R_Outer^N)

并获取半径(如 this 2D case )

R = R_NPow^1/N

然后计算结果点坐标:

res_x1 = R * ux1 
res_x2 = R * ux2 
... 
res_xn = R * uxn 

(1) Muller, M. E. “关于在 - 维球体上均匀生成点的方法的注释。”通信。副教授。电脑。马赫。 2, 19-20, 1959 年 4 月。