基于MCE-CA-Markov模型的森林景观格局演变和模拟预测——以宁乡市为例

陈佳楠，唐代生，贾剑波

（中南林业科技大学 a.林学院；b.水土保持与荒漠化防治湖南省高等学校重点实验室；c.南方林业生态应用技术国家工程实验室，湖南 长沙 410004）

森林景观是以森林生态系统为主体构成的景观，是景观生态学研究的一个重要分支。森林景观格局及其时空演变研究已成为森林景观空间分析的重要内容[1-2]。研究森林景观格局的基础是森林景观类型的划分，这也是选择森林景观模型的关键标准[3]。国外学者对景观生态学的关注起步较早，常用的景观分类系统有大尺度欧洲景观分类系统[4]、城市林业分类系统[5]、景观生态分类系统（湿地、森林和农业）[6]。此外，也有一些国外学者根据研究目的的不同，建立了其他的景观分类系统。Dumas等[7]根据森林与房屋的关系，把研究区内的森林景观分成5 种城市森林景观类型。Gond 等[8]使用多变量分析，将圭亚那盾的森林景观类型划分为33 类。针对我国森林景观类型的划分，国内大部分学者根据森林资源的特点，基于土地类型、植被型、优势树种等因子构建了森林景观分类系统，如周亚东等[9]根据地类、林种、优势树种和起源划分森林景观类型，胡美娟等[10]结合植被型划分森林景观类型，吴见等[11]则根据土地利用现状划分森林景观类型；或从森林景观发挥功能的角度，对森林景观类型进行分类，如何友军等[12]依据森林景观在城市中不同功能划分森林景观类型。

研究森林景观格局多运用景观指数法和景观动态模拟模型分析法[13-15]。景观指数法可从数量上反映森林景观格局的现状，借助景观动态模拟模型，能在更大空间尺度和时间维度对景观格局变化进行研究。常用的景观动态模型有Markov 模型、CA 模型、FORET 模型、LANDIS 模型等[16]。CA-Markov 模型将CA 模型模拟复杂系统空间变化的能力和Markov 模型长期预测的优势相结合，能够有效地模拟预测不同情景下景观格局的空间变化过程[17-19]。加入多准则评价（MCE）方法的MCE-CA-Markov 模型，能够量化自然和社会经济等影响因子，补充CA 模型的转换规则，使模型更加优化[20-21]，而基于MCE-CA-Markov 模型对森林景观格局的研究较少。

宁乡市作为湖南省会次中心城市、长株潭城市群的重要组成部分，是湖南省林业重点县，且于2018年启动省级森林城市创建工作，其自然地理区位、经济发展水平与森林资源特点典型性强。近年来城市建设及土地开发利用对其森林景观格局的演变产生了严重的影响，研究森林景观格局演变和模拟预测具有现实意义。本研究以湖南宁乡市为研究对象，以RS 和GIS 为技术手段，在确定森林景观类型的基础上，借助转移矩阵和景观指数，分析森林景观结构特点、景观格局动态变化特征，进而引入高程、坡度等因子，构建MCECA-Markov 模型，模拟预测未来的森林景观格局，揭示区域空间结构演变规律，目的是优化景观资源配置，为宁乡市森林城市创建提供参考。

1 材料与方法

1.1 研究区概况

宁乡市位于湘中北、洞庭湖南缘，介于111°53′～112°46′E、27°55′～28°29′N 之间，东靠长沙市望城区、岳麓区，南接韶山市、湘乡市，西与娄底市、安化县接壤，北与桃江县、益阳市毗邻，总面积2 906 km2，下辖4 个街道21个镇4 个乡，地处雪峰山余脉东缘，属湘中丘陵向滨湖平原的过渡地区，地势西高东低，以丘陵为主，兼有山地、岗地和平原。该区属亚热带大陆性季风湿润性气候区，四季分明，年平均气温16.8℃。境内水系发达，主要有沩江、靳江两大水系。

根据《湖南植被》区划，该区属中亚热带常绿阔叶林带的北部亚热带湘赣植被区，主要植被类型为针叶林、常绿阔叶林、针阔混交林、竹林等，主要树种有杉木Cunninghamia lanceolata、马尾松Pinus massoniana、湿地松Pinus elliottii、毛竹Phyllostachys edulis、枫香Liquidambar formosana等。据统计，2017年宁乡市林地总面积为1 286.83 km2，森林覆盖率为42.35%。其中：有林地为1 143.47 km2，疏林地为3.18 km2，灌木林地为87.11 km2，未成林地为11.82 km2，苗圃地为0.85 km2，无立木林地为29.68 km2，宜林地为10.72 km2。地理位置如图1所示。

图1 研究区地理位置Fig.1 Geographical location of study area

1.2 数据源与处理

林地与非林地在较长时空范围内相互转化，并统一于森林景观，参照《湖南省森林资源规划设计调查技术规定》地类划分，将研究区域的森林景观类型分为8 类（表1）。选取研究区2009、2013、2017年3 期Landsat 遥感 影像（表2）和DEM 数据，数据均来源于中国科学院计算机网络信息中心地理空间数据云平台(http://www.gscloud.cn)，分辨率均为30 m。对3 期遥感影像分别进行辐射定标、大气校正、影像镶嵌、影像裁剪预处理后，结合监督分类和目视解译，得到研究区不同研究年份的森林景观类型分类结果。采用随机抽样，创建随机点295 个，用同期的Google Earth在线地图验证样点森林景观分类结果，精度均为80%以上（表3），满足研究需要。坡度数据由DEM 数据在GIS 软件中的表面分析得到；GDP 数据来源于中国科学院地理科学与资源研究所(http://www.resdc.cn/DOI)；河流、道路、行政中心的矢量数据由Google Earth 上提取得到，并将其转换为栅格数据。所有栅格数据统一坐标系，栅格大小统一为30 m×30 m。

表1 森林景观类型的划分Table 1 Classification of forest landscape types

表2 遥感影像数据列Table 2 List of remote sensing image date

表3 宁乡市森林景观分类的精度Table 3 Accuracy of forest landscape classification in Ningxiang city

1.3 研究方法

1.3.1 转移矩阵

利用转移矩阵来分析景观格局动态变化能够有效地表现出某一时段内景观组分之间数量变化的关系，分析区域景观变化的结构特征和各类型变化的方向[22-23]，能直观反映研究初期阶段和末期阶段的景观分类结构，并体现研究期间内各景观类型的转变状况。

不同景观面积占比计算方法如下：

其中

式中：Pij为景观分类系统中基于各级分类的第i类景观在第j年的面积比例；Sij为景观分类系统中基于各级分类的第i类景观在第j年的面积；Ts为评价区域总面积；n为景观类型数；i、j表示研究初期与研究末期的景观。

1.3.2 景观指数

景观指数可以包括景观格局的有效信息，能够定量描述景观的差异。计算景观指数简单快捷，因此景观指数法格外受研究者青睐[24-26]。根据研究需求，从类型水平选取斑块数量（NP）、斑块密度（PD）、平均斑块面积（AREA_MN），从景观水平选取香农多样性指数（SHDI）、香农均匀度指数（SHEI)、聚集度指数（AI）、蔓延度（CONTAG）和连接度（CONNECT）共8 个景观指数，对宁乡市森林景观格局进行研究。选取的指数及生态学意义见表4～5，选取指数的计算公式和详细信息见文献[27]。

表4 类型水平指数的选取Table 4 Selection of the patch-level landscape index

表5 景观水平指数的选取Table 5 Selection of the landscape-level landscape index

1.3.3 MCE-CA-Markov 模型

1.3.3.1 CA 模型

CA 模型即元胞自动机（Cellular automata）模型，由元胞及其状态、元胞空间、元胞邻域、转换规则4 个部分组成。CA 模型具有时间、空间、状态离散型，任意元胞变量只存在有限且离散状态，并依据相同的转变规定进行同步修正，状态改变的规则在时间、空间上均为局部特征，普通CA 模型为：

式中：S为元胞有限、离散状态集合；N为元胞的邻居；f为局部映射元胞的转换规则。

1.3.3.2 Markov 模型

Markov 模型也被称为Markov 链，通过分析系统里每个状态的转移概率，可以预测对象的未来状态。由于该模型的无后效性，Markov 链广泛用于模拟土地利用动态变化[28]。在土地利用覆盖变化研究中，利用Markov 模型，预测土地覆盖变化的计算式为：

其中

式中：S(T)为T时刻土地利用结构状态；S(T0)为T0时刻土地利用结构状态；Pij为转移概率矩阵，研究初期到末期由类型i转移为j的概率，且0＜Pij＜1。

CA 模型和Markov 模型均为时间离散、状态离散的动力学模型，但是Markov 模型预测没有引入空间变量，而CA 模型的状态变量则与空间位置紧密相连，二者具有一定的局限性。将两个模型有效地结合，一方面可以更加容易地制定转换规则，另一方面可以受人为因素的影响变小，能同时从数量和空间上较好地对森林景观的分布进行预测，进而模拟景观格局的演变，具有一定的科学性和实用性。

1.3.3.3 MCE 模型

转换规则是CA 模型的核心，利用IDRISI 软件中的MCE 模块来生成各景观类型转移适宜性图集，进而确定元胞在下一时刻的状态，从而来定义演化规则或标准。MCE 模型包括两部分：约束条件和影响因子。约束条件是适宜性图像构建的特定限制规则，设置部分景观类型禁止向其他景观类型转化。影响因子则通过FUZZY 模块标准化，统一将连续的或离散的数据转化为0～255 值域内的标准化数据。把约束条件和标准化因子图层导入MCE 模块，使用Collection Editor 工具将各森林景观类型适宜性图像合并成一个适宜性图集。

水域景观作为稀缺资源，往往禁止向其他景观类型转化。建设用地景观作为已建成区域，往往很少转化为其他景观类型。因此，本研究将水域和建设用地作为约束条件。

影响因子的选择取决于影响景观类型的因子，本研究根据宁乡市森林资源特点和获得的数据，结合前人的研究，选取了高程、坡度、河流、道路、GDP 和行政中心6 个影响因子[29]。根据研究区内不同景观类型对各因子条件的适宜性不同，对不同因子条件的适宜性分级。其中，耕地在坡度小于5°时最适宜，5°～25°适宜性逐渐减小，大于25°为不适宜；在距河流1 000 m 内为最适宜，大于1 000 m 适宜性逐渐减弱；在距行政中心1 000 m内最适宜，1 000～3 000 m 内适宜，大于3 000 m适宜性逐渐减弱；建设用地在坡度大于25°时为不适宜；在距道路500 m 内最为适宜，500～1 000 m为适宜，大于1 000 m 适宜性逐渐减弱；乔木林地、竹林地、灌木林地、其他林地受限制因子较少，高程一般是其主要的影响因子；其他非林地可以任意转化。

2 结果与分析

2.1 宁乡市森林景观格局

遥感解译后，得到宁乡市2009、2013、2017年3 期的森林景观类型分类结果（表6）。经过统计分析，得到森林景观空间分布特征。从森林景观类型构成来看，在林地中以乔木林地为主；在非林地中，耕地在每期数据中占比均最高，高达51%。各类森林景观类型的构成见表6。

表6 宁乡市森林景观类型的构成Table 6 Composition of forest landscape types in Ningxiang city

结合图2可知，林地主要分布在宁乡市西部和东南部的山地，占44%～45%。其中乔木林地占比最高，占35%～37%；竹林地、灌木林地、其他林地占比较小。在非林地中，耕地主要分布在宁乡市东北部的平坦地区，为研究区内占比最高的森林景观类型，占49%～51%；主要河流从西部向东部流经宁乡市；而建筑用地大多分布在宁乡市东北部平坦地区，其余的建设用地零星分布在耕地之中；其他非林地多为裸土地，占比极小。

图2 2009—2017年宁乡市森林景观分类图Fig.2 Forest landscape classification map in Ningxiang city from 2009 to 2017

2009—2013年，各森林景观类型均有变化，耕地面积变化量最大，减少了45.01 km2；其他林地面积减少量次之，为33.71 km2；灌木林地则减少了15.65 km2。与此同时，乔木林地、竹林地、建设用地、水域、其他非林地分别增加了39.50、18.00、32.18、2.84、1.85 km2。总体变化量绝对值较前时间段相比有所减少。建设用地、灌木林地、水域面积分别增加了29.81、4.78、0.79 km2；其余森林景观类型面积都有所减少，其中其他林地面积减少了25.59 km2，乔木林地、竹林地、耕地和其他非林地的减少量都较少。从2009—2017年总体来看，其他林地和耕地急剧减少，建设用地和乔木林地大量增加，说明宁乡市在城市的扩张中占用了大量的耕地和林地，但同时也加强了对森林的经营与管理，使得乔木林地面积增加。

2.2 转移矩阵

为了使各森林景观之间相互转换过程得到更好的体现，整理得到不同研究时段的森林景观转移矩阵，表7～8 揭示了各森林景观之间的流转方向以及空间过程。由表7可知，在2009—2013年间，各森林景观类型面积转出和转入量各不相同。乔木林地主要流向了灌木林地、其他林地、耕地，其中流向耕地的面积最大，为112.17 km2；竹林地、灌木林地、其他林地主要流向了乔木林地和耕地。结合乔木林地面积增加的结果可知，实施封山育林政策和对荒山荒地植树造林后，其他林地大量转化为乔木林地，有效地增加了乔木林的面积。建设用地的增加主要来源于耕地和乔木林地，其中耕地的转化量更大。耕地除了转化为建设用地外，也向乔木林地转化了部分面积，说明退耕还林政策得到了有力实施。另外，水域和其他非林地变化量不大。

表7 2009—2013年宁乡市森林景观类型面积转移矩阵Table 7 Transfer matrix of forest landscape type area in Ningxiang city from 2009 to 2013 km2

由表8可知，2013—2017年各森林景观类型面积的转出和转入量也有所不同。乔木林地主要流向耕地、建设用地和灌木林地，其中乔木林地流向耕地最为明显；灌木林地则主要流向了耕地；其他林地主要流向乔木林地和耕地，其中其他林地转化为乔木林地的面积占2017年其他林地总面积的30.96%。耕地主要流向乔木林地和建设用地，与2009—2013年相比，流向建设用地的面积显著减少，说明对耕地的保护更为严格。与此同时，由乔木林地转化为建设用地的面积与前一研究时段相比也有所减少。竹林地、水域和其他非林地变化量不大。

表8 2013—2017年宁乡市森林景观类型面积转移矩阵Table 8 Transfer matrix of forest landscape type area in Ningxiang from 2013 to 2017 km2

2009—2017年，森林景观转出以林地和非林地中的耕地为主。建设用地的增加主要来源于耕地和林地，在宁乡市城镇化的扩张过程中不可避免地导致路网框架的扩大和居民住宅用地的增加，但是随着人们越来越注重保护生态环境，耕地和林地的减少速度得到了缓解。与此同时，随着封山育林政策的有效实施，森林群落演替良好，原有的其他林地等转化为乔木林地，且造林绿化工作的有效开展，减少了荒山荒地的面积，增加了乔木林地的面积。

2.3 景观格局变化特征

2.3.1 类型水平

在类型水平上，宁乡市的森林景观格局变化特征结果如表9所示。由表9可知，乔木林地的斑块数量最大，其余按从大到小排序为竹林地、耕地、水域、建设用地、其他林地、灌木林地、其他非林地；各类森林景观的斑块密度大小排序为乔木林地＞竹林地＞耕地＞水域＞建设用地＞其他林地＞灌木林地＞其他非林地；耕地的平均斑块面积远大于其它森林景观类型，乔木林地的平均斑块面积次之，其它森林景观类型的平均斑块面积相差不大。

表9 宁乡市类型水平景观指数值Table 9 The patch-level landscape index value in Ningxiang city

从2009—2017年宁乡市森林景观的斑块数量、斑块密度和平均斑块面积变化分析可知，乔木林地和耕地的斑块数量和斑块密度均增加，且平均斑块面积减少，表明其破碎化程度增大，多年的变化较为稳定；竹林地、灌木林地、其他林地的斑块数量和斑块密度均减少，表明其景观破碎化程度减弱，更加聚集；建设用地的斑块数量和斑块密度增加较为明显，表明在2009—2017年间，宁乡市的城市化建设较快，使得建设用地的数量增加；水域和其他非林地的斑块数量和斑块密度少量增加。

2.3.2 景观水平

在景观水平上，宁乡市的森林景观格局变化特征结果如表10所示。由表10可知，2009—2013年，宁乡市景观格局整体变化波动较大。蔓延度指数从2009年的57.692 4 下降到2013年的56.669 6，表明研究区各类景观连接度降低，分布分散。而2013—2017年各类指数变化不大，宁乡市整体景观格局变化相对稳定。从整体来看，香农多样性指数和香农均匀度指数在2009—2017年间先增加后减少，总的来看为增加；聚集度指数和蔓延度指数则表现为先减少再略微增加的趋势，总的来看为减少；连接度指数则先从0.229 3 增加到0.240 4，再略微减少到0.232 3，总的来看为减少。总体而言，研究区各景观类型在景观中呈均衡化趋势分布，景观聚集度下降，各类森林景观异质性朝着多样性发展。

表10 宁乡市景观水平景观指数值Table 10 The landscape-level landscape index value in Ningxiang city

2.4 宁乡市森林景观格局演变模拟

以宁乡市2013年森林景观分类图为初始年，利用适宜性图集以及2009—2013年森林景观转移概率矩阵，得到2017年宁乡市森林景观分类模拟结果（图3）。

图3 2017年宁乡市森林景观模拟分类结果Fig.3 Simulated forest landscape classification map of Ningxiang city in 2017

将2017年宁乡市森林景观模拟分类结果和2017年宁乡市森林景观遥感解译分类结果对比分析，对模拟结果精度进行评估，Kappa 系数为0.837 4，说明MCE-CA-Markov 模型的模拟效果较好，具有可行性。

2.5 宁乡市森林景观格局演变预测

以宁乡市2017年森林景观分类图为初始年，在CA-Markov 模块下导入2013—2017年马尔科夫转移概率矩阵和适宜性图集，间隔年份为4 a。循环次数设置为4，采用默认的摩尔型邻域5×5 滤波器。将以上参数设置完成后，获得2021年森林景观模拟分类图，模拟结果见图4。

图4 2021年宁乡市森林景观预测分类结果Fig.4 Predicted forest landscape classification map of Ningxiang city in 2021

由表11可知，与2017年宁乡市各森林景观类型相比，乔木林地、竹林地、灌木林地、其他林地、其他非林地面积均有减少，分别减少了36.07、6.70、9.84、7.55、0.27 km2，乔木林地减少的面积最多；与之相反，耕地和建设用地分别增加了26.8 km2和32.66 km2，水域面积也有少量增加。

表11 2021年宁乡市森林景观类型面积Table 11 Forest landscape area of Ningxiang city in 2021

由表12可知，2017—2021年期间，宁乡市增加的建设用地主要由耕地转化而来，与之相应地，增加的耕地主要由乔木林地和灌木林地转化得到，说明在城市建设的发展进程中，保护耕地的力度在不断提高。部分其他林地转化为乔木林地，说明城市的造林绿化工作在进一步推进，其他森林景观类型的转化量不大。

表12 2017—2021年森林景观类型面积转移矩阵Table 12 Transfer matrix forest landscape type area from 2017 to 2021 km2

3 结论与讨论

林地与非林地在较长时空范围内可相互转化，森林景观格局演变与非林地紧密相关。参照《湖南省森林资源规划设计调查技术规定》地类划分标准将宁乡市森林景观类型分为8 类，解译3 期Landsat 遥感影像得到2009—2017年森林景观类型分类结果。采用转移矩阵和景观指数方法，并结合高程、坡度等6 个影响因子构建MCECA-Markov模型，研究森林景观格局演变和模拟预测。

基于森林景观空间分布特征及转换过程，在2009—2017年期间，其他林地和耕地面积显著减少，乔木林地和建设用地面积增加，表明宁乡市城市化进程加快，且加强了森林经营管理，乔木林地面积增加。2013年后，耕地转出速度明显下降，表明耕地保护力度加强。同时得益于封山育林和荒山荒地造林，使得其他林地大量转为乔木林地。

从景观指数变化角度分析，在2009—2017年期间，乔木林地破碎化程度最大，8年间呈现先减弱后加剧的趋势；竹林地破碎化程度减弱；灌木林地破碎化程度也呈现出先减弱后加剧的趋势；耕地具有高聚集性；建设用地破碎化程度显著加剧；其他景观类型变化不显著。总之，研究区各景观在景观格局中呈均衡化趋势分布，景观聚集度下降，各类景观异质性朝着多样性方向发展。

结合高程、坡度等6 个影响因子，构建MCECA-Markov 模型，将2017年宁乡市森林景观类型的遥感解译结果和模型模拟结果进行对比，检验模拟结果的Kappa 系数为0.837 4，说明该模型对宁乡市森林景观格局模拟效果较好。2021年宁乡市森林景观格局预测结果与2017年各森林景观类型对比：乔木林地、竹林地、灌木林地、其他林地、其他非林地均减少，主要转化为耕地和建设用地，而建设用地随着城市发展持续增加，耕地受到严格保护。

本研究基于2009、2013、2017年3 期遥感影像数据，对其森林景观格局现状和变化过程进行分析；构建MCE-CA-Markov 模型并验证模型的可行性，模拟分析2021年宁乡市森林景观格局的变化情况，研究结果对于宁乡市森林城市创建具有一定的参考价值。由于研究区域范围及研究目的的不同，在数据来源和景观类型划分上有所侧重[12,14,30]。因本研究基于Landsat 遥感影像在县级空间尺度上将森林景观类型划分为8 类，未来研究可考虑选用更高分辨率的遥感影像数据，对各类森林景观类型划分进一步细化。在影响景观类型因子选择方面，结合前人的研究，只选取了高程、坡度、河流、道路、GDP 和行政中心6 个影响因子，可增加如人口[29]等其它影响因子进一步完善MCE-CA-Markov 模型，更好地预测未来森林景观格局变化。有关森林景观格局演变规律与驱动因素的相关关系等方面，还有待深入研究。