Классифицировать массив Python ближайшим "семенным" регионом?

У меня есть растр экологических мест обитания, который я преобразовал в двумерный массив numpy Python (example_array ниже). У меня также есть массив, содержащий "семенные" регионы с уникальными значениями (seed_array ниже), которые я хотел бы использовать для классификации областей среды обитания. Я бы хотел, чтобы "вырасти" мои семенные регионы в "области моей среды обитания" таким образом, чтобы местам обитания присваивался идентификатор ближайшего семенного региона, измеряемый "через" области ареала обитания. Например:

Мой лучший подход использовал функцию ndimage.distance_transform_edt для создания массива, изображающего ближайшую "семенную" область для каждой ячейки в наборе данных, которая затем была заменена обратно в массив среды обитания. Однако это не очень хорошо работает, поскольку функция не измеряет расстояния "через" мои области обитания, например внизу, где красный круг представляет собой неправильно классифицированную ячейку:

Ниже приведены образцы массивов для моей среды обитания и данных семян, а также пример такого вида вывода, который я ищу. Мои фактические данные намного больше - более миллиона мест обитания/семян. Любая помощь будет высоко оценена!

import numpy as np
import scipy.ndimage as ndimage
import matplotlib.pyplot as plt

# Sample study area array
example_array = np.array([[0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0],
                          [0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 1, 1],
                          [0, 0, 0, 0, 1, 1, 0, 0, 0, 1, 1, 1],
                          [0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0],
                          [1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 1, 1],
                          [1, 1, 0, 1, 0, 0, 0, 0, 1, 1, 1, 1],
                          [1, 1, 1, 1, 0, 0, 1, 1, 1, 0, 0, 1],
                          [1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0],
                          [1, 1, 1, 1, 0, 0, 1, 1, 1, 0, 0, 0],
                          [1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0],
                          [1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0],
                          [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]])

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

seed_array = np.array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
                       [0, 0, 0, 0, 4, 0, 0, 0, 0, 0, 0, 0],
                       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 0],
                       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
                       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
                       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
                       [0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
                       [0, 1, 1, 1, 0, 0, 2, 2, 0, 0, 0, 0],
                       [0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0],
                       [0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0],
                       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
                       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]])

# Plot seeds
plt.imshow(seed_array, cmap="spectral", interpolation='nearest')

desired_output = np.array([[0, 0, 0, 0, 4, 0, 0, 0, 0, 0, 0, 0],
                           [0, 0, 0, 4, 4, 4, 0, 0, 0, 3, 3, 3],
                           [0, 0, 0, 0, 4, 4, 0, 0, 0, 3, 3, 3],
                           [0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 3, 0],
                           [1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 3, 3],
                           [1, 1, 0, 1, 0, 0, 0, 0, 2, 2, 3, 3],
                           [1, 1, 1, 1, 0, 0, 2, 2, 2, 0, 0, 3],
                           [1, 1, 1, 1, 1, 2, 2, 2, 2, 0, 0, 0],
                           [1, 1, 1, 1, 0, 0, 2, 2, 2, 0, 0, 0],
                           [1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0],
                           [1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0],
                           [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]])

# Plot desired output
plt.imshow(desired_output, cmap="spectral", interpolation='nearest')

Ответ 1

Вы можете использовать сегментацию водораздела от scikits-image:

Дистанционное преобразование

from scipy import ndimage as nd
distance = nd.distance_transform_edt(example_array)

Сегментация водораздела

from skimage.morphology import watershed, square
result = watershed(-distance, seed_array, mask=example_array, \
                   connectivity=square(3))

Результат

subplot(1,2,1)
imshow(-distance, 'spectral', interpolation='none')
subplot(1,2,2)
imshow(result, 'spectral', interpolation='none')

В качестве другого варианта и после вашего первоначального подхода вы можете использовать водораздел, чтобы найти связанных соседей с ближайшими семенами. Как вы упомянули в вопросе:

Рассчитайте расстояние до семян:
```
distance = nd.distance_transform_edt(seed_array == 0)
```

Вычислите водораздел на расстоянии:

result = watershed(distance, seed_array, mask=example_array, \
                   connectivity=square(3))

Результат графика:

figure(figsize=(9,3))
subplot(1,3,1)
imshow(distance, 'jet', interpolation='none')
subplot(1,3,2)
imshow(np.ma.masked_where(example_array==0, distance), 'jet', interpolation='none')
subplot(1,3,3)
imshow(result, 'spectral', interpolation='none')

Дальнейшее обсуждение: Метод Watershed пытается вырастить регионы из посевных пиков, протекающих через градиент изображения. Поскольку ваше изображение является двоичным, регионы будут одинаково расширяться по всем направлениям от посевных точек и, таким образом, дают вам точку между двумя регионами. Для получения дополнительной информации о водоразделе см. wikipedia.

В первом примере преобразование расстояния вычисляется в исходном изображении, и, таким образом, регионы расширяются равномерно от семян, пока они не достигнут точки разделения в середине.

Во втором примере преобразование расстояния вычисляется из всех пикселей в любую из посевных точек, а затем применяет водораздел в этом пространстве. Watershed в основном назначит каждый пиксель ближайшему семену, но он добавит ограничение связи.

ПРИМЕЧАНИЕ разность значков в картах расстояний как в плане, так и в воде.

ПРИМЕЧАНИЕ В картах расстояний (левое изображение на обоих графиках) синий означает закрытие, где красный означает далеко.