사면 단위와 격자 단위 기반의 산사태 위험도와 위험등급 분석 및 비교

분석 단위 추출 및 중첩분석 단계는 검증 자료로 사용될 산사태 피해지 자료의 구축과 공간정보의 명확성을 확보하기 위해 수행한다. 먼저, 산사태 위험성 평가 검증을 위한 산사태 인벤토리는 행정안전부의 생활안전지도에서 제공하는 산사태 발생 이력 정보(’16~’23)를 활용하였다. 다만, 산사태 발생 이력 정보 내에는 발생일시, 발생원인, 피해지 주소만이 존재함으로 웹사이트를 통해 확인되는 엄정면의 자료를 별도로 확인하였으며, 확인된 192개소의 발생 이력을 토대로 Google Earth Pro(Ver 7.3.6., Google, California, U.S.A.)의 Time Slider 기능을 활용하여 시계열 분석을 수행하였으며 피해 원형이 확인되는 131개소의 피해지를 직접 구축하였다(Figure 3).

Figure 3. Distribution of landslide damage areas.

Download Original Figure

구축한 각 산사태 피해지 공간정보는 피해 범위 내 최대 표고로 확인되는 지점을 발생원으로 가정하여 격자 단위와 사면 단위 분석 결과에 대한 검증 자료로 이용하였다. 산사태 발생원의 추출은 ArcGIS Pro(Ver 3.4.2., ESRI, California, U.S.A.)를 이용하여 산사태 피해지 공간정보와 공간해상도 10 m급의 수치표고모형을 중첩 분석하여 산사태 피해지별 최대 표고점을 추출하였다(Figure 4).

Figure 4. extract landslide source point.

Download Original Figure

격자 단위는 산사태 위험지도 산정기법을 통해 획득하였으며, 사면 단위의 추출은 GRASS GIS(Ver 8.3.1., GRASS development Team, Open Source Geospatial Foundation, U.S.A.)환경 하에 Alvioli et al.(2016)이 개발한 r.slopeunits 모듈을 이용하였다. 구동에 필요한 매개변수 설정은 표준화된 값의 제공이 되지 않아, 계수를 조정하며 0.25 m 항공정사영상을 토대로 실제 지형과의 일치성이 가장 높은 경험적 계수를 적용하였다(Table 2). t는 흐름 축적의 초기 임계치(m²), a는 사면의 최소면적(m²), c는 경사 방향으로의 원형 분산 정도, r은 t에 대한 반복 감소 계수, m은 최대 반복 횟수이다.

이후, Table 1의 열거된 자료를 중첩 분석하여, 분석 단위별 산림환경특성을 검토하였다.

다만, r.slopeunits는 수문학적 분류를 통한 사면 단위 추출을 기반함으로 공간해상도가 적절하지 않거나 수계학적 연결 관계의 탐색이 어려운 경우 등에 의해 단일 셀이 분석 범위로 인식되거나 국소 사면으로 분할되는 등의 기하학적 문제가 발생한다. 이에 따라, 본 연구에서는 식 1에서 식 3을 통해 기하학적 문제를 자동 검출하고 보정 처리하였다(Figure 5).

Figure 5. Removal of geometric errors in slopes.

Download Original Figure

식에서, SU는 전체 사면의 집합, s는 특정사면, A_t는 면적 최솟값(m²), SU_small은 A_t 이하 사면의 부분집합, ${\partial }_{s_i}$는 특정사면 i의 경계, L은 ${\partial }_{s_i}$와 ${\partial }_{s_j}$간의 길이, N(s_i)는 특정사면 s_i와 경계를 공유하는 이웃 사면의 집합, A_q는 이웃 사면 중 면적이 가장 큰 사면, s_merged는 처리 사면이다.

해당 식에서 설정할 수 있는 값은 A_t이며, 이 연구에서는 10 m급 공간해상도에서 산출될 수 있는 가장 작은 기하학적 오류 면적인 1,000 m²를 면적 최솟값으로 설정하였다.

분석 단위별 산사태 위험도 평가는 산사태 위험지도의 위험도 산정 방안과 같고, 이는 식 4와 식 5를 통해 산정된다(Woo et al., 2014; NIFoS, 2023). 산사태 위험등급은 식 5로부터 산출되는 H에 대해 Table 3의 기준을 적용하여 구분한다.

식에서, LS는 산사태 발생확률, S는 사면 경사(°), L은 사면 길이(m), C는 사면 곡률, D는 경급, SD는 토심(cm), A는 사면 방위, B는 모암, F는 임상, H는 산사태 위험도, a는 상부사면 기여면적, β는 경사(Radian)이다. 범주형 자료에 대한 인자값 설정은 NIFoS(2023)를 참고하여 적용하였다.

다만, 사면 단위의 산사태 위험도 평가는 기하학적 문제가 처리된 사면에 대해 중첩분석을 수행하되 연속형 자료는 사면 내 자료의 평균값을 이용하고, 범주형 자료는 사면 내 자료의 최빈값을 이용하여 산사태 위험도와 위험등급을 산정하고, 기술통계 분석을 수행하였다.

분석 단위별 위험도 및 위험등급 검토는 격자 단위와 사면 단위로 산출된 위험도와 위험등급을 검증 자료와 대조하고 각 분석 단위별 평가 결과를 비교 및 통계 분석하기 위해 R(Ver 4.2.2., R Foundation, Vienna, Austria)을 이용해 연속형 자료(위험도)와 범주형 자료(위험등급)를 구분하여 통계 분석하였다.

위험도의 비교는 분석 단위별 일치도를 평가하기 위해 윌콕슨 부호 순위 검정(Wilcoxon signed-rank test)을 수행하였으며, 검정 결과에 대한 실질적 크기의 산정을 위해 Cohen’s d를 이용하여 효과 크기를 분석하였다.

위험등급의 비교는 두 범주형 자료의 등급 분포양상과 차이를 확인하기 위해 카이제곱 검정(Chi-Square Test)을 이용하여 동질성(Homogeneity)과 독립성(Independence)을 분석하였다(Wilcoxon, 1945; Cohen, 1988).

발생원을 기반한 격자 단위의 중첩분석 결과는 Figure 6과 Table 4와 같다.

Figure 6. Overlay analysis results of target areas including landslide source points(n=131). (6a: Slope, 6b: Slope length, 6c: Curvature, 6d: TWI, 6e: Soil depth, 6f: Forest type, 6g: DBH class, 6h: Geology, 6i: Aspect).

Download Original Figure

발생원과 중첩된 격자 단위의 평균 경사는 약 23.9±7.4°로, 급경사지(20~25°)에서 다수 분포한다. 평균 사면 길이는 약 51.6±53.8 m로 변동성이 크게 나타나며 이는 사면 하부 또는 소규모 산사태가 포함되어 나타나는 결과로 판단된다. 평균 곡률은 약 −0.89±2.78로 Suzuki(1977)에 의한 분류 결과, 오목한 지형에서 다수 분포하나 자료의 변동성이 크게 나타난다. 평균 TWI는 약 3.1±1.3으로 전체 대상지의 평균 TWI 값인 약 2.9±1.4보다 소폭 높은 경향이 확인되어 상대적으로 습윤한 지형으로 나타났다. 평균 토심은 약 52±22.3 cm로 얕은 토심에서 발생한 것으로 나타나나 이는 입력자료의 한계로 전토심이 아닌 유효토심 자료를 사용한 결과로 판단된다. 임상은 침엽수림이 최빈값으로 확인되며, 경급은 치수, 모암은 흑운모 화강암, 방위는 남동이 최빈값으로 나타났다.

사면 단위 추출 결과, 대상지는 총 3,496개의 산지 사면으로 구분된다(Figure 7). 추출된 산지 사면 내 산사태 피해 발생 사면의 중첩분석 결과는 Table 5와 같다.

Figure 7. Slope unit extraction results and representative sites.

Download Original Figure

평균 면적은 약 10,302±9,856 m²이며, 5,000 m² 이하의 사면이 2,263개로 약 65% 수준인 것으로 나타났다. 사면의 평균 경사는 약 22.6±6.3°로 격자 단위 중첩분석 결과에 비해 다소 낮게 나타났다. 이는 분석 단위의 확대로 인해 산지 하부, 경사가 감소하는 구간이 포함된 결과로 판단된다. 평균 사면 길이는 약 60.5±98.8 m로 격자 단위 중첩분석 결과에 비해 높게 나타나나 변동성 또한 상승하는 것으로 확인된다. 소규모 사면이 포함되어 나타나는 결과로 판단된다.

평균 사면 곡률은 약 0.16±1.13으로 볼록한 유형의 사면이 다수 추출되었다. 볼록한 사면의 추출이 다수 확인되는 것은 사면 단위 추출 기법에 기인하는 것으로, 사면 단위는 계곡선과 능선을 기준으로 사면을 분류함으로, 계곡선이 경계면으로 작용하는 사유로 판단된다.

사면 평균 TWI는 약 3.0±2.9로 격자 단위 중첩분석 결과에 비해 소폭 낮게 나타나나 변동성은 크게 나타난다. 이는 경사의 경우와 같이 분석 범위의 확대에 따른 결과로 판단된다. 사면 내 평균 토심은 약 48.7±16.6 cm로 격자 단위 추출에 의한 발생원 중첩분석 결과에 비해 낮게 나타났다. 임상은 활엽수림이 최빈값으로 나타났으며, 경급은 중경목, 모암은 흑운모 화강암, 사면 방위는 북서로 확인된다.

분석 단위별 산사태 위험도 결과에 대한 Shapiro-Wilk 정규성 검정 결과, 사면 단위의 위험도 결과에서 p<0.05의 수준으로 정규성 가정을 만족하지 않는 것으로 나타났다(Table 8).

이에 따른, 윌콕슨 부호 순위 검정 결과, p<0.05의 수준으로 분석 단위별 위험도 분석 결과 간의 유의성이 존재하며, 이는 Zhuang et al.(2016)과 Huang et al.(2024)이 시사한 바와 같이, 격자 단위와 사면 단위의 분석 방법의 차이에 기인한 것으로 판단된다[Figure 9(a)].

Figure 9. Visualization results of the statistical analysis on landslide hazard. (9a: Violin plot including the results of the Wilcoxson signed-rank test, 9b: Cohen’s d results for individual observations).

Download Original Figure

다만, Eboh et al.(2021)이 시사한 바와 같이, 윌콕슨 부호 순위 검정 결과의 유의성이 존재하는 것은 윌콕슨 부호 순위 검정이 데이터의 분포를 고려하지 않고 중앙값을 고려하기 때문으로, 비모수적 효과 크기 지표를 고려하여 결과의 실질적인 차이를 분석해야 한다.

이에 따라, Cohen’s d를 이용한 효과 크기 분석 결과, d(효과 크기)는 0.161의 수준으로 격자 단위와 사면 단위의 산사태 위험도 평가 결과는 방법론적 차이에 의한 통계적 유의성은 존재하나 실질적인 차이가 존재하지 않는 것으로 나타났다[Figure 9(b)]. 즉, 격자 단위와 사면 단위 간 산사태 위험도 평가의 통계적 차이에도 불구하고 실질적인 차이가 없으므로, 예측 모형의 정확도와 현장 적용에 큰 영향이 없는 것으로 판단된다.

카이제곱 검정을 이용한 분석 단위별 산사태 위험등급 검토 결과, 동질성의 유의 확률은 p>0.05의 수준이며 χ² = 4.04로, 통계적 유의성이 나타나지 않았다[Figure 10(a)]. 즉, 두 분석 단위에서 산출된 위험등급 분포에는 통계적으로 유의한 차이가 존재하지 않아 두 분석 단위별 위험등급 평가 체계의 일관성이 존재한다고 볼 수 있다.

Figure 10. Chi-square test results for landslide hazard grade. (10a: Homogeneity test results using bar graph, 10b: Independence test results using heatmap plot).

Download Original Figure

위험등급 분포가 유사한 경향으로 나타난 것은, 분석 단위별 중첩분석 결과와 같이 지형적 동질성에 기인한 것으로 판단되며, 효과 크기 분석 결과에 근거하여 실질적인 차이가 부재한 결과로 판단된다.

독립성 검정 결과, 유의 확률은 p<0.05의 수준이며 χ² 은 83.95로, 분석 단위별 위험등급의 분포는 통계적 유의성이 존재하는 것으로 나타났다. Heatmap plot으로 결과를 확인할 경우, 특히 두 분석 단위에서 1등급과 5등급으로 분류되는 조합의 표준화 잔차(Standardized residual)가 각각 약 5.40과 5.04로 확인되어, 동일 산사태 피해지에 대해 두 분석 단위 모두 고위험(1등급)과 저위험(5등급)으로 분류되는 빈도가 기대치보다 크게 나타났다[Figure 10(b)].

다만, 중간 위험등급(2~4등급)의 경우는 음의 표준화 잔차가 확인되어 관측 빈도가 기대치보다 낮은 경향을 보였다. 이는 두 분석 단위가 고위험과 저위험과 같은 극단적 평가에 대해서는 일관된 분류 경향을 보이나 중간 위험등급에 대한 평가 분류 결과는 상대적으로 불일치가 존재함을 의미한다.

이러한 두 분석 단위 간의 비독립성은 산사태 위험등급 산정에 선행되는 산사태 위험도 평가 과정의 차이와 평균화로부터 기인하는 것으로 판단되어 분석 목적에 따른 분석 단위의 설정이 필요할 것으로 판단된다.

특히, Woodard et al.(2023)의 선행 연구를 참고할 경우, 공간해상도에 따라 공간적인 추세와 통계적 특성이 상이해질 가능성이 존재한다. 이는, 격자 단위 분석 방법인 산사태 위험지도 제작 방법을 사면 단위에 적용할 경우, 공간적 세분성의 저하와 사면 단위별 불균질성으로 정밀도가 저하될 가능성이 존재하고, 입력 인자의 평균화 또는 단순화(최빈값 이용)에 대한 추가적 정밀도 문제가 존재함으로 한계점이 존재하여, 추후 사면 단위를 기준으로 한 산사태 위험도 평가 모형 개발 연구가 필요할 것으로 판단된다.

또한, 사면 단위 분석의 경우, 분석 단위의 면적이 상이함으로 이에 따른 통계적 왜곡이 발생할 가능성이 존재한다. 다만, 관련 연구사례가 미비한 실정으로 이 연구에서는 적용하지 못한 한계점이 존재한다. 향후 면적-가중 평균(Area-|Weighted Mean) 또는 면적 기반의 보정 통계량을 산출 및 적용하여 사면 단위의 상이한 면적 차이를 고려한 연구가 필요할 것으로 판단된다.

다만, Chen and Wu(2016)은 위험도 평가 및 재난 관리 정책 수립 측면에서, 로지스틱 회귀 분석 모형과 같은 통계적 모형에 사면 단위 기반 산사태 위험도 분석을 적용할 시 효율성과 실용성을 보유함을 시사한 바 있으며, Woodard et al.(2024)은 사면 단위를 통한 산사태 위험 구역 산출 시, 사면 단위가 격자 단위에 비해 직관적이고 명확한 산사태 위험 구역의 구분을 가능하게 하여, 산사태 피해의 사전 예방과 조기 대응 차원에서의 보수적 접근이 가능함을 밝힌 바 있다.

특히, 산림청은 산사태 발생 우려지역에 대한 기초 및 실태조사를 수행하고 있으며, 피해 유형(산사태·토석류)에 따라 평가를 달리한다.

다만, 기초조사 중 토석류 유형이 유역을 기초 공간 단위로 하는 것에 반해 산사태 유형은 공간 단위 설정이 미비한 실정이다.

이에 따라, 토석류 유형에서 산사태 위험지도를 통해 평가하는 황폐 발생원에 대한 검토가 사면 단위에서는 배제되는 실정이고, 사면을 대상으로 조사하는 산사태 유형은 사면의 설정부터 평가까지 작업자의 경험과 지식에 좌우된다.

따라서, 공간적 기초 단위의 설정을 통해 조사와 평가 방법을 일원화할 필요성이 존재한다. 또한, 사면 단위를 도입할 경우, 공간 단위별 통합 산림 관리체계의 마련이 가능할 것으로 판단된다. 이는 기존에 목적별로 상이하게 분리되어 관리되는 공간적 단위를 (1단계)유역 단위→(2단계)사면 단위→(3단계)격자 단위로 연결하여 공간적 기준을 세분화하고 기초 공간 단위를 유역으로 설정하여 사면 단위와 격자 단위를 유역의 하부 요소로 관리하는 등의 체계적이며 효율적인 산림 관리방안을 수립할 수 있을 것으로 기대된다.

종합하여, 산사태 위험도 및 위험등급 평가에서 사면 단위 기반의 분석을 수행할 경우, 정책·행정적 관리 측면에서의 효율성과 실용성을 보유하나 격자 단위와 달리 일관된 분석 체계의 유지 및 국지적 세부 특성 산출의 한계점을 나타낸다. 따라서, 재해 관리적 측면 또는 심층적 분석과 같은 재해 평가적 측면 등을 고려하여 활용 목적에 따른 분석 단위의 설정이 필요하며 분석 단위의 변경에 따른 통계적 차이와 분류 특성 등을 고려해야 할 것으로 판단된다.