Obliczanie odległości między unikatowymi regionami tablic Pythona?

Question

Mam rastrowe z zestawem unikalnych łatek/regionów ID, które przekonwertowałem na dwuwymiarową tablicę Pythona. Chciałbym, aby obliczyć parami odległości euklidesowych między wszystkimi regionami, aby uzyskać minimalną odległość oddzielającą najbliższe krawędzie każdej poprawki rastrowej. Ponieważ tablica była pierwotnie rastrem, rozwiązanie musi uwzględniać przekątne odległości między komórkami (zawsze mogę zamienić odległości zmierzone w komórkach z powrotem na metry, mnożąc przez rozdzielczość rastra).

Eksperymentowałem z funkcją cdist z scipy.spatial.distance zgodnie z sugestią w this answer to a related question, ale jak dotąd nie udało mi się rozwiązać problemu, korzystając z dostępnej dokumentacji. Jako wynik końcowy najlepiej byłoby mieć tablicę 3 x X w postaci "od ID, do ID, odległości", w tym odległości między wszystkimi możliwymi kombinacjami regionów.

Oto zbiór danych próbki przypominający mój danych wejściowych:

import numpy as np 
import matplotlib.pyplot as plt 

# Sample study area array 
example_array = np.array([[0, 0, 0, 2, 2, 0, 0, 0, 0, 0, 0, 0], 
          [0, 0, 2, 0, 2, 2, 0, 6, 0, 3, 3, 3], 
          [0, 0, 0, 0, 2, 2, 0, 0, 0, 3, 3, 3], 
          [0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 3, 0], 
          [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 3], 
          [1, 1, 0, 0, 0, 0, 0, 0, 3, 3, 3, 3], 
          [1, 1, 1, 0, 0, 0, 3, 3, 3, 0, 0, 3], 
          [1, 1, 1, 0, 0, 0, 3, 3, 3, 0, 0, 0], 
          [1, 1, 1, 0, 0, 0, 3, 3, 3, 0, 0, 0], 
          [1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0], 
          [1, 0, 1, 0, 0, 0, 0, 5, 5, 0, 0, 0], 
          [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 4]]) 

# Plot array 
plt.imshow(example_array, cmap="spectral", interpolation='nearest') 

Answer 1

Odległości między oznaczonymi regionach obrazu może być obliczony z następującego kodu,

import itertools 
from scipy.spatial.distance import cdist 

# making sure that IDs are integer 
example_array = np.asarray(example_array, dtype=np.int) 
# we assume that IDs start from 1, so we have n-1 unique IDs between 1 and n 
n = example_array.max() 

indexes = [] 
for k in range(1, n): 
    tmp = np.nonzero(example_array == k) 
    tmp = np.asarray(tmp).T 
    indexes.append(tmp) 

# calculating the distance matrix 
distance_matrix = np.zeros((n-1, n-1), dtype=np.float) 
for i, j in itertools.combinations(range(n-1), 2): 
    # use squared Euclidean distance (more efficient), and take the square root only of the single element we are interested in. 
    d2 = cdist(indexes[i], indexes[j], metric='sqeuclidean') 
    distance_matrix[i, j] = distance_matrix[j, i] = d2.min()**0.5 

# mapping the distance matrix to labeled IDs (could be improved/extended) 
labels_i, labels_j = np.meshgrid(range(1, n), range(1, n)) 
results = np.dstack((labels_i, labels_j, distance_matrix)).reshape((-1, 3)) 

print(distance_matrix) 
print(results) 

ta zakłada całkowitą Identyfikatory i w razie potrzeby należy je przedłużyć. Na przykład, na podstawie danych z badań nad obliczona odległość jest matryca,

# From 1    2   3   4    5   # To 
[[ 0.   4.12310563 4.   9.05538514 5.  ] # 1 
[ 4.12310563 0.   3.16227766 10.81665383 8.24621125] # 2 
[ 4.   3.16227766 0.   4.24264069 2.  ] # 3 
[ 9.05538514 10.81665383 4.24264069 0.   3.16227766] # 4 
[ 5.   8.24621125 2.   3.16227766 0.  ]] # 5 

a pełna moc można znaleźć here. Zauważ, że zajmuje to odległość euklidesową od środka każdego piksela. Na przykład odległość między strefami 1 i 3 wynosi 2,0, a są one oddzielone 1 pikselem.

To podejście typu brute-force, w którym obliczamy wszystkie odległości parami między pikselami w różnych regionach. Powinno to wystarczyć dla większości aplikacji. Jeśli jednak potrzebujesz lepszej wydajności, spójrz na scipy.spatial.cKDTree, która byłaby bardziej wydajna w obliczaniu minimalnej odległości między dwoma regionami, w porównaniu do cdist.