[GeoPandas] 사용 시 알아두면 좋은 점들

개요

이 포스팅의 원문이 되는 글은 사실 3년 전에 쓰여졌습니다. 당시에는 원래 쓰던 언어도 바꾸고(C++ -> Python) 일하는 환경도 바꾸고 다루는 데이터도 바뀌는 등 큰 변화를 겪고 있었어요. 근 1년 간은 Raster를 집중적으로 다뤘습니다만 원래는 Vector를 중점적으로 다뤘습니다. 그 때의 추억도 소환할 겸, 복기도 하고, 다시 공부하는 느낌으로 포스팅을 쓰게 되었습니다.

오늘 포스팅은 좀 어려운 내용들을 다룹니다. 여러분들은 포스팅을 보면서 아래의 표정을 짓게 될 겁니다.

네, 이해합니다. 석사 과정 당시에 공부한 내용들이 이 포스팅을 통해 빛을 발하는군요.

(+) 2024.08.05 기준 원문을 기반으로 작성한 글입니다. 추후 python 코드로 좀 더 쉽게 이해할 수 있게 예시를 덧붙이겠습니다.

개념

이 포스팅에서 설명할 키워드들입니다.

• join 연산 관련 주의할 점들
• Area 연산 관련 주의할 점들
• Distance 연산 관련 주의할 점들
• Dissolve 연산 관련 주의할 점들
• Pandas concat vs append vs merge
• sindex

join

overlay

해당 연산의 종류는 아래와 같습니다.

• intersection
• union
• identity
• symmetric_difference
• difference
• identity (오! 새로 생겼나?)
  이 연산들은 집합 연산을 생각하면 이해하기에 편합니다. 저는 identity 연산을 처음 보는 거 같아서 뭔지 살펴봤는데요, 다음과 같은 연산을 하는 옵션이군요.
  

  

두 GeoDataFrame df1, df2에 대해 순서대로 입력하고 'identity'옵션을 주면 결과는 {df1인데 df2에도 속하는 부분의 geometry} + {나머지 df1} 가 됩니다. 그림에서 확인할 수 있듯 출력 결과에 geometry가 5개 있음을 확인할 수 있습니다.

해당 연산에서 주의할 점은 연산의 결과에 따라 도형이 변화되어 새 도형으로 반환된다는 점입니다.
union결과만 봐도 그렇습니다. 저는 이 연산의 결과가 단순 두 GeoDataFrame 내의 geometry들의 concat이 아닐까 생각했습니다만 아니었습니다!

아래 글을 참고해서 썼습니다.
Set operations with overlay

overlaps

이 함수에서 고려해야 할 옵션이 여러가지 있긴 한데 저는 이 연산을 다음과 같은 상황에서 사용합니다.

• 용례 : 관심 지역(Area Of Interest) 내에 어떤 객체들이 포함되는지 알아보기 위해서 주로 사용합니다.
• 결과물 : list of Boolean

같이 비교할 연산들은 아래와 같습니다.

• contains : A가 B를 포함하는가?
• crosses : A가 B를 가로지르는가? (contains은 아님)
• disjoint : A와 B가 겹치지 않는 부분 구하기
• touches : A의 경계선과 B의 경계선이 한 점에서 만나는가?
• within : A는 B에 포함되는가? (contain과의 차이는 주어와 목적어의 순서입니다.)

각 연산을 명확하게 구분하려면 아래의 생소한 개념이 필요합니다.
최소한 도형을 구성하는 요소에 대해 알아두신다면 GeoPandas, PostGIS의 공간 연산 파트를 보실 때 이해가 수월하실 겁니다.

도형을 이루는 요소에는 interior, boundary, exterior가 있습니다.
도형을 규정하는 외곽을 boundary, boundary에 둘러싸인 내부를 interior, boundary 밖을 exterior라고 합니다.

예시를 볼까요?

Polygon에 대해서 어떻게 이 세 요소를 규정하는지 보여주고 있습니다.
그리고 Line에 대해서도 이 세 요소를 규정할 수 있습니다.

출처 - PostGIS 기술 문서

그리고 이 개념을 바탕으로 다음의 더 생소한 개념을 설명하겠습니다.
!어려운 개념이기 때문에 이런게 있다고만 알아두셔도 좋습니다!

9-intersection-model

두 도형이 교차하는 경우의 수를 9가지 경우로 구분할 수 있다는 모델입니다.
위키백과
도형의 3가지 구성 요소를 바탕으로 어떻게 도형이 '겹친다'고 말할 수 있는지를 구분하고 있습니다.
아래 그림은 9-intersection-model의 도식화입니다.
이렇게 보니까 굉장히 생소하군요.

우리는 이제 도형이 겹친다 라는 문장을 단순히 생각할 수 없게 되었습니다.
interior까지 겹치는 걸까? boundary만 겹치는 걸까? 아니면 포함되는건가?
이런 질문들이 꼬리에 꼬리를 물 겁니다.

예시를 하나 들겠습니다.
아래 상황에서의 9-intersection-model은 어떻게 작용할까요?

답은 아래와 같습니다.
표에서 숫자는 intersection 여부를 판단했을 때, 연산의 결과에 대한 차원 수 입니다.
0차원은 점, 1차원은 선, 2차원은 면입니다. 그리고 F는 해당 없음 입니다.

sjoin

이 연산은 이름 답게 spatial join의 줄임말입니다. 공간적으로 겹치는지 아닌지 확인하고, 겹친다면 그 다음은 일반적인 join처럼 속성의 합으로 결과를 보여줍니다.
용례는 다음과 같습니다.

• 한국 시, 도의 경계 shp 파일과 전국 중증외상센터 위치 POI 파일을 각각 읽어들여서 sjoin을 수행한 결과 : 각 시,도에 위치한 중증외상센터 위치 정보 및 속성 정보

Area

• 자주 쓰는 좌표계 중 EPSG:4326은 unit이 degree이고 EPSG:3857은 unit이 meter입니다.(metre라고 부를 때가 있습니다.)
• EPSG:3857은 공간을 m단위로 일정하게 나누고, EPSG:4326은 공간을 degree단위로 일정하게 나눴습니다.
• m^2으로 결과를 얻기 위해서는 반드시 데이터를 EPSG:3857로 reprojection을 해줘야 합니다.
• GeoPandas는 to_crs()를 통해서 reprojection이 가능합니다.
• 도형의 경우 shapely 라이브러리를 많이 쓸 텐데요, 자체적으로 projection 관련 함수가 없어서 pyproj를 이용해서 reprojection을 합니다.
• 혹시나 모르고 EPSG:4326과 같은 degree가 unit인 좌표계를 사용하는 중에 모르고 area 함수를 쓰면 경고문이 뜹니다. 다시 한번 강조합니다. 에러가 아닌 경고문이 뜹니다.
• 그리고 GeoPandas에서의 모든 공간 연산은 3차원 기하를 지원하지 않습니다. GeoPandas는 shapely에 의존하는데요, shapely에서 3차원 기하 연산을 지원하지 않습니다.
• 혹시나 PostGIS에서 데이터를 뽑았을 때 결과물이 PolygonZ같이 끝에 Z가 붙는 도형이 나왔다면 이는 3차원 도형을 기술하는 경우입니다. 다만 Z축 값이 모두 0이라면 Polygon과 같은 2차원 도형으로 변환하면 됩니다.

Distance

• 50m 이내의 거리를 재는 거라면 EPSG:4326에서 두 좌표의 차를 구해 대략적으로 쓰는 111.12km(1 degree to km)를 곱해줍니다. 하지만! 이건 대략적으로 확인하려고 쓰는 값입니다. 대략적으로 이정도일거 같다~ 하고 얘기할 때 저 값을 자주 곱합니다. 하지만 정확한 방법은 아닙니다.
• 왜 50m냐구요? 교수님에게 들은 팁입니다. (feat. 교수님의 노하우)
• 50m 이상의 거리를 측정하는 거라면 EPSG:3857로 옮겨서 계산하세요.

Dissolve

• 일종의 aggregation 입니다. (공통 속성을 가진 geometry끼리 합침)
• aggregation 결과를 하나의 도형으로 합치고, 속성도 하나의 row로 합칩니다.
• 결과값은 multi~타입으로 나오게 됩니다. ex) polygon -> multipolygon
• 그렇기에 count를 해보면 1로 나옵니다.
• 용례 : Merge geometry into single geometry라고 검색할 때는 보통 join의 개념에서 단순 row를 합치는 방법들이 보통 나옵니다. 그렇기 때문에 합치고 싶은 대상을 명확히 정의하고 검색하세요.

Pandas concat vs append vs merge

여러 파일을 읽어와서 하나의 GeoDataFrame으로 처리하고 싶을 때 고려할 사항입니다.
GeoPandas는 (이름에서도 알 수 있지만) Pandas를 상속받고, 이에 GIS 특화 함수와 속성을 구현합니다. 그래서 왠만한 Pandas 함수가 먹힙니다.
(그리고 이런 말 할 때 마다 강조하는 거지만, 왠만한 이란 말은 예외도 존재한다는 말입니다.)

# shp 데이터의 경로를 받아와서 이를 한 데이터프레임으로 만들 때, 반복적으로 append를 부를 수 있지만 버퍼를 지속적으로 사용하기 때문에
# 많은 수의 데이터를 한꺼번에 읽어와 한 데이터프레임으로 저장하려면 아래와 같이 concat를 써보세요
gdf = gpd.GeoDataFrame(pd.concat([gpd.read_file(i) for i in target_file_list], 
                            ignore_index=True), crs=gpd.read_file(target_file_list[0]).crs)

이 항목을 작성할 때 참고한 글이 있습니다.
pandas-dataframe-concat-vs-append
여기서 갖가지 연산들이 어떤 차이점을 가지는 지 알아볼 수 있습니다.

• Concat gives the flexibility to join based on the axis(all rows or all columns)
• Append is the specific case(axis=0, join='outer') of concat
• Join is based on the indexes (set by set_index) on how variable =['left','right','inner','outer']
• Merge is based on any particular column each of the two dataframes, this columns are variables on like 'left_on', 'right_on', 'on'

정리하면 concat을 기본으로, append는 concat의 한 경우입니다.
Join은 index 기반으로 이루어지고 Merge는 대상이 되는 두 DataFrame의 컬럼을 기반으로 이루어집니다.

sindex

DB에는 index가 있습니다. GIS는 spatial index가 있습니다. 2차원을 기반으로 어떻게 데이터를 인덱싱할까요?
가장 대표적인 예시는 공간을 직사각형을 기준으로 나누는 R-tree가 있습니다.
이런 식으로 공간을 직사각형으로 나누고 + 계층적으로 나눕니다.

또 다른 기준도 있습니다. h3 index입니다. h3는 우버에서 제안한 공간 인덱싱 방법입니다.
h3는 공간을 정육각형으로 나누고 + 계층적으로 나눕니다. 캘리포니아를 어떻게 h3로 나눴는지 그 예시를 보시죠.

H3: Uber’s Hexagonal Hierarchical Spatial Index
왜 정육각형일까요? 축구공을 생각해보시면 이해가 수월하실 겁니다.

결론

오랫만에 공간 연산 관련 내용들을 복습하니 감회가 새롭습니다.
각 문제마다 적절한 공간 연산이 존재합니다. 어떤 상황에서 어떤 공간 연산을 사용할 지 방법을 찾는 게 GIS 엔지니어의 역할 중 하나입니다.
Raster도 다루기에 매력적인데, 실생활의 문제를 해결하기 위해 Vector 데이터를 고르고, 어떤 공간 연산을 적용해서 문제를 풀 지 설계해서 그 결과를 보는 일은 즐겁습니다. 이번 글은 내용이 좀 어렵기 때문에 언제든 댓글 등으로 궁금한 점이 있으시다면 문의 주세요.
읽어주셔서 감사합니다.