산림조사 자료의 표본점 할당 모델링을 통한 고해상도 산림공간정보 구축

참조자료로 제7차 국가산림자원조사(NFI, national forest inventory), 2019년 지리산국립공원 공원자원조사(NPRS, national park resource survey) 및 서울대학교 남부학술림(SNUF, Seoul National University Nambu Forest) 조사 자료를 이용하였다(Table 1, Table 2). 총 111개 조사구에서 101종, 12,900개체에 관한 자료가 수집되었다. 조사 위치, 면적과 임목의 수종, 수고, 흉고직경을 공통적으로 활용하였고, 이를 바탕으로 임상, 임분고, 경급을 정의하였다.

임상은 5개 유형 즉, 신갈나무림, 참나무림, 활엽수림, 침엽수림, 침활혼효림으로 나누었으며, 지리산의 최우점종인 신갈나무(Quercus mongolica)를 비롯한 참나무류의 분포를 그 외의 임상과 구분하였다. 각 분류는 순서대로 신갈나무, 참나무류, 기타 활엽수, 침엽수의 수종별 흉고단면적 비율이 75% 이상인 경우에 해당하며, 나머지는 침활혼효림으로 분류하였다. 임분고는 각 표본점의 상위 3개 임목의 평균 수고를 기준으로 0~15 m, 15~21 m, 21 m 이상의 3개 범주로 분류하였다. 경급은 임목의 흉고직경을 소경목(6~18 cm), 중경목(18~30 cm), 대경목(30 cm 이상)으로 구분한 다음, 흉고단면적 비율이 가장 높은 등급을 할당하였다. 각 조사구의 좌표를 기준으로 기후, 지형, 위치 등 환경 정보를 추가하여 최종 참조자료를 구축하였다(Table 2).

식생변수를 정의하기 위해서는 조사방법, 조사항목, 조사면적 및 자료구조 등 여러 측면에서 상이한 세 기관의 자료를 일관성 있게 통합하는 자료 처리 과정이 중요하였다. NFI 자료는 전국에 분포한 표본점 중 지리산국립공원 경계로부터 500 m 이내에 있는 34개를 선택하였다. NFI의 4개의 부표본점으로 구성된 집락표본점을 중점에 위치한 하나의 표본점으로 취급하였으며, 일부 부표본점이 산림 이외의 토지에 위한 경우 나머지의 평균적인 정보를 이용하였다. 각 표본점에서 수종별로 단위면적당 흉고단면적을 산출하기 위하여 NFI의 기본조사원(0.04 ha)과 흉고직경 30 cm 이상의 대경목을 조사하는 대경목조사원(0.08 ha)을 구분하여 계산하였다. 또한, 임분 내 대표성이 높은 수종만 포함하기 위해서 수관급이 우세목, 준우세목, 중층목에 해당하는 임목만 선택하고, 피압목, 폭목은 분석에서 제외하였다. NPRS 자료는 지리산 국립공원 내에 40개의 표본점에서 교목층(수고 >8 m), 아교목층(2~8 m), 관목층(<2 m), 초본층으로 구분된 정보를 제공하며, 이 중 교목층과 아교목층만 이용하였다. 20 m × 20 m의 교목층 조사구와 10 m × 10 m 사분면의 아교목층 조사구를 구분하여 단위면적당 흉고단면적을 계산하였다. SNUF의 자료는 지리산 반야봉 이하 남서쪽 영역에 집중적으로 분포한 37개의 영구표본구에서 층위 구분 없이 조사되었다. NPRS와 SNUF 자료는 흉고직경 6 cm 이하의 임목을 포함하는데, NFI의 기준에 맞추어 흉고직경 6 cm 이상의 임목만 선택하였다.

목표자료를 구축하고자 다양한 출처의 래스터 자료를 수집하였다. 임상도는 벡터 기반 자료이므로 격자 형태의 래스터 자료로 변환하였다. 식생변수 항목이 참조자료와 일치하도록 대축척 임상도(1:5,000, 이하 임상도)의 임상, 임분고, 경급 속성을 재분류하였다(Table 2). 연구지역의 26개 임상을 신갈나무림, 참나무림, 활엽수림, 침엽수림, 침활혼효림의 5개 유형으로 단순화하였다. 임분고는 2 m 단위로 구분된 기존 자료를 참조자료와 동일한 3개 범주로 재분류하였다. 경급은 소경목, 중경목, 대경목으로 구분된 기존자료의 속성을 그대로 이용하였다.

지형변수는 NASA SRTM 30 m 수치표고모형 자료(Farr et al., 2007)를 100 m 해상도로 리샘플링하여 고도, 경사, 사면향을 계산하였다. 위경도는 격자 중심점의 좌표를 추출하였다. 최고기온, 최저기온, 강수량은 기상청의 2000~ 2019년 기후 자료를(Kim et al., 2022), 광합성유효광량(PAR)은 기상청 태양광자원지도(Kim et al., 2020)의 2016년 7월~2021년 6월의 평균 경사면 일사량을 각각 격자별로 평균 내어 이용하였다.

랜덤 포레스트는 무작위로 구성된 의사결정나무들을 결합해 예측 성능을 높이는 앙상블 기법이며(Breiman, 2001), 이 연구에서는 식생, 기후, 지형, 위치에 관한 총 12개 변수를 모형 학습에 활용하였다(Table 2). 이 연구는 다중 반응변수 예측 기능이 추가된 R yaImpute 패키지(Crookston and Finley, 2008)의 ‘변형 랜덤 포레스트’를 활용했으나, 간략히 ‘랜덤 포레스트’로 지칭하고 있다. 참조자료의 조합에 따라 총 7개의 랜덤 포레스트 모형을 구성하고, 표본점 할당을 통해 7개의 표본점 고유번호의 래스터 자료를 획득하였다(Table 3). 모형 1(NFI), 모형 2(NPRS), 모형 3(SNUF)은 각각 34, 40, 37개의 표본점을 가지는 단일 자료를 기반으로 하였으며, 다기관 자료를 통합한 모형인 모형 4(NFI+NPRS), 모형 5(NFI+SNUF), 모형 6(NPRS+SNUF), 모형 7(NFI+NPRS+SNUF)은 각각 73, 71, 77, 109개의 표본점을 활용하였다. 계산 부담과 정확도를 최적화하는 사전 분석을 통해 각 모형이 식생변수마다 100개씩, 총 300개의 의사결정나무를 학습하도록 구성하였다(Rigatti, 2017).

이후, 학습된 참조자료의 식생변수를 반응변수로, 목표자료의 12개 변수를 예측변수로 설정해 모든 격자점에 가장 유사한 표본점을 예측하였다. 개별 의사결정나무가 산출한 표본점 후보 중 가장 여러 번 선택된 고유번호가 최종적으로 할당되며, 이를 통해 생성된 고유번호 래스터는 산림조사 자료와 연계되어 각 격자점에서 수종, 흉고직경, 수고, 밀도 등 구체적인 임목정보를 활용할 수 있다(Riley et al., 2016).

산림조사 자료와 임상도의 일치율은 전반적으로 낮았으며(Table 5), 이는 해당 지역의 높은 모델링 난이도를 나타낸다. 점대점 비교에서는 임상, 임분고, 경급의 일치율이 각각 34.2%, 42.3%, 46.8%이었으며, 점대면 비교에서는 각각 61.8%, 57.7%, 55.9%로 더 높았다. 임상도는 천연림 구획의 최소기준을 0.5 ha 이상으로 설정하지만, 산림조사는 0.04~0.16 ha에서 수행되었으므로 공간적 규모의 차이가 낮은 일치율의 원인 중 하나일 것이다(Choi et al., 2024). 이때, 두 평가방식의 큰 일치율 차이는 지리산의 산림 환경이 복잡하고 임상도의 공간적 패턴의 이질성이 높다는 것을 보여준다.

모형의 임상도에 대한 예측 성능을 F1 점수로 평가한 결과, 세 기관의 자료를 모두 결합한 모형 7이 임상(0.66), 임분고(0.81), 경급(0.65) 모두에서 가장 우수한 성능을 보였다(Table 6). 그러나 최저점의 경우 임상은 SNUF 모형의 0.34, 임분고는 NFI 모형의 0.55, 경급은 NPRS 모형의 0.41로 나타나, 단일 자료 모형 간에 평가항목마다 예측 성능의 차이가 있음을 보여준다. 전반적으로, 다기관 자료를 통합한 모형이 단일 자료 모형에 비해 각 변수의 예측 성능이 우수했다.

이러한 경향은 모형 예측 결과의 공간적 분포 패턴에서도 확인된다(Figure 2). 단일 자료 모형에서는 기관별 조사 자료의 특성이 크게 반영되었으나 여러 자료를 활용한 모형에서는 조사 자료의 특성이 완화되어 임상도와 유사한 패턴이 나타났다(Yoccoz et al., 2001; Zipkin and Saunders, 2018). 단일자료가 예측한 임상의 공간적 패턴을 보면, NFI는 침활혼효림, NPRS는 참나무림, SNUF는 활엽수림을 지리산의 주요 임상으로 예측하였다. NFI와 NPRS 모형은 지리산국립공원 경계 주위에 분포하는 소나무림을 나타내는 저지대 침엽수림 분포를 잘 예측하였으나, SNUF 모형은 저지대 식생을 잘 예측하지 못하였다. 임상도는 지리산 천왕봉과 반야봉 부근의 아고산침엽수림 분포를 잘 드러내고 있으나, 표본점 할당 모델링은 고지대 침엽수를 예측하지 못하였다. 이는 모든 자료에 고지대 조사가 부족하기 때문에 나타난 현상이다.

Figure 2. Forest attribute maps derived from the target image (Forest type map) and outputs from models 1 to 7: (a) Vegetation class, (b) height class, and (c) diameter class.

Download Original Figure

모형 예측 결과와 임상도의 임상, 임분고, 경급의 비율을 비교한 결과, 입력자료에 따른 차이가 정량적으로 드러났다(Figure 3). 단일 자료 모형을 비교할 때, NFI 모형은 침활혼효림을, SNUF 모형은 활엽수림을 과대 추정하였으며, NPRS 모형이 임상도와 가장 유사한 분포를 도출하였다. 이러한 차이는 NPRS 조사가 398~1,680 m의 가장 넓은 고도 범위에서 수행된 반면, SNUF는 622~1,350 m의 제한된 고도 범위의 산림을 집중적으로 조사한 결과로 해석된다. NFI는 481~1,477 m의 고도 범위에서 수행되었으나, 연구에 활용된 34개 조사지점 중 6개가 지리산국립공원 경계 외부에 위치하여 저지대의 소나무 및 낙엽송 군락지를 다수 포함한 점이 영향을 미친 것으로 파악된다(Figure 1).

Figure 3. (a) Vegetation, (b) height, and (c) diameter distributions: Forest type map (target) vs. predictions for models 1 to 7 with percentages (left axis) and counts (right axis). Missing classes were assigned a count of 1.

Download Original Figure

모든 모형은 임상도에 비해서 신갈나무, 참나무림, 침엽수림의 비율을 과소 추정하고, 활엽수림과 침활혼효림은 과대 추정하였다(Figure 3). 이 연구에서는 임상을 정의할 때 참조자료의 흉고단면적 비율을 기준으로 삼았으나, 임상도는 주로 항공영상의 수관 면적 비율을 이용한다. 수관 상층을 차지하지 못하는 비우점 활엽수종의 비율은 흉고단면적을 기준으로 할 때 더 높게 계산된다. 임상을 정의하는 기준의 차이로 인한 오차를 개선하기 위해서는, 향후 산림조사에서 수관 면적에 대한 정보를 파악하거나 기존 자료로부터 수관 면적을 추정할 수 있어야 한다.

임분고의 경우, NPRS 모형은 15~21 m 범위, SNUF 모형은 21 m 이상 범위의 식생 비율을 과대 추정하였으며, NFI 모형은 임상도와 비교적 유사한 비율을 도출했다. 다른 모형이 고도가 높을수록 임분고를 낮게 예측하는 것과 달리, NFI 모형은 저지대에서 더 높은 임분고를 예측한 것이 예외적이었다(Figure 2). 경급에서는 NFI 모형이 소경목 비율은 낮고 대경목 비율은 높은 것에 반해, NPRS 모형은 소경목과 대경목 비율이 모두 높고 중경목 비율은 낮게 나타났다. 임분고의 오차는 각 기관의 수고 측정 방식의 차이에서, 경급의 오차는 조사지점 선정의 차이에서 주로 비롯되었을 수 있다.

모든 111개 표본점에서 모형 예측 결과와 산림조사 자료 간의 일치율을 평가한 결과, 세 기관의 자료를 모두 결합한 모형 7이 전반적인 일치율이 가장 높았다(Table 7). 임상, 임분고, 경급 각각에서 점대점 61.3%, 62.2%, 63.1% 및 점대면 82.9%, 79.3%, 75.7%의 일치율을 기록하였다. 단일 자료 모형 중에서는 NPRS 모형이 모든 변수에서 가장 우수한 성능을 보였다. 그러나 입력자료의 수가 많을수록 모형에 미치는 영향이 커지므로, 표본점 수가 많은 NPRS 자료가 평가에서 유리함을 인지해야 한다. SNUF 모형은 일부 지역에 조사 자료가 한정되어 다른 기관 조사 자료의 임상을 상호 예측하는데 어려움이 있었다.

학습에 사용된 표본점에서의 일치율을 별도로 평가한 결과, 전반적으로 61.3% 이상의 수준을 보였으나, 경향은 전체 표본점에 대한 재현 능력과 상반되었다(Table 8). 단일 자료 모형은 일치율이 76.5% 이상으로 대부분 높았으나, 여러 자료를 결합할수록 특정 지점에 대한 국지적 예측 정확도가 낮아졌다. 즉, 다기관 자료를 활용하면 전체 지역에 대한 예측 능력은 향상되지만, 실제 산림조사 지점에 정확한 값을 할당하는 능력은 떨어질 수 있음을 의미한다. 이러한 오차가 발생할 가능성은 임상도와 산림조사 자료가 불일치하거나, 자료가 일부 지역에 편중되어 있을 경우, 기관별 측정 방식에 차이가 있을 경우에 증가할 수 있다(Roberts et al., 2017). 따라서, 산림조사 단계에서부터 공간적 표본 분포를 균형 있게 확보하는 전략이 중요하며, 모델링 과정에서도 단순한 자료 결합에만 의존하기보다 기관별 측정 편차를 보정하는 절차를 거쳐야 한다(Morrison, 2016). 표본점 할당 모델링에서는 전체 예측 성능뿐 아니라 입력자료의 다양성과 지점별 공간 분포가 국지적 오차에 미치는 영향을 종합적으로 고려할 필요가 있다.

최종적으로 선정된 모형은 세 기관의 조사 자료를 모두 결합한 모형 7로서, 임상도와 산림조사 자료 전체 지점에 대해 모든 식생변수에서 상당한 예측 성능을 보였다(Table 6, Table 7). 다기관 자료 통합은 기관별 자료의 특성을 상쇄하여 임상도의 식생 분포와 산림 구조를 안정적으로 예측하는 데 기여하였다. 단, 학습에 이용된 표본점에서의 국지적 예측에서 오차가 증가할 수 있다는 점에서 모형 해석에 유의해야 한다.