sf : R과 지도(2)

이번 글은 저번에 이어서, sf 패키지를 활용한 지도 다루기에 대해 알아보자 합니다. 저번 글은 지리정보(GIS)를 다루기 위한 정말 기초적인 부분이였다면, 이번 글은 보다 실용적인 부분에서 어떻게 파일을 불러오고, 어떤 활용법을 쓸 수 있는지에 대해 소개해 볼 것입니다.

 

1. 데이터프레임과 SF

앞서 다뤘듯이 sf패키지는 geometry 형식으로 폴리곤(면)과 포인트(점)을 다룹니다. 그러나 앞서 글에서 봤듯이 이를 선언하기엔 너무 많은 작업을 요구하는 것이 현실이기도 합니다. 그렇다면 sf용 객체를 쉽게 선언할 방법은 없을까요? 다행히도 st_as_sf라는 함수로 데이터프레임을 sf형식으로 변환할 수 있습니다. 

st_as_sf 함수

예를 들어서, x좌표와 y좌표 두 가지의 정보를 칼럼으로 가진 데이터프레임이 있는데 이를 point 형식의 객체로 변환하고 싶다면, 다음과 같이 해주면 됩니다.

5 * 2 데이터프레임을 Point로 변환한 모습

여기서 coords는 변환하는 데이터의 어떤 칼럼을 좌표로 쓸지입니다. 예를 들어 df에 id라는 칼럼이 추가되있다면, 위 코드를 실행했을때 x와 y는 geometry로 변환되고 id라는 칼럼은 속성값으로 남게됩니다.

id가 칼럼으로 들어간 sf객체

그런데 여기서 한 가지 의문이 더 생깁니다. 이와 같이 변형된 값이 위경도를 표현하는 건지, 실제 지도값도 이렇게 들어가는지 입니다. 실제로 지도의 좌표를 표현하는 방법은 여러가지 이유(예 : 직선거리에 대한 직관적 표현 및 간결한 계산)로 인해 다양하며, 이를 투영법(projection)이라고 합니다. 

 

2. 지도 불러오기와 투영법 지정

그렇다면 투영은 어떻게 하는 걸까요? sf 패키지에서는 crs(coordinatecoordinate reference system)를 이용하여 조정하며, st_transform 혹은 위의 st_as_sf에서의 crs파라미터 등으로 조정하게 됩니다. 이를 위해 실제 데이터를 먼저 불러와보겠습니다.

구글 지도나 클릭센스 등 일부 프로그램의 경우 kml 형식을 채용하지만, 역시 제일 보편적인 형식은 .shp형식일 것입니다. sf 패키지에서는 데이터를 읽어주기 위해서 read_sf 함수를 쓰게 되는데, 이번 글에서는 shp에 대해서만 해보겠습니다. 사용하는 shp 파일은 http://www.gisdeveloper.co.kr/?p=2332에서 다운받으시면 됩니다(2019년 법정동, 시군구 지도)

불러온 두 데이터 중 시군구 데이터

read_sf의 경우 tibble형식으로 받게 되며(파라미터가 디폴트값일 경우), options라는 파라미터를 통해 인코딩을 수정할 수 있습니다.

밑의 class에서 볼 수 있듯이, tibble dataframe이면서 sf객체의 속성을 가지고 있는 것을 알 수 있습니다. 재밌는 점은 이 때문에 filter, mutate 같은 일반적인 dplyr문법이 먹힘은 물론, groupby 연산이 geometry를 자동으로 포함해서 계산한다는 점입니다. 예를 들어 SIG_CD의 앞 두 글자(=시도 코드)를 기준으로 그룹을 묶어 summarise로 요약하면, 자연스래 해당 그룹의 geometry가 합쳐집니다.

위 코드의 결과에 따라 시/도 단위로 합쳐진 3개의 폴리곤

이 sf객체는 ggplot의 geom_sf를 통해서도 표현될 수 있는데, geom_sf가 받는 aesthetic은 geometry로, sf객체의 geometry칼럼을 대응시켜주면 됩니다. 그 외에도 color(폴리곤의 선 색깔), fill(폴리곤 안 색깔) 등 옵션을 이용해 다양한 방법으로 지도를 표현할 수 있습니다.

예제코드에서 그려진 서울 시군구 지도

그런데 위의 지도나 객체 head에서 보이듯, 지도의 좌표값이 기존에 아는 위경도값과 매우 다름을 알 수 있습니다. 이는 shp데이터가 위경도가 아닌 다른 좌표계를 쓰고 있는 것으로, 이제 crs를 수정해줘야 할 차례입니다. 기존의 좌표계는 객체 호출시 콘솔에 보이는 proj4string 혹은 st_crs 함수를 통해 확인 가능하며, st_transform 함수를 통해서 변경이 가능합니다. 여기서 crs = 4326은 미리 저장되어 있는 crs값 중 epgs:4326(위경도) 값을 사용하겠다는 의미입니다.

기존의 crs와 위경도로 변경된 폴리곤

 

3. 지도 데이터 쓰기

지도 데이터를 불러왔다면 당연히 다시 쓸 방법도 알아야 할 것입니다. read.csv의 반대가 write.csv이듯, read_sf는 write_sf와 쌍을 이루고 있습니다. 여기서 주의할 점은 read_sf에서는 options파라미터로 인코딩 등을 조정했다면, write_sf에서는 layer_options를 사용합니다. 또한, 쓸 때 인코딩을 바꿔서 내보냈다면, 불러올 때 내보냈을 때의 인코딩을 기준으로 불러와야합니다. 만약 인코딩을 맞추지 않는다면 아래와 같이 글자가 깨지게 될 것입니다.

서울시 시군구 지도를 write_sf한 후 다시 불러온 모습

 

4. 응용

이제 지금까지 알아본 바를 이용해 간단한 응용을 해볼 것입니다. 첫번째는 지도에서 겹치는 구간만 빼내는 것입니다. 저번 글에서 확인했었던 st_intersects와 st_intersection을 활용하게 됩니다. 예를 들어 종로구와 겹치는 법정동의 부분만 찾아보고 싶다면 st_intersection을 활용해 겹치는 부분만 생성된 폴리곤을 만들고, 겹치는 법정동이 누구인지에 대해 가져오고 싶다면 st_intersects를 이용해 해당 법정동을 가져오게 합니다. 이를 응용하면 아래와 같이 종로구에 겹치는 법정동과 종로구의 경계선 표시, 그리고 겹치는 부분에 대한 색칠을 한 그래프를 만들 수 있게 됩니다.

종로구와 겹치는 법정동 지도

다른 예제로는 거리 계산에 대한 문제입니다. KATEC좌표계 등 유클리디안 거리를 기반으로 하는 좌표계가 있는가하면 위경도처럼 Great circle distance를 기준으로 하는 경우도 있기에 좌표간의 거리를 계산할 경우 주의를 해야 하는데, 이 경우 st_distance 함수는 좌표계에 맞는 방법을 자동으로 찾아서 계산하는 함수로 매우 유용합니다.

유클리디안 거리와 great circle distance 계산

 

마치며

이번에는 지도 관련 프로젝트를 꽤 장기간 진행하게 되면서 경험한 것 중 기초적이지만 R에서 지도를 다룰때 알면 편리한 것들에 대해서만 정리해보았습니다. 저도 처음 지리정보 데이터를 python이나 R에서 접했을 때 생소하고 다루기도 어려워 꽤 많은 시간을 잡아먹었던 기억이 있어서, 이 글로 정보를 대략적이나마 접하고 지도를 다루게 되면 저보다 쉽게 해결해나갈 수 있을 것이라고 생각됩니다.

 

# 예제코드는 모두 git repository로 옮겨서 관리하게 되었습니다. https://github.com/JunmoNam/applebox_blog

이번 예제 : https://github.com/JunmoNam/applebox_blog/blob/master/R/17)%20geometry(2).R