서 론
소나무속은 중생대 백악기에 한반도에 처음 등장했다(Gong, 2008). 다양한 수종과 공존해왔던 소나무는 사람들의 산지이용과 훼손이 늘어나면서 산림 내 비중이 점차 증가했다. 특히 조선시대에는 국가적으로 필요한 목재자원 확보를 위해 소나무 중심 보호 정책이 적극적으로 시행되었고 과도한 산림이용으로 인한 황폐화가 복합적으로 영향을 주어, 1910년 당시 현재 대한민국에 해당하는 지역의 산림 상당 부분(나무가 자라는 산지의 78%)이 소나무숲이었고 그 중 80% 이상이 어린나무로 구성된 미성숙한 산림이었다(Bae and Kim, 2019). 교란이 많았던 한반도 산지에서 소나무가 넓게 자리잡은 것은 훼손 지역의 천이 초기 수종으로서 소나무가 산림을 재구성하는데 큰 역할을 하였기 때문이다. 교란 이후 산림의 2차 천이 과정을 통해 수종변화 과정이 지속되었고 현재 소나무림 면적 비율은 25%로 이전보다 상당히 감소했지만, 우리나라 산림생태계 내에서 단일 수종 기준으로 소나무가 여전히 가장 넓은 면적을 점유하고 있다.
일반적으로 천이 초기에 조성되는 선구수종은 소나무와 같은 내음성이 약한 양수이고, 그 다음으로 참나무류와 같이 내음성이 있는 수종에 의해 구성되는 절차를 밟는다(Spurr and Cline, 1942). 소나무는 입지・토양 환경의 좋고 나쁨을 가리지 않고 한반도 전역에서 잘 자라는데, 이러한 능력 때문에 산림의 회복과 유지, 생태적 안정화를 도모하는데 큰 역할을 해왔다. 유기물 함량이 적고 모래질 입자가 많으며 산도가 높고 건조한 지역의 경우, 소나무는 잘 생육하지만 활엽수는 잘 생육하지 못해 내음성 활엽수종이 소나무림 입지에 침입하여 경쟁에서 이겨 번성하기는 쉽지 않았다(Lee et al., 1999). 그러나 우리나라 산림의 토양 상태가 점차 양호해짐에 따라, 소나무는 비옥한 환경에서 잘 자라는 참나무류 등 활엽수와의 경쟁에 밀려 분포 면적이 지속적으로 감소하고 있다. Kim et al.(2000)은 월악산 소나무림 조사 결과를 통해 굴참나무로 바뀌다가 신갈나무와 졸참나무로 천이가 진행될 것으로 판단했다. 또한 건조한 지역의 숲틈에서 소나무 치수가 유지되어 소나무림 유지 가능성이 있다고 보았다. Lee et al.(2006)이 가야산 국립공원 홍류동 계곡 소나무림의 15년간 변화를 분석한 결과, 대부분의 지역이 아교목층을 중심으로 서어나무, 굴참나무, 졸참나무로 천이가 진행되고 있음을 확인하였다. Lee et al.(2009)이 속리산국립공원 소나무림을 대상으로 17년간의 식생구조 변화를 분석한 결과, 천이잠재성이 낮은 소나무군집(소나무 순림), 천이 잠재성이 높은 소나무군집(아교목층에 낙엽활엽수가 우점하는 소나무림), 천이진행중인 소나무 군집(소나무-산벚나무군집, 소나무-졸참나무군집), 도태중인 소나무군집(졸참나무-소나무군집, 갈참나무-소나무군집), 낙엽활엽수혼효군집으로 발달하는 것으로 파악되었다. Lee et al.(2012)의 남한산성 소나무림의 19년간 식생구조 변화 분석 결과 소나무림→소나무 및 참나무류 혼효림→참나무류림→서어나무림으로 천이가 진행되었고, 시간이 흐를수록 참나무류의 세력이 확대되었다. Jang et al.(2013)은 용문사계곡 소나무림의 22년간 변화를 분석을 통해 소나무는 세력 유지 또는 약화되었고, 경쟁관계에 있는 졸참나무, 서어나무, 까치박달의 세력이 증가하였다고 보고하였다. 용문사 계곡 소나무림의 천이 경향은 소나무림→졸참나무→까치박달, 서어나무로 예측되었다. Lee et al.(2014)의 한라산국립공원 영실 소나무림의 22년간 식생구조 변화 관찰에서는 소나무의 우점도가 지속적으로 유지되고 있고 아교목층에서 소나무와 서어나무류 세력이 증가하고 있어 천이가 진행되지 않고 있음을 확인했다. 소나무림의 자연적 쇠퇴 현상은 다양한 지역에서 확인되고 있으나, 지역에 따라 편차가 있음이 확인되었다. 국가 통계에 따르면 소나무림을 포함한 침엽수림은 1980년대 전국 산림의 약 50%를 차지하고 있었지만 2020년 현재 37%까지 감소했는데, 과거 10년 동안 침엽수림이 침엽수림으로 그대로 유지되는 비율이 경기(서울・인천), 충남(대전・세종) 지역이 가장 낮고 경북(대구), 강원 지역이 높았던 것으로 분석되었다(Kim et al., 2022).
산림생태계의 천이과정이 진행될수록 산림군집의 구조상의 복잡성, 외부교란으로부터의 안정성, 산림군집의 성숙도가 증가하며 이는 종다양성과 정비례하는 경향이 높은 것으로 인식되고 있다(Odum, 1969; Loucks 1970; Bazzaz, 1979). 하지만 천이 과정에 따른 종다양성의 변화 경향은 임상, 천이 계열, 교란, 관측 기간 등에 따라 차이를 보인다. 미국 동부 뉴잉글랜드 지역 온대림은 100여년 간 영급이 높을수록 Simpson 다양성 지수가 높아지는 경향을 보인 반면(Pan et al., 2018), 독일 남동부 온대림에서 300여년 동안의 장기적인 천이 과정에서 영양 단계별 생물다양성을 비교한 Hilmers et al. (2018)의 연구에서는 천이 초기부터 중반 200년까지 식물종의 풍부도와 밀도가 감소하는 경향을 보이다 산림이 쇠퇴하기 시작하면서 다시 증가하는 U자형 패턴을 갖는다고 하였다. 국내에서 강원도 평창군 일대의 묵밭이 산림으로 천이되는 과정을 관찰한 Lee (2006)은 천이 초기부터 종다양성이 증가하는 경향을 보이며 10-20년 사이에 식물의 종풍부도 및 Shannon 다양성 지수가 가장 높고 이후 다소 감소하는 경향을 보인다고 하였다. 이와 달리 광양 지역 낙엽활엽수림을 개벌한 곳에서 이차 천이가 진행되는 과정을 10여년 간 관찰한 Oh et al. (2004)의 연구에서는 남사면과 북사면 모두 목본 수종의 종 수는 미약하게 증가하지만 Shannon 다양성 지수는 지속적으로 감소한다고 하여 교란의 성질에 따른 종다양성 변화 경향의 차이를 보였다. 한편 천이가 충분히 진행된 광릉 및 계방산 지역의 천연활엽수림에서는 지난 10-15년 동안 하층 식생의 종풍부도가 지속적으로 줄어든 것으로 나타나(Cheon et al., 2014; Kim et al., 2021) 주로 2차 천이가 진행된 국내 임분들은 천이 초기에 종풍부도가 증가하다가 성숙림에 도달하게 되면 종풍부도가 감소하는 경향을 보인다는 것을 알 수 있다. 이를 볼 때 천이 초기에 해당하는 소나무림도 천이가 진행되면 단기적으로 종다양성이 증가하는 경향을 보일 것으로 예상되는데 설악산 계곡 지역에 발달한 소나무림에서 각 천이 계열을 조사한 Chun et al. (2007)은 천이가 진행될수록 산림 교목성 수종이 하층에 유입되면서 종풍부도와 균등도가 점차 증가한다고 하였다. 이와 유사하게 Lee et al.(2006)이 가야산 소나무림을 15년 간격으로 조사한 결과에서도 아교목층과 관목층의 출현 종이 늘어나면서 종다양성이 증가하는 경향을 보였다.
소나무림의 자연적 쇠퇴가 가시화되면서 과거 소나무가 우리 국민들에게 제공해왔던 다양한 가치들이 함께 쇠퇴할 위험에 처해 있다. 천연기념물, 시도기념물, 보호수, 왕릉숲, 사찰숲, 마을숲, 유적지, 휴양・교육공간 등의 중요 소나무와 소나무림이 전국 곳곳에 분포되어 있고, 문화재용 목재생산림, 소나무 유전자원보호림 등 보호구역, 소나무 종자 육성・관리 채종림 등 소나무 자원 육성과 관련된 소나무림도 매우 많다. 이러한 역사・문화적, 자원・생태적 측면의 중요 소나무림의 경우 지속적인 관리를 통해 소나무림 가치를 유지시킬 필요가 있으나, 자연적인 변화 압력의 현황과 특성을 이해해야 산림 관리의 정확한 진단과 처방이 가능하다. 또한 자연적인 소나무림 쇠퇴 정도에 따라 소나무림 유지 관리지역 또는 신규 수종 육성 지역 등의 새로운 미래 산림관리 방향을 마련하는 것 또한 필요하다. 이를 위해서는 소나무림의 자연적 변화 과정이 어떤 특성에 따라 진행되고 있는지, 지역적으로는 어떻게 차이가 나는지에 대한 명확한 이해가 필요하다.
따라서, 본 연구에서는 우리나라 전역의 소나무림이 어떠한 변화 과정을 겪고 있는지를 파악하고, 소나무림 변화 특성에 기반하여 어느 지역이 미래 소나무림 변화 압력이 상대적으로 클지 예측하는 것을 목적으로 하였다. 이를 위해, 시계열 국가산림자원조사 자료를 이용하여 전국 산림의 10년 동안의 소나무림 변화 정보를 확보하였고, 소나무림 변화에 영향을 줄 수 있는 다양한 변수(지형, 임분, 교란, 기후)를 구축하고, 이를 기반으로 소나무림 변화와 유지에 영향을 주는 변수를 평가하였다. 또한 유효변수들을 이용하여 모형을 구축하고 10년 경과에 따른 미래의 소나무림의 변화 모습을 예측했다.
재료 및 방법
소나무림 변화의 정보는 국가산림자원조사(National Forest Inventory, NFI) 시계열 조사 자료를 통해 취득하였다. 국가산림자원조사는 국가산림통계산출을 위하여 1970년대부터 수행되고 있는 산림청의 조사 사업으로, 2006년 제5차 조사부터는 전국적으로 고정된 표본점을 설치하여 5년 주기로 반복적인 변화모니터링 기반의 현장조사를 수행하고 있다. 따라서 해당 자료는 전국 산림의 임분 단위 변화를 파악할 수 있는 국가 규모의 유일한 현장조사 자료라고 할 수 있다. 본 연구에서는 과거 10년 동안의 자연적 임분 동태 변화에 따른 소나무림 변화 상태를 확인하기 위해 제5차 NFI 자료(2006~2010)와 제7차 NFI(2016~2020) 자료를 비교 대상 자료로 활용하였다.
소나무림의 자연적 변화 특성을 파악하기 위해, 표본점 위치가 훼손되지 않아 산림의 변화 모니터링이 가능한 표본점, 인위적인 산림관리 활동 이력이 없고 병해충이나 산불 등 교란이력이 없는 표본점을 선정하였으며, 조사정보에 기록되지 않았으나 산림관리 활동이 있었을 것으로 추정되는 치수림 및 1영급 임분 표본점은 분석에서 제외하였다. NFI6과 NFI7 모두에서 모니터링이 가능하다고 평가되어 최종적으로 NFI5와 NFI7 자료간의 비교가 가능한 표본점은 총 12,833개소였으며, 인위적 활동과 교란이 있을 것으로 추정된 표본점을 제외한 총 9,220개 표본점 자료를 이용하여 분석을 수행했다. 이후 상세 분석 과정에서는 해당 표본점 중에서 소나무류 중요치 등급 등에 따라 추가 선별하여 모델 구축에 활용하였다.
표본점의 대표수종을 파악하기 위해 임분 상층(우세목, 준우세목)을 구성하는 교목을 대상으로 중요치를 계산하였다. 중요치는 상대적 본수 비율과 상대적 흉고단면적합 비율의 평균값으로 계산되었으며, 소나무류[소나무(PD), 곰솔(PT)]의 중요치에 따라 Class1(IV>0.75), Class2(0.5< IV≤0.75), Class3(0.25<IV≤0.5), Class4(IV≤0.25)로 구분하였다. 이 중 Class1에 해당하는 표본점이 실제 임상도에서 구분하고 있는 수종 구분 기준(임분 상층 75% 이상을 점유하고 있는 수종)에 따른 소나무림과 가장 유사하다고 볼 수 있다.
중요치 그룹을 기준으로 NFI5와 NFI7 사이 기간 동안의 중요치 변화 추세를 파악하였다. 또한 소나무림 중요치 변화와 함께 임분 전체(교목, 관목 포함)의 생물다양성 변화가 어떻게 변화하였는지 함께 파악하였다. 생물다양성은 종 다양성 지수인 Shannon 다양성 지수로 다음의 식 1과 같이 산출하였다(Shannon, 1948).
여기서 pi는 i종의 상대 피도(흉고단면적)를 의미하며, S은 표본점 내 출현한 전체 목본성 수종의 종 수를 의미한다.
임분 내 소나무림의 변화 요인을 파악하여 이를 설명할 수 있는 모형을 개발하기 위해, NFI5 전체 표본점 중에서 소나무 임분으로 규정되는 Class1 표본점을 대상으로 분석을 수행하였다. 과거(NFI5) Class1이었던 표본점이 10년이 경과한 후(NFI7)에도 Class1로 유지가 되는 표본점과 중요치의 감소에 따라 Class2~4로 변화하는 표본점을 구분하여 이 두 가지 경우에 영향을 미치는 요인을 파악하였다. 모형 적용을 통해 최종적으로 소나무림의 변화 또는 미변화 여부를 결정해야 하므로, 소나무림 중요치 변화에 대한 연속형 모형이 아닌 변화 여부를 통계적으로 결정할 수 있는 분류형 모형을 개발하는 것을 목표로 하였다.
모형의 종속변수는 소나무 중요치 등급 유지를 0(소나무림 유지), 중요치 등급 변화를 1(소나무림 변화)로 설정했으며, 소나무림 변화를 설명할 변수로는 지형변수, 임분변수, 교란변수, 기후변수를 구축했다(Table 1). 소나무림 변화 모형 구축 자료는 표본점 조사 정보를 활용하였고, NFI 조사를 통해 확보되지 않은 자료들은 공간정보로 제작하여 해당 위치에 해당하는 자료를 추출하여 활용하였다.
지형의 경우 고도, 사면향 자료는 표본점 조사정보를 이용하였다. 사면 방향에 따른 건조도 영향을 정량적으로 비교하기 위해 남서사면이 건조도가 상대적으로 높고 북동사면이 상대적으로 습윤한 특성을 고려하여 남서사면(SW)을 180, 북동사면(NE)을 0으로 지정하고 전체 사면향 값을 이에 따라 변환하여 활용했다. 즉, 사면향의 값이 클수록 건조, 작을수록 습윤한 상태를 말한다. 지형학적수분지수(TWI, Topographical Wetness Index)는 사면 형태에 따른 수분조건을 반영한 지수(Beven and Kirkby, 1979)로서, 전국 30 m 수치표고모델을 이용하여 산출했다. 임분 변수는 표본점의 임분조사 정보의 영급, 수관밀도, 상층의 소나무 상대중요치, 생물다양성 지수(Shannon index) 자료를 이용하였다. 임분의 유지와 변화에 영향을 줄 수 있는 외부 교란 요인으로 표본점 주변지역(반경 100 m)의 산림 내 소나무림 비율과 비산림 비율을 활용하였다. 해당 변수 선정을 위해 전국 임상변화 요인 분석에 대한 선행 연구에서 기도출된 주요 변수들을 포함하고(Kim et al., 2022), 여기에 생물다양성 지수를 추가로 고려하여 소나무림 변화/미변화에 영향을 주는 변수를 평가하였다. 소나무림 및 비산림 면적 비율 자료는 NFI5 조사 시기를 고려하여 대축척임상도가 처음 제작되어 배포된 2012년 제작 임상도를 이용하여 자료를 제작했다.
기후 부문에서는 소나무림 분포지의 일반적인 기후특성을 반영하기 위해 기상청의 격자기상자료(MK-PRISM ver1.2)가 구축되어 있는 2000~2019년의 평균기온과 강수량 자료를 이용하여 20년 평년 연평균 자료를 구축해 지역별 기후특성을 반영하는 자료로 활용했다. 격자기상자료는 기상청 기후정보포털(http://www.climate.go.kr)에서 다운로드하여 이용했다.
소나무림 유지 또는 변화에 영향을 주는 요인을 파악하고 이를 기반으로 현재 기준으로 이후 소나무림 변화 가능성을 파악하기 위해 표본점 조사정보 기반의 소나무림 변화 모형을 구축하였다. 우선, 소나무림 우점도 변화에 영향을 줄 수 있는 지형, 임분, 교란, 기후 변수들에 대해 독립변수 t-test을 통해 통계적 유의성을 검증하고, 모형 구축에 활용할 수 있는 변수를 우선 선정했다(유의확률이 0.05보다 작은 변수).
관련 변수를 복합적으로 고려하여 변화, 미변화를 구분할 수 있는 후보모형으로 Random Forest(RF)(Breiman, 2001)와 Support Vector Machine(SVM)(Cortes and Vapnik, 1995)을 이용했다. SVM의 경우 데이터들을 구분하는 알고리즘에 따라 Linear SVM과 Radial SVM 모형을 고려하였다.
최적모형 선정을 위해 소나무림 변화에 대한 동일한 입력자료를 구축하고 세 가지 기계학습 모형을 적용했다. 소나무 중요치 그룹이 변화하지 않은 표본점 개수는 변화한 표본점 개수보다 훨씬 개수가 많기 때문에, 소나무림 변화/미변화 모형 산출 시에 개수가 많은 미변화 표본점의 특징이 상대적으로 영향력이 커지는 문제가 발생했다. 변화 표본점의 특징을 모형 구축에 반영하기 위해, 심플 개수가 적은 경우 비교 대상이 되는 샘플의 개수와 동일하게 맞춰주는 Upsampling 기법을 적용했다. 또한 세 모형 간의 비교를 위해 K-fold Cross Validation 방법을 적용하고, 분류정확도를 평가하는 전체정확도(Overall Accuracy), 카파계수(Kappa Coefficient), 민감도(Sensitivity), 특이도(Specificity)를 이용하여 모형 성능을 평가했다. 통계분석과 모형구동은 R 프로그램으로 수행하였다. 소나무림 변화 예측 최종 모형 구축을 위해 모형 구축과 검증을 위한 샘플을 구분했던 것을 모두 포함해 전체 표본점 자료를 이용해 RF, SVM (Linear), SVM (Radial) 모형을 도출했다.
소나무림 변화를 설명할 수 있는 최종 모형을 이용하여 전국 소나무림에 대한 미래 소나무림 변화를 전망하기 위해 모든 입력자료(지형, 교란, 기후)를 30 m 격자해상도의 공간정보로 구축했다. 30 m 수치표고모델 자료를 이용하여 사면향, 지형학적수분지수 자료를 제작하였다. 임분 외부 교란특성 정보 구축을 위해 2022년 현재 기준 임상도를 이용하여 100 m 반경 내 산림 중 소나무림 면적 비율, 100 m 반경 내 비산림 면적 비율 지도를 제작하였다. 또한 임분 내부 교란 특성을 반영하는 소나무림 내 임분상층의 소나무류 중요치 정보는 NFI7 자료로부터 구했다. NFI7의 소나무림 표본점을 대상으로 소나무류(소나무, 곰솔)의 중요치를 계산하고 해당 점 자료를 역거리가중(IDW, Inverse Distance Weighted) 기법을 이용하여 공간적으로 확장했다(Figure 1). 공간보간법 적용을 위해 소나무림 중요치의 공간적 자기상관성을 분석한 결과, 전역적 모란지수(Global Moran’s I)가 0.22로 산출되어 거리에 따른 군집화 경향이 있는 것으로 파악되었다(p<0.05). 기후정보는 기존에 구축되어 있는 1 km 해상도 자료를 30 m 자료로 Resampling하여 공간해상도를 통일해서 준비했다.
준비된 현재 기준 공간 입력자료와 최종 도출된 3가지 모형을 이용해서 현재 기준으로부터 10년 후의 미래 전국 소나무림 변화를 예측하고, 세 모형에서 모두 소나무림이 변화된 것으로 추정된 경우를 소나무림 변화지역으로 설정하여 최종 결과를 도출하였다.
결과 및 고찰
NFI5와 NFI7 간의 변화분석이 가능하고 인위적 활동과 교란이 없을 것으로 추정된 9,220개 표본점 자료를 기반으로 소나무림의 변화 경향을 파악하였다. 임분상층 소나무 중요치에 따른 등급의 변화를 파악한 결과, 소나무가 대표수종인 Class1 표본점이 10년 후에 등급 유지되는 비율은 75%이며, 나머지 25%는 중요치가 감소하여 등급이 변경되었다. 소나무류가 생육하는 혼효림에 속하는 Class2, Class3 역시 10년 경과 후에 67.3%, 72.0%만이 등급이 유지되었다. 과거 10년 동안 전체 산림 내에서 소나무류의 중요치가 전반적으로 감소하는 방향으로 변화하였다(Table 2).
※ Class1: IV > 0.75, Class2: 0.5 < IV ≤ 0.75, Class3: 0.25 < IV ≤ 0.5, Class4 : IV ≤ 0.25
임분 상층의 소나무류 중요치 감소 현상은 지역적으로 차이를 보였다. 전남, 제주, 전북, 충남 등이 중요치 감소량이 크게 나타났고, 경북, 충북, 경기, 경남 지역의 감소량이 상대적으로 적었다. 주로 서해 남부지역 및 제주 지역의 변화가 크게 나타난 것으로 요약될 수 있다. 소나무류(소나무, 곰솔)와 참나무류(신갈나무, 졸참나무, 굴참나무)와의 경쟁 현황 파악을 위해, 전체 수종이 아닌 소나무류와 참나무류 만을 대상으로 한 소나무류 상대중요치 변화를 살펴보았을 때, 경기, 전남, 전북, 충남이 가장 변화가 컸고, 제주, 강원, 경북, 경남이 변화가 상대적으로 적었다. 경기도의 경우 타 지역보다 참나무류와의 경쟁이 상대적으로 매우 큰 것으로 추정되며, 제주의 경우 소나무류의 전반적인 중요치 감소는 매우 크게 발생하였으나 참나무류와의 경쟁은 크지 않아 다른 요인이 관여하고 있는 것으로 파악되었다(Figure 2).
임분상층의 소나무 중요치 등급이 변화하는 것과 동시에 임분 전체의 생물종 구성도 변화하였다. NFI5 Class1 표본점이 NFI7에서도 Class1로 그대로 유지(소나무 순림이 유지)되는 표본점의 경우에도 평균적인 생물다양성 지수는 0.11로 약간 증가하였다. 또한 중요치 감소에 따라 등급이 변경(소나무 순림이 혼효림으로 변화)되는 경우에는 생물다양성 지수의 변화량이 더 크게 나타났다. 소나무가 분포하는 표본점 Class1~3이 주로 변화하는 구간(Table 2의 회색부분)은 전반적으로 생물다양성 지수가 증가한 것으로 파악되었다(Table 3).
| NFI7 | Class1 | Class2 | Class3 | Class4 | Average |
|---|---|---|---|---|---|
| NFI5 | |||||
| Class1 | 0.11 | 0.27 | 0.69 | 0.73 | 0.16 |
| Class2 | −0.08 | 0.09 | 0.22 | 0.35 | 0.12 |
| Class3 | −0.47 | −0.01 | 0.06 | 0.09 | 0.06 |
| Class4 | - | 0.07 | 0.14 | 0.04 | 0.04 |
전체 소나무류 등급 중에서 소나무 순림에 해당하는 Class1의 변화를 집중적으로 분석하였다. Class1 표본점에 대한 소나무류 중요치 등급 변화(이하 “소나무림 변화”)는 지형, 교란, 기후 등 여러 요인과의 관련성을 확인할 수 있었다. 소나무림 변화와 미변화 사이의 통계적 차이(p<0.05)가 확인된 변수는 사면향(Aspect), 지형학적수분지수(TWI), 소나무면적비율(PDPT100), 비산림면적비율(NF100), 소나무중요치(IV), 생물다양성지수(SI), 연강수량(PRCP) 이었다.
소나무림 변화 지역은 미변화 지역보다 사면향 지수가 낮고, 지형학적수분지수(TWI)가 높으며, 강수량(PRCP)이 큰 경향을 보였다. 즉, 지형적, 기후적으로 수분조건이 상대적으로 양호한 지역에서 변화가 많이 나타났다. 또한 소나무림 변화 지역은 주변의 소나무면적비율(PDPT100)이 낮고, 비산림비율(NF100)은 높으며, 해당 임분 내 소나무류 중요치(IV)가 낮고 생물다양성지수(SI)는 높았다. 즉, 임분 내・외부의 소나무의 중요치가 상대적으로 낮고 비산림 비율이 높아 외부교란 가능성이 많은 지역에서 변화가 많았다. 이를 요약하면, 교란요인이 많으며 수분조건이 좋은 지역에서 소나무림의 중요치 감소가 많이 발생했다고 볼 수 있다(Figure 3).
영급과 수관밀도 등 임분 변수에서는 소나무림 변화지역이 상대적으로 영급과 수관밀도가 낮았으나 통계적인 유의성은 낮은 것으로 분석되었다. 기존의 Kim et al.(2022)의 임상 변화 요인 분석에서도 영급과 수관밀도가 상대적으로 낮은 지역이 침엽수림에서 활엽수림으로 전환되는 것과 관련성이 높다고 분석된 바 있으나, 본 연구에서는 외부 교란 가능성이 높은 표본점을 제외하는 과정에서 영급이 낮은 표본점이 상당부분 제외됨에 따라 임분 특성은 세밀하게 고려하지 못하게 된 것으로 판단된다.
소나무림 변화와 미변화를 분류하기 위한 모형을 구축하기 위해 통계분석을 통해 유의성이 확인된 7개 변수 중에서 변수 간 상관성이 높은 변수를 제외하고 최종적으로 6개 변수를 모형 입력변수(Aspect, TWI, PDPT100, NF100, IV, PRCP)로 선정하였다. 이 변수들을 이용하여 RF, SVM (Linear), SVM(Radial) 모형을 적용한 결과, 소나무림 변화/미변화 분류 결과 추정치에 대한 K-fold Validation의 전체 정확도 평균은 0.72~0.74, Kappa계수 평균은 0.44~0.48로 유사한 것으로 나타났다. Sensitivity와 Specificity에서는 SVM(Linear)가 가장 균등하게 추정하는 것으로 확인되었으며 RF와 SVM(Radial)의 경우 약간의 편향성을 보였다(Table 4). 개별 모형의 불확도 영향을 줄이기 위해 세 모형의 결과를 조합하여 활용하도록 하였다.
| Model | Overall Accuracy | Kappa Coeff. | Sensitivity | Specificity |
|---|---|---|---|---|
| RF | 0.72 | 0.44 | 0.88 | 0.55 |
| SVM (Linear) | 0.74 | 0.48 | 0.74 | 0.74 |
| SVM (Radial) | 0.73 | 0.47 | 0.71 | 0.76 |
최종 도출된 RF, SVM (Linear), SVM (Radial) 모형을 이용하여 현재 기준으로부터 10년 후의 미래 소나무림 변화를 예측하고 세 모형에서 모두 변화로 추정된 경우를 소나무림 변화지역으로 판정하고 결과를 도출하였다. 그 결과, 전국적으로 현재 기준 소나무림 순림의 약 14.2%가 향후 10년 후에는 자연적인 임분동태 변화에 따라 혼효림으로 전환될 가능성이 높을 것으로 추정되었다. 지역적으로 보면 제주와 경기가 42.8%, 26.9%로 가장 높았고 경북과 강원이 8.8%, 13.8%로 가장 낮았다(Table 5).
광역시와 시군 기준으로 구획을 세분화해서 보았을 때, 전반적인 지역적 패턴이 나타났으며 경기, 충청, 전남 지역을 따라서 변화의 민감도가 큰 것으로 확인되었다. 경북, 경남지역의 변화 민감도는 가장 낮았다(Figure 5).
결 론
본 연구에서는 국가 규모의 현장조사 자료를 이용하여 우리나라 소나무림의 전반적인 변화 경향을 파악하고 지역별 차이를 비교하였다. 그 결과, 과거 10년 동안 우리나라 소나무림은 전반적으로 중요치가 감소하는 방향으로 변화해왔으며, 이와 동시에 임분 내 생물다양성이 점차 증가하고 있는 것으로 나타났다.
소나무림의 변화는 지형적, 기후적으로 수분조건이 상대적으로 양호한 지역에서 많이 나타났는데, 이러한 결과는 습윤한 지역에서 참나무류의 생육이 더욱 좋아지면서 기존에 점유하고 있던 소나무류와의 경쟁에서 이기게 되는 수종경쟁 메커니즘(Bae et al., 2012; Lee et al. 2015; Park et al., 2017)과 관련성 있는 결과로 볼 수 있다. 수종의 생육 경쟁에는 양분, 수분, 빛이 관여하는데, 양분과 수분 조건이 좋아져 빛 경쟁으로 수종 경쟁을 할 경우 음수인 참나무류가 양수인 소나무류의 수관을 결국 점유하여 도태시키는 현상이 발생하게 된다. 따라서 소나무가 경쟁력을 보이는 지역은 양분과 수분 조건이 열악하여 소나무가 다른 활엽수보다 생육 활동의 상대우위를 점유하거나, 수고가 매우 높이 자라 활엽수종의 임분 상층 점유가 어려운 지역이 해당될 수 있다.
또한, 소나무림 변화는 임분 내・외부의 교란 요인이 상대적으로 많은 곳에서 발생했다. 내부적으로는 임분 내 소나무의 중요치, 즉 다른 수종이 어느 정도 섞여 있는지에 따라 이후 변화의 동력이 발생한다는 의미이다. 즉, 소나무림일지라도 상층부가 모두 소나무로 구성되어 있는 임분 보다는 여러 가지 수종이 혼재되어 있는 경우가 변화의 가능성이 더 높으며, 외부적으로는 주변지역에 소나무림이 많은 지역보다는 다른 수종의 산림이 혼재되어 있는 경우, 순수한 산림지역보다는 비산림지역 등 외부 교란이 있는 지역에서 변화 가능성이 더 높다는 것이다. 이는 변화의 잠재력 측면에서 볼 수 있는 부분이다. 우리 산림은 전반적으로 생물다양성이 증가되는 방향으로 변화하고 있기 때문에, 이러한 소나무림의 변화의 잠재력 또는 점차 커지고 있다고 볼 수 있다. 즉, 과거 소나무림이 변화했던 속도보다 이후에는 더 빠른 속도로 변화할 수 있는 가능성이 있다.
이러한 특성에 따라 현재 전국적으로 분포해 있는 소나무림의 변화 취약성은 지역에 따라 다르게 예상되었다. 경기, 충청, 전남 지역을 따라서 변화의 민감도가 크고, 경북, 강원지역의 변화 민감도는 상대적으로 낮았다. 따라서 소나무림 변화 민감도가 큰 지역 내에 위치한 보호와 육성 대상 중요 소나무림의 경우, 좀더 강력한 소나무 보호와 인위적 관리 대책이 시행될 필요가 있다. 한편, 변화 민감도가 적은 지역의 경우 소나무 생육이 안정적이므로 기존의 소나무 중심의 산림 유지 및 활용 방안이 효율적이며, 변화 민감도가 높은 지역의 경우 타 수종으로 육성 방향을 전환하는 대책 마련이 필요하다.
본 연구의 소나무림 변화 분석 및 미래 전망 결과는 몇 가지 고려할 점들이 있다. 우선 본 연구결과는 자연적인 임분동태 변화에 따른 소나무림 변화 압력만을 추출하여 분석하였으므로, 병해충과 산불 등 교란 요인은 고려하지 않았다. 따라서 본 결과를 통해 산출된 변화 전망 결과와 현실에서의 소나무림 변화 특성은 차이가 날 수 있음을 유의해야 한다. 소나무재선충병에 따른 경남지역의 집중적인 소나무림 피해, 강원・영동지방의 산불로 인한 피해 등을 고려하면 소나무림의 피해와 감소량은 더 클 것이며, 지역적 변화 패턴도 다르게 나타날 수 있다.
또한 NFI5와 NFI7 표본점 자료 기준으로 소나무림이 과거 10년 사이 약 25% 정도가 혼효림으로 전환된 것과 비교했을 때 미래 추정치는 14.2%로 다소 낮게 평가 된 측면이 있다. 이는 표본점 자료 분석과 전국 대상의 공간자료 분석의 차이가 있을 수 있음을 고려해야 한다.
과거 산림통계와 비교하였을 때, 2000년 기준 침엽수림 면적이 약 2,711천ha이었고 2020년에는 2,320천ha로 감소하였는데 이는 평균적으로 10년에 약 7.5%의 침엽수림 면적 감소율이 있었던 것으로 계산된다. 물론 이는 전체 침엽수림에 대한 변화율이며 소나무림에 대한 변화율은 아니지만 한반도 산림의 상당한 부분의 침엽수림이 소나무로 구성되어 있다는 점을 고려하였을 때, 본 연구의 추정 결과가 과거 실제 침엽수림 변화율보다는 다소 크게 추정되었다고 평가할 수 있다. 그러나 실제 산림통계와의 비교 시에는 침엽수종의 수종 구성, 신규 조림 면적 반영 등 다양한 요인이 통계 안에 포함되어 있음을 고려해야 하며, 또한 우리나라 산림의 생물다양성이 점차 증가하고 있어 과거보다 미래에 소나무림의 변화가 더 가속화 될 가능성도 검토해보아야 한다. 한편으로, 표본점 변화 분석 결과에서 확인된 혼효림에서 소나무림으로 역방향으로 변화하는 지역에 대한 부분이 고려되지 못한 부분은 본 연구의 한계로 지적될 수 있다. 혼효림에서 소나무림으로 변화하는 표본점은 소나무림이 혼효림으로 변화하는 표본점 개수의 약 16% 정도로 적은 비중을 차지하며 이러한 변화는 산림생태계의 자연적 변화 과정에서 주로 나타나는 프로세스는 아니지만, 해당 지역에 대한 특징을 파악하고 적용하면 소나무림과 혼효림이 변화하는 양방향의 변화를 보다 명확하게 고려할 수 있을 것이다.
본 연구에서 제시한 자연적 임분동태 변화에 따른 소나무림 변화 민감도 평가 결과는 소나무림 감소에 대한 지역적인 변화압력의 차이를 도출함으로써 전국 소나무림에 대한 지역 기반의 입체적 관리 방안을 마련하는데 기여할 수 있다. 더 나아가 소나무림 변화에 대한 보다 세밀한 분석을 위해 병해충, 산불 등 다양한 교란 요인과의 통합적 분석, 소나무림과 혼효림의 양방향 변화에 대한 정밀 분석, 생물다양성 변화와 산림구성 변화의 관계에 대한 연구, 수종별 중요치 감소량에 대한 정량적 모형 구축 연구가 추가적으로 수행되어야 한다.
