日长变化、大气角动量和ENSO年际信号的相关分析及2020—2021年拉尼娜事件

孔昭洋，周永宏,3，许雪晴，安显然

(1.中国科学院 上海天文台，上海200030； 2.中国科学院大学，北京100049； 3.中国科学院 行星科学重点实验室，上海200030； 4.上海市空间导航与定位技术重点实验室，上海200030)

1 引言

厄尔尼诺-南方涛动(El Nino-Southern Oscillation,ENSO)是指大约每2～7 a出现在热带太平洋地区的气候系统年际振荡。作为一种海—气耦合现象，ENSO在海洋中表现为每隔几年出现在赤道中—东太平洋表层海水的暖或冷异常，厄尔尼诺(El Nino)对应暖异常，拉尼娜(La Nina)对应冷异常。在大气中，ENSO表现为南方涛动，即为赤道东、西太平洋之间海平面气压存在的此消彼长现象。

地球自转变化不仅表征了太阳、月球和其他行星引力等外力矩的影响，还反映了固体地球内部地核与地幔的相互作用，以及地球表面海洋和大气等流体圈层的耦合过程。在年际以及更高频时间尺度上，大气和海洋等地球表层流体活动是影响地球自转速率变化的主要原因[1−5]。其中，作为年际尺度上最强的气候变率，ENSO是导致日长年际变化的重要因子。

在不考虑外力矩作用的情况下，固体地球和大气、海洋等流体圈层可近似看成保守的动力学系统。根据总角动量守恒原理，大气和海洋的动力学过程所引起的角动量变化，会使固体地球的角动量产生相反的变化，从而使地球自转速率或日长(length of day,LOD)发生相应改变。因而，日长变化(ΔLO D)与ENSO事件之间有密切的联系。

针对ΔLO D、大气角动量(atmospheric angular momentum,AAM)与ENSO的联系，许多学者开展了相关研究：Langley等人[6]最早给出了AAM与ΔLO D之间的定量关系。Stefanick[7]基于1963―1973年的数据指出年际尺度上LO D的变化、AAM与南方涛动指数(southern oscillation index,SOI)具有相似性。在经历了1982―1983年的超强厄尔尼诺事件后，Rosen和Salstein[8]在AAM和ΔLO D中检测到了1982―1983年的超强暖事件信号。Eubanks等人[9]、Chao[10]通过超前滞后相关得到南方涛动指数超前日长年际变化1～2月，表明年际分量的ΔLO D主要由ENSO引起。此后，Dickey等人[11−13]认为年际ΔLO D主要由与AAM轴向分量有关的χ3激发，同时指出大气角动量存在起源于赤道并向中高纬度传播的趋势，且与ENSO密切相关；之后还提出海表面温度(sea surface temperature,SST)的低频变化可能是年际和年代际尺度AAM发生变化的激发源；以及极向的热带温度梯度(tropical temperature gradient,TTG)在太平洋中部海区(Nino3.4)平均SST出现最大异常后的1～2月达到峰值，因此认为是TTG所引起的热成风异常驱动了大气角动量中的风项，从而解释了ENSO对年际ΔLO D的驱动。近年来，de Viron和Dickey[14]对两种不同类型的厄尔尼诺事件(东部型和中部型)与ΔLO D的关系进行了AAM与力矩方面的研究，认为东部型ENSO的AAM信号比中部型大2倍以上，并用力矩解释了这一差异，同时认为东太平洋型ENSO具有更大的山脉力矩和摩擦力矩；山脉力矩越大说明AAM异常越强，但与此异常相关的风速在地表产生较强的负摩擦力矩又对AAM异常起到抵消作用。Lambert等人[15]对1980年以来3次超强厄尔尼诺事件鼎盛时期的AAM和力矩与ΔLO D的关系进行了分析，认为在1982―1983年和1997―1998年的东部型事件期间，ΔLO D主要受山脉力矩驱动，而对2015―2016年的东部—中部混合型事件，摩擦力矩补偿了较弱的山脉力矩而驱动ΔLO D异常。

基于以上研究，本文采用1962年1月―2021年1月内国际地球自转与参考架服务(International Earth Rotation and Reference Systems Service,IERS)的日长变化观测值序列、基于美国国家环境预报中心/大气研究中心(National Centers for Environmental Prediction/National Center for Atmospheric Research,NCEP/NCAR)再分析资料集R1归算的大气角动量函数，与表征ENSO的海洋尼诺指数(ocean Nino index,ONI)和经过标准化的负南方涛动指数(minus southern oscillation index,MSOI)进行对比，综合分析地球自转变化、AAM和ENSO在年际尺度上的关联性，尝试给出过程中对应的物理解释，同时在ΔLO D年际分量中检测到2020―2021年拉尼娜事件的信号。

2 资料与数据处理

2.1 厄尔尼诺-南方涛动

由于ENSO本身所具有的多样性和复杂性[16,17]，给出准确的定义并测定其持续时间和强度是一项有相当难度的工作。常用的判定标准一般来自美国国家海洋与大气管理局(National Oceanic and Atmospheric Administration,NOAA)，将热带太平洋分成了Nino1+2,Nino3,Nino3.4,Nino4海区(如图1所示)，其中Nino3.4海区是业务上定义的热带太平洋中部区域[18]。NOAA对ENSO的定义为：热带太平洋中部(Nino3.4)的季节性(任意连续3月的平均值)海表温度(SST)比对应区域的气候态平均值高0.5℃(厄尔尼诺)或低0.5℃(拉尼娜)。而对厄尔尼诺事件和拉尼娜事件起止时间和强度的判定主要以海洋尼诺指数(ONI)为标准，O N I指数的每个数据点均使用连续3月Nino3.4的海表面温度异常(sea surface temperature anomaly,SSTA)滑动平均值，如果连续出现至少5月的O N I|S S T A|≥0.5℃便会被定义为一次厄尔尼诺或者拉尼娜事件，这些冷暖异常事件通常具有广泛的影响。值得注意的是，目前O N I指数的SSTA基于第5个版本的NOAA扩展重构海温资料(extended reconstructed sea surface temperature version 5,ERSSTV.5)，该版本以每5年滚动更新一次气候态的方式，以尽可能地消除Nino3.4海区SST长期变暖趋势对ENSO监测的影响[19]。

图1 Nino指数分区与塔西提(Tahiti)和达尔文(Darwin)位置示意图

南方涛动指数(S O I)用来衡量南方涛动在大气中的表征，为塔西提岛(热带东太平洋)与达尔文港(热带西太平洋)海平面气压差的月平均标准化指标，见图1。有时S O I也是衡量赤道太平洋东风强度的指标，其正值对应拉尼娜相关的赤道东风异常，负值对应厄尔尼诺相关的赤道西风异常[20]。图2为1962年1月―2021年1月标准化的S O I与Nino3.4 SSTA的月平均序列对比(数据来自NOAA)，可见两者呈负相关且气压的改变相比于海温异常的变化更具高频噪声波动。值得注意的是：两者在对ENSO强度的描述方面并不完全相似，如对1982―1983年、1997―1998年和2015―2016年这3次超强厄尔尼诺事件中最大强度的反映，S O I有递减趋势，而Nino3.4 SSTA则显示递增趋势。但另一方面，这两个指数通常具有良好的时间同步性，图3显示了ENSO在海洋与大气中的不同表征指数在不同修正方法下的相关性，时间区间选取基本同上，可见月平均的S O I与Nino3.4 SSTA(绿色点划线)、3月滑动平均的S O I与O N I(Nino3.4 SSTA月数据的3月滑动平均值)(红色点划线)、S O I与O N I的年际(1～10 a)分量(蓝色点划线)均显示最大负相关同步，且相关系数由于序列平滑性的增加依次增强为-0.73,-0.84和-0.92。

图2 1962年1月—2021年1月的南方涛动指数(S O I)与Nino3.4海表温度异常

图3 ENSO在海洋与大气中的不同表征指数在不同修正方法下的相关性

为了与O N I所示的ENSO冷暖相位有相同的对应，本文将南方涛动指数取负得到负南方涛动指数(MSOI)。为了去除滤波过程中端部效应的影响，本文先将时间序列用自回归(auto-regressive,AR)模型进行左右端部各延拓30月，再通过带通1～10 a的2阶巴特沃斯滤波器，最后去除端部效应影响较大的AR模型延拓部分，得到年际分量的时间序列[21]。

2.2 日长变化

绕地球自转轴旋转速率的变化通常用日长变化(ΔLO D)表示：

式中，ΔLO D表示地球自转日长的变化，Δω表示地球自转角速度的变化，ω0表示地球自转的平均角速度。由式(1)可知地球自转的日长变化(ΔLO D)与地球自转角速度的变化(Δω)成反比[22]。

国际地球自转与参考架服务(IERS)提供了自1962年1月1日以来包含日长变化(ΔLO D)在内地球自转参数(earth orientation parameters,EOP)的日观测数据。这些数据最早由光学天体测量方法得到，自20世纪70年代起逐渐由激光测月(lunar laser ranging,LLR)、激光测卫(satellite laser ranging,SLR)、甚长基线干涉(very long baseline interferometry,VLBI)、全球卫星导航定位系统(global navigation satellite system,GNSS)以及星基多普勒轨道和无线电定位组合系统(Doppler orbitography and radio positioning integrated by satellite,DORIS)等现代空间大地测量技术测量得到[23]。目前最新的EOP 14 C04序列，对日长变化的测量精度达到微秒水平。

本文选取1962年1月1日―2021年1月31日的ΔLO D时间序列进行分析。考虑地球系统内角动量守恒，需要去除潮汐外力矩，根据IERS规范计算并扣除固体潮中带谐项对LO D的贡献，得到非潮汐的ΔLO D日数据，然后转为月平均数据[24,25]。扣除潮汐后的日长变化序列及其分解序列显示于图4，其中，图4a)的浅蓝色实线和粉色实线分别表示非潮汐的多尺度日长变化序列和通过多项式拟合得到的日长变化“十年”尺度项；图4b)红色实线表示：对非潮汐日长变化序列扣除“十年”尺度项和季节项后，通过前述的AR模型扩展、巴特沃斯滤波器滤波，以及端部效应扣除得到的与ENSO对应的年际尺度ΔLO D时间序列；图4c)中绿色实线表示通过叠加周年项、半年项和1/3年项拟合得到的季节项；图4d)深蓝色实线为ΔLO D的亚季节等高频项。

图4 1962年1月—2021年1月的月平均ΔLO D多时间尺度分解示意图

2.3 大气角动量

为了研究大气角动量对日长变化的激发情况，本文使用Eubanks给出了轴向的大气角动量函数(atmospheric angular momentum function,AAMF)χ3分别对于压力项和风项的表达式[26]：其中，u为纬向西风，ps为地表大气压力，R为地球半径，Ω为地球平均自转速率，g为重力加速度，取9.8 m·s−2，Cm为地幔轴向主转动惯量，λ为经度，φ为纬度。

上式中需要用到的地表大气压力和水平风场等气象要素资料来自美国国家环境预报中心/美国大气研究中心(NCEP/NCAR)的再分析资料集R1，数据输出间隔为6 h(每日以GMT 00:00点起始)，网格分辨率为2.5°×2.5°，风场垂直方向从1 000～10 hPa，共17层[27,28]。值得注意的是：压力项的计算中，由于选取的是1962―2021年的长期数据，所以应该考虑海平面高度随气压的变化。因此，本文使用反变气压(inverted barometer,IB)近似，即假设海洋对大气压力变化存在响应[29,30]。在风项计算时，本文依据前人工作考虑了地形因素的影响，即从不同地形起伏所对应的高度垂直积分至10 hPa[31]。轴向大气角动量函数(AAMFχ3)为压力项与风项之和。为了与日长变化的观测单位(ms)相一致，还需将χ3乘上转换因子8.64×107ms·rad−1。最后，将所得序列扣除与ΔLO D类似的拟合季节项后，经AR模型延拓和1～10 a频段滤波以及端部效应扣除，得到与ENSO指标和ΔLO D相对应的AAM年际分量。

3 对比分析结果

图5显示的是ΔLO D,AAM与MSOI的年际分量对比，可见在年际尺度上三者有相似的起伏趋势，典型的1982―1983年、1997―1998年和2015―2016年超强厄尔尼诺事件分别体现在这三者的振幅上。图6给出了以上年际分量时间序列的超前滞后相关情况，其中相关系数检验均超过99%的显著性水平[32]。可见MSOI分别超前于AAM(蓝色实线)和ΔLO D(红色实线)约1月，且MSOI与AAM和ΔLO D的最大相关系数分别达到0.69和0.55。另一方面，AAM与ΔLO D同步(绿色实线)，相关系数达到0.70。这反映了在年际尺度上，固体地球与表层大气海洋系统之间的强烈耦合作用，同时也验证了在外力矩为0时大气与固体地球角动量守恒下的瞬时转换过程。图6所得结论与前人研究基本一致[7,11−13,33]。

图5 1962年1月—2021年1月间ΔLO D,AAM,MSOI的年际分量对比

图6 1962年1月—2021年1月间的年际分量时间序列的超前滞后相关情况

为了进一步分析年际尺度上日长变化、大气角动量激发和南方涛动的时频特征，本文引入正则化Morlet小波[34]。设时间序列为f(t)，定义小波函数为：

其中，Ψ(t)为基本小波，a为确定频率特征的伸缩尺度因子，b为时间域上的平移因子，小波变化可以同时用于分析准周期时变序列的时间域和频率域特征。图7给出ΔLO D,AAM和MSOI在年际分量上的小波谱结果，其中红色对应正相位，蓝色对应负相位，颜色越深表示强度越大。可见ΔLO D,AAM和MSOI在年际分量上具有相似的时频分布结构。

图7 1962年1月—2021年1月间ΔLO D,AAM和MSOI的年际分量小波分析谱

以上统计分析结果反映了ENSO通过影响大气角动量进而导致地球自转速率发生改变的过程。总体来说，当厄尔尼诺出现时，MSOI指数为正(S O I指数为负)，表征赤道东风带的信风减弱，出现西风异常，这将导致纬向沃克(Walker)环流减弱；同时，由于厄尔尼诺期间的热带太平洋暖海水异常加热大气，使得经向哈德来(Hadley)环流增强，进一步使副热带西风急流增强。图6显示年际分量的MSOI超前于AAM与ΔLO D1月左右，可能是由于热带海气相互作用通过大气的经向环流输送使副热带出现异常响应所需要的时间。在厄尔尼诺期间热带太平洋上空出现的西风异常以及副热带西风急流的增强等过程共同导致了AAM因西风分量的增加而增强。由于固体地球与大气的系统总角动量守恒，AAM的增加使得地球自转减速，ΔLO D增加。当厄尔尼诺通过皮耶克尼斯正反馈(Bjerknes feedback)[35]达到峰值后开始逐渐消亡时，其产生的大气环流异常也会逐渐恢复，在此期间AAM和ΔLO D会相对降低，因此通常出现AAM和ΔLO D极大值。前人在研究中也提及了类似的反馈过程[8,36,37]。

反之，在拉尼娜期间，由于赤道东太平洋出现冷海水异常，赤道东风增强，纬向沃克环流增强，经向哈德来环流减弱，副热带西风急流减弱。热带和副热带的东风分量增强，导致AAM减弱，根据地气系统角动量守恒，地球自转速率增加，ΔLO D降低。同理，当拉尼娜开始消亡时，AAM和ΔLO D往往会相对增加，而通常导致AAM和ΔLO D极小值的出现。

值得注意的是，尽管ENSO是气候系统年际尺度的最强变率，也是影响日长变化年际分量的重要原因，但由于导致地球自转速率变化的原因非常复杂，所以ENSO并不能完全解释年际尺度的日长变化。平流层准两年振荡(quasi-biennial oscillation,QBO)等气候模态和一些地球内部物理过程也可能对日长变化的年际分量产生影响[26,38]。

图8显示了大致时间段内ΔLO D年际分量与O N I的对比情况。其中，如O N I部分所示，截至2021年1月末，2020年7月至2021年1月期间Nino3.4海区SSTA连续3月滑动平均值所得的O N I指数已经满足连续5月S S T A≤-0.5℃的拉尼娜事件判定条件。O N I指数的峰值出现在2020年11月(Nino3.4 OND SSTA平均)达到-1.3℃，根据NCAR的标准[18]，是一次中等强度的拉尼娜事件。同时，将图中的ΔLO D年际分量与O N I序列(厄尔尼诺事件(红色填充)和拉尼娜事件(蓝色填充))进行对比，可以发现，年际分量的ΔLO D并不完全与以O N I指数为判定标准的ENSO冷暖事件在起止时间和强度上严格对应，这可能是由于两种数据分辨率不同、数据处理方法导致的误差以及影响ΔLO D因素的复杂性所致。按照基于O N I的判定标准，从1962年1月至2020年1月，共出现了19次厄尔尼诺事件以及16次拉尼娜事件，如果忽略上述导致O N I和ΔLO D对应上出现误差因素的影响，几乎所有的厄尔尼诺和拉尼娜事件均使年际分量的ΔLO D出现对应的极大值和极小值。对于目前仍在持续的2020年夏季至2021年春季中等强度拉尼娜事件，年际ΔLO D在对应时期出现幅度约-0.18 ms的变化(图8a)右侧)，我们将持续关注本次事件的后续进展。

图8 ΔLO D年际分量与O N I的对比图

4 总结与讨论

本文基于天文观测的日长变化序列、大气角动量和表征ENSO的S O I与O N I指数，分析了三者在年际时间尺度上的关联性，并尝试给出物理过程描述，同时检测到年际日长变化中的2020―2021年拉尼娜事件信号。

由ENSO、大气角动量及日长变化得到的互相关关系，除了验证O N I与MSOI这两个分别表征EN与SO的指数具有很好的相关性和同步性以外，还显示了ΔLO D,AAM和MSOI在统计意义上的超前滞后关系。可知ΔLO D与AAM基本同步，反映了大气与固体地球角动量在系统总角动量守恒条件下的瞬时转换过程。而MSOI分别超前于AAM和ΔLO D约1月，则从统计角度说明，ENSO是ΔLO D在年际尺度上的重要激发源。

基于现有的观测资料，从1962年1月至2021年1月，共出现了19次厄尔尼诺事件和17次拉尼娜事件。由于厄尔尼诺和拉尼娜会导致热带及副热带大气环流出现异常(反映为AAM的变化)，因此也会使LO D出现相应改变。对于2020年夏季至2021年春季发生的中等强度拉尼娜事件，反映在地球自转上是转动速率加快，且年际日长变化序列在对应时期出现了幅度约―0.18 ms的变化。截至2021年2月底，本次拉尼娜事件仍在延续，我们将继续关注相关的天文和气象发展过程。

致谢

感谢国际地球自转与参考架服务(IERS)和美国国家环境预报中心/大气研究中心(NCEP/NCAR)提供的数据资料。感谢审稿专家及编辑部老师对本文提出的宝贵意见。