从随机样本（python）构建一个近似均匀的网格

Posted 2023-03-12

【中文标题】从随机样本（python）构建一个近似均匀的网格

【英文标题】：Build an approximately uniform grid from random sample (python)

【发布时间】：2013-11-09 19:31:17

【问题描述】：

我想从采样数据构建一个网格。我可以使用机器学习 - 聚类算法，例如 k-means，但我想限制中心大致均匀分布。

我想出了一种使用 scikit-learn 最近邻搜索的方法：随机选择一个点，删除半径 r 内的所有点，然后重复。这很好用，但想知道是否有人有更好（更快）的方法。

针对 cme​​ts 我尝试了两种替代方法，一种结果慢得多，另一种几乎相同......

方法0（我的第一次尝试）：

def get_centers0(X, r): 

    N = X.shape[0]
    D = X.shape[1]
    grid = np.zeros([0,D])
    nearest = near.NearestNeighbors(radius = r, algorithm = 'auto')

    while N > 0:
        nearest.fit(X)
        x = X[int(random()*N), :]
        _, del_x = nearest.radius_neighbors(x)
        X = np.delete(X, del_x[0], axis = 0)
        grid = np.vstack([grid, x])
        N = X.shape[0]

    return grid

方法一（使用预计算图）：

def get_centers1(X, r): 

    N = X.shape[0]
    D = X.shape[1]
    grid = np.zeros([0,D])
    nearest = near.NearestNeighbors(radius = r, algorithm = 'auto')
    nearest.fit(X)
    graph = nearest.radius_neighbors_graph(X)

    #This method is very slow even before doing any 'pruning'

方法二：

def get_centers2(X, r, k): 

    N = X.shape[0]
    D = X.shape[1]
    k = k
    grid = np.zeros([0,D])
    nearest = near.NearestNeighbors(radius = r, algorithm = 'auto')

    while N > 0:
        nearest.fit(X)
        x = X[np.random.randint(0,N,k), :]

        #min_dist = near.NearestNeighbors().fit(x).kneighbors(x, n_neighbors = 1, return_distance = True)
        min_dist = dist(x, k, 2, np.ones(k)) # where dist is a cython compiled function
        x = x[min_dist < 0.1,:]

        _, del_x = nearest.radius_neighbors(x)
        X = np.delete(X, del_x[0], axis = 0)
        grid = np.vstack([grid, x])
        N = X.shape[0]

    return grid

按如下方式运行：

N = 50000
r = 0.1
x1 = np.random.rand(N)
x2 = np.random.rand(N)
X = np.vstack([x1, x2]).T

tic = time.time()
grid0 = get_centers0(X, r)
toc = time.time()
print 'Method 0: ' + str(toc - tic)

tic = time.time()
get_centers1(X, r)
toc = time.time()
print 'Method 1: ' + str(toc - tic)

tic = time.time()
grid2 = get_centers2(X, r)
toc = time.time()
print 'Method 1: ' + str(toc - tic)

方法0和方法2差不多……

Method 0: 0.840130090714
Method 1: 2.23365592957
Method 2: 0.774812936783

【参考方案1】：

我不确定您到底要做什么。您提到想要创建“近似网格”或“均匀分布”，而您提供的代码选择点的子集，使得成对距离不大于r。

几个可能的建议：

如果你想要的是一个近似网格，我会构建你想要近似的网格，然后查询每个网格点的最近邻。根据您的应用程序，您可能会进一步将这些结果修剪为与网格点的距离大于对您有用的剪切点。

如果您想要的是从点中提取的近似均匀分布，我会在每个点进行核密度估计 (sklearn.neighbors.KernelDensity)，然后进行随机子从每个点的局部密度的倒数加权的数据集中选择。

如果你想要的是一个点的子集，使得没有成对距离大于r，我将首先构建一个半径为r的radius_neighbors_graph，它将一次性为您提供所有太靠近的点的列表。然后，您可以使用类似于您上面编写的剪枝算法来根据这些稀疏图距离删除点。

希望对你有帮助！

【讨论】：

我在你的观点 3 之后。不知道 radius_neighbors_graph 会试一试并报告。

对于我想到的样本量，图形方法似乎要慢得多....

【参考方案2】：

我想出了一个非常简单的方法，它比我之前的尝试效率更高。

这只是简单地循环数据集并将当前点添加到网格点列表中，前提是它与所有现有中心的距离大于 r。这种方法比我以前的尝试快了大约 20 倍。因为没有涉及外部库，所以我可以在 cython 中运行这一切......

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
def get_centers_fast(np.ndarray[DTYPE_t, ndim = 2] x, double radius):

    cdef int N = x.shape[0]
    cdef int D = x.shape[1]
    cdef int m = 1
    cdef np.ndarray[DTYPE_t, ndim = 2] xc = np.zeros([10000, D])
    cdef double r = 0
    cdef double r_min = 10
    cdef int i, j, k

    for k in range(D):
        xc[0,k] = x[0,k]

    for i in range(1, N):
        r_min = 10
        for j in range(m):
            r = 0
            for k in range(D):
                r += (x[i, k] - xc[j, k])**2
            r = r**0.5
            if r < r_min:
                r_min = r
        if r_min > radius:
            m = m + 1
            for k in range(D):
                xc[m - 1,k] = x[i,k]

    nonzero = np.nonzero(xc[:,0])[0]
    xc = xc[nonzero,:]

    return xc

按如下方式运行这些方法：

N = 40000
r = 0.1
x1 = np.random.normal(size = N)
x1 = (x1 - min(x1)) / (max(x1)-min(x1))
x2 = np.random.normal(size = N)
x2 = (x2 - min(x2)) / (max(x2)-min(x2))
X = np.vstack([x1, x2]).T

tic = time.time()
grid0 = gt.get_centers0(X, r)
toc = time.time()
print 'Method 0: ' + str(toc - tic)

tic = time.time()
grid2 = gt.get_centers2(X, r, 10)
toc = time.time()
print 'Method 2: ' + str(toc - tic)

tic = time.time()
grid3 = gt.get_centers_fast(X, r)
toc = time.time()
print 'Method 3: ' + str(toc - tic)

新方法的速度提高了大约 20 倍。如果我提前停止循环（例如，如果 k 次连续迭代未能产生新中心），它可以变得更快。

Method 0: 0.219595909119
Method 2: 0.191949129105
Method 3: 0.0127329826355

【参考方案3】：

也许您只能在每 k nearest 对象以加快处理速度。大多数时候邻里结构应该不会有太大变化。

【讨论】：

好点。我有一个替代版本，我只在开头放置 nearest 对象，然后跟踪我到目前为止删除的点。虽然它实际上速度较慢，但​​我认为问题在于，当您改装时，随着剩余样品的缩小，您的速度会加快。你的想法可能会解决这个问题。我会试试的。

【参考方案4】：

听起来您正在尝试重新发明以下内容之一：

集群特征（参见 BIRCH） 数据气泡（请参阅“数据气泡：层次聚类的质量保持性能提升”） 树冠预聚类

即这个概念已经被发明了至少 3 次，只是略有不同。

从技术上讲，它不是集群。 K-means 也不是真正的聚类。

将其描述为矢量量化更为恰当。

【讨论】：

谢谢，我想会是这样。我不认为您可以将我指向执行此操作的特定 python 库？