基于移动窗口分析法的矿区景观格局时空演变研究

梅昭容，李云驹，康 翔，魏善宝,潘剑君

(1.南京农业大学公共管理学院，南京 210095；2.国家磷资源开发利用工程技术研究中心，昆明 650500；3.南京农业大学资源与环境学院，南京 210095)

0 引言

矿产资源开采在给地区带来巨大经济利益的同时，也造成了土地不断被采场、排土场、堆占区等工矿用地取代，使得矿区成为景观要素在空间上变化最剧烈的地区之一。景观生态学将地表景观要素在空间上的排列描述为景观格局，认为其在很大程度上影响了生态学过程，可以用来推测生态学过程的特征[1-3]。我国矿产资源储量丰富，近年来随着经济发展引发的工业原料需求增加造成了采矿规模的进一步扩大，对地区生态环境及后续土地利用规划均造成了明显影响。在此背景下，深入研究景观格局的时空变化规律，对发掘自然与人类活动影响下景观的变化走向，动态调控景观过程的演变模式，以及为矿区的环境评估及政府部门制定土地利用与生态保护计划提供依据均有重要意义[4-6]。

随着地理信息系统(geographic information system，GIS)等数字化技术在景观生态学领域中的应用，利用特征指数分析区域景观格局变化与驱动机制的研究得以迅速发展，在不同时间与空间尺度均取得了许多成果[7-8]。其中，Whittaker[9]提出的移动窗口分析法有效实现了特征指数的可视化，在河谷沿岸的景观格局空间梯度[10]、水电站开发对景观梯度的影响[11]、城镇景观的异质性与空间格局变化[12-15]等研究中取得了一系列成果。由此可见将移动窗口分析法与时间序列上的景观指数相结合，不仅有助于直观了解景观格局的变化过程，也有利于进一步分析采矿等人类密集活动对周边景观格局的影响机制。但是目前针对矿区景观的研究仍以依靠景观指数进行宏观分析为主[4-6，16-17]，景观格局时空演变的分析十分缺乏。

本研究以位于昆明市附近的某露天磷矿区为研究区，基于2007—2017年3期遥感影像，结合移动窗口分析法与GIS技术，在确定窗口合理分析尺度的基础上，获取10 a来研究区景观破碎化与多样性指数的空间分布图。随后分析区域在露天采矿背景下的景观格局时空变化模式，及不同强度的采矿活动对景观格局改变的影响机制，以期为矿区的景观格局时空分析提供参考，并为环境监测与制定后续工作计划提供基础数据支撑。

1 研究区概况及数据源

1.1 研究区概况

研究区位于昆明市西南部，紧邻滇池南岸，地理位置在N24°40′～24°47′，E102°28′～102°35′之间，总面积约98.74 km2。区域内山地与高原盆地交错分布，地势中间高而四周低，中部山地保存着较为完整的大面积林地；城镇集中分布在南北两侧海拔较低地区，耕地较为破碎地分布其间。研究区地处亚热带高原季风气候带，年平均气温在16.5 ℃左右，年降雨量为1 450 mm。区域内磷矿资源丰富，开采条件优越，现有昆阳磷矿、海口磷矿与尖山磷矿3个主要矿区。其中规模最大的昆阳磷矿于1965年建矿，自滇池西岸向西，被划分为6个采区，目前主要开采第3—5采区，排土场与尾矿库复垦工作是后续工作的重点；海口磷矿始建于1966年，初建规模为4个采区，第3采区是目前主要的采矿段；规模最小的尖山磷矿建于2003年，自滇池西岸向西至白塔山，共由数十家大小不一的国营与民营矿山组成。自2007年以来的10 a间，在区位、政策与自然条件的影响下，各矿区均处在开采作业持续稳定的中期阶段，到2017年，已实现年原矿开采820万t，采剥总量3 500万m3的规模。矿山规模的扩大在带动了采矿及其衍生产业的高速发展的同时，也造成了地区生态压力的日趋加大。

1.2 数据源及其预处理

本研究选择的基础数据为2007年11月6日的CBERS-02B、2013年4月4日的资源三号(ZY-3)与2017年3月19日的高分二号(GF-2)卫星影像(下载自中国资源卫星应用中心网站http://www.cresda.com.cn)。将3景影像的全色波段与多光谱波段分别进行精确地理配准后，采用Gram-Schmidt法进行融合处理[18-19]，融合后影像的空间分辨率分别为2.36 m，2.1 m与1 m。参考国家土地利用分类标准，并利用面向对象分类法[20]，将研究区分为灌草地、林地、工矿用地、耕地、居民地以及水体6类主要景观，得到景观类型矢量图(图1(a)—(c))。为确保分类精度，选取650个有代表性的样本点(图1(d))与分类结果进行混淆矩阵分析(表1)，以评估分类结果。

(a)2007年景观类型 (b)2013年景观类型

(c)2017年景观类型 (d)样本点分布情况

图1 研究区2007年、2013年和2017年景观类型图、剖面线设置与样本点分布情况
Fig.1Landscapetype,transectsettingandsamplepointdistributioninstudyareain2007,2013and2017

表1 分类精度评价Tab.1 Precision evaluation of classified results

样本点由Google Earth影像目视解译得到，尽量保证其分布均匀且覆盖全区。将目视解译、全球定位系统(global positioning system，GPS)实地定位采样结果与企业提供的土地规划图等数据相结合，确定样本点的实际景观类型。由于灌草地与耕地、工矿用地与居民地容易被相互错分，因此选取样本点较多，650个样本点中，共选取150个灌草地样本点、120个耕地样本点、100个林地样本点、120个工矿用地样本点、120个居民地样本点及40个水体样本点。验证结果显示，面向对象分类法对各景观类型的分类精度均高于75%，3期景观类型图Kappa系数均大于0.85，显示出分类结果与实际景观的良好一致性，能够满足研究需要。随后在景观类型矢量图的基础上，利用景观指数粒度效应分析法[3，21-22]，确定最适宜分析粒度，将3幅景观类型矢量图转化为栅格大小为10 m×10 m的栅格图，进行后续分析。

2 研究方法

2.1 景观指数选取

根据该矿区的景观类型及分布特点，本文选择从景观破碎化与多样性2个角度衡量区域景观格局时空演变情况。参考前人研究[3，12-14]，选取斑块密度(patch density，PD)和最大斑块指数(largest patch index，LPI)衡量景观破碎化程度。其中：PD表示单位面积内的斑块总数，是景观异质性的直接反映；LPI用景观中最大斑块的面积除以景观总面积得出，其数值反映出景观中是否存在优势斑块。计算公式分别为

(1)

(2)

式中：A为景观总面积；i为斑块数量；an为景观中每个斑块的面积；Ni为景观中的斑块数量。

选取香农多样性指数(Shannon’s diversity index，SHDI)反映景观多样性程度，其值反映了景观要素的多少及各要素占景观的比例，具体公式为

(3)

式中：Pi为第i种景观类型占景观总面积的比例；m为景观类型总数。

2.2 移动窗口分析法

本文采用移动窗口分析法生成景观指数栅格图。先设定一个窗口半径，从左上角开始移动，每移动一格，计算窗口范围内的景观指数数值，并赋予中心栅格，形成景观指数栅格图。因此窗口半径，即分析尺度十分重要。此前有关研究表明，当窗口幅度过小时，相邻像元间的景观指数差距较大，使部分地区结果失真[12]。因此本研究分别以550 m，750 m，950 m，1 200 m和1 500 m为半径，在Fragstats4.2软件中进行筛选。

2.3 剖面线分析

为定量反映不同规模与强度的采矿活动对景观格局变化的影响机制，本研究两两连接3个露天矿区的采矿中心，建立3个剖面线，利用ArcGIS10.3软件中的Spatial Analysis模块将景观指数图中的信息提取至剖面线中进行分析。

3 结果与分析

3.1 景观格局基本特征

分析研究区不同时期各景观类型的斑块面积(total area，TA)与斑块所占景观总面积比例(percentage of landscape，PLAND)，有助于反映区域整体景观格局的分布现状及演变规律[10,12]，具体如表2所示。由表2可知，研究区内林地的面积较大，且其面积占比在10 a间始终保持在30%以上，是区域优势景观；灌草地则为面积变化最剧烈的景观类型，其面积呈明显的先减小后增加趋势；工矿用地面积在10 a间迅速增加，面积占比也在2017年超过了20%；其余景观类型面积较小且相对稳定。结合图1，本研究认为研究区内灌草地与林地的相互转化，以及工矿用地持续侵占周边灌草地，是景观构成变化的主要原因。

表2 2007—2017年间各类景观面积及所占比例Tab.2 Area of each landscape types and their percentages during 2007—2017

3.2 景观格局时空变化模式

3.2.1 特征尺度的选取

根据景观指数曲线的波动情况来选取特征尺度。由于景观指数在各剖面线中具有相似变化趋势，本研究选择其中一条剖面线进行分析。从图2中看出，曲线的波动幅度随窗口的减小而增大，当半径减少至550 m时，曲线出现明显的锯齿状波动，当半径增加至950 m时曲线趋于平滑。由于随着窗口半径的增加，数据量迅速加大，且窗口半径在1 500 m时景观指数曲线的波动幅度与950 m时差距不大，因此本研究选择将950 m作为移动窗口特征尺度。

(a)PD (b)LPI (c)SHDI

图2 不同尺度下的景观指数曲线

Fig.2Landscapeindexunderthedifferentscales

3.2.2 特征尺度下的景观破碎化时空演变

研究区2007—2017年间景观破碎化空间分布及变化状况如图3和图4所示。

(a)2007年PD分布 (b)2013年PD分布 (c)2017年PD分布

(d)2007—2013年PD变化情况 (e)2013—2017年PD变化情况

图3 2007—2017年间研究区PD分布及变化

Fig.3DistributionandchangeofPDinstudyareaduring2007—2017

(a)2007年LPI分布 (b)2013年LPI分布 (c)2017年LPI分布

(d)2007—2013年LPI变化情况 (e)2017—2013年LPI变化情况

图4 2007—2017年间研究区LPI分布及变化

Fig.4DistributionandchangeofLPIinstudyareaduring2007—2017

由图3—4可知，10 a间区域景观破碎化的高值区从四周向中心迁移，破碎化程度逐年显著，在前6 a间的增幅较为剧烈，后4 a的分布模式则趋于稳定。2007—2013年间，PD增加最迅速的地区为西部草地与海口磷矿、尖山磷矿周边，受陡坡种植、城镇扩张及磷矿开采侵占周边土地的影响，这2个区域的斑块被分割得十分破碎；PD降低最迅速的地区位于长期开展同类建设活动的昆阳磷矿与南部城镇。2013—2017年间，矿区的排土与运输工作造成的森林资源破坏导致位于3个矿区间的中部林地PD迅速增加，海口磷矿与尖山磷矿由于仍处于持续开采扩张中，PD也呈小幅上升趋势；昆阳磷矿则由于用地类型趋于稳定，PD持续降低。LPI计算结果表明，采矿活动强度加大导致3个工矿用地中心逐渐形成LPI高值区，斑块团聚效应逐年显著。区域LPI整体显著下降的趋势则可能与该地除采矿活动外，缺乏其他明显的经济建设活动有关。

3.2.3 特征尺度下的景观多样性时空演变

SHDI空间分布及变化情况如图5所示，总体呈现出中间低而四周高的分布格局，整个研究区的景观多样性逐渐变得十分显著。

(a)2007年SHDI分布 (b)2013年SHDI分布 (c)2017年SHDI分布

(d)2013—2007年SHDI变化情况 (e)2017—2013年SHDI变化情况

图5 2007—2017年内研究区景观多样性空间分布图

Fig.5Spatialdistributionmapoflandscapediversityinstudyareaduring2007—2017

2007—2013年间受农业与采矿活动侵占土地的影响，西部的耕地与尖山、昆阳磷矿边缘处的景观复杂程度迅速增加；2013—2017年间SHDI增加最剧烈的地区则主要分布于昆阳磷矿北部林地。而长期开展同类经济活动的露天矿区内部、南部建城区以及人类干扰少的山地最深处，景观多样性程度在10 a间持续降低，生态状况稳定。

3.3 采矿活动对景观格局变化的影响

经统计，10 a间3个主要矿山的采矿规模均有明显增加，其中昆阳磷矿的面积由约5.8 km2增加至约8.5 km2，涨幅为46 %；海口磷矿的面积由约3.7 km2增加至约5.8 km2，涨幅为56 %；尖山磷矿的面积由3.5 km2增加至约5.7 km2，涨幅为62 %。到2017年，3个矿山的年磷矿石产量分别为260万t，250万t与100万t。目前开采规模最大的昆阳磷矿开始进入开采稳定期，扩张速度相对稳定；而规模最小，正处在作业中前期的尖山磷矿，采矿强度增加迅速。为进一步探究不同强度采矿活动对景观格局影响的异同，本文采用剖面线分析法进行后续分析。

3.3.1 采矿活动对景观破碎化的影响

如图6所示，PD在各剖面线上的数值与波动幅度均逐年大幅提高，表明采矿活动对景观异质性程度的影响大幅上升。昆阳—海口剖面线与尖山—昆阳剖面线上，PD在采场与林地交界处达到最大值后开始减小，在约2～3 km处的林地中心达到最小值。海口—尖山剖面线上的PD走势在各年份间有较大变化，2007年PD的最高值出现在距海口磷矿2 km处的灌草地，与采矿活动无明显关系，2013与2017年，PD的最高值则出现在尖山磷矿周边1 km左右地区。与此同时，LPI的走势与PD大致相反。海口—尖山剖面线上破碎化程度显著，无明显的斑块团聚中心，昆阳—海口剖面线与尖山—昆阳剖面线上，则在2013年后在中心林地形成大面积LPI高值区。本研究认为，采矿活动是造成矿区周边地区景观破碎化的重要原因，从范围上看，扩张速度缓慢的昆阳磷矿的影响范围仅在周边1.5 km，而正处于采矿活动中、前期的海口与尖山磷矿的影响范围在2 km以上；从强度上看，昆阳磷矿的影响强度最大，海口磷矿则最小。各矿区对2.5 km以外区域的影响较小，区域内景观异质性最低，斑块面积较大，生态质量较好。

(a)海口—尖山剖面线PD (b)海口—尖山剖面线LPI

(c)昆阳—海口剖面线PD (d)昆阳—海口剖面线LPI

(e)尖山—昆阳剖面线PD (f)尖山—昆阳剖面线LPI

图6 景观破碎化指数沿剖面线的变化趋势

Fig.6Changeoflandscapefragmentationindexalongtranscets

3.3.2 采矿活动对景观多样性的影响

如图7所示，10 a间SHDI的上下波动程度在2007—2013年间增加明显，在2013—2017年间保持稳定，总体而言采矿活动对多样性的影响半径小于破碎化。在昆阳—海口剖面线与尖山—昆阳剖面线上，SHDI均在距矿区中心约1～1.5 km左右与林地的交界处出现第一个最高值拐点，随后迅速降低，在距各矿区约3 km处的中心林地达到最小值，最高值与矿区中心的距离随时间增加而增加。该走势表明采矿活动对景观组成的影响主要集中在矿区边缘，且影响范围随采矿规模的扩大而增加。在海口—尖山剖面线上无明显SHDI低值区，最高值在2013年后出现在距离尖山磷矿1.8 km左右的居民地，整条剖面线的多样性程度十分显著。总体而言，昆阳磷矿与海口磷矿的开采活动对区域景观构成的影响半径在1 km内，尖山磷矿的持续开采则对其周边2 km左右区域的景观组成均造成影响。

(a)海口—尖山剖面线 (b)昆阳—海口剖面线 (c)尖山—昆阳剖面线

图7 SHDI沿剖面线的变化趋势

Fig.7ChangeofSHDIalongtranscets

4 结论与讨论

本研究在GIS技术与景观生态学理论的支持下，结合移动窗口法空间分析技术，以昆明市某露天磷矿区为例，探讨了采矿背景下，区域景观格局在2007—2017年间的时空演变趋势，及采矿活动对景观格局变化的影响。与此前针对矿区景观的研究相比[4-6,16-17]，该方法有助于直观了解区域景观格局的演变模式并定量分析人类活动对其的干扰。

1)从景观整体构成上看，林地是该矿区最主要的景观类型，灌草地的面积变化剧烈，工矿用地面积持续增加，而其余景观类型的面积则相对稳定。由此认为灌草地与林地间的相互转化，以及工矿用地持续侵占周边灌草地是区域景观组成变化的主要原因。

2)从景观格局指数空间分布情况来看，研究区在2013年前后形成稳定的景观格局，受缺乏其他明显经济活动影响，10 a间除露天矿区与山体林地内部斑块团聚效应显著，景观类型单一外，整个研究区的景观破碎化与多样性程度增加明显。

3)分析采矿活动对景观格局的影响发现，采矿活动的影响取决于距离矿区中心的距离与采矿强度的大小。昆阳磷矿在1 km半径内对景观格局造成剧烈影响，海口磷矿与尖山磷矿的影响半径与影响程度则逐年增加。因此在后续工作中，尖山磷矿应注意对扩张速度的监控，海口磷矿应以实现矿山扩张与生态稳定为目标，昆阳磷矿则应注意保护外排土场边缘的林地，各矿均应减少对2.5 km外林地的侵扰。

研究结果的精度通常受诸多因素的影响。考虑到研究区大小，本研究仅依照一级分类标准分类，更为细致的分类标准是否会影响景观格局的精度仍值得进一步研究。移动窗口的尺度对研究结果也十分重要，随着窗口半径的增加，景观曲线趋于平滑，但是各斑块间的差异性有所降低，是否存在一个能使曲线更平滑，并精确反映区域景观异质性的窗口幅度，也有待进一步发掘。