火成岩中长石含量与其特征光谱间响应关系研究
杨长保1, 高文博1,*, 侯光宇2, 李星喆1, 高曼婷1
1. 吉林大学地球探测科学与技术学院, 吉林 长春 130026
2. 中国人民解放军32023部队, 辽宁 大连 116000
*通讯联系人 e-mail: 905920143@qq.com

作者简介: 杨长保, 1972年生, 吉林大学地球探测科学与技术学院副教授

摘要

长石是地表岩石最重要的造岩矿物, 在地壳中的比例高达60%, 几乎是所有火成岩的主要矿物成分。 随着高光谱技术的发展, 国内外众多学者研究主要造岩矿物含量与其特征光谱的响应关系, 对遥感岩矿识别以及矿化蚀变信息提取提供了多种可能性。 该研究以USGS光谱库里18个火成岩样本为基础数据, 研究长石的特征光谱及其与含量之间的定量关系。 通过原始光谱反射率及其变换(包括小波三层分解高频分量、 小波二层分解、 去包络线后光谱、 去包络线后小波三层分解高频分量及去包络线后小波二层分解), 研究其与长石的含量之间的相关关系, 结果表明: (1)分析六种光谱反射率的变换, 去包络线后小波三层分解高频分量的光谱反射率与长石含量的相关关系最好, 且相关系数正负不断变化, 根据相关系数极值获得长石的特征谱带为431, 570, 972, 1 456, 1 856, 2 292.9和2 481 nm; (2)原始光谱反射率与长石含量的相关性曲线趋势较为平缓, 而经小波分解得到的高频分量后, 趋势明显, 经去包络线及小波分解得到高频分量后, 相关性曲线的变化趋势愈加明显, 由此可见, 自变量的微小变化就会引起因变量变化, 当岩石中长石的含量极小时, 小波分解处理能够提高模型的精度。 将长石含量与特征光谱的关系量化, 运用多元逐步线性回归分析以及最小二乘法建模, 建立6个线性回归模型和6个最小二乘法回归模型, 结果表明: (1)去包络线后的光谱比原始光谱建立的回归模型精度更高, 经过小波二层分解后的低频分量建模的回归模型精度优于未进行小波分解的光谱, 其中去包络线后小波二层分解低频分量建立的回归模型效果最佳。 (2)多元线性回归建立的模型精度优于最小二乘法, 同时筛选对因变量影响较大的自变量, 972, 1 456, 1 856, 2 292.9以及2 481 nm。 因此选择去包络线后的光谱进行多元线性回归法进行分析长石含量与光谱反射率之间关系, 考虑到不同的特征吸收波段对长石含量的影响因子不同, 可以利用长石的特征光谱定量反演某一区域内的长石的含量, 对识别矿物具有重要意义。

关键词: 火成岩; 长石; 特征光谱; 相关性分析; 多元回归分析; 小波分析
中图分类号:P585 文献标志码:A
Response Relationship between Feldspar Content and Characteristic Spectra in Igneous Rocks
YANG Chang-bao1, GAO Wen-bo1,*, HOU Guang-yu2, LI Xing-zhe1, GAO Man-ting1
1. College of Earth Exploration Science and Technology, Jilin University, Changchun 130026, China
2. Unit 32023 of the People’s Liberation Army of China, Dalian 116000, China
Abstract

Feldspar is the most important rock-forming mineral in surface rocks, accounting for up to 60% in the earth's crust. With the development of high-spectral technology, many scholars at home and abroad have been studying the response of the main building of rock and mineral content and the spectrum of the features, and it offers a variety of possibilities for remote sensing, mineralization, and mineralization. It is based on data of 18 igneous samples in the USGS spectral library to study the quantitative relation between the characteristic spectra and the content of the feldspar. Through the original spectral reflectivity and the transformation (including the three layers of the small wave to break down the high frequencies, the little wave layers, the spectrum of the wave to the back of the line, and then the little wave layer after the cable, and then the little wave layer after the cable, and then the small wave layer of the wave, and then the second layer of the cable) study the correlation between the high and the long rocks, and it turns out: (1) to analyze the transformation of six spectral reflectivity, the relation between the spectral reflectivity of the high frequency and the long rock content of the small wave three layers of the envelope, and the correlation coefficients are the best, and the correlation coefficients are constantly changing, and based on the relative coefficients of the relative coefficients, the high value of the long rock is 4, 31, 570, 972, 1 456, 1 856, 2 292.9, 2 481 nm; (2) the correlation between the original spectral reflectance and feldspar content curve trend is relatively flat, and after the wavelet decomposed high frequency component, through to the envelope and wavelet decomposition for the high frequency component, the change trend of correlation curve is becoming ever more obvious, therefore, the independent variable of the small changes will cause changes in the dependent variable, very hour when the content of the rock feldspar, wavelet decomposition processing can improve the accuracy of the model. The relation between the content of feldspar and the characteristic spectrum is quantified, using a multi-element stepwise linear regression analysis and a least-square method of modeling, establishing six linear regression models and six least-squares regression models, and the results show that: (1) the spectrum of the spectrum behind the envelope is more accurate than the original spectrum, and the low-frequency portion of the lower wave part of the wave is better than the one that's not done with the small wave decomposition, which is the best way to get to the back model of the low-frequency portion of the small wave after the envelope. (2) the multilinear regression model is better than the least-squares, and the variables that affect the variables that affect the larger variables, 972, 1 456, 1 856, 2 292.9 and 2 481 nm. In view of that relation between the content of feldspar and the spectral reflectance of the long stone and the influence factors on the content of feldspar in different absorption band, it is important to use the characteristic spectrum of feldspar to quantitatively invert the content of feldspar in a region.

Keyword: Igneous rock; feldspar; Characteristic spectrum; Correlation analysis; Multiple regression analysis; wavelet analysis
引 言

高光谱遥感技术具有丰富的光谱信息, 在可见光范围内, 不同的矿物具有不同的光谱特征谱带, 利用这些特征谱带, 可以有效的识别岩石中的矿物成分, 反演矿物在岩石中的质量百分比。 2012年, 樊磊等通过高光谱数据对黑龙江省北部地区进行岩石特征波段的提取, 通过对应分析方法找出与岩石相关的特征波段, 解释其物理意义并将样品中岩石类型进行分类, 提供了岩石与成像光谱数据精细分类的依据[1]; 2015年杨长保等通过高光谱数据研究岩石中黑云母的含量与光谱的吸收深度之间的关系, 并用数学模型的方式定量分析岩石的光谱吸收深度与矿物含量之间的关系[2]。 借鉴前人的研究, 发现岩类随着造岩矿物的不同, 光谱会出现一些系统的区域性变异变化, 但也存在局部的稳定性和规律性[3, 4, 5, 6, 7, 8], 通过研究确定了岩石中主要造岩矿物长石的特征光谱, 并通过建立长石含量与光谱之间的定量关系来为遥感岩矿识别以及矿化蚀变信息提取提供一种崭新的思路。

1 实验部分
1.1 样本数据

本次研究所用数据是以美国地质调查局提供的标准数字光谱库为基础数据, 选取火成岩样本18块, 其中火山岩3块(andesite1, andesite2, basalt), 侵入岩15块。 岩石样本包括各类矿物的百分含量以及用ASD地物光谱仪测量的0.35~2.5 μ m波长范围内的反射率。 表1是18块岩石主要造岩矿物的百分含量。

表1 18个岩石样本主要矿物含量说明 Table 1 Main mineral content description of 18 rock samples
1.2 光谱数据去包络线处理

岩石的光谱测量值包含有背景信息以及所需要的特征光谱信息, 光谱数据的去包络线处理目的是为了抑制背景信息, 从而放大光谱上的吸收特征。 “ 包络线” 定义为逐点直线连接局部波谱的最大值, 并使出现在峰值点上的外角大于180° , 再以原始光谱上的值除以包络线上对应的值, 即得到去除包络线的光谱值, 如式(1)[9]

R'(λi)=R(λi)Rc(λi)(1)

其中, R'(λ i)为去包络线后的光谱反射率, R(λ i)为某波段在吸收谷处的反射率, Rc(λ i)为对应波长处包络线上的值。

1.3 光谱数据的小波分解重构

在离散小波变换(DWT)中, 在空间Vj=VJ-1+WJ-1上表示信号, 也就是对于每一个在Vj上表示的信号x(t)能用上面提到的两个空间中的基函数来表示。

x(t)=kcA0(k)ϕj, k(t)=kcA1(k)ϕj-1, k(t)+kcD1(k)ωj-1, k(t)(2)

其中, A1(k)和D1(k)是尺度度量空间j-1的系数, j-1空间是由尺度度量空间j的系数A0(k)分解而来。 类似的, 可以用A1(k)和D1(k)重构A0(k)。

在尺度度量空间j对系数A0(k)进行分解得到在尺度度量空间j-1的两个系数A1(k)和D1(k)通过重构得到系数A0(k)。

如上图的分解与重构, 可以通过小波变换算法计算, 当小波与尺度在空间内是正交的, 可以用内积计算公式得到系数cA1(k)和D1(k)

cA1(k)=< x(t), ϕj-1, k(t)> =< ncA0(n)j, n(t), ϕj-1, k(t)> =ncA0(n)< ϕj, n(t), ϕj-1, k(t)> (3)

下面是内基计算方法的具体计算公式

具体系数计算过程如下

cA1(k)=nh0(n-2k)cA0(n)cD1(k)=< x(t), ωj-1, k(t)> =< ncA0(n)ϕj, n(t), ωj-1, k(t)> =ncA0(n)< ϕj, n(t), ωj-1, k(t)> (5)

cD1(k)=nh1(n-2k)cA0(n)(6)

对于上面的小波分解过程, 通过分别设计高通滤波器以及低通滤波器两组滤波器的系数可实现, 重构只是分解的逆过程。

在MATLAB中实现小波分解重构的算法处理。 具体操作为利用小波分解进行多层分解, 得到细节特征有效的各光谱曲线。 依据光谱曲线, 可得原始光谱中的噪声以及微小的吸收特征是影响高频系数的主要因素。 利用sym5小波基对原始光谱曲线进行三层小波分解得出其低频分量以及高频分量。 图1和图2分别是原始光谱曲线和去包络线后进行小波分解重构的结果。 图1是某个样本的原始光谱反射率进行小波分解重构后的结果。 (a)是样本的原始光谱曲线, (b)为样本光谱曲线小波分解低频分量光谱, (c)为样本光谱曲线小波分解高频分量光谱, (d)为样本光谱曲线小波分解重构后保留的特征光谱。 图2为样本反射率光谱进行包络线去除后进行小波分解重构的结果。 其中(a)是某样本去包络线后的光谱曲线, (b)是某样本包络线去除后小波分解低频分量光谱, (c)是样本包络线去除后小波分解高频分量光谱, (d)是样本包络线去除后小波分解重构后保留的特征光谱。

图1 原始光谱曲线小波分解结果Fig.1 The result of the wavelet decomposition of the original spectrum curve

图2 去包络线后光谱曲线小波分解结果Fig.2 De-envelope curve de-spectrum curve wavelet decomposition results

2 结果与讨论

在光谱库中筛选出18个样本, 对样本光谱进行去包络线, 小波分解等预处理, 并运用多元逐步回归法进行建模。

2.1 光谱数据的相关性分析

将岩石中长石的含量与光谱的反射率, 小波分解高频波段, 小波二层分解, 去包络线, 去包络线后小波分解高频分量以及去包络线后小波二层分解后的光谱反射率分别做双变量相关性分析, 随着波长的变化, 长石的含量与各种算法下的光谱反射率的相关性系数也不同, 相关系数曲线趋势图如图3和图4所示。 从表2可以看出, 岩石中长石的含量与岩石的反射率的相关系数为0.407, 呈正相关关系。 但对原始光谱进行小波二层分解后, 提高了相关系数为0.595, 同时岩石中长石含量与岩石的反射率相关系数呈正负不断变化, 这能有效的区分不同波段对长石含量的影响。 这为后面建立回归方程所需的特征波段提供理论依据。

图3 原始光谱与其小波分解后的光谱相关系数曲线Fig.3 Correlation coefficient of the original spectrum and its wavelet decomposition

图4 去包络线与其小波分解后的光谱相关系数曲线Fig.4 De-envelope the spectral correlation coefficient following its wavelet decomposition

表2 长石含量与最大相关波段分析 Table 2 Analysis of feldspar content and maximum correlation band

综合图3和图4可以看出, 对原始光谱去包络线后小波分解高频的光谱与长石含量的相关性最好, 相关系数范围为[-0.7, 0.7], 其次是原始光谱经小波分解得到高频分量后的光谱, 其与长石含量的相关性范围是[-0.7, 0.6]。 研究发现经小波分解得到高频分量后, 光谱反射率与长石含量相关性系数曲线趋势有较为明显的变化, 经去包络线及小波分解得到高频分量后, 相关性系数曲线的变化趋势愈加明显, 这说明: 自变量的微小变化就会引起因变量变化, 当岩石中长石的含量极小时, 通过小波分解处理能够提高模型的精度。 由表2以及图3和图4等相关系数趋势图可知, 当波长处于431, 570, 972, 1 456, 1 856, 2 292.9以及2 481 nm处相关系数达到极显著水平, 说明在这7个波段处, 长石的含量与光谱反射率之间的相关关系最好, 故岩石光谱在这7处存在特征吸收波段, 可见火成岩中长石的含量与岩石光谱之间的关系复杂, 具有多个特征吸收谱带。

2.2 长石含量与光谱反射率的回归方程建立

将以上分析结果中的波段作为自变量, 并设为X1X6, 长石的含量作为因变量进行回归分析, 选用enter法和两阶最小二乘法建立回归方程, 并使用判断系数R2、 调整后的判断系数R2以及F统计量来判断模型的稳定性。 分别用原始的反射率光谱, 小波分解后高频分量, 小波二层分解, 去包络线, 小波分解后的高频分量以及去包络线后的小波二层分解作为自变量, 长石的含量作为因变量, 进行多元回归分析, 得到12个方程。 通过对表3表4所列的12个方程分析可得, (1)线性回归分析模型的精度要优于最小二乘法回归模型, 在线性回归分析模型中, 去包络线小波重构分解得到方程的F统计量最大, 为43.827, 且调整后的R2最大为0.922; 在最小二乘法回归模型中, 原始光谱小波分解得到的方程中, F统计量最大为16.623, 调整后的R2统计量最大为0.859; (2)在线性回归分析中, 经过小波二层分解得到的方程的精度比小波三层分解高频得到的方程精度高, 可见在小波分解中, 三层分解的高频信号携带大部分噪声, 不建议使用小波分三层分解进行反演; (3)在线性回归分析中, 去包络线光谱与去包络线光谱小波二层分解得到的方程能够自动筛选对因变量作用明显的波段, 并将对因变量影响微小的波段予以剔除, 其中972, 1 456, 1 856, 2 292.9以及2 481 nm这些波段是对长石含量影响最大; (4)最小二乘法回归分析中, 原始光谱的小波三级分解、 去包络线光谱、 去包络线光谱的小波三级分解以及去包络线光谱的小波二级分解的F统计量均小于2.9, 调整后的R2统计量均小于0.5, 不具有明显的统计意义, 故在反演模型中将其剔除。

表3 长石含量多元线性回归分析方程 Table 3 Multivariate linear regression analysis equations for feldspar content
表4 长石含量最小二乘法回归方程 Table 5 The least square regression equation for feldspar content
3 结 论

在研究岩石中长石的含量与光谱反射率之间的关系时发现, 去包络线处理后的光谱反射率与长石的相关性比未经过反射率处理的相关性强; 经过小波分解后, 相关性均有显著提高, 但去包络线后再小波分解得到的相关系数明显更高, 范围为[-0.7, 0.7], 据此结合相关系数趋势图, 得到在 431, 570, 972, 1 456, 1 856, 2 292.9以及2 481 nm处长石的含量与光谱反射率的相关性系数达到极值, 岩石在这七处有特征光谱。

研究长石含量与长石特征光谱之间的回归分析时, 运用多元逐步线性回归分析以及最小二乘法建模, 建模精度通过调整后的判定系数(Adjusted R2)和F统计量的大小判断, 采用enter线性回归分析以及最小二乘法得到12个回归方程, 分析可得如下结论: (1)在线性回归分析中小波二级分解比小波三级分解高频信号建立的模型精度高, F统计量大, 且调整后的R2大, 可见小波三级分解高频信号携带大量噪声, 不建议使用; (2)在线性回归分析中, 去包络线光谱与去包络线光谱小波二层分解得到的方程能够自动筛选对因变量作用明显的波段, 并将对因变量影响微小的波段予以剔除; (3)线性回归分析所建立的模型精度要优于最小二乘法所建立的模型精度; 因此选择线性回归分析的方法, 对去包络线后在小波二级分解的光谱建立的模型, 是精度最优解。 将长石的含量与光谱的反射率进行线性回归, 结合不同特征波段的影响权重可以利用长石的特征光谱定量反演某区域内长石的含量并能有效识别矿物成分。

参考文献
[1] FAN Lei, ZHAO Wen-ji, GONG Zhao-ning, et al(樊磊, 赵文吉, 宫兆宁, ). Journal of Jilin University·Earth Science Edition(吉林大学学报·地球科学版), 2012, 42(2): 575. [本文引用:1]
[2] YANG Chang-bao, ZHANG Chen-xi, LIU Fang, et al(杨长保, 张晨曦, 刘舫, ). Spectroscopy and Spectral Analysis(光谱学与光谱分析), 2015, 35(9): 2583. [本文引用:1]
[3] XU Jing, JIANG Ping-an, WU Hong-qi(胥静, 蒋平安, 武红旗). Hubei Agricultural Science(湖北农业科学), 2015, 54(1): 43. [本文引用:1]
[4] SUN Ya-qin, TIAN Shu-fang, WANG Xing-zhen, et al(孙娅琴, 田淑芳, 王兴振, ). Geoscience(现代地质), 2016, 30(1): 239. [本文引用:1]
[5] QIN Kai, CHEN Jian-ping, ZHAO Ying-jun, et al(秦凯, 陈建平, 赵英俊, ). Uranium Geology(铀矿地质), 2018, 34(3): 174. [本文引用:1]
[6] GUO Bang-jie, ZHANG Jie-lin, WU Ding(郭帮杰, 张杰林, 武鼎). Science Technology and Engineering(科学技术与工程), 2018, 18(17): 125. [本文引用:1]
[7] YIN Zi-yao, LIU Tang, WANG Zhen(殷子瑶, 刘唐, 王震). Beijing Surveying and Mapping(北京测绘), 2018, 32(7): 788. [本文引用:1]
[8] Gomez C, Adeline K, Bacha S, et al. Remote Sensing of Environment, 2018, 204: 18. [本文引用:1]
[9] Elsy Ibrahim, Pierre Barnabé, Erick Ramanaidou, et al. International Journal of Applied Earth Observations and Geoinformation, 2018, 73: 653. [本文引用:1]