西山煤田生态变化的遥感监测与评价

周玄德，窦文章,,，李国栋

(1.新疆大学资源与环境科学学院，新疆 乌鲁木齐 830046；2.北京博雅方略旅游研究院，北京 100044；3.北京大学软件与微电子学院，北京100871；4.北京大学战略研究所，北京 100871)

能源是人类历史发展过程中不可或缺的要素之一，其中煤炭资源发挥着极其重要的作用。在我国发展过程中，煤炭资源助推了社会经济等各方面的发展，然而伴随着煤炭能源的枯竭以及在使用过程中带来的生态环境等问题，日益引起人们的重视。大量的煤炭开采区，均面临着大气污染严重、空气质量差等问题，给人们的生活带来极大的影响。很多学者展开了对煤炭区域的生态环境治理对策等研究工作。刘佳骏等[1]对中国主要煤炭基地生态环境脆弱度进行了判别，分析了生态环境的约束因子，给出了不同煤炭绿色开采的措施;张海娜等[2]从煤炭资源型城市入手，对引发生态环境问题的原因进行了分析，提出了改善煤资源型城市生态环境问题的对策。在对煤炭开发区生态环境问题定性分析的同时，很多学者也展开了大量的定量评价性的研究，其中主要以煤炭开发区的生态风险评价体系的构建及实证分析为主[3]。万伦来等[4]通过系统动力学原理构建了煤炭开发对生态系统的胁迫模型，并对淮南市进行分析；张思锋等[5]通过半定量研究给出榆林煤炭开采区的风险源区域；路苹[6]、王广成等[7]构建了煤炭区生态环境的评价体系并展开研究。关于煤炭开发区的工作还包括对土壤修复、生态修复等方面的研究[8]，生态修复技术体系为我国煤炭基地的可持续发展提供了建议。

随着科学技术的进步，多源遥感数据可以实现实时生态监测，操作相对简单、成本较低，可以很好地用于煤炭开发区生态环境的监测[9]。基于此，本文将选取山西省太原市西山煤田作为研究区，通过遥感方法对生态环境变化进行监测，分析时空变化特征，揭示生态指数的变动规律。

1 研究区概况

西山煤田位于山西省中部，太原市西侧，主要以炼焦煤为主，是山西省六大煤田之一。本文选取了西山煤田的东部主要区域，范围为东经111°58′4″～112°30′0″，北纬37°24′0″～38°2′24″，面积3 250.71 km2，主要涵盖了古交矿区，还包括大量的人类活动较频繁的居住区，其中北部主要是古交矿区，南部为居住区，中间区域主要是地形相对复杂的过渡区域，人类活动较少。

2 数据来源与研究方法

数据选取2004年9月22日、2010年9月23日的Landsat TM影像，2016年9月23日Landsat8 OLI影像，行列号为125/034，分辨率30 m，日间数据，云量较少。数据来源于中国科学院计算机网络信息中心地理空间数据云平台 (http:∥www.gscloud.cn)。

生态指数的评价选用徐涵秋等人提取的遥感生态指数(RSEI)的计算方法[10-12]。该方法主要是由绿度、湿度、热度、干度四个指标构成，即通过这四个指标的合成构建遥感生态指数。其中，绿度指标是指植被指数，湿度指标是指湿度分量，热度是通过热红外波段反演的地表温度，干度是由建筑指数和裸土指数构成，详见公式(1)～(11)。

NDVI=(ρNIR-ρR)/(ρNIR+ρR)

(1)

WET=C1ρB+C2ρG+C3ρR+C4ρNIR+

C5ρSWIR1+C6ρSWIR2

(2)

NDBSI=(IBI+SI)/2

(3)

IBI=

(4)

(5)

式中：NDVI、WET、NDBSI分别为绿度、湿度、干度，主要通过波段运算计算得到；ρB、ρG、ρR、ρNIR、ρSWIR1、ρSWIR2分别为Blue、Green、Red、NIR、SWIR1、SWIR2波段的反射率；C1、C2、C3、C4、C5、C6依次为计算湿度的系数，在TM影像中，C1=0.0315、C2=0.2021、C3=0.3102、C4=0.1594、C5=-0.6806、C6=-0.6109，在OLI影像中，C1=0.1511、C2=0.1973、C3=0.3283、C4=0.3407、C5=-0.7117、C6=-0.4559。

LST=T/[1+(λT/ρ)lnε]

(6)

(7)

L6/10=gain×DN+bias

(8)

式中：LST为地表温度;T为传感器温度值;K1、K2为定标系数；L6/10为TM/TIRS的热红外波段的辐射值，对应B6、B10波段；DN为像元灰度值;gain、bias依次为波段增益值、偏置值；λ为热红外波段的中心波长；ρ为波尔兹曼常数；ε为地表比辐射率，取值参见文献[11]。在TM影像中，K1=607.76 W/(m2·sr·μm)，K2=1 260.56 K，gain=0.055,bias=1.18243，λ=11.4350 μm；在TIRS影像中，K1=774.89 W/(m2·sr·um)，K2=1 321.08 K，gain=3.342，bias=0.1,λ=10.9000 μm；ρ=1.438×10-2m·K。

RSEI=f(G,W,T,D)

(9)

RSEI0=

(1-{PCA{[f(NDVI,WET,LST,NDBSI)]}})

(10)

(11)

式中：RSEI为遥感生态指数，是通过四个指标综合得到；PCA为主成分分析，PC1、PC2、PC3、PC4依次为主成分分析后得到的第一、第二、第三、第四主成分；RSEI0为遥感生态指数初始值，即对PCA进行正规化处理后的值；RSEI为RSEI0归一化后的值；RSEI0-min、RSEI0-max分别为RSEI0的最小值、最大值。

表1显示，在对绿度、湿度、热度、干度四个指标进行主成分分析发现，PC1、PC2的累计贡献率都在85%以上，能够包含反应遥感生态指数的主要信息，因此选用PC1、PC2根据其贡献率，加权求和得到研究年份的遥感生态指数，见图1。

表1 各主成分特征值Table 1 Characteristic values of principal components

图1 生态指数空间分布图(2004年、2010年和2016年)Fig.1 Spatial distribution map of ecological indices in 2004,2010,2016

3 结果分析

3.1 生态指数的基本统计特征

根据式(1)～(11)，分别计算了2004年、2010年、2016年绿度、湿度、热度、干度、遥感生态指数的相关统计值，见表2。研究发现，各类指标均表现为一定程度的波动变化。绿度指标2004年为0.345，2010年为0.476，增长了0.131，2016年为0.435，相比2010年下降了0.041，整体上表现为先上升后下降的倒“V”型；湿度指标2004年为-0.095，2010年为-0.09，2016年为-0.03，相比2004年增长了67.87%，保持着持续增长的趋势；地表温度指标在2004年、2010年、2016年的值依次为25.068、19.903、26.868，整体表现为先下降后上升的“V”型变化趋势；干度指数也由2004年的0.096下降到2010年的0.075，2016年相比2010年上升了0.0128，达到0.088，与地表温度的变化趋势相一致，也表现为“V”型变化趋势。

2004年遥感生态指数为0.639，2010年该指数为0.694，增长了0.056，2016年为0.655，相比2010年下降了0.039，总体表现为倒“V”型态势。通过四个指标共同作用产生的体现的生态环境状况的遥感生态指数也表现一定程度的波动性。研究区遥感生态指数的变化特征也刚好是四个指标变化特征的体现。绿度、湿度对遥感生态指数起到正向的促进作用，它们两者的上升一定程度的促进遥感生态指数的上升，反之，绿度、湿度的下降也会带动遥感生态指数的下降；热度、干度指标对遥感生态指数起到反向的抑制作用，即它们两者的变动方向与遥感生态指数的变动方向相反，因此它们两者的“V”型变化特征会对遥感生态指数起到倒“V”型的作用。因此，遥感生态指数的变化是四个指标的综合体现，只是相互作用的强度有所差别。

3.2 生态指数的空间分布特征

在生态指数的空间分布中，将研究区遥感生态指数按照0.2等间隔进行分割，依次构成了差、较差、中、良、优，共五个等级，也构成了生态指数的不同类型，为深一步的分析构建基础。下面的分析所指的不同类型的遥感生态指数，就是划分的五类生态指数。为了引用方便，分别给予编号，依次为A、B、C、D、E。分析中给出了5种不同类型遥感生态指数的空间分布图，然后对各类型的区域面积构成及级别变化情况进行分析。

从图2和表3可以发现，研究区不同等级遥感生态指数变动较明显。E级生态指数区域面积增长较明显， 2004年该区域面积仅89.23 km2， 2010年为295.40 km2，2016年达到了341.32 km2，由总面积的2.76%增长到10.54%，该类生态指数由相对分散向集聚式条带分布，主要集中在古交主矿区与南部人类密集区之间的过渡区。D级生态指数在整个区域中一直占主导，占总面积的一半以上，2004年、2010年、2016年依次分别占总面积的63.74%、75.94%、58.77%，该类区域东南部区域减少明显。C级生态指数的区域面积中，2004年1 080.02 km2，主要集中分布在西北部的古交矿区；2010年下降到481.31 km2，主要分布在京昆高速沿线区域；2016年达到了919.64 km2，主要分布在东南区域，该区域属于密集居住区。A级、B级生态指数所占区域面积较少，比重较低，整体变化不大。研究区生态指数等级的变化过程，主要体现在西北区域的主矿区生态的好转、东南居住区域生态恶化的状态。

表2 生态指数的基本统计值Table 2 Basic statistical values of the ecological index

图2 生态指数的空间等级分布图Fig.2 Spatial level distribution map of the ecological index

为了进一步细化不同类型遥感生态指数的转化情况，分别统计不同时间段的级差变化情况。表4中，分别给出了2004～2010年、2010～2016年、2004～2016年三个时间段的遥感生态指数的等级变化情况。研究遥感生态指数等级变化中，主要以保持原状态为主，2004～2010年不变区域占总面积的66.86%，2010～2016年不变区域占总面积的73.27%，2004～2016年不变区域占总面积的56.27%。2004～2010年，遥感生态指数升级的区域面积达941.2974 km2，降级区域的面积为131.7852 km2，生态指数改善明显；2010～2016年，遥感生态指数升级的区域面积达186.3162 km2，降级区域的面积为679.0941 km2,生态指数下降明显；2004～2016年，遥感生态指数升级的区域面积达871.6842 km2，降级区域的面积为544.4631 km2,生态指数升级趋势较明显。同时，不同类型生态指数升级、降级过程中，主要一个级差变化为主，跳级的区域极少。

3.3 不同类型遥感生态指数的重心偏移

为较准确地把握不同类型遥感生态指数的移动方向，主要采用ArcGIS10分别测算各类遥感生态指数的重心坐标、偏移距离及方位角，见表5和表6。

表3 遥感生态指数各等级面积及占比Table 3 Area and proportion of each level of remote sensing ecological index

表4 遥感生态指数等级变化统计Table 4 Statistics on the level change of remote sensing ecological index

表5 不同类型遥感生态指数重心坐标Table 5 Center of gravity coordinates of different types of remote sensing ecological indices

表6 不同类型遥感生态指数重心偏移量及变动方位角Table 6 Center of gravity offset and variable azimuth of different types of remote sensing ecological index

2004～2010年，各类型生态指数的偏移距离差异较大。A级、B级生态指数的偏移距离较大，分别为37.172 km、31.255 km，其次是C级生态指数，偏移距离为13.594 km，D级、E级生态指数的偏移距离相对较小，分别为4.293 km和3.492 km。2010～2016年，各类型生态指数的偏移量相比上一阶段均表现为缩小的趋势。其中，A级、B级生态指数的中心偏移量依然最高，分别为15.270 km、18.814 km,C级以上的偏移量均较小，E级的重心偏移量仅2.495 km。

从两个不同时段各类生态指数偏移距离的大小来看，总体上体现了生态指数越差偏移距离越大的规律性。主要是由于低级别类型的生态指数相对较小，很小的变化都会带来重心较大的偏移量。同时也说明了研究区的生态持有一定程度的稳定性。特别是2010年以来，生态系统的稳定性有所加强。

为了比较直观的展示不同等级遥感生态指数的具体重心位置及转移方向，绘制了重心转移变动图，见图3。研究发现，各级生态指数的重心的集聚性明显。其中，A级生态指数的重心先由北部的矿区移动到南部的居住区，2016年又开始向北移动，但仍在人类居住区范围内；B级生态指数的重心，最初位于古交市中心，2010年该重心向南移动，2016年继续向南移动，达到人类居住相对密集的区域；C级生态指数重心由北向南一直移动，2016年位于京昆高速边缘；D级生态指数的重心，表现为由南向北移动，2016年位于古交矿区与南部村庄密集区域的过渡地带；E级生态指数的重心，也表现为由南向北的移动，与D级重心区域的位置较接近。A级、B级、C级生态指数的重心总体表现为由北向南移动，2016年的重心基本位于居民相对密集的村庄，说明了2004年矿区的生态较差，后来逐渐向人类居民密集分布的南部移动，该区域生态环境变差明显。D级、E级为生态较好的区域，均表现为由南向北的移动方向，而且始终位于主矿区和居民密集分布的过渡地带，说明了生态优质区域在远离人类活动区，向北转移，该区域主要由于地形较复杂，人类活动干扰相对较少。

图3 不同类型遥感生态指数重心偏移Fig.3 Center of gravity shift of different types of remote sensing ecological index

3.4 生态指数与四个指标的相关性

在城市遥感生态指数的评价体系中，为了量化四个指标对遥感生态指数的作用，通过计算两两之间的相关系数进行分析，见表7。

结果显示，绿度、湿度与遥感生态指数的相关系数大于0， 表现为正相关； 热度、干度指标与遥感生态指数的相关系数小于0，表现为负相关。在正相关的指标中，绿度指标的相关系数总体高于湿度指标的相关性，说明在促进遥感生态指数升级过程中，绿度指标的贡献度更大。在负相关指标中，热度指标的相关性绝对值一致高于干度指标的相关性，说明了在阻碍遥感生态指数的升级过程中，热度的影响更明显。通过计算，得到四个指标在三年中与遥感生态指数相关系数的平均值。

通过对平均作用强度取绝对值发现，对遥感生态指数的作用强度中，热度的作用最显著，相关度高达0.859，植被的作用次之，相关系数为0.799，干度的作用强度0.707，湿度的作用最小，相关系数为0.596。由此说明，在优化升级生态指数过程中，降低温度的作用比较明显，其次提高植被覆盖度也可以很好的调整生态环境，同时也降低了裸土的面积，弱化了裸土的负面作用。

表7 遥感生态指数与各指标的相关系数Table 7 Correlation coefficient between remote sensing ecological index and each indicator

4 结 论

1) 2004～2016年西山煤田区域遥感生态指数及其相关指标均表现为波动变化趋势。遥感生态指数、绿度、热度、干度均表现为先上升后下降的倒“V”变化趋势，湿度则表现为持续的增长。

2) 西山煤田生态指数以良为主，其中优级生态指数区域增长明显。在不同类型的遥感生态指数中，良级区域的面积2004年为2 063.87 km2，占总面积的63.74%，2010年面积为2 459.03 km2，占总面积的75.94%，2016年面积为1 903.16 km2，占总面积的58.77%，一致处于主导地位。E级生态指数区域面积增长较明显，由总面积的2.76%增长到10.54%。

3) 2004～2016年西山煤田生态指数的变化体现在西北区域主矿区生态好转、东南居住区域生态恶化。研究区不同类型的生态指数变化中，一半以上的区域仍保留原来的状态，2004～2010年不变区域占总面积的66.86%，2010～2016年不变区域占总面积的73.27%，2004～2016年不变区域占总面积的56.27%。

4) 2004～2016年不同类型的遥感生态指数重心均发生偏移。重心偏移距离较大的为A级、B级，其中，2004～2010年偏移距离分别为37.172 km、31.255 km，2010～2016年偏移距离分别为15.270 km、18.814 km。其他类型的生态指数重心偏移较少。

5) 西山煤田生态指数评价系统中，热度的作用最强，绿度次之。遥感生态指数与构成指标的相关系数分析中，热度与遥感生态指数的平均相关强度为0.859，植被与遥感生态指数的相关强度为0.799，干度、湿度的作用强度依次分别是0.707、0.596。