如何以矢量化方式平均給定距離內的所有座標

Question

我確實找到了計算點羣集的中心座標的方法。然而，當初始座標的數量增加時（我有大約100 000個座標），我的方法非常慢。

瓶頸是代碼中的for循環。我試圖通過使用np.apply_along_axis來刪除它，但發現這只不過是一個隱藏的Python循環。

是否有可能以矢量化的方式檢測並平均出各種大小的過於接近點的聚類？

import numpy as np 
from scipy.spatial import cKDTree 
np.random.seed(7) 
max_distance=1 

#Create random points 
points = np.array([[1,1],[1,2],[2,1],[3,3],[3,4],[5,5],[8,8],[10,10],[8,6],[6,5]]) 

#Create trees and detect the points and neighbours which needs to be fused 
tree = cKDTree(points) 
rows_to_fuse = np.array(list(tree.query_pairs(r=max_distance))).astype('uint64') 

#Split the points and neighbours into two groups 
points_to_fuse = points[rows_to_fuse[:,0], :2] 
neighbours = points[rows_to_fuse[:,1], :2] 

#get unique points_to_fuse 
nonduplicate_points = np.ascontiguousarray(points_to_fuse) 
unique_points = np.unique(nonduplicate_points.view([('', nonduplicate_points.dtype)]\ 
               *nonduplicate_points.shape[1])) 
unique_points = unique_points.view(nonduplicate_points.dtype).reshape(\ 
              (unique_points.shape[0],\ 
              nonduplicate_points.shape[1])) 
#Empty array to store fused points 
fused_points = np.empty((len(unique_points), 2)) 

####BOTTLENECK LOOP#### 
for i, point in enumerate(unique_points): 
    #Detect all locations where a unique point occurs 
    locs=np.where(np.logical_and((points_to_fuse[:,0] == point[0]), (points_to_fuse[:,1]==point[1]))) 
    #Select all neighbours on these locations take the average 
    fused_points[i,:] = (np.average(np.hstack((point[0],neighbours[locs,0][0]))),np.average(np.hstack((point[1],neighbours[locs,1][0])))) 

#Get original points that didn't need to be fused 
points_without_fuse = np.delete(points, np.unique(rows_to_fuse.reshape((1, -1))), axis=0) 

#Stack result 
points = np.row_stack((points_without_fuse, fused_points)) 

預期輸出

>>> points 
array([[ 8.  , 8.  ], 
     [ 10.  , 10.  ], 
     [ 8.  , 6.  ], 
     [ 1.33333333, 1.33333333], 
     [ 3.  , 3.5  ], 
     [ 5.5  , 5.  ]]) 

EDIT 1：爲循環創建變量

#outside loop 
points_to_fuse = np.array([[100,100],[101,101],[100,100]]) 
neighbours = np.array([[103,105],[109,701],[99,100]]) 
unique_points = np.array([[100,100],[101,101]]) 

#inside loop 
point = np.array([100,100]) 
i = 0 

：1環與期望的結果

步驟1的實施例

步驟2：檢測其中一個獨特的點的points_to_fuse陣列中出現的所有位置

locs=np.where(np.logical_and((points_to_fuse[:,0] == point[0]), (points_to_fuse[:,1]==point[1]))) 
>>> (array([0, 2], dtype=int64),) 

步驟3：創建點的陣列，並且在這些位置處的相鄰點並計算平均

一個完整的運行後

array_of_points = np.column_stack((np.hstack((point[0],neighbours[locs,0][0])),np.hstack((point[1],neighbours[locs,1][0])))) 
>>> array([[100, 100], 
      [103, 105], 
      [ 99, 100]]) 
fused_points[i, :] = np.average(array_of_points, 0) 
>>> array([ 100.66666667, 101.66666667]) 

環路輸出：

>>> print(fused_points) 
>>> array([[ 100.66666667, 101.66666667], 
      [ 105.  , 401.  ]]) 

Answer 1

瓶頸不是必需的循環，因爲所有的街區都不一樣大小。

陷阱是points_to_fuse[:,0] == point[0]在循環中觸發二次複雜性。您可以通過按索引排序點來避免這種情況。

爲例做，即使它並沒有解決整個問題（的rows_to_fuse產生後）：

sorter=np.lexsort(rows_to_fuse.T) 
sorted_points=rows_to_fuse[sorter] 
uniques,counts=np.unique(sorted_points[:,1],return_counts=True) 
indices=counts.cumsum() 
neighbourhood=np.split(sorted_points,indices)[:-1] 
means=[(points[ne[:,0]].sum(axis=0)+points[ne[0,1]])/(len(ne)+1) \ 
for ne in neighbourhood] # a simple python loop. 
# + manage unfused points. 

另一改進是，如果你想加快代碼來計算與numba手段，但我認爲現在的複雜性是最佳的。