Двоичные / унарные функции предиката, перекрестно сравнивающие ВСЕ объекты со ВСЕМИ другими объектами в python

Я задал очень похожий вопрос раньше. Поскольку решение arcpy очень громоздко, теперь я ищу в основном ту же функцию в geopandas. Возникает вопрос: каков самый быстрый / лучший способ применить функцию бинарного предиката (например touches ), где каждая функция x сравнивается с всеми остальными функциями либо x, либо другого набора данных y. Я ожидал, что результат будет аналогичен поведению по умолчанию в R:

Если y отсутствует, эффективно вызывается st_predicate(x, x), и возвращается квадратная матрица с диагональными элементами st_predicate(x[i], x[i]).

Чтобы проиллюстрировать это с помощью фиктивных данных и функции st_overlaps():

library(sf)

b0 = st_polygon(list(rbind(c(-1,-1), c(1,-1), c(1,1), c(-1,1), c(-1,-1))))
a0 = b0 * 0.8
a1 = a0 * 0.5 + c(2, 0.7)
a2 = a0 + 1
a3 = b0 * 0.5 + c(2, -0.5)
x = st_sfc(a0,a1,a2,a3)

plot(x)


st_overlaps(x)
#> Sparse geometry binary predicate list of length 4, where the predicate was `overlaps'
#>  1: 3
#>  2: 3
#>  3: 1, 2
#>  4: (empty)

Как добиться аналогичного поведения в python / geopandas? Очевидно, geopandas автоматически выравнивает x и _16 _ / _ 17_, и выполняется поэлементное сравнение (см. этот вопрос SO и этот выпуск на github). В python запуск x.overlaps(x) просто возвращает серию панд с четырьмя True значениями.

import geopandas as gpd

x.overlaps(x)
0      True
1      True
2      True
3      True

python r sf gis geopandas
person Ratnanil    schedule 08.07.2019    source источник

Ответы (2)


arrow_upward
0
arrow_downward

Идиоматический способ Python выразить это с помощью понимания списка, например чтобы создать список, состоящий из кортежа (index: (overlapping-index)), вы должны написать

[ ( ind, 
    [ind2 for ind2, g2 in enumerate(series) if g.overlaps(g2)] 
  ) 
  for ind, g in enumerate(series) ]

Результат:

[(0, [2]), (1, [2]), (2, [0, 1]), (3, [])]

Но, как отмечает Мартинфлайс, это не суперэффективный способ сделать это, поскольку он не использует никакого пространственного индексирования.

Вы можете повысить производительность, используя операцию наложения, см. http://geopandas.org/set_operations.html

person Michael Entin    schedule 09.07.2019

arrow_upward
1
arrow_downward

Это определенно не самый быстрый способ, поскольку это простой итератор, но если ваши данные не огромны, он может сработать.

import geopandas as gpd
from shapely.geometry import Polygon

b0 = Polygon([(-1,-1), (1,-1), (1,1), (-1,1)])
a1 = Polygon([(1.5,0.2), (2.5,0.2), (2.5,1.2), (1.5,1.2)])
a2 = Polygon([(0,0), (2,0), (2,2), (0,2)])
a3 = Polygon([(1.5,-1), (2.5,-1), (2.5,-0.2), (1.5,-0.2)])

series = gpd.GeoSeries([b0, a1, a2, a3])

results = {}
for poly in series.iteritems():
    results[poly[0]] = []
    for poly2 in series.drop(poly[0]).iteritems():
        if poly[1].overlaps(poly2[1]):
            results[poly[0]].append(poly2[0])

Это даст вам подсказку с вашими ценностями.

{0: [2], 1: [2], 2: [0, 1], 3: []}

Однако имейте в виду, что он проверяет A-> B, а затем B-> A и что он также проверяет полигоны, даже если они явно далеко. Чтобы ускорить это, вы можете использовать пространственный индекс rtree, чтобы проверять только те, которые могут перекрываться, вместо того, чтобы проверять каждый многоугольник друг с другом (дважды).

person martinfleis    schedule 08.07.2019