numpy配列内のラベル付きコンポーネント間の最小エッジ - エッジ間ユークリッド距離

Question

大きなnumpy配列には多数の異なる形式があり、numpyとscipyを使用してそれらの間のエッジ - エッジ間のユークリッド距離を計算したいとします。

ノート：私は検索を行いました。これは他の質問と同様に、配列内のラベル付きパッチ間の最小距離を取得するためのもので、ポイント間または別個の配列間ではありません。

私の現在のアプローチはKDTreeを使って動作しますが、大規模な配列ではひどく非効率です。基本的には、各ラベル付きコンポーネントの座標を検索し、他のすべてのコンポーネント間の距離を計算しています。最後に、平均最小距離が一例として計算される。

私は、よりスマートなアプローチを、Pythonを使用して、好ましくは余分なモジュールなしで探しています。

import numpy 
from scipy import spatial 
from scipy import ndimage 

# Testing array 
a = numpy.zeros((8,8), dtype=numpy.int) 
a[2,2] = a[3,1] = a[3,2] = 1 
a[2,6] = a[2,7] = a[1,6] = 1 
a[5,5] = a[5,6] = a[6,5] = a[6,6] = a[7,5] = a[7,6] = 1  

# label it 
labeled_array,numpatches = ndimage.label(a) 

# For number of patches 
closest_points = [] 
for patch in [x+1 for x in range(numpatches)]: 
# Get coordinates of first patch 
    x,y = numpy.where(labeled_array==patch) 
    coords = numpy.vstack((x,y)).T # transform into array 
    # Built a KDtree of the coords of the first patch 
    mt = spatial.cKDTree(coords) 

    for patch2 in [i+1 for i in range(numpatches)]: 
     if patch == patch2: # If patch is the same as the first, skip 
      continue 
     # Get coordinates of second patch 
     x2,y2 = numpy.where(labeled_array==patch2) 
     coords2 = numpy.vstack((x2,y2)).T 

     # Now loop through points 
     min_res = [] 
     for pi in range(len(coords2)): 
      dist, indexes = mt.query(coords2[pi]) # query the distance and index 
      min_res.append([dist,pi]) 
     m = numpy.vstack(min_res) 
     # Find minimum as closed point and get index of coordinates 
     closest_points.append(coords2[m[numpy.argmin(m,axis=0)[0]][1]]) 


# The average euclidean distance can then be calculated like this: 
spatial.distance.pdist(closest_points,metric = "euclidean").mean() 

---

EDIT はただ、提案されたソリューション@morningsunテストし、それは巨大な速度向上です。しかし、戻り値はわずかに異なります：

# Consider for instance the following array 
a = numpy.zeros((8,8), dtype=numpy.int) 
a[2,2] = a[2,6] = a[5,5] = 1 

labeled_array, numpatches = ndimage.label(cl_array,s) 

# Previous approach using KDtrees and pdist 
b = kd(labeled_array,numpatches) 
spatial.distance.pdist(b,metric = "euclidean").mean() 
#> 3.0413115592767102 

# New approach using the lower matrix and selecting only lower distances 
b = numpy.tril(feature_dist(labeled_array)) 
b[b == 0 ] = numpy.nan 
numpy.nanmean(b) 
#> 3.8016394490958878 

EDIT 2

ああ、それを考え出しました。 spatial.distance.pdistは適切な距離行列を返さないため、値が間違っていました。ここ

Answer 1

標識オブジェクトの距離行列を見つけるための完全ベクトル化方法は次のとおり

import numpy as np 
from scipy.spatial.distance import cdist 

def feature_dist(input): 
    """ 
    Takes a labeled array as returned by scipy.ndimage.label and 
    returns an intra-feature distance matrix. 
    """ 
    I, J = np.nonzero(input) 
    labels = input[I,J] 
    coords = np.column_stack((I,J)) 

    sorter = np.argsort(labels) 
    labels = labels[sorter] 
    coords = coords[sorter] 

    sq_dists = cdist(coords, coords, 'sqeuclidean') 

    start_idx = np.flatnonzero(np.r_[1, np.diff(labels)]) 
    nonzero_vs_feat = np.minimum.reduceat(sq_dists, start_idx, axis=1) 
    feat_vs_feat = np.minimum.reduceat(nonzero_vs_feat, start_idx, axis=0) 

    return np.sqrt(feat_vs_feat) 

このアプローチは、Nがゼロピクセルの数であり、O（N 2 ）メモリを必要とします。これがあまりにも厳しい場合は、1つの軸に沿って「逆ベクトル化」することができます（forループを追加します）。