地球自转运动可以用地球定向参数(Earth Orientation Parameters, EOP)表征。地球定向参数包括岁差、章动、UT1-UTC及极移的两个分量x_p和y_p。地球定向参数是实现地球参考坐标系与天球参考坐标系相互转换的必需参数, 在天文地球动力学研究、深空探测及卫星导航等领域有重要应用。甚长干涉基线测量(Very Long Baseline Interferometry, VLBI)、多普勒卫星定轨和无线电定位(Doppler Orbitography and Radiopositioning Integrated by Satellite, DORIS)以及全球卫星导航系统(Global Navigation Satellite System, GNSS)等空间测地技术是测量地球定向参数的主要手段, 其中UT1-UTC和极移的测量精度可以分别达到10 μs和100 μas^[1], 但复杂的资料处理过程致使地球定向参数的获取存在一定的延迟。由于深空探测等领域对地球定向参数实时测量值有重要需求, 而卫星自主导航又对地球定向参数中长期预报值有重要需求, 因此对地球定向参数进行高精度的短期与中长期预报是非常必要的。在地球定向参数的几个分量中, UT1-UTC是变化最快、最难预报的一个分量。

目前有多种UT1-UTC预报方法, 其中多数方法结合最小二乘(Least Squares, LS)外推和其它模型进行UT1-UTC预报^[2-4], 这些方法的思路为, 首先采用最小二乘外推模型提取UT1-UTC序列中的周期性成分, 并对其进行外推, 然后利用神经网络(Neural Network, NN)或者自回归(Auto Regressive, AR)模型等对最小二乘拟合残差进行建模、预报, 最后再将周期项和残差项的外推值相加得到UT1-UTC预报值。实际应用中发现, 在利用最小二乘外推模型对UT1-UTC观测资料进行拟合时, 在拟合序列的两端存在发散畸变现象, 这种现象在数据处理中称为端部效应^[5-7]。端部效应使残差项与周期项的预报值出现偏差, 最终导致UT1-UTC预报值不准确。本文针对UT1-UTC预报中最小二乘拟合出现的端部畸变现象, 在采用最小二乘外推模型对UT1-UTC序列进行拟合之前, 先利用端部延拓方法对UT1-UTC观测资料进行数据延拓, 即在UT1-UTC序列的两端增加应用时序分析方法延拓出的若干数据点, 形成一个新序列, 然后用新序列求解最小二乘外推模型系数, 最后再基于最小二乘外推模型对UT1-UTC序列中的周期项进行外推。数值分析表明, 通过在UT1-UTC观测序列端部增加统计延拓数据, 可以有效地抑制端部效应, 从而改进UT1-UTC的预报效果。

1 预报方法 1.1 数据预处理

本文UT1-UTC观测资料来源于国际地球自转与参考系服务(International Earth Rotation and Reference Systems Service, IERS)发布的EOP-08-C04序列, 其中UT1-UTC数据的采样间隔为1d。UT1-UTC观测资料中含有闰秒(Leap Second)及多种周期项、准周期项, 其中, 对于周期为5 d~18.6 a的62个固体地球带谐潮汐项应用IERS协议(IERS Conventions)给出的经验公式予以扣除^[8], 扣除62个固体地球带谐潮汐项后的UT1-UTC称作UT1R-UTC, 然后再去掉UT1R-UTC序列中的闰秒, 获得UT1R-TAI序列。本文对UT1-UTC的预测实质上是针对UT1R-TAI的预测。

1.2 最小二乘外推模型

UT1-UTC序列扣除62个固体地球带谐潮汐项后, 还含有长期趋势项、周年项、半周年项等周期性变化成分, 长期趋势项与周年项、半周年项利用如下模型进行拟合、外推:

$ {f_{{\rm{UT1R- TAI}}}} = a + bt + \sum\limits_{i = 1}^2 {\left[{{c_i}{\rm{cos}}(2{\rm{ \mathsf{ π} }}t/{T_i}) + {d_i}{\rm{sin}}(2{\rm{ \mathsf{ π} }}t/{T_i})} \right], } $

(1)

其中, T₁, T₂分别表示半周年项和周年项的振荡周期, 取T₁=182.62 d, T₂=365.24 d^[2]; a, b为趋势项参数; c₁, d₁为半周年项参数; c₂, d₂为周年项参数。这些未知参数可以根据最小二乘法求解。

1.3 神经网络模型

反向传播(Back Propagation, BP)神经网络是当前应用较为广泛的神经网络算法之一, 但这种算法存在收敛速度慢、过拟合及泛化能力弱等不足。鉴于此, 本文选用近年来发展起来的极限学习机(Extreme Learning Machine, ELM)神经网络算法对最小二乘拟合残差序列进行建模、预报。与反向传播神经网络算法相比, 极限学习机算法训练速度非常快, 不需要人为设置网络参数, 且具有很强的泛化能力^[9], 极限学习机神经网络在地球定向参数预报中已经成功应用^[10], 其原理如下所述。

给定N个训练样本(x_i, y_i), 其中输入向量x_i=[x_i1, x_i2, …, x_im]^T∈R^m, 输出向量y_i=[y_i1, y_i2, …, y_in]^T∈Rⁿ, m, n分别表示输入向量维数与输出向量维数, 则具有h个隐含层节点的单隐层前馈神经网络能够以零误差逼近这N个样本, 即有下式成立:

$ {\mathit{\boldsymbol{y}}_i} = \sum\nolimits_{j = 1}^h {{\mathit{\boldsymbol{\beta }}_j}} g({\mathit{\boldsymbol{w}}_j}\cdot{\mathit{\boldsymbol{x}}_i} + {b_j}), i = 1, 2, \ldots, N, $

(2)

其中, w_j=[w_j1, w_j2, …, w_jm]^T为第j个隐含层节点和输入节点之间的连接权值; β_j=[β_j1, β_j2, …, β_jn]^T为第j个隐含层节点和输出节点之间的连接权值; b_j为第j个隐含层节点的偏置; g(x)为激活函数。(2)式可以表示为矩阵相乘的形式:

$ \mathit{\boldsymbol{Y}} = \mathit{\boldsymbol{H\beta }}, $

(3)

其中, H为隐含层输出矩阵:

$ \mathit{\boldsymbol{H}} = {\left[{\begin{array}{*{20}{c}} {g\left( {{\mathit{\boldsymbol{w}}_1}\cdot{\mathit{\boldsymbol{x}}_1} + {b_1}} \right)}& \ldots &{g({\mathit{\boldsymbol{w}}_h}\cdot{\mathit{\boldsymbol{x}}_1} + {b_h})}\\ \vdots&\ldots&\vdots \\ {g\left( {{\mathit{\boldsymbol{w}}_1}\cdot{\mathit{\boldsymbol{x}}_N} + {b_1}} \right)}& \ldots &{g({\mathit{\boldsymbol{w}}_h}\cdot{\mathit{\boldsymbol{x}}_N} + {b_h})} \end{array}} \right]_{N \times h}}, $

(4)

$ \mathit{\boldsymbol{\beta }} = {\left[{\begin{array}{*{20}{c}} {{\mathit{\boldsymbol{\beta }}_1}^{\rm{T}}}\\ {{\mathit{\boldsymbol{\beta }}_2}^{\rm{T}}}\\ \vdots \\ {{\mathit{\boldsymbol{\beta }}_h}^{\rm{T}}} \end{array}} \right]_{h \times n}}, \mathit{\boldsymbol{Y}} = {\left[{\begin{array}{*{20}{c}} {{\mathit{\boldsymbol{y}}_1}^{\rm{T}}}\\ {{\mathit{\boldsymbol{y}}_2}^{\rm{T}}}\\ \vdots \\ {{\mathit{\boldsymbol{y}}_N}^{\rm{T}}} \end{array}} \right]_{N \times n}}. $

(5)

极限学习机神经网络的训练步骤:(1)随机选取输入层权值w_j、偏置b_j, j=1, 2, …, h; (2)计算隐含层输出矩阵H; (3)计算输出层权值β=H^†Y, 其中H^†为矩阵H的Moore-Penrose广义逆。

1.4 预报过程

由于最小二乘外推和神经网络组合模型(LS + NN)不仅具有较高的预报精度, 还具有较好的长期预测稳定性, 故基于最小二乘外推和神经网络组合模型的UT1-UTC短期预测值精度较高, 可以很好地反映UT1-UTC的变化规律^{[1, 3]}。作者曾将端部数据延拓的思想和方法应用于日长变化资料分析中, 在改善最小二乘拟合的端部效应方面进行了研究, 取得了较为明显的效果^[6-7]。本文尝试将端部数据延拓的思想和方法应用于UT1-UTC预报, 以限制UT1-UTC资料拟合中最小二乘拟合出现的端部畸变现象。

最小二乘拟合端部效应的UT1-UTC预测方法和常规方法的差别之处在于, 该方法在对UT1-UTC观测序列建立趋势项及周期项最小二乘外推模型之前, 先应用统计学方法在UT1R-TAI序列首部和尾部进行数据延拓, 以抑制最小二乘拟合端部畸变, 预报过程如下:

(1) 首先通过(1)式对UT1R-TAI序列作最小二乘拟合, 建立周期项及趋势项外推模型, 然后利用神经网络模型对最小二乘拟合残差序列进行建模、预报, 最后联合最小二乘外推模型和神经网络模型在UT1R-TAI序列首部和尾部分别外推若干数据点, 这样UT1R-TAI序列加上首部和尾部外推的若干数据点就构成了一个新序列。

(2) 利用新序列求解最小二乘外推模型系数, 即用新序列重新建立趋势项及周期项最小二乘外推模型, 然后再结合最小二乘外推模型和神经网络模型对UT1R-TAI序列作预报。

(3) 将周期为5 d~18.6 a的62个固体地球带谐潮汐项及闰秒恢复到UT1R-TAI序列的预测值中即可得到UT1-UTC的预测值。

2 数值分析

选取2001.1.1~2016.9.20期间的UT1-UTC观测序列进行数值分析, 其中2010.1.1~2016.9.20为预测期, 建模数据长度为10年, 每隔7 d预报1次, 总共进行了300次预测。图 1(a)、(b)分别给出了2001.1.1~2016.9.20时段的UT1-UTC测量序列及其扣除闰秒和62个固体地球带谐潮汐项的UT1R-TAI序列。

为了验证端部延拓方法对最小二乘拟合端部效应的改善效果, 首先比较端部延拓前后最小二乘拟合效果。图 2绘出了2001年1月1日至2010年12月31日时段UT1R-TAI序列及其端部数据延拓点, 其中, 端部延拓数据点数为360, 即在首尾两端各延拓180个数据点。从图 2可以看到, UT1R-TAI端部延拓数据点与原始序列符合得较好, 可以很好地反映原始UT1R-TAI序列的变化规律。图 3绘出了端部数据延拓前后2001年1月1日至2010年12月31日时段UT1R-TAI序列的最小二乘拟合序列, 基于端部效应改善的最小二乘(Edge-effect Corrected Least Squares, ECLS)表示端部效应修正的最小二乘拟合序列。为了更加直观、清楚地展示端部数据延拓方法对最小二乘拟合端部畸变的改善效果, 对图 3作局部放大处理, 图 4(a)绘出了端部数据延拓前后最小二乘拟合序列首部前180个历元的最小二乘拟合残差数据点, 图 4(b)绘出了端部数据延拓前后最小二乘拟合序列尾部最后180个历元的最小二乘拟合残差数据点。

从图 3和图 4可以发现, 与直接对原始UT1R-TAI序列进行拟合相比, 端部数据延拓后最小二乘拟合残差序列在首部和尾部更加接近于0, 换言之, 端点延拓后最小二乘拟合的UT1R-TAI序列在首部、尾部和原始序列吻合得更好, 这说明端点延拓方法能够有效地抑制最小二乘拟合出现的端部畸变现象。

为了检验端部效应改善的最小二乘外推模型对UT1-UTC预报的改善效果, 分别利用最小二乘外推和神经网络组合模型和基于端部效应改善的最小二乘和神经网络组合模型对UT1-UTC作1~360 d时长预报, 采用平均绝对误差(Mean Absolute Error, MAE)作为预测性能评价指标, 其计算公式可以写为

$ MA{E_i} = \frac{1}{M}\sum\limits_{j = 1}^M {\left| {P_j^i-O_j^i} \right|}, $

(6)

其中, P_j, O_j分别为j点的UT1-UTC预测值及观测值; i为预报长度; M为预报期数, 本文总共作了300期的预测, 即M=300。

图 5给出了最小二乘外推和神经网络组合模型与基于端部效应改善的最小二乘和神经网络组合的UT1-UTC预测平均绝对误差对比图, 表 1给出了最小二乘外推和神经网络组合模型与基于端部效应改善的最小二乘和神经网络组合模型在不同预报长度下的平均绝对误差统计结果。从图 5及表 1可以发现, 对于1~30 d的短期预报, 基于端部效应改善的最小二乘和神经网络组合模型的预报精度相对于常规最小二乘外推和神经网络组合模型的预报精度并无改善, 但从第30 d开始, 基于端部效应改善的最小二乘和神经网络组合模型的预报精度明显优于常规最小二乘外推和神经网络组合模型, 预报精度最大提高了27%, 且一直保持在10%以上, 这说明与常规最小二乘外推模型相比较而言, 端部效应改善的最小二乘外推模型对于UT1-UTC中长期预报具有更明显的优势, 同时也从侧面反映出最小二乘拟合的端部畸变现象对UT1-UTC中长期预报的影响更大。

3 结论

本文提出了一种顾及最小二乘拟合端部效应的UT1-UTC预报方法, 这种方法与基于常规最小二乘外推模型的UT1-UTC预报方法区别在于, 该方法在进行最小二乘拟合之前, 首先在原始序列两端增加统计延拓数据, 然后再对数据延拓后的新序列进行最小二乘拟合, 目的是将最小二乘拟合存在的端部畸变搬移到新序列的两端, 从而抑制最小二乘拟合序列的端部畸变。数值分析表明, 最小二乘拟合的端部效应主要影响UT1-UTC中长期预报结果, 通过在观测资料的两端增加用统计学方法延拓的外推数据点, 然后再进行最小二乘拟合, 能够有效地改善最小二乘拟合出现的端部畸变现象。与常规最小二乘外推模型相比, 基于端部效应改善的最小二乘外推模型的UT1-UTC中长期预报精度提升明显。