广西防城金花茶国家级自然保护区建立对区域土地利用动态的影响

黄秋霞，廖南燕，杨元征,3,4，俎佳星,3,4，王嘉丽,3,4，蔡文华,3,4，杨 健,3,4，苏宏新,3,4，杨 婷

（1. 南宁师范大学地理科学与规划学院，南 宁530001；2. 广西壮族自治区防城金花茶国家级自然保护区管理处，广西 防城港 538021；3. 北部湾环境演变与资源利用教育部重点实验室，南宁 530001；4. 广西地表过程与智能模拟重点实验室，南宁 530001）

自然保护区的建立是保护珍稀动植物及其生境的重要举措（Maiorano et al., 2008; Lung et al.,2010; Radeloff et al., 2010），其在缓解土地利用变化、保护生态多样性方面起积极作用（Stein et al.,2008;Xu et al.,2017）。自1956年中国建立第一个自然保护区（鼎湖山国家级自然保护区，面积为11.33 km2）以来，自然保护区数量和面积不断增加，建设体系日益成熟（张镱锂等，2015），了解自然保护区是否有效实现保护目的逐渐成为热点话题，但其建设成效评估体系仍不完善（祝萍等，2018）。目前，评价保护区建设成效的方法很多，至今尚未形成统一标准。土地利用变化能够直观表现自然保护区保护成效（郭子良等，2021），但简单分析保护区成立以来的土地利用类型面积变化并不能达到有效评价的目的。从长时间序列深入了解保护区及其邻近区土地利用类型动态演变，可以综合评估保护区建设成效，为土地资源规划提供决策支撑，对维持生态系统平衡和资源环境的可持续发展具有重要意义（Lambert et al.,2017）。

土地利用是链接人类社会经济活动与自然生态过程的纽带（Liu et al.,2014），其利用方式转变对生物多样性（García et al., 2007）、生态系统服务（田义超等，2015；刘晓娟等，2019；Chen et al.,2020; Cao et al., 2021）、资源的可持续利用（姚静韬等，2020）以及气候变化都有深远影响（Bu et al.,2014;Huang et al.,2015）。近年来，国内外学者对自然保护区土地利用变化做了大量研究，主要利用地理信息系统技术，结合数理统计与景观生态学方法分析保护区建立后土地利用和景观格局动态变化（Zhang et al.,2017；张宏静等，2021），探讨造成该动态变化的驱动因素（Su et al., 2015；满苏尔·沙比提等，2016，2017），并预测未来土地利用变化趋势及其对生态系统服务的影响（井云清等，2016；Zhang et al.,2019;Li et al.,2021）。但目前研究较多集中在保护区内部，较少针对保护区内部与外部土地利用动态差异。

广西防城金花茶国家级自然保护区（简称“金花茶保护区”）位于广西壮族自治区防城港市境内，于1986年经广西壮族自治区人民政府批准，建立广西防城上岳金花茶自治区级自然保护区，1994年晋升为国家级自然保护区，主要保护对象为珍稀濒危金花茶组植物及其生境，是中国乃至世界金花茶植物的分布中心。目前对金花茶保护区的研究主要集中在植物多样性（庞洁等，2008）、植物资源保护及开发（谢钊等，2019）、昆虫区系分析（莫璐瑕等，2014）等方面，对其土地利用变化的研究比较少。保护区的建立使人类活动逐渐减少，对金花茶生境质量提升起一定的积极作用，但在当地经济发展的需求下，部分人类活动对金花茶保护区土地利用格局仍然造成影响（黄晓园等，2019）。保护区及邻近区域土地利用变化会威胁保护区生态安全与生态系统服务功能。因而，有必要探索保护区的建立对区域土地利用产生的动态影响。鉴于此，以1986—2019 年7 期Landsat5 TM/Landsat8 OLI 遥感影像为数据源，采用决策树分类方法进行土地利用分类，探讨金花茶保护区非核心区和保护区外部2 km邻近区的土地利用类型动态演变规律。以期为保护区的建设及监管、资源开发及生态安全保护提供一定的数据支撑。

1 研究区概况与方法

1.1 研究区概况

金花茶保护区位于广西壮族自治区防城港市防城区境内（李广超等，2021），总面积为9 098.62 hm2，其核心区面积为1 479.10 hm2，非核心区面积为7 619.52 hm2，非核心区包括保护区内实验区和缓冲区两部分，其中缓冲区面积为3 459.21 hm2，实验区面积为4 160.31 hm2，跨越那梭、华石、大菉和扶隆4个乡镇，地理位置为21°44′-21°48′N、108°03′-108°12′E（图1），保护区内人口约1 900人（实地走访调查所得），粮食作物以玉米和水稻为主，经济作物主要有红薯、蔬菜、八角、肉桂。保护区地处十万大山山脉，属热带季风气候，全年气候温和湿润，冬短夏长（李广超等，2021），年均温约21.8℃，年均降雨量高达2 900 mm。植被类型以阔叶混交林、亚热带常绿阔叶林为主。为更好地评估保护区建设与发展对当地生态环境的保护力度，将研究区划分为保护区内的核心区（简称“核心区”）、保护区非核心区（简称“非核心区”），保护区外部邻近区（简称“邻近区”）三部分（见图1）。其中，邻近区为在保护区外部建立2 km 的缓冲区，外部缓冲区面积与保护区面积大致相同。通过对比保护区外部邻近区与非核心区土地利用变化趋势，可为综合评价保护区建设成效提供数据支撑。

图1 研究区示意Fig.1 Spatial distribution of the study area

1.2 数据来源

Landsat 系列数据具有时间跨度长、重访周期短、覆盖范围广等优点，被广泛用于土地利用时空格局演变研究中。基于Google Earth Engine（GEE）云平台①https://earthengine.google.com/获取金花茶保护区及其邻近区1986、1994、1999、2004、2009、2014、2019 年7 期Landsat5 TM/Landsat8 OLI 地表反射率产品，影像行列号为125/45，空间分辨率为30 m。选取当年云量最小影像，消除云雨等因素的影响。该数据产品已进行大气校正与正射校正处理。

1.3 研究方法

因核心区人为干扰较少，土地利用状况稳定，主要对比非核心区和邻近区的土地利用变化。采用决策树分类方法进行土地利用遥感解译，并基于分类结果对比非核心区和邻近区的土地利用演变规律。为定量评价土地利用时空演变规律，选取单一土地动态度、综合土地利用动态度和土地利用程度综合指数表征土地利用变化动态度，利用数量转移轨迹与质心迁移轨迹定量展示土地利用时空演变特征。

1.3.1 土地利用类型解译 基于金花茶保护区土地利用实际情况，并结合Landsat5 TM/Landsat8 OLI遥感数据的特点，参考中国科学院土地利用/土地覆盖分类系统，利用ENVI 5.3 软件进行决策树分类，将研究区土地利用划分为林地、草地、耕地、水域、建设用地和未利用地6类（未利用地包括裸地、裸岩等），并依据野外调查数据和Google Earth高分辨率遥感影像进行精度验证。

1.3.2 土地利用类型动态演变 不同土地利用类型的变化速度可以通过土地利用动态度定量描述（邬亚娟等，2020），单一土地动态度K可用于分析某一用地类型在一定时间范围内的数量变化情况（张浚茂等，2019），综合动态度LC可以表征研究区不同时间节点全部用地类型的整体变化速率（张浚茂等，2019），计算公式为：

式中：Uai、Ubi分别表示研究初期和研究末期某单一土地利用类型i的面积；T为研究时长；LUi为研究期初期第i类土地利用类型的面积；ΔLUi-j为研究期内第i类土地利用类型转为其余土地利用类型面积的绝对值（张浚茂等，2019）。

土地利用程度可以反映人类活动对土地利用变化的影响程度（樊嘉琦等，2018），计算公式为：

式中：La为研究区土地利用程度综合指数，数值越大，表明该区域人类活动越强烈；Ri为地类i的土地利用程度的分级指数，参考自文献（邬亚娟等，2020）（表1）；Ci为地类i占研究区总面积的百分比；n为土地利用类型的数量。

表1 土地利用程度分级赋值Table 1 Grades of land use degree classification

1.3.3 土地利用数量变化轨迹 为认识非核心区和邻近区不同时段内各土地利用类型之间转移的方向和数量，计算不同时段的土地利用变化转移矩阵（张佰发等，2020）。转移矩阵的公式为：

式中：Yij表示某时期内第i种土地利用类型转化为第j种土地类型的面积百分比；n为土地利用类型的数量。

通过R语言中的NetworkD3 包，采用桑基图（Sankey diagram）将土地利用转移矩阵进行可视化，构建土地利用数量变化轨迹模型。桑基图是一种显示物质或者能量流动的可视化方法（姜婷婷等，2016），在材料分析、金融、生态学等多个领域被广泛运用。

1.3.4 土地利用迁移轨迹 为有效反映土地利用类型的时空迁移过程，采用标准差椭圆法（The Standard Deviational Ellipse,SDE）计算各地类的质心坐标。其质心坐标计算公式为（刘涛等，2013）：

式中：SDE(x，y)表示某一地类的算数平均中心坐标，即重心；（xi，yi）表示某地类第i个图斑的空间区位坐标；n为土地利用类型的数量。

利用ArcGIS 10.5 叠加各时段不同地类的重心坐标，得到邻近区和非核心区土地利用重心迁移轨迹，运用欧式距离方法计算得到不同土地利用类型重心迁移距离。

2 结果与分析

采用决策树分类法进行遥感解译得到金花茶保护区及邻近区1986、1994、1999、2004、2009、2014和2019年7期土地利用数据（图2），根据野外调查数据和Google Earth 高分辨率遥感影像对土地利用分类结果进行精度验证，7 期遥感图像的分类结果总体精度均＞85%，Kappa 系数均＞0.81，分类精度能满足研究需要。

图2 1986—2019年研究区土地利用Fig.2 Land use pattern in the study area during 1986-2019

2.1 土地利用总体结构变化

由7个时期土地利用时空分布情况可知，保护区（核心区和非核心区）土地利用类型主要为林地；未利用地和草地呈片状分布于中部、中偏西部海拔较高区域，呈现一定的集聚性；耕地呈点带状环保护区边缘较平缓地势分布（见图2）。邻近区土地利用类型主要为林地、耕地；林地成片状集中分布，东部、西部分布较多，南部分布较为破碎；南部有一条明显的耕地带，呈西南—东北走向，建设用地成片状分布，与耕地的分布具有一定的空间协同性。整体而言，金花茶保护区及其邻近区土地利用类型主要为林地、耕地，二者面积之和占总面积的93%以上（图3），呈现集聚性分布；林地面积分布最广，占比为76.73%～79.46%，主要分布在研究区中部及西南部；其次是耕地，面积比例为15.60%～18.36%，呈连片状分布，南部分布较多；水域、未利用地、建设用地、草地占总面积比例较少，面积比例均小3%，其中建设用地面积分布最少。近35 年间1986-2019 年土地利用类型发生较大变化，建设用地面积占比增长幅度最大，且呈现逐年上升趋势，面积占比从1986年的0.19%增长至2019 年的0.93%，2004—2009 年增长了0.29%，较为迅速；林地呈现波动性上升趋势，从1986 年的78.09%上升到了2019 年的78.82%；草地面积呈现略微增加趋势，面积占比从1986年的1.30%增长至2019年的1.68%；未利用地面积缩减幅度最大，面积占比从1986 年的2.02%减少至2019 年的1.07%，2004-2009 年面积占比减少了1.87%；耕地和水域呈现先增加后减少趋势，总体上面积占比分别减少了0.76%和0.13%。

图3 金花茶保护区及邻近区土地利用面积变化趋势Fig.3 Land use changes at various stages in Golden Camellia Reserve and its adjacent areas

2.2 土地利用动态变化特征

根据非核心区和邻近区单一土地利用类型动态度可以看出（表2），非核心区建设用地单一土地利用类型动态度居于第一，其中，建设用地1986-2019年平均动态度为161.60%，总体呈现显著扩张趋势，1994—1999年动态度最高为56%，2009年之后动态度为负值且呈现下降趋势；草地的动态度相对于其余用地类型一直处于较高水平，年平均动态度为1.84%，总体位居第二且呈现扩张趋势；林地呈现微增趋势；未利用地1986—2019 年平均动态度为-1.58%，呈现波动减少趋势；水域和耕地的动态度波动较平稳，呈微小降低趋势。邻近区建设用地的动态度一直处于较大值（见表2），1986-2019年平均动态度为11.42%，呈现扩张趋势；其余5 种土地利用类型的动态度波动较为平缓，林地呈现微小增加趋势，草地、未利用地、水域和耕地呈微小缩小趋势。综上，非核心区和邻近区建设用地的动态度最高，但非核心区的建设用地动态度更大，2个地区近35年的单一土地利用类型动态变化除草地外其他类型趋于一致，非核心区草地为正值且排名第二，邻近区草地则为负值。

表2 邻近区和非核心区单一土地利用动态度Table 2 The rate of land use change at various stages in adjacent areas and non-core areas%

1986—2019年，邻近区综合土地利用动态度曲线波动较大（图4-a），呈现先急剧减少再上升到最高值，之后波动下降的趋势。1994—1999年处于较低值，1999—2004 年上升到最高值0.63%，表明该时期人类活动对邻近区土地利用变化造成较大影响；2004年以后呈波动下降趋势，土地利用变化趋于平缓。非核心区1986—2019年综合土地利用动态度呈现先增加后递减的趋势，2004—2009年处于最高值0.55%，该时期土地利用变化剧烈。整体而言，邻近区综合土地利用动态度比非核心区变化更剧烈，表明保护区内部非核心区土地利用格局更为稳定，受人为活动影响较小，邻近区人地矛盾更为剧烈，2 个区域2009—2019 年综合动态度都趋于下降，表示土地利用变化都趋于平缓。

图4 邻近区和非核心区1986—2019年综合土地利用动态度（a）和土地利用程度综合指数（b）Fig.4 Comprehensive land use dynamic degree(a)and comprehensive index of land use degree(b)in the adjacent and non-core areas from 1986 to 2019

邻近区和非核心区土地利用程度综合指数变化趋势相似，均呈现先上升再下降，之后较为平缓波动趋势（图4-b）。邻近区在1999 年达到最大值229.58后略微下降，2004年以后平缓上升，这与该时期土地利用方式转换相一致，表示土地利用程度高的建设用地在扩张，未利用地等土地利用程度较低的地类面积占比减少，土地利用强度显著提高。非核心区土地利用程度综合指数在1999年达到最大值，随后呈波动下降趋势，该区域土地利用程度综合指数始终比邻近区低，表明非核心的土地利用受人为干扰较少，建设用地、耕地等高土地利用强度所占比例得到控制并逐渐减少。

2.3 土地利用数量转移轨迹

从金花茶保护区土地利用数量转换轨迹看（图5） （空间上的土地利用类型转移情况见文末附图1），邻近区1986—2019 年土地利用类型转换频繁（图5-a）。其中，林地转换量最为显著，主要流向耕地、草地，1986—1994 年流向耕地面积高达1 423.33 hm2，1994—1999 年变化相对较少，其他时期转移面积波动减缓；耕地以流出为主，主要流向为林地、建设用地、草地，1986—1994年流向林地面积最大，达到1 149.14 hm2；建设用地流入行为明显，主要由耕地、未利用地、林地流入，由耕地流入面积逐年增加。就土地利用类型转换的轨迹线而言，1986—1994 年不同地类之间转化最为活跃，转换轨迹多样化；1994—1999年不同地类之间转换相对较少；1999—2019年4个时段不同土地利用类型间的转变幅度较为一致。整体上，1986—2019年邻近区地类转换轨迹线比较稠密，各地类相互转换比较多样，不同类型间变化较活跃。

附图1 金花茶国家级自然保护区邻近区、非核心区1986—2019年土地利用变化情况Accompanying Fig.1 Land use change in adjacent and non-core areas of the Fangcheng Camellia National Nature Reserve from 1986 to 2019

非核心区1986—2019年土地利用类型之间转换较少（图5-b）。林地流入行为最为明显，主要由耕地、草地、未利用地流入，耕地转入面积最大，1986—1994 年流入面积达最大值212.75 hm2，其余年间转换面积相对较少；草地以流入为主，主要由林地、未利用地、耕地转入，2004—2009年由未利用地转入52.73 hm2；未利用地流出明显，主要转化为草地、林地，2004—2009 年转换为草地104.22 hm2；耕地、水域和建设用地与其他土地利用类型之间的转移变化相对较少。与邻近区对比，非核心区轨迹线较为稀疏，说明该区域土地利用类型变化较少，且转变类型更简单。综上，保护区的建设与管理显著效果，当地生态环境得到有效保护。

图5 1986—2019年土地利用数量变化轨迹Fig.5 Land use change trajectory from 1986 to 2019

2.4 土地利用时空迁移分析

土地利用重心迁移轨迹分析结果表明（图6），邻近区林地重心呈现三角形状迁移特征，总体向东北方向迁移，1999—2004年迁移幅度最大；耕地重心呈现“V”字型迁移特征，总体向西南方向迁移，1999—2004年迁移幅度最大；建设用地重心向西北方向迁移，1986—1994年迁移幅度最大；草地重心向东北方向迁移，2004—2009年迁移幅度最大；未利用地重心向东北方向迁移，1986—1994年迁移幅度最大；水域重心迁移幅度较小，向东北方向迁移。耕地与建设耕地总体迁移方向都为向西迁移。

图6 邻近区、非核心区土地利用重心迁移轨迹Fig.6 Migration trajectory of land use gravity center in adjacent areas and non-core areas

非核心区6种土地利用类型重心整体上都是向东北方向迁移，建设用地重心迁移幅度最大，在1986—1994年迁移距离最远；林地在1999—2004年重心迁移幅度最大；耕地重心呈现“W”字型迁移特征，2004—2009年迁移幅度最大；草地在1986—1994年重心迁移幅度最大；未利用地在2014—2019年重心迁移幅度最大；水域重心呈现“Z”字型迁移特征。

3 讨论

综上，保护区周边邻近区土地利用变化剧烈，土地利用类型转换频繁，非核心区土地利用变化较平缓，土地利用类型相互转换较少，表明保护区建立后管护政策的贯彻实施，管护力度不断加强，减少了人为活动的干扰，保护区的生态环境质量得到改善。1986—2019年非核心区和邻近区建设用地动态度变化最大，都为扩张趋势，主要由耕地流入，保护区为加强开发管理对自然保护区内的公路进行扩建，这是非核心区用地动态度变化的主要原因。邻近区土地利用变化主要受那梭农场和冲皇沟生态旅游区的发展建设影响，西南方向原有居民地向西北方向不断扩张，说明西部人类活动逐渐增多，社会经济活动活跃，需对保护区西部加强管理；退耕还林、沼气池建设工程等政策的实施，使非核心区和邻近区林地面积都有所增加，非核心区的林地因得到保护增幅更为明显，2 个区域林地都往东北方向迁移，与森林资源消耗的减少、东北方向耕地退耕有一定关联；两个区域的耕地均呈减少趋势，与建设用地、林地的动态演变具有一定协同性；非核心区未利用地主要分布在海拔较高山顶区域；随着保护区生态环境的改善，西南方向较低海拔的未利用地转换为林地和草地，邻近区的未利用地主要转换为林地与耕地。

从时间梯度看，邻近区综合动态度在1994—1999 年处于较低值，土地类型转换较少，这与90年代以后外出务工人员大量增加，土地利用受人为干扰下降有一定关系，当地居民通过务工收入水平显著提高，扩大住宅面积，土地利用程度高的建设用地在扩张，使得土地利用程度上升；1999—2004年，那湾生态园的建设、那梭农场特色产业的开发，建设用地在西部显著扩张，人类活动剧烈，土地利用类型之间转换剧烈，以致土地综合动态度较高，2002年退耕还林政策的实施，使得耕地面积减少，土地利用程度下降；2003年国家重点公益林补助试行之后，当地居民建设沼气池的积极性提高，减少了森林资源消耗，退耕还林政策的继续实施，有效控制了耕地的面积，但建设用地不断扩张，土地利用动态度趋于波动增长，土地利用程度平缓提升。非核心区受道路建设、退耕还林等政策的影响，土地利用变化较强烈，土地利用类型之间转换较频繁，1986—2009年土地利用动态度一直处于增长趋势，2009年之后呈现下降趋势，表明人类活动逐渐较少，土地利用格局趋于稳定。

整体而言，保护区的设立，退耕还林还草等工程的实施，管护工作力度的不断加强，使得保护区内部土地利用受人为干扰下降，土地利用格局更为稳定，建设用地、耕地等高土地利用强度所占比例得到控制并逐渐降低。同时，保护区的建立也促进了周边生态旅游、矿产资源、水资源的开发以及经济作物产业的发展，促使周边人类活动增加，土地利用变化剧烈，土地利用动态度一直处于较高水平。总体上，保护区的管理取得较为显著的成效，内部生态环境得到改善，但内部建设用地的缓慢增加与资源的不合理生产开发，仍然对保护区的生境质量产生负面影响。保护区在未来的管理建设中，应合理控制区内人类活动，强化土地集约化发展，有效实施生态移民搬迁、异地安置等政策。

金花茶保护区位于海拔较高区域，其地形特点在一定程度上也驱动着土地利用类型及其动态演变（崔步礼等，2011）。为更精准地评估保护区的建设成效，提取1986—2019年非核心区坡度＜6°的土地利用数据，并进行土地利用动态度的计算，1986-2019 年7 个年份土地利用强度分别为224.41、221.54、223.96、223.02、224.15、222.01、221.20，均小于邻近区的226.99、229.58、229.59、226.51、228.72、227.68、228.86，证明保护区的建设对生态环境的提高有促进作用。

4 结论

利用1986—2019年遥感影像对金花茶保护区及其邻近区进行土地利用遥感解译，并基于GIS空间分析与R语言统计分析方法，探讨了研究区近35年来土地利用动态演变规律，得到的主要结论有：

1）研究区土地利用类型整体分布具有一定的集聚性，保护区土地利用类型以林地为主，外部邻近区以林地和耕地为主，近35年间土地利用类型发生较大变化，建设用地面积占比增长幅度最大，未利用地面积缩减幅度最大，其余地类呈现波动增减趋势。

2）自保护区建立以来，邻近区土地利用变化强烈，非核心区土地利用格局更为稳定；非核心区和邻近区建设用地的动态度一直处于较大值，其余5 种土地利用类型的动态度波动较为平缓；邻近区综合土地利用动态度比非核心区变化更剧烈；邻近区和非核心区土地利用程度综合指数变化趋势相似，均呈现先上升再下降，之后较为平缓波动趋势。

3）邻近区近35年间土地利用类型之间转换频繁，相对而言，非核心区土地利用类型相互转换较少，2 个区域均是耕地与林地之间的转移交换最为剧烈。邻近区转换轨迹线比较稠密，各地类相互转换呈多样化，非核心区轨迹线较为稀疏。

4）土地利用重心迁移具有一定的方向性，非核心区6 种土地利用类型总体上都是向东北迁移，邻近区建设用地与耕地重心向西迁移，其余地类重心向东北迁移，表明建设用地、耕地在已有基础上不断扩张。

5）保护区的管理取得较为显著的成效，保护区内部土地利用受人为干扰下降，土地利用格局更为稳定，生态环境得到改善，但内部建设用地的缓慢增加与资源的不合理生产开发，仍然对保护区的生境质量产生负面影响。

本研究定量对比分析了金花茶保护区非核心区和外部邻近区在保护区建设与发展历程中土地利用动态演变，探究保护区的建立对区域土地利用动态的影响，可为保护区的建设及监管、资源开发及生态安全保护提供一定的数据支撑。