土地利用景观格局演变作为土地利用/土地覆被变化（LUCC）在土地利用数量、结构和布局上的综合体现，是LUCC研究的重点内容之一。矿区景观是以人类采矿及其他生产、生活活动为主要特征的景观，煤炭资源的大规模开发利用在推动区域社会经济发展的同时，也造成了较严重的地面塌陷、矸石堆积等生态环境问题，加剧了矿区生态系统的脆弱程度；矿区塌陷地复垦治理工程的实施，对矿区生态系统的结构、组成和功能也会产生一定的影响，这些影响相互联系、错综复杂，既有直接的关系也有间接的影响。进行土地利用景观格局演变的驱动力分析，可揭示区域土地利用的演变规律及驱动机制，以调整土地利用政策，改善区域生态环境。

传统的土地利用变化驱动力分析研究主要运用因子分析、主成分分析、相关关系分析等方法，通过选择影响土地利用变化的各种自然与人文因子来建立多元线性回归模型，从而定量地进行土地利用覆被变化的驱动力分析。但当因变量是一个二分类或多分类变量而不是连续变量时，线性回归方法就不再适用，这时通常采用Logistic回归模型，当因变量是二分类变量时，采用二元Logistic回归模型，当因变量是多分类变量时，采用多元Logistic回归模型^[1]。多元Logistic回归模型主要应用于农牧交错带^[2]、矿区^[3]等生态脆弱区的土地利用变化驱动因素分析；二元Logistic回归模型主要用于耕地变化^[4]、居民点变化^[5]的驱动力分析、城市景观格局演变的驱动力分析^[6-9]，以及干旱区^[10]、矿区^[11-12]、流域^[13-18]等生态环境脆弱区的土地利用景观格局演变的驱动力分析；基于二元Logistic回归模型的结果常进行景观格局情景模拟。前人的研究结果表明Logistic回归模型能较好地揭示景观格局演变的主要驱动力及其作用机理^{[13, 16-17]}。

本研究以地面塌陷较严重的微山矿区作为研究区，基于某一特定时间段内矿区土地利用变化的分析结果^[19-22]，利用CLUE-S模型及二元Logistic逐步回归模型，进行各土地利用类型与驱动因子的回归分析，分析各驱动因子对解释变量的贡献大小，以预测未来土地利用变化趋势^{[7, 10-17, 23-25]}，从而为景观格局优化及当地政府制定科学合理的土地利用策略提供理论依据，对土地资源的合理开发、利用与保护也具有重要的现实意义。

1 材料与方法 1.1 研究区概况

微山县位于山东省济宁市南部，地处东经116° 34'~117°24'、北纬34°27'~35°20'，南北相距120 km，东西相距8~30 km，总面积1 779.8 km²；其中湖面面积1266 km²，占全县总面积的三分之二，由北到南依次为南阳湖、独山湖、昭阳湖、微山湖，统称南四湖。区内辖3个街道、10个镇、2个乡和1个县经济开发区（2014年行政区划）。在选取研究范围时，主要选择煤矿较集中、地表塌陷较严重的地区，同时考虑研究区的连贯性，最终选取微山县的11个乡镇作为研究区，总面积1 176.86 km²。研究区地理位置及范围见图 1。

1.2 数据来源与处理

选取2000年8月21日、2005年9月4日、2010年9月18日的Landsat5 TM遥感影像和2016年9月2日的Landsat8 OLI遥感影像作为数据源，其云量均小于2%，空间分辨率为30 m；其他辅助数据包括：微山县土地利用现状图、微山县矿区分布图、微山县土地利用总体规划（2006—2020年）、微山县矿产资源总体规划（2008—2015年）、济宁市采煤塌陷地治理规划（2016—2030年）等。

参照《土地利用现状分类和编码》（GB/T 21010— 2007），考虑遥感影像的分辨率，结合研究区地面塌陷积水较严重的实际情况，将研究区景观类型分为耕地、其他农用地、城乡建设用地、塌陷积水区、水域、滩涂沼泽6类。利用ENVI 5.3，采用监督分类中的支持向量机分类，人工目视判读与决策树分类，NDVI、MNDWI等指数模型进行影像分层分类提取，各种分类方法综合运用，以提高分类精度，最终得到2000、2005、2010、2016年4期土地利用类型分布图（图 2）。在4期分类图上分别随机选取200个检查点，通过实地调查并结合相近年份的土地利用现状图和天地图网站上的高清卫星图，获取检查点的实际土地利用类型，利用ENVI软件计算混淆矩阵和Kappa系数，经验证Kappa系数均在0.8以上，能满足后续研究的需要。

解译后的30 m×30 m的栅格大小超出了CLUE-S模型运行对数据量的上限要求，故本研究将解译后的30 m的栅格数据分别重采样为60、90、120、150 m共4种尺度，以进行各土地利用类型与驱动因子的二元Logistic回归分析。

1.3 驱动因子的选择及数据处理

为了解释矿区土地利用与其驱动因子之间的关系，需从影响土地利用变化的诸多因素中准确筛选出能够反映研究区土地利用变化的驱动因子，并进行驱动因子与各土地利用类型之间关系的定量分析。

基于2005年微山县土地利用现状图，提取主要铁路公路、主要河流、城乡居民点、工矿、沟渠等矢量数据，在ArcGIS中转化为栅格数据，并分别重采样为60、90、120、150 m，利用空间分析工具中的欧式距离计算每个栅格像元距主要铁路公路、距主要河流、距城乡居民点、距工矿、距沟渠的最近距离，以得到5个距离驱动因子。

高程、坡度等地形因素对土地利用景观格局演变的影响也较大，本研究高程数据下载的ASTER GDEM V2数据，其空间分辨率为30 m，需在ArcGIS中分别重采样至60、90、120、150 m，并利用栅格表面分析工具得到不同尺度下的坡度图。

将不同尺度下的各土地利用类型数据及7个驱动因子数据进行二元Logistic回归分析，采用受试者工作特征曲线（Receiver operating characteristic curve，简称ROC曲线）对各尺度下回归分析的结果进行精度验证。经过对比分析，60 m下重采样过程中，各土地利用类型数据及驱动因子数据信息损失量较小，同时回归分析结果精度较高，故本研究采用60 m×60 m的栅格大小进行驱动力分析。

60 m下各驱动因子栅格图如图 3所示。

1.4 各土地利用类型与驱动因子的回归分析模型

回归分析作为标准的统计学分析方法，在诸多行业和领域的数据分析中都发挥着非常重要的作用。在回归分析中，要求被解释变量为数值型变量，解释变量可以为分类变量，也可以为连续变量；当被解释变量是二分类变量时，通常采用二元Logistic回归模型，即被解释变量的取值只包含两个类别，如好或坏，发生或不发生，常用Y=1或Y=0表示，X表示解释变量，则P（Y=1|X）表示在X的条件下Y=1的概率。某种土地利用类型在特定空间位置上出现时数值为1，不出现时为0，就是一个典型的二分类变量。故本研究采用二元Logistic逐步回归模型分析计算各种土地利用类型在特定空间位置的发生概率P。

Logistic回归模型的数学表达式^[26]：

(1)

式中：P（/1-P）称为事件的“发生比”（Odds），用Ω表示，即某事件发生概率与不发生概率之比；β₀为常量；β₁为解释变量与被解释变量之间的相关系数。将上式进行数学变换，得：Ω=P（/1-P）=e^{β₀ + β₁X}

在二元Logistic回归分析中，将研究区划分为许多栅格单元，各土地利用类型作为被解释变量，引起土地利用变化的驱动因子作为解释变量，运用Logistic逐步回归计算各种土地利用类型在特定空间位置的发生概率，以分析各土地利用类型与驱动因子之间的关系，其数学表达式为：

(2)

经过数学变换，可得到如下函数关系式：

(3)

式中：P_i为每个栅格单元可能出现某种土地利用类型i的概率；X_1，i~X_n，i为与该土地利用类型i相关的各驱动因子；β₀为常量；β₁~β_n为各驱动因子与土地利用类型i之间的相关系数。

将所有驱动因子的栅格数据文件利用ArcGIS软件转化为ASCⅡ文件格式，再转化为CLUE-S模型能用的sc1gr^*.fil文件格式；利用解译后的2005年60 m的土地利用类型分布栅格数据，将6种土地利用类型分别提取出来，有这种地类的空间位置设置为“1”，没有的设置为“0”，存为6种地类的栅格格式文件，再转为ASCⅡ格式，进而制作CLUE-S模型用的cov0.0、cov1.0、cov2.0、cov3.0、cov4.0、cov5.0共6种单一地类文件。将6种cov^*.0文件和驱动因子sc1gr^*.fil文件导入CLUE-S模型中，利用File Convert模块进行数据转换，生成stat.txt文件，以进行下一步的回归分析。

在进行二元Logistic逐步回归时，自变量筛选方法采用基于最大似然估计的向前逐步回归法（Forward：LR），变量进入方程的依据是比分检验统计量，剔除出方程的依据是最大似然估计原则下的似然比卡方^[1]。步进概率中Entry和Removal分别为0.05和0.10，表示回归系数比分检验的概率P值小于0.05时，相应变量是与土地利用分布格局相关性较高的因子，可进入回归方程；概率P值大于0.10时，表示相应变量对解释土地利用格局影响不显著，应剔除出回归方程。

2 结果与分析 2.1 土地利用类型数量变化特征分析

从表 1可知，2000—2016年耕地、其他农用地和水域是研究区的主要土地利用类型，3种土地利用类型约占区域总面积的80%。

从表 2可知耕地2000—2005年、2005—2010年、2010—2016年面积呈现减少、增加、减少的趋势；城乡建设用地2000—2005年单一动态度数值最大，面积增幅最大，2005—2010年面积呈缓慢增加趋势，2010—2016年间城乡建设用地单一动态度为负值，面积有所减少；塌陷积水区2000—2005年面积增幅最大，2005—2010年面积增幅变缓，2010—2016年面积减少，整个研究期内塌陷积水区面积共增加4.73 km²；其他农用地和水域面积总体呈增加趋势，滩涂沼泽呈减少趋势。

2.2 各驱动因子与土地利用类型的回归结果分析

回归系数的显著性检验采用Wald统计量，如果某驱动因子Wald统计量越大，且对应的概率P值小于给定的显著性水平0.05时，则认为该驱动因子对土地利用类型的变化影响越显著^[5]。各驱动因子与6种土地利用类型的回归结果见表 3~表 6。

2.2.1 回归结果分析

由表 3~表 6可以看出，距居民点的距离是影响耕地及城乡建设用地变化的最主要驱动因子；影响塌陷积水区变化的驱动因子中，距工矿的距离和高程的Wald统计量较大，故这两个因子是影响塌陷积水区变化的最主要驱动因子；距居民点的距离也是影响水域变化的最主要驱动因子，究其原因，微山湖作为境内的湖泊淡水资源，人们利用距离居民点较近的湖泊淡水资源进行水产养殖，增加收入来源，故水产养殖等其他农用地面积的增加会导致湖泊水域面积的减少。

表 3~表 6中的回归系数β是由Logistic回归方程得出的相关系数。回归系数β有正有负，正负值表示各个驱动因子X对Logit（P）的作用方向。β值为正，说明解释变量X_i增加一个单位时，土地利用类型发生比随X_i的增加而增加；β值为负，说明解释变量Xi增加一个单位时，土地利用类型发生比随X_i的增加而减小^[26]。

从表 3可知，耕地与距居民点及河流的距离、高程、坡度呈负相关，说明随着距居民点及河流的距离增大，高程、坡度的增加，耕地面积呈减少趋势；耕地与距道路、沟渠的距离呈正相关，说明主要道路周边及离沟渠较近的地区是耕地面积增加的主要区域。

由表 4可知，其他农用地与距主要道路、河流的距离呈负相关，说明离主干道及主要河流越远，其他农用地面积呈减少的趋势，特别是其他农用地中的水产养殖用地，离河流越远，面积会大幅减少。

由表 5可知，城乡建设用地与距居民点、工矿用地及主干道路的距离呈显著负相关，说明居民点、工矿用地及主干道路周边是城乡建设用地面积扩张的主要区域范围；塌陷积水区与距工矿用地、主要河流、道路的距离呈显著的负相关，说明工矿周边及河流、主干道路周边是塌陷积水区面积增加的主要区域。

表 6的结果表明，水域及滩涂沼泽与距主要沟渠、道路的距离呈显著负相关。

2.2.2 回归结果检验

采用ROC曲线对驱动因素的解释能力进行检验。ROC曲线是以真阳性率（灵敏度）为纵坐标，假阳性率（1-特异性）为横坐标绘制的曲线。通过计算ROC曲线下的面积可以衡量模型的精度，曲线下的面积值介于0.5~1之间，面积值愈接近1，说明诊断或预测的精度愈高。一般情况下，ROC曲线下的面积值大于0.7，表明所选因子有较好的解释能力。图 4为各地类与驱动因子回归分析结果的ROC检验图。

从图 4可以看出，除其他农用地略低（ROC曲线下面积值为0.639），各土地利用类型的ROC曲线下的面积均超过0.7，其中城乡建设用地、水域、塌陷积水区、滩涂沼泽和耕地ROC曲线下的面积值分别为0.940、0.903、0.815、0.795和0.793，表明进入回归方程的因子对城乡建设用地和水域这两种土地利用类型的空间分布格局具有很好的解释能力，各地类的概率分布与实际的地类分布较一致，回归模型模拟效果较好；其他农用地ROC曲线下的面积值略低的原因可能是其他农用地包括的地类较多，所以其模拟效果稍差。

3 讨论

利用遥感及地理信息系统技术进行数据源及驱动因子的数据处理，以某种土地利用类型在特定空间位置上是否出现作为因变量，采用二元Logistic逐步回归模型，可有效地对引起各土地利用类型空间变化的驱动因子进行定量分析，以揭示各土地利用类型与驱动因子之间的相互关系及驱动机制。

矿区土地利用景观格局是人类在长期的生产生活过程中逐步演化积累形成的，其受人口、社会经济发展的影响较大，所以不同的社会经济形势下，土地利用景观格局的形成和演变机理随着社会化进程的加快也在逐渐发生改变，其主导因子也相应地发生变化。本研究在选取驱动因子时，由于数据获取的限制，仅选取了7个驱动因子进行回归分析，未来在开展类似研究时，应尽可能多地选取自然、社会经济等驱动因子进行土地利用景观格局演变的驱动力分析。

4 结论

（1）2000—2016年研究区耕地面积呈波动式减少，城乡建设用地及塌陷积水区面积2000—2005年增幅较大、增速较快。

（2）耕地的变化受距居民点的距离、距河流的距离等社会经济、区位可达性因子的影响较大，同时也受高程、坡度等自然驱动因子的影响，两者综合作用使得耕地发生变化。城乡建设用地的变化主要受社会经济驱动因子的影响较大。塌陷积水区的变化是区位因子、自然、社会经济驱动因子综合作用的结果。

（3）利用SPSS软件对采煤塌陷区2000—2016年各地类与驱动因子进行Logistic回归建模分析，并利用ROC曲线对回归结果进行检验，除其他农用地略低，各土地利用类型的ROC曲线下的面积均超过0.7，表明进入回归方程的因子能较好地解释土地利用景观格局的变化。

（4）基于CLUE-S模型及二元Logistic逐步回归模型进行矿区土地利用变化的驱动力分析，从引起各土地利用类型变化的影响因子中筛选出最主要驱动因子，并分析各驱动因子对各地类变化的影响程度，可为当地政府制定科学合理的土地整治、保护等相关政策提供科学依据。