集合经验模态分解和交叉小波在赣江流量与太阳黑子活动关系中的应用

万智巍， 贾玉连， 洪祎君， 蒋梅鑫

(1.江西师范大学 地理与环境学院 鄱阳湖湿地与流域研究教育部重点实验室，江西 南昌 330022； 2.中国科学院地理科学与资源研究所 陆地表层格局与模拟院重点实验室， 北京 100101)

1 研究背景

河流的水文要素变化过程是气候变化和下垫面要素交互作用的结果，具有复杂的时空变化演变规律[1]。近年来，随着全球极端水文、气候事件频发，了解河川径流时空变化规律成为社会的迫切需求[2-3]。地球作为一个复杂的开放系统，太阳辐射是其重要的外部能量输入和影响因素。研究表明，全球和区域降水在长时间尺度上受到太阳活动的影响，但其具体的影响过程和相互关系在不同区域有着不同的表现形式和规律特征[4]。Neff等[5]的研究论证了全新世以来，太阳活动与印度季风之间复杂的耦合关系。Zhang等[6]的研究也表明，太阳活动的变化对东亚区域降水和河流流量有着复杂的影响机制。窦睿音等[7]的研究表明，太阳黑子数与关中平原的整体干湿状况的关系存在转换现象，20世纪70年代之前二者为反向相关，之后二者之间的关系转换为正相关。

河流流量的变化是区域气象水文活动的综合体现，反映了流域范围内整体的干湿状况，因此流量变化成为气候变化区域响应的指示器[8]。相关研究表明，尽管太阳黑子活动与河流流量之间存在相关性，但是由于气象水文过程属于非线性、非稳定的复杂系统，如利用传统统计方法分析其相关性可能会受到原始序列中噪音的干扰[9]，难以探究二者之间在不同时域和频域上的复杂相关关系。因此，本文结合集合经验模态分解(Ensemble Empirical Mode Decomposition, EEMD)[10]和交叉小波变换(Cross Wavelet Transform, CWT)[11]，提出一种改进的相关分析方法——EEMD-CWT综合分析，并将其应用于赣江长时间流量序列和同期太阳黑子数序列之间的相关性分析。文中结果一方面可以为赣江流域水资源和干湿状况长期变化规律和影响因素的探讨提供基础资料，另一方面可以为提高气象水文系统与外部强迫之间的非线性关系研究提供参考和思路。

2 数据与方法

2.1 流量数据与太阳黑子数据

赣江的流量数据来自1950-2015年南昌外洲水文站的实测资料。外洲水文站位于江西省南昌市范围内的赣江下游右岸(115.84°E、28.63°N)，控制面积为80 948 km2，占整个赣江流域范围的99%以上[12]，因此外洲水文站的流量数据可以反映整个赣江流域的产流情况。外洲水文站设立于1949年10月，从1950年开始有连续的流量观测数据，是目前赣江流域最长的流量数据序列(图1(a))。

太阳黑子数据来自比利时皇家天文台数据中心，该研究机构是世界权威的太阳活动观测中心，因此本文选取其发布的1950-2015年太阳黑子相对数(图1(b))作为反映太阳活动状况的序列(数据下载网址为http://www.sidc.be/silso/datafiles)。

图1 1950-2015年赣江年均流量和太阳黑子数

由图1可以观察到，流量序列与太阳黑子序列之间的相关关系并不明显。因此，利用SPSS 20软件，对流量数据与太阳黑子数据进行初步的相关分析。结果表明，二者之间的Pearson相关系数为-0.118(P=0.347)。

2.2 研究方法

2.2.1EEMD方法 经验模态分解(EmpiricalModeDecomposition,EMD)是一种利用各信号周期不同的振幅频率特征，将原始序列分解为一系列具有不同时间尺度本征模函数(IntrinsicModeFunction,IMF)的方法，被广泛应用于水文、气象、地球物理等领域的多尺度时频分析[13-15]。集合经验模态分解(EEMD)则是EMD方法的改进算法，通过在原始信号中加入噪声减弱了模态之间的混叠，从而提高了分解的有效性[16-17]。

2.2.2CWT方法 交叉小波变换(CrossWaveletTransform,CWT)是一种结合交叉谱分析和小波分析的新方法，可以同时在时域和频域上分析两个信号之间的相关性[11]。交叉小波能量谱可以反映出两个信号经过小波分析之后具有相同能量谱的区域，因此可以揭示出两个信号在不同时频域上相互作用的程度。其具体计算方法为：

(1)

2.2.3EEMD-CWT综合分析 本文提出的EEMD-CWT综合分析的思路是：首先利用EEMD方法将原始信号分解为不同尺度的互相正交的本征模函数(IMF)分量，以实现信号的平稳化，降低原始信号中的噪音影响。然后利用多窗谱分析(Multi-TaperMethod,MTM)[18]确定原始信号中所存在的各种周期。最后，选择具有相同物理意义的周期尺度下的本征模函数(IMF)分量进行交叉小波变换，得到二者的交叉小波能量谱和小波凝聚谱，并以此来分析两个序列在不同时频域上的相关性。

3 结果与分析

3.1 流量的多尺度特征和周期

对1950-2015年赣江流量序列进行EEMD分解(图2)，共得到6个IMF分量，其中IMF1～5为不同时间尺度下流量序列的周期项，IMF6为整个序列的长期变化趋势项。由图2可以看出，IMF1分量的波动振幅最大，频率也最高，但是周期性较差。这一信号代表了原始序列中的主要变化因素，是水文序列非线性和非平稳特征的体现。EMD方法的提出者——Huang等[18]的研究指出，如果趋势项为单调上升或单调下降函数则说明原始信号具有明显的非平稳性。由图2可以清楚的观察到IMF6为单调上升函数，这也论证了赣江近66年来的流量序列具有明显的非平稳性。IMF2～5这4个分量的周期性和平稳性已经得到了很大程度上的改善，表现出了比较明显的周期变化规律。并且随着分解尺度的递进，其振动幅度和频率逐渐下降，信号的平稳性逐渐增强。如IMF4和IMF5已经比较接近于正弦波。

(图中纵坐标为相对值)

通过对1950-2015年赣江流量序列的多窗谱(MTM)分析表明(图3)，近66年来赣江流量序列具有明显的22、8、4.5和3.7a周期(通过α=0.05显著性检验)，其中11a周期也比较明显，但是未能通过相关检验。谱值最大的22a周期即为太阳活动的Hale周期[19]，也就是常见的太阳黑子11a周期的2倍。另外，准4～8a的周期很可能与ENSO的活动周期有关[20]。

(虚线表示95%置信区间)

3.2 太阳黑子的多尺度特征和周期

对1950-2015年太阳黑子数序列进行EEMD分解(图4)，共得到6个IMF分量，其中IMF1～5为不同时间尺度下太阳黑子数序列的周期项，IMF6为整个序列的长期变化趋势项。由图4可以看出，太阳黑子数的变化过程较为平稳，具有明显的周期性特征。在IMF1分量中，除1950-1960年这段时间外，其余时间其振幅都较为平稳。IMF2～5这4个分量具有非常明显的周期变化特征，各个时间尺度下的变化过程都具有类似正弦波的结构。IMF6趋势项的变化过程显示，1950-1970年为太阳黑子数上升阶段，1970-2015年为下降阶段。

(图中纵坐标为相对值)

进一步的多窗谱分析(MTM)表明(图5)，太阳黑子数序列具有显著的36、11、5.4和3.2a周期(通过α=0.05显著性检验)。其中太阳黑子数的11a周期表现最为明显，其谱值明显大于其他几个周期。

(虚线表示95%置信区间)

3.3 流量与太阳黑子交叉小波分析

现代研究表明[21-22]，太阳活动异常对地球气候系统、降水、流域流量的影响主要是通过影响环流和季风强度实现的。如过去120年以来南亚和东亚地区太阳辐射的增强通过影响沃克环流和哈德来环流，最终导致区域降水和流域流量的增加。因此可以通过太阳活动和流量之间的关系，进一步揭示二者在不同尺度下的相互影响过程。

通过前文的分析可知，赣江流量序列与太阳黑子数序列的相关系数仅为-0.118，且未能通过α=0.1的显著性检验，不能直接利用未经处理的原始序列进行二者之间的相关分析。但是，经过EEMD分解后，各IMF分量的平稳性有了很大程度的改善，周期性特征也更加明显。与此同时，流量和太阳黑子的多窗谱分析(MTM)表明，二者都具有明显的太阳活动周期，因此可以选取具有与太阳活动周期相关的IMF分量进行对比分析。以此为基础，经过交叉小波分析得到赣江流量和太阳黑子数之间的交叉小波能量谱和小波凝聚谱。

图6为1950-2015年赣江流量序列IMF2分量与太阳黑子数序列的IMF2分量对比。由图6可知，赣江流量变化与太阳黑子数序列之间的相关关系在不同时间段内是不一样的。其中，1950-1960年呈负相关、1961-1973年呈正相关、1974-1984年呈负相关、1985-1995年呈正相关、1996-2007年为负相关(部分时段相差1/4个周期)、2008-2015年为正相关(部分时段相差1/4个周期)。近66年来，大体上二者以11a为周期，在负相关和正相关之间交替变换。

(图中纵坐标为相对值)

图7为1950-2015年赣江流量与太阳黑子数序列的交叉小波谱和小波凝聚谱，交叉小波能量谱中颜色越偏红表示能量谱密度值越大。由图7(a)可以观察到在11a周期附近赣江流量与太阳黑子数的相关关系最为显著，其中黑色实线所围合区域表示通过了95%置信度检验，灰色阴影区则表示可能会受到边界效应的影响。图中→表示二者之间变化的位相一致(正相关)，←表示位相相反(负相关)，↑和↓则表示二者之间提前或落后1/4个周期。在图7(a)中的11a周期尺度上，1950-2015年二者之间由负相关逐渐顺时针转为相差1/4个周期，再转为正相关。小波凝聚谱(图7(b))中相关系数较高的区域与交叉小波能量谱(图7(a))中基本一致，也是在11a尺度上达到最大值。其中1950-2000年，二者的相关系数在0.6以上，2001-2015年二者的相关系数也在0.4以上。小波凝聚谱可以发现在32a周期尺度上，二者的相关系数值在整个研究时段都在0.8以上，但是考虑到本研究的长度只有66a，而且由于边界效应的存在，这一相关关系可能并不具有实际意义。

图7 1950-2015年赣江流量与太阳黑子数序列的交叉小波谱和小波凝聚谱

4 结 论

本文基于赣江外洲水文站1950-2015年的长时间流量序列和同期太阳黑子数序列，结合EEMD及CWT综合分析方法对流量与太阳活动之间的多尺度相关关系进行了分析。研究结果表明：

(1)EEMD-CWT综合分析可以有效利用EEMD方法将原始信号进行分解，得到具有平稳性的IMF分量，并借助多窗谱分析(MTM)识别出信号中的不同周期，选择具有相同物理意义的赣江流量序列IMF分量和太阳黑子IMF分量，再利用交叉小波变换(CWT)得到二者在时域和频域上的相关关系变化过程。

(2)1950-2015年赣江流量序列的EEMD分解结果表明，该序列具有明显的非平稳性，其中IMF1分量的波动振幅最大，频率最高，代表了原始序列中的主要变化因素，是水文序列非线性和非平稳特征的体现。

(3)1950-2015年太阳黑子数序列的EEMD分解结果表明，太阳黑子数的变化过程较为平稳，具有明显的周期性特征。趋势项的变化过程显示，1950-1970年为太阳黑子数上升阶段，1970-2015年为下降阶段。

(4)多窗谱分析(MTM)表明，近66年来赣江流量序列具有22、8、4.5和3.7a周期，太阳黑子数序列具有36、11、5.4和3.2a周期。

(5)交叉小波分析(CWT)表明，在11a周期尺度上，1950-2015年赣江流量与太阳黑子数之间由负相关逐渐顺时针转为相差1/4个周期，再转为正相关。其中1950-2000年，二者的相关系数在0.6以上，2001-2015年相关系数在0.4以上。

[1]GATESWL.Anoverviewoftheresultsoftheatmosphericmodelintercomparisonproject(AMIPI)[J].BulletinoftheAmericanMeteorologicalSociety,1999，80(2)：29-55.

[2] 夏 军，谈 戈. 全球变化与水文科学新的进展与挑战[J]. 资源科学,2002，24(3)：1-7.

[3] 范垂仁，夏 军，张利平，等. 中国水旱灾害长期预报[M]. 北京：中国水利水电出版社，2008.

[4] 张兰生，方修琦，任国玉. 全球变化[M]. 北京：高等教育出版社，2000.

[5]NEFFU，BURNSSJ，MANGINIA，etal.Strongcoherencebetweensolarvariabilityandthemonsooninomanbetween9and6KyrAgo[J].Nature,2001，411(6835)：290-293.

[6]ZHANGP，CHENGH，EDWARDSRL，etal.Atestofclimate，sun，andculturerelationshipsfroman1810-Yearchinesecaverecord[J].Science,2008，322(5903)：940-942.

[7] 窦睿音，延军平. 关中平原太阳黑子活动周期与旱涝灾害的相关性分析[J]. 干旱区资源与环境,2013，27(8)：76-82.

[8]CHENF，YUANYJ，WEIWS，etal.Tree-ringrecordedhydroclimaticchangeinTienshanMountainsduringthepast500years[J].QuaternaryInternational,2015,358：35-41.

[9]YEHJ，SHIEHJ，HUANGNE.Complementaryensembleempiricalmodedecomposition：anovelnoiseenhanceddataanalysismethod[J].AdvancesinAdaptiveDataAnalysis,2010，2(2)：135-156.

[10]WUZhaohua，HUANGNE.EnsembleEmpiricalModeDecomposition：Anoise-assisteddataanalysismethod[J].AdvancesinAdaptiveDataAnalysis,2009,1(1)：1-41.

[11]GRINSTEDA，MOOREJC，JEVREJEVAS.Applicationofthecrosswavelettransformandwaveletcoherencetogeophysicaltimeseries[J].NonlinearProcessesinGeophysics,2004，11(5/6)：561-566.

[12] 顾朝军，穆兴民，高 鹏，等. 赣江流域径流量和输沙量的变化过程及其对人类活动的响应[J]. 泥沙研究,2016,41(3)：38-44.

[13]JIF，WUZ，HUANGJ，etal.Evolutionoflandsurfaceairtemperaturetrend[J].NatureClimateChange,2014，4(6)：462-466.

[14]CHENYangkang，MAJitao.RandomnoiseattenuationbyF-Xempiricalmodedecompositionpredictivefiltering[J].Geophysics,2013，79(3)：81-91.

[15]BREAKERLC，LOORHR，CARROLLD.Trendsinseasurfacetemperatureoffthecoastofecuadorandthemajorprocessesthatcontributetothem[J].JournalofMarineSystems,2016，164：151-164.

[16]WUZ，HUANGNE，CHENX.Themulti-dimensionalensembleempiricalmodedecompositionmethod[J].AdvancesinAdaptiveDataAnalysi,2009，1(3)：339-372.

[17]MANNME，LEESJM.Robustestimationofbackgroundnoiseandsignaldetectioninclimatictimeseries[J].ClimaticChange,1996，33(3)：409-445.

[18]HUANGNE，SHENZ，LONGSR，etal.Theempiricalmodedecompositionandthehilbertspectrumfornonlinearandnon-stationarytimeseriesanalysis[J].ProceedingsoftheRoyalSocietyofLondonA(Mathematical，PhysicalandEngineeringSciences),1998，454(1971)：903-995.

[19]RINDD.Thesun'sroleinclimatevariations[J].Science,2002，296(5568)：673-677.

[20] 王钟睿，钱永甫. 江淮梅雨的多尺度特征及其与厄尔尼诺和大气环流的联系[J]. 大气科学学报. 2004，27(3)：317-325.

[21] 葛全胜，刘路路，郑景云，等. 过去千年太阳活动异常期的中国东部旱涝格局[J]. 地理学报，2016，71(5):707-717.

[22]BHATTACHARYYAS，NARASIMHAR.RegionaldifferentiationinmultidecadalconnectionsbetweenIndianmonsoonrainfallandsolaractivity[J].JournalofGeophysicalResearchAtmospheres，2007，112(24):177-180.