考虑最小二乘拟合端部效应的UT1-UTC预报

雷雨,蔡宏兵,高玉平

(1.中国科学院 国家授时中心,西安 710600；2.中国科学院 时间频率基准重点实验室,西安 710600)

0 引言

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

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

1 预报方法

本节首先介绍最小二乘和AR组合模型预报UT1-UTC的原理,在此基础之上,建立考虑最小二乘拟合端部效应的UT1-UTC预报算法。

1.1 数据预处理

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

1.2 最小二乘外推模型

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

fUT1R-TAI=a+bt+c1cos(2πt/T1)+d1sin(2πt/T1)+c2cos(2πt/T2)+d2sin(2πt/T2),

(1)

式(1)中,T1、T2分别表示半周年项和周年项的振荡周期,取T1=182.62 d、T2=365.24 d[4],a、b代表趋势项参数,c1、d1代表半周年项参数,c2、d2代表周年项参数,这些未知参数可以根据最小二乘法求解。

1.3 自回归模型

AR模型是对平稳时间序列{zt,t=1,2,…,n}建立的一个概率统计模型,它根据变量自身的历史变化规律来建立统计模型,其数学公式为

(2)

式(2)中,p为模型阶数,εt为白噪声,φ1,φ2,…,φp为模型参数,可以通过求解Yule-Walker方程来获得[18-19]。

运用AR模型的关键之处在于选取模型阶数p。已有的定阶准则包括传递函数准则、预测误差准则及信息论准则[11]。本文选用赤池信息量准则(Akaike information criterion,AIC)来确定模型阶数p,AIC准则的目标函数可以表示为

(3)

在模型阶数p与模型参数φ1,φ2,…,φp确定以后,可以通过如下方式对时间序列作多步外推：

(4)

1.4 预报过程

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

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

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

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

图1为考虑最小二乘拟合端部效应的UT1-UTC预测方法的流程图。

图1 考虑最小二乘拟合端部效应的UT1-UTC预测流程

2 数值分析

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

图2 UT1-UTC观测序列及其UT1R-TAI序列

为了验证端部延拓方法对最小二乘拟合端部效应的改善效果,首先比较端部延拓前后最小二乘拟合效果。图3绘出了端部数据延拓前后2001-01-01至2010-12-31时段UT1R-TAI序列的最小二乘拟合序列,其中,端部延拓数据点数为200,即在首尾两端各延拓100个数据点,ECLS(edge-effect corrected least squares)表示端部效应修正的最小二乘拟合序列。为了更加直观地、清楚地展示端部数据延拓方法对LS拟合端部畸变的改善效果,对图3作局部放大处理,图4(a)绘出了端部数据延拓前后最小二乘拟合序列首部前100个历元的最小二乘拟合残差数据点,图4(b)绘出了端部数据延拓前后最小二乘拟合序列尾部最后100个历元的最小二乘拟合残差数据点。

图3 UT1R-TAI的LS与ECLS拟合序列

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

图4 UT1R-TAI的LS与ECLS拟合的首尾两端残差序列

为了进一步检验端部数据延拓方法对最小二乘拟合端部效应的改善效果,分别利用LS+AR模型和ECLS+AR模型对UT1-UTC作1～360 d时长预报,采用平均绝对误差(mean absolute error,MAE)作为预测性能评价指标,其计算公式为[20]

(5)

式(5)中,Pj和Oj分别为j点的UT1-UTC预测值及观测值,i表示预报长度,σi表示预报长度为i时的UT1-UTC预测平均绝对误差,M表示预报期数,本文总共作了286期的预测,即M=286。

图5给出了LS+AR方法和ECLS+AR方法的UT1-UTC预测平均绝对误差对比图,表1给出了LS+AR方法和ECLS+AR方法在不同预报跨度下的平均绝对误差统计结果,其中将1～15 d跨度的预报称为短期预报,将大于15 d跨度的预报称为中长期预报。从图5及表1可以发现,对于1～15 d的短期预报,ECLS+AR模型的预报精度相对于常规LS+AR模型的预报精度并无改善,但从第15 d开始,ECLS+AR模型的预报精度明显优于常规LS+AR模型,精度最大提高了34%,且一直保持在15%以上,这说明与常规LS+AR模型比较而言,ECLS+AR模型对于UT1-UTC中长期预报具有更明显的优势,同时也从侧面反映出最小二乘拟合的端部畸变现象对UT1-UTC中长期预报的影响更大。

图5 LS+AR模型与ECLS+AR模型的UT1-UTC预报MAE对比

表1 LS+AR模型与ECLS+AR模型的UT1-UTC预报MAE统计

3 结语

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