第 56 卷 第 1 期 天 文 学 报 Vol.56 No.1 2015 年 1 月 ACTA ASTRONOMICA SINICA Jan., 2015 doi: 10.15940/j.cnki.0001-5245.2015.01.007 基 于 高 斯 过 程 的 日 长 变 化 预 报 雷 雨 1,2,3 赵 丹 宁 1,3 高 玉 平 1,2 蔡 宏 兵 1,2 (1 中 国 科 学 院 国 家 授 时 中 心 西 安 710600) (2 中 国 科 学 院 时 间 频 率 基 准 重 点 实 验 室 西 安 710600) (3 中 国 科 学 院 大 学 北 京 100049) 摘 要 由 于 日 长 (length-of-day, LOD) 变 化 具 有 复 杂 的 时 变 特 性, 传 统 线 性 模 型 如 最 小 二 乘 外 推 模 型 时 间 序 列 分 析 模 型 等 的 预 报 效 果 往 往 不 甚 理 想, 所 以 将 一 种 新 型 的 机 器 学 习 算 法 高 斯 过 程 (Gaussian processes, GP) 方 法 用 于 LOD 变 化 预 报, 并 将 预 报 结 果 同 利 用 反 向 传 播 神 经 网 络 (back propagation neural networks, BPNN) 和 广 义 回 归 神 经 网 络 (general regression neural networks, GRNN) 的 预 报 结 果 以 及 地 球 定 向 参 数 预 报 比 较 竞 赛 (Earth Orientation Parameters Prediction Comparison Campaign, EOP PCC) 的 预 报 结 果 进 行 对 比. 结 果 表 明, GP 用 于 LOD 变 化 预 报 是 高 效 可 行 的. 关 键 词 天 体 测 量, 时 间, 方 法 : 数 据 分 析 中 图 分 类 号 : P127; 文 献 标 识 码 : A 1 引 言 LOD 变 化 是 表 征 地 球 自 转 变 化 的 一 个 重 要 参 数, 它 是 指 天 文 意 义 上 的 1 d 和 标 准 日 长 86 400 s 之 间 的 差 异, 反 映 了 地 球 自 转 速 率 的 变 化. LOD 和 极 移 (polar motion, PM) 统 称 为 地 球 自 转 参 数 (Earth rotation parameters, ERP). ERP 是 实 现 天 球 参 考 系 和 地 球 参 考 系 之 间 相 互 转 换 的 必 需 参 数, 在 深 空 探 测 卫 星 精 密 定 轨 和 天 文 地 球 动 力 研 究 等 领 域 都 有 重 要 应 用 [1]. 现 代 测 地 技 术 ( 甚 长 干 涉 基 线 (Very Long Baseline Interferometry) 全 球 卫 星 导 航 系 统 (Global Navigation Satellite Systems, GNSS) 和 卫 星 激 光 测 距 (Satellite Laser Ranging, SLR) 等 ) 被 广 泛 应 用 于 地 球 自 转 变 化 的 常 规 监 测 中, 提 供 了 高 时 空 分 辨 率 和 高 精 度 的 观 测 资 料. 然 而, 由 于 复 杂 的 数 据 处 理 过 程, 由 现 代 测 地 技 术 获 取 的 ERP 往 往 需 要 延 迟 几 天 甚 至 2 个 星 期, 所 以 对 ERP 进 行 实 时 快 速 的 预 报 成 为 一 项 值 得 深 入 研 究 的 课 题. LOD 变 化 的 精 确 预 报 是 ERP 预 报 中 的 难 点 之 一, 特 别 是 在 厄 尔 尼 诺 (El Niño) 事 件 发 生 期 间, 热 带 季 风 的 变 化 导 致 LOD 变 化 出 现 大 幅 振 荡. LOD 变 化 的 高 精 度 实 时 快 速 预 报 引 起 了 越 来 越 多 学 者 的 关 注. LOD 变 化 主 要 由 潮 汐 项 和 非 潮 汐 项 2 部 分 组 成, 2014-07-02 收 到 原 稿, 2014-08-04 收 到 修 改 稿 国 家 自 然 科 学 基 金 项 目 (10573019) 资 助 leiyu@ntsc.ac.cn
54 天 文 学 报 56 卷 潮 汐 项 可 以 由 国 际 地 球 自 转 与 参 考 系 服 务 (International Earth Rotation and Reference Systems Service, IERS) 协 议 给 出 的 模 型 精 确 确 定 [2], 而 非 潮 汐 项 中 的 半 年 项 和 周 年 项 等 季 节 性 变 化 主 要 是 由 固 体 地 球 和 全 球 大 气 海 洋 以 及 地 下 水 之 间 的 角 动 量 交 换 引 起 的 [3]. 学 者 们 在 ERP 预 报 方 面 已 经 做 了 许 多 研 究, 提 出 了 各 种 预 报 模 型, 包 括 最 小 二 乘 (least squares, LS) 外 推 模 型 [4] LS 外 推 模 型 和 自 回 归 (autoregressive, AR) 模 型 的 组 合 (LS+AR) [4 5] 卡 尔 曼 滤 波 联 合 大 气 角 动 量 (Kalman Filter+OAM) [6 7] LS 外 推 模 型 和 人 工 神 经 网 络 (artificial neural networks, ANN) 模 型 的 组 合 (LS+ANN) [3,8 11] 模 糊 推 理 系 统 (fuzzy-inference systems, FIS) [12] 以 及 离 散 小 波 变 换 (discrete wavelet transform, DWT) 和 自 协 方 差 (autocovariance, AC) 模 型 的 组 合 (DWT+AC) [13] 等. 为 了 对 比 不 同 模 型 的 预 报 效 果, 维 也 纳 理 工 大 学 大 地 测 量 与 地 球 物 理 研 究 所 从 2005 年 10 月 1 日 至 2008 年 2 月 28 日 组 织 了 全 球 性 的 地 球 定 向 参 数 预 报 比 较 竞 赛, 2 yr 多 的 预 报 结 果 表 明, 没 有 一 种 模 型 既 适 合 于 ERP 所 有 分 量 的 预 报 又 适 合 于 所 有 跨 度 的 预 报 [14]. 受 多 种 激 发 因 素 的 影 响, 地 球 自 转 变 化 呈 现 复 杂 的 非 线 性 不 规 则 变 化 特 性, 因 此 采 用 非 线 性 的 预 报 方 法 对 其 进 行 预 报 在 理 论 上 更 为 合 理 [8 11]. ANN 是 逼 近 复 杂 非 线 性 函 数 的 一 种 有 效 工 具, 所 以 有 许 多 学 者 将 其 应 用 于 ERP 预 报 中, 并 取 得 了 显 著 的 预 报 效 果 [3,5,8 11]. 但 是 ANN 存 在 一 些 缺 点, 例 如 网 络 拓 扑 结 构 难 以 确 定 训 练 过 程 存 在 过 学 习 现 象 迭 代 过 程 易 陷 入 局 部 最 优 收 敛 速 度 较 慢. 此 外, ANN 的 优 化 目 标 是 基 于 经 验 的 风 险 最 小 化, 无 法 保 证 网 络 的 泛 化 能 力 [15]. 针 对 上 述 缺 点 很 多 学 者 提 出 了 改 进 措 施, 例 如 将 ANN 和 其 他 人 工 智 能 (artificial intelligence, AI) 算 法 结 合 进 行 网 络 优 化, 如 遗 传 算 法 与 ANN 的 组 合 粒 子 群 算 法 与 ANN 的 组 合 等, 然 而 这 些 算 法 均 存 在 一 定 的 不 足, 仍 处 于 不 断 的 尝 试 研 究 阶 段. GP 是 近 年 来 发 展 起 来 的 一 种 新 型 机 器 学 习 算 法, 它 有 着 严 格 的 统 计 学 习 理 论 基 础, 对 处 理 高 维 数 小 样 本 非 线 性 等 复 杂 问 题 具 有 很 好 的 适 应 性, 且 泛 化 能 力 强 [16 17]. 与 ANN 和 支 持 向 量 机 相 比, GP 具 有 容 易 实 现 超 参 数 自 适 应 获 取 非 参 数 推 断 灵 活 以 及 输 出 具 有 概 率 意 义 等 优 点 [16 17]. 现 已 成 为 机 器 学 习 领 域 的 研 究 热 点, 并 在 许 多 领 域 得 到 了 成 功 的 应 用 [18 19]. 本 文 将 GP 算 法 用 于 LOD 变 化 预 报, 研 究 表 明 : GP 算 法 用 于 LOD 变 化 预 报 是 可 行 的, 且 预 报 效 率 和 精 度 较 高. 2 GP 的 基 本 原 理 及 其 用 于 LOD 变 化 预 报 的 方 法 2.1 GP 回 归 的 基 本 原 理 高 斯 过 程 又 称 正 态 随 机 过 程, 其 任 意 有 限 变 量 集 合 都 有 着 联 合 高 斯 分 布 的 特 性, 即 对 于 任 意 的 变 量 x 1, x 2,, x n 与 其 对 应 的 函 数 f(x 1 ), f(x 2 ),, f(x n ) 的 联 合 概 率 分 布 服 从 n 维 高 斯 分 布. 高 斯 分 布 的 全 部 统 计 特 征 完 全 由 它 的 均 值 函 数 m(x) 和 协 方 差 函 数 C(x, x ) 来 确 定, 一 般 记 为 f(x) GP(m(x), C(x, x )). 若 给 定 训 练 样 本 集 D = {(x i, y i ) x i R d, y i R, i = 1, 2,, n}, 其 中 d 为 向 量 x i 的 维 数, 则 对 于 测 试 样 本 输 入 x, GP 模 型 的 预 测 值 为 y = k T (x )K 1 y, (1)
1 期 雷 雨 等 : 基 于 高 斯 过 程 的 日 长 变 化 预 报 55 其 中 k(x ) = [C(x, x 1 ), C(x, x 2 ),, C(x, x n )] 为 测 试 样 本 输 入 和 训 练 样 本 输 入 值 之 间 的 1 n 阶 协 方 差 矩 阵, K 是 训 练 样 本 输 入 值 之 间 的 n n 阶 协 方 差 矩 阵, K ij = C(x i, x j ). (1) 式 表 明, GP 模 型 可 以 根 据 协 方 差 函 数 和 测 试 样 本 进 行 预 测. 协 方 差 函 数 在 GP 回 归 模 型 中 起 到 关 键 作 用, 它 表 达 了 一 种 样 本 间 的 相 似 性, 对 所 要 学 习 的 函 数 提 供 了 假 设 信 息, 协 方 差 函 数 必 须 是 半 正 函 数. 常 用 的 协 方 差 函 数 为 平 方 指 数 函 数 [17], 即 [ C(x i, x j ) = σf 2 exp 1 d ( ) (m) x i x (m) 2 ] j + σ 2 l nδ 2 ij, (2) m m=1 其 中 x (m) i x (m) j 分 别 是 向 量 x i x j 的 第 m 个 分 量, l 1, l 2,, l d 为 协 方 差 的 尺 度 因 子, σf 2 为 信 号 方 差, σn 2 为 噪 声 方 差, 称 参 数 集 θ = {l 1, l 2,, l d, σ f, σ n } 为 超 参 数 ; δ ij 为 Kronecker 函 数, 当 i = j 时 δ ij = 1, 当 i j 时 δ ij = 0. 2.2 基 于 GP 回 归 模 型 的 LOD 变 化 预 报 2.2.1 数 据 预 处 理 本 文 所 用 的 LOD 变 化 数 据 来 自 IERS 发 布 的 EOP C04 序 列, 采 样 间 隔 为 1 d. LOD 变 化 序 列 中 周 期 为 5 d 18.6 yr 的 固 体 地 球 带 谐 潮 汐 项 可 以 通 过 IERS 协 议 给 出 的 经 验 模 型 精 确 确 定 [2], 近 周 日 和 半 周 日 海 洋 潮 汐 项 不 作 修 正, LOD 变 化 的 长 期 趋 势 项 季 节 性 变 化 的 周 年 和 半 年 项 等 根 据 下 述 线 性 模 型 确 定 [4] : 4 f LODR (t) = a + bt + [c i,1 cos(2πt/p i ) + c i,2 sin(2πt/p i )], (3) i=1 其 中 LODR 表 示 经 过 固 体 地 球 带 谐 潮 修 正 后 的 LOD 变 化 序 列, p 1 p 2 p 3 和 p 4 分 别 表 示 半 年 项 周 年 项 9.3 yr 项 和 18.6 yr 项 的 周 期, 取 p 1 = 182.62 d p 2 = 365.24 d p 3 = 3396.732 d p 4 = 6793.464 d, t 为 协 调 世 界 时 (Coordinated Universal Time, UTC), 在 拟 合 时 单 位 转 换 为 d. a b 表 示 长 期 趋 势 项 的 参 数, c 1,1 c 1,2 表 示 半 年 项 的 参 数, c 2,1 c 2,2 表 示 周 年 项 的 参 数, c 3,1 c 3,2 表 示 9.3 yr 项 的 参 数, c 4,1 c 4,2 表 示 18.6 yr 项 的 参 数, 这 10 个 未 知 参 数 通 过 最 小 二 乘 法 求 得. 经 过 上 述 数 据 预 处 理 后 的 剩 余 部 分 为 含 有 非 线 性 成 分 的 残 差 序 列, 主 要 包 括 海 洋 近 周 日 半 日 潮 项 以 及 不 规 则 的 短 周 期 成 分. 图 1 从 上 至 下 依 次 绘 出 了 1990 2010 年 期 间 LOD 变 化 的 原 始 序 列 带 谐 潮 序 列 线 性 模 型 拟 合 序 列 以 及 残 差 序 列. 本 文 采 用 GP 回 归 模 型 对 残 差 序 列 进 行 建 模 和 预 测, 将 线 性 模 型 的 预 报 结 果 和 残 差 序 列 的 预 报 结 果 相 加 即 可 获 得 最 终 的 LOD 变 化 预 报 值. 2.2.2 建 模 和 预 报 GP 模 型 的 建 模 过 程 就 是 通 过 对 样 本 数 据 D 的 训 练, 确 定 协 方 差 函 数 的 超 参 数. 超 参 数 的 选 取 方 法 主 要 有 交 叉 检 验 法 贝 叶 斯 推 理 法 和 最 大 似 然 法 [13]. 本 文 采 用 最 大 似 然 法 选 取 超 参 数, 即 任 意 给 定 超 参 数 的 初 值, 采 用 共 轭 梯 度 优 化 算 法 求 取 训 练 样 本 对 数 似 然 函 数 的 最 大 值, 从 而 得 到 似 然 函 数 最 大 值 所 对 应 的 超 参 数 即 为 最 优 超 参 数. 对 数 似 然 函 数 的 形 式 为 ln p(y x, θ) = 1 2 yt K 1 y 1 2 ln K n ln(2π). (4) 2
56 天 文 学 报 56 卷 Raw series/ ms Tidal term/ ms Fitting of linear model/ ms Residual/ ms 5 0 5 1 0 1 5 0 5 1 0 (a) (b) (c) (d) 1 1990 1992 1994 1996 1998 2000 2002 2004 2006 2008 2010 year 图 1 LOD 变 化 的 原 始 序 列 (a); 带 谐 潮 项 (b); 线 性 模 型 拟 合 项 (c); 残 差 项 (d) Fig. 1 The raw series (a); tidal terms (b); fitting terms of linear model (c); and residual terms (d) of the LOD variations 除 了 协 方 差 函 数 及 其 超 参 数, 样 本 的 输 入 和 输 出 方 式 也 非 常 重 要. 按 以 下 方 式 构 建 样 本 的 输 入 和 输 出 : 在 训 练 阶 段, 样 本 的 输 入 和 输 出 方 式 为 {ξ(t d), ξ(t d + 1),, ξ(t 1)} ξ(t), t = d + 1, d + 2,, n, (5) 其 中 {ξ(i), i = 1, 2,, n}, 表 示 LOD 变 化 序 列 经 数 据 预 处 理 后 的 残 差 序 列, 根 据 经 验 确 定, 本 文 取 d = 4. 在 预 报 阶 段, 预 报 跨 度 为 k = 1, 2,, d, d + 1, 时 样 本 的 输 入 和 输 出 方 式 分 别 为 {ξ(n d + 1), ξ(n d + 2),, ξ(n)} ξ(n + 1) k = 1 {ξ(n d + 2), ξ(n d + 3),, ξ(n + 1)} ξ(n + 2) k = 2 {ξ(n), ξ(n + 1),, ξ(n + d 1)} ξ(n + d). k = d { ξ(n + 1), ξ(n + 2),, ξ(n + d)} ξ(n + d + 1) k = d + 1., (6) 其 中 ξ(n + k) 表 示 残 差 序 列 为 k 的 预 报 值. 可 以 看 出 : 当 k 2 时, 每 一 次 预 报 都 需 要 前 几 步 的 预 报 值 作 为 输 入, 所 以 这 种 样 本 输 入 和 输 出 模 式 实 质 上 是 一 种 递 推 迭 代 预 报 方 式, 其 优 点 是 对 于 连 续 的 多 天 预 报, 只 需 进 行 一 次 模 型 计 算, 而 不 像 跨 度 模 式 等 需 要 进 行 多 次 模 型 计 算 [20 21], 预 报 效 率 大 大 提 高. 此 外, 这 种 模 式 的 数 据 利 用 率 高, 对 于 长 度 为 n
1 期 雷 雨 等 : 基 于 高 斯 过 程 的 日 长 变 化 预 报 57 的 基 础 序 列, 可 构 建 n d 个 训 练 样 本. 2.2.3 精 度 评 定 指 标 采 用 均 方 根 误 差 (root mean square error, RMSE) 和 平 均 绝 对 误 差 (mean absolute error, MAE) 作 为 预 报 结 果 的 精 度 评 定 指 标, 其 计 算 公 式 分 别 为 RMSE i = 1 N (P i j Oj i N )2, (7) j=1 MAE i = 1 N N P i j Oj i. (8) j=1 其 中 i 为 预 报 跨 度, N 为 预 报 期 数, P i j O i j 分 别 表 示 第 j 期 的 第 i d LOD 变 化 的 预 测 值 和 实 际 值. 3 实 验 分 析 首 先 将 1990 年 1 月 1 日 至 1999 年 12 月 31 日 的 LOD 变 化 残 差 序 列 用 于 GP 模 型 的 训 练, 然 后 用 训 练 好 的 GP 回 归 模 型 对 2000 2001 年 的 LOD 变 化 残 差 序 列 进 行 1 10 d 15 d 20 d 25 d 30 d 60 d 90 d 120 d 150 d 360 d 跨 度 的 预 报 ( 与 Schuh 等 [3] 和 张 晓 红 等 [8 9] 的 预 报 时 间 段 相 同 ). 图 2 给 出 了 基 于 GP 方 法 的 预 报 跨 度 为 1 d 的 LOD 变 化 残 差 的 预 报 曲 线 (a) 和 预 报 误 差 (b) 图, 图 2(a) 中 虚 线 和 实 线 分 别 代 表 残 差 预 报 值 和 观 测 值. 0.6 (a) Residual/ ms Prediction error/ ms 0.3 0.0 0.3 Prediction Observation 0.6 2000.1 2000.2 2000.3 2000.4 2000.5 2000.6 2000.7 2000.8 2000.9 year 0.05 0.05 (b) 2000.1 2000.2 2000.3 2000.4 2000.5 2000.6 2000.7 2000.8 2000.9 year 图 2 跨 度 为 1 d 的 LOD 变 化 残 差 的 预 报 结 果 (a) 和 预 报 误 差 (b) Fig. 2 The prediction results of the residual of LOD variations (a) and the predicted errors (b) at the prediction horizon of 1 d
58 天 文 学 报 56 卷 同 时 本 文 将 基 于 GP 方 法 的 LOD 变 化 预 报 结 果 同 Schuh 等 [3] 使 用 的 反 向 传 播 神 经 网 络 (back propagation neural networks, BPNN) 和 张 晓 红 等 [8 9] 使 用 的 改 进 的 BPNN 以 及 广 义 回 归 神 经 网 络 (general regression neural networks, GRNN) 的 预 报 结 果 进 行 对 比, 预 报 RMSE 值 见 表 1. 从 图 2(a) (b) 和 表 1 可 以 看 出, GP 用 于 LOD 变 化 预 报 是 可 行 和 有 效 的, 随 着 预 报 跨 度 的 增 大, 预 报 精 度 有 所 降 低. [3] 表 1 GP 预 报 结 果 与 BPNN 预 报 结 果 (Schuh 等 的 预 报 结 果 ) 改 进 的 BPNN 及 GRNN 预 报 结 果 [8 9] ( 张 晓 红 等 的 预 报 结 果 ) 的 比 较 ( 单 位 : ms) Table 1 The comparison of the prediction results of GP with those of the BPNN (Schuh et al. [3] ), the modified BPNN, and the GRNN (Zhang et al. [8 9] ) (unit: ms) Prediction day GP BPNN [3] Modified BPNN [8 9] GRNN [8 9] 1 0.027 0.019 0.027 0.037 2 0.057 0.049 0.073 0.074 3 0.078 0.074 0.093 0.097 4 0.097 0.097 0.110 0.117 5 0.112 0.121 0.131 0.134 6 0.126 0.142 0.148 0.151 7 0.138 0.159 0.162 0.164 8 0.149 0.174 0.170 0.174 9 0.158 0.184 0.176 0.179 10 0.166 0.193 0.185 0.187 15 0.190 0.246 0.211 0.204 20 0.201 0.251 0.217 0.210 25 0.207 0.249 0.215 0.211 30 0.210 0.245 0.219 0.217 60 0.232 0.292 0.219 0.222 90 0.266 0.306 0.231 0.226 120 0.260 0.314 0.229 0.226 150 0.252 0.330 0.237 0.233 180 0.257 0.361 0.234 0.234 210 0.268 0.397 0.241 0.236 240 0.257 0.377 0.236 0.236 270 0.251 0.386 0.231 0.240 300 0.222 0.402 0.249 0.247 330 0.280 0.372 0.262 0.254 360 0.273 0.347 0.245 0.250
1 期 雷 雨 等 : 基 于 高 斯 过 程 的 日 长 变 化 预 报 59 为 了 更 加 直 观 地 比 较 4 种 方 法 的 预 报 精 度, 图 3 绘 出 了 不 同 跨 度 的 预 报 精 度. 从 图 中 可 以 看 出, 在 短 期 (1 30 d) 预 报 中, 除 了 当 跨 度 为 1 3 d 时 GP 的 预 报 精 度 略 低 于 BPNN 预 报 精 度 外, 其 它 跨 度 的 预 报 精 度 均 高 于 另 外 3 种 方 法 的 预 报 精 度. 对 于 中 期 (1 360 d) 预 报, GP 的 预 报 精 度 仍 优 于 BPNN 的 预 报 精 度, 但 低 于 改 进 BPNN 和 GRNN 的 预 报 精 度. 在 预 报 效 率 上, 由 于 GP 模 型 参 数 ( 协 方 差 参 数 ) 可 以 自 适 应 获 取, 而 不 必 像 ANN 技 术 需 要 对 训 练 样 本 反 复 训 练 才 能 得 到 最 优 网 络 参 数, 故 训 练 速 度 较 ANN 技 术 要 快. 此 外, 因 为 本 文 采 用 的 GP 预 报 模 式 只 需 建 模 一 次 便 可 实 现 LOD 变 化 的 多 天 连 续 预 报, 因 此 预 报 所 用 时 间 较 少, 对 于 1 360 d 的 连 续 预 报, 训 练 时 间 和 预 报 时 间 之 和 一 般 在 20 min 以 内, 保 证 了 算 法 的 实 时 性. 0.40 0.35 0.30 GRNN BPNN Modified BPNN GP (a) 0.50 0.40 GRNN BPNN Modified BPNN GP (b) RMSE/ ms 0.25 0.20 0.15 RMSE/ ms 0.30 0.20 0.05 0 5 10 15 20 25 30 day in future 0 60 120 180 240 300 360 day in future 图 3 GP 预 报 精 度 与 BPNN 改 进 BPNN 和 GRNN 预 报 精 度 对 比. (a) 短 期 (1 30 d) 预 报, (b) 中 期 (1 360 d) Fig. 3 预 报 The comparison of the prediction accuracies of GP with those of the BPNN, the modified BPNN, and the GRNN. (a) The short-term (1 30 d) prediction, and (b) the medium-term (1 360 d) prediction 为 了 与 EOP PCC 的 预 报 结 果 进 行 比 较, 选 取 了 1990 年 1 月 1 日 至 2005 年 9 月 30 日 的 LOD 变 化 数 据 作 为 基 础 序 列, 预 报 2005 年 10 月 1 日 到 2008 年 2 月 28 日 ( 与 EOP PCC 预 报 时 间 段 相 同 ) 1 360 d 跨 度 的 LOD 变 化 值, 统 计 了 预 报 结 果 的 MAE, 并 与 EOP PCC 的 结 果 进 行 了 对 比, 对 比 结 果 见 图 4 6. 图 4 6 中 不 同 颜 色 和 不 同 形 状 线 条 分 别 代 表 参 与 EOP PCC 的 不 同 团 队 所 得 的 预 报 误 差, 参 与 此 项 竞 赛 的 团 队 详 细 情 况 参 见 文 献 [14]. 在 图 4 5 中 红 色 实 线 代 表 Gross 团 队 的 预 报 误 差, 粉 色 实 线 代 表 Kalarus 团 队 的 预 报 误 差, 蓝 色 虚 线 代 表 Akyilmaz 团 队 的 预 报 误 差, 蓝 色 点 划 线 代 表 Kosek 团 队 的 预 报 误 差, 绿 色 实 线 绿 色 虚 线 和 绿 色 点 划 线 代 表 Zotov 团 队 的 预 报 误 差, 黄 色 实 线 代 表 Pasynok 团 队 的 预 报 误 差, 蓝 色 实 线 代 表 Mendes Cerveira 团 队 的 预 报 误 差, 黑 色 实 线 代 表 本 文 预 报 误 差 ; 在 图 6 中 黑 色 方 形 线 条 代 表 Mendes Cerveira 团 队 的 预 报 误 差, 黑 色 三 角 形 线 条 代 表 Kosek 团 队 的 预 报 误 差, 黑 色 五 角 星 线 条 代 表 Gross 等 的 预 报 误 差, 黑 色 圆 形 线 条 代 表 本 文 预 报 误 差, 其 中 超
60 天 文 学 报 56 卷 短 期 (1 10 d) 预 报 精 度 较 高 的 团 队 是 Gross 和 Kalarus 团 队, 短 期 (1 30 d) 预 报 精 度 较 高 的 团 队 是 Gross Kalarus 和 Kosek 团 队, 而 参 与 中 期 (1 360 d) 预 报 竞 赛 的 只 有 Gross Kosek 以 及 Mendes Cerveira 3 个 团 队. 从 图 4 6 的 比 较 中 可 以 看 出, 对 于 1 4 d 的 预 报, GP 方 法 的 预 报 精 度 低 于 排 在 第 1 位 的 Gross 等 和 排 在 第 2 位 的 Kalarus 等 的 预 报 精 度, 从 第 5 d 开 始, GP 的 预 报 精 度 优 于 Kalarus 等 的 预 报 精 度, 但 仍 低 于 排 在 第 1 位 的 Gross 等 的 预 报 精 度 ; 对 于 短 期 (1 30 d) 预 报, GP 的 预 报 精 度 仅 次 于 排 在 第 1 位 的 Gross 等 的 预 报 精 度, 与 并 列 排 在 第 2 位 的 Kalarus 等 和 Kosek 等 的 预 报 精 度 大 致 相 当 ; 对 于 中 期 预 报, GP 的 预 报 效 果 则 不 如 EOP PCC. 0.40 0.35 0.30 MAE/ ms 0.25 0.20 0.15 0.05 1 2 3 4 5 6 7 8 9 10 day in future 图 4 超 短 期 (1 10 d) MAE 对 比 Fig. 4 The comparison of the predicted MAE for the ultra short-term (1 10 d) MAE/ ms 0.40 0.35 0.30 0.25 0.20 0.15 0.05 0 5 10 15 20 25 30 day in future 图 5 短 期 (1 30 d) MAE 对 比 Fig. 5 The comparison of the predicted MAE for the short-term (1 30 d)
1 期 雷 雨 等 : 基 于 高 斯 过 程 的 日 长 变 化 预 报 61 0.50 0.40 MAE/ ms 0.30 0.20 0 60 120 180 240 300 360 day in future 图 6 中 期 (1 360 d) MAE 对 比 4 讨 论 与 总 结 Fig. 6 The comparison of the predicted MAE for the medium-term (1 360 d) 本 文 根 据 实 验 验 证 了 GP 用 于 LOD 变 化 的 可 行 性 和 有 效 性. 与 ANN 方 法 相 比, GP 方 法 较 容 易 实 现, 并 且 它 不 需 要 太 多 的 先 验 信 息, 只 需 事 先 选 择 适 当 的 协 方 差 函 数, 其 超 参 数 在 训 练 过 程 中 便 可 以 自 适 应 地 确 定, 从 而 可 以 避 免 预 报 的 人 为 主 观 性, 提 高 预 报 结 果 的 可 信 度. 通 过 实 例 发 现, GP 方 法 用 于 LOD 变 化 预 报 可 以 取 得 较 好 的 预 报 效 果. 通 过 与 ANN 预 报 结 果 以 及 与 EOP PCC 预 报 结 果 的 比 较 发 现, GP 方 法 的 短 期 (1 30 d) 预 报 精 度 较 高, 但 中 期 (1 360 d) 预 报 精 度 则 不 如 ANN 预 报 精 度 和 EOP PCC 预 报 精 度, 这 可 能 是 由 本 文 所 使 用 的 递 推 预 报 模 式 的 误 差 累 积 效 应 引 起 的, 对 此 可 以 尝 试 以 下 两 种 方 法 对 GP 中 期 预 报 精 度 进 行 改 进, 一 种 是 改 进 样 本 输 入 方 式, 如 采 用 连 续 输 入 方 式 或 者 跨 度 输 入 方 式 [20 21], 另 外 一 种 方 法 是 对 GP 模 型 进 行 在 线 训 练, 本 文 暂 不 对 此 进 行 讨 论. 在 预 报 效 率 方 面, 由 于 本 文 采 用 递 推 预 报 模 式 进 行 预 报, 故 只 需 一 次 模 型 计 算 便 可 以 实 现 多 步 预 报, 极 大 提 高 了 预 报 效 率. 对 于 跨 度 为 1 360 d 的 预 报, 应 用 GP 方 法 预 报 LOD 变 化 只 需 要 20 min 左 右 的 时 间, 而 应 用 ANN 则 需 要 数 小 时, 预 报 效 率 大 大 提 高, 这 对 于 ERP 的 实 时 快 速 预 报 具 有 重 要 的 现 实 意 义. 此 外, 基 于 GP 方 法 预 报 LOD 变 化 所 需 训 练 样 本 数 量 远 远 小 于 ANN 所 需 样 本 量, 这 在 历 史 数 据 较 少 的 情 况 下 进 行 LOD 变 化 的 预 报 具 有 更 高 的 现 实 意 义. 由 于 协 方 差 函 数 类 型 超 参 数 选 取 方 法 以 及 样 本 输 入 方 式 等 对 GP 方 法 的 预 报 效 果 都 有 一 定 的 影 响, 因 此, 如 何 从 上 述 角 度 来 优 化 GP 模 型 以 进 一 步 提 高 LOD 变 化 的 预 报 精 度 是 我 们 下 一 步 的 研 究 重 点. 对 此 我 们 将 另 行 文 讨 论. 致 谢 感 谢 IERS 提 供 的 LOD 变 化 资 料, 对 中 国 科 学 院 上 海 天 文 台 郑 大 伟 研 究 员 提 供 的 帮 助 表 示 由 衷 的 感 谢! [1] Gambis D, Luzum B. Metro, 2011, 48: 165 [2] McCarthy D D, Petit G. ITN, 2003, 13: 605 参 考 文 献
62 天 文 学 报 56 卷 [3] Schuh H, Ulrich M, Egger D, et al. JGeod, 2002, 76: 247 [4] Tomasz N, Kose W. JGeod, 2008, 82: 83 [5] 许 雪 晴, 周 永 宏. 飞 行 器 测 控 学 报, 2010, 29: 70 [6] Freedman A P, Steppe J A, Dickey J O, et al. JGR, 1994, 99: 6981 [7] Gross R S, Eubanks T M, Steppe J A, et al. JGeod, 1998, 72: 215 [8] 张 晓 红, 王 琪 洁, 朱 建 军, 等. 天 文 学 报, 2011, 52: 322 [9] Zhang X H, Wang Q J, Zhu J J, et al. ChA&A, 2012, 36: 86 [10] 王 琪 洁, 廖 德 春, 周 永 宏. 科 学 通 报, 2007, 52: 1728 [11] Wang Q J, Liao D C, Zhou Y H. ChSBu, 2008, 53: 969 [12] Akyilmaz O, Kutterer H. JGeod, 2004, 78: 82 [13] Kosek W, Kalarus M, Johnson T J. ArtSa, 2005, 40: 119 [14] Kalarus M, Schuh H, Kosek W, et al. JGeod, 2010, 84: 587 [15] Samanta B, Al-Balushi K R, Al-Araimi S A. JAGI, 2003, 16: 657 [16] Seeger M. IJNS, 2004, 14: 69 [17] 何 志 昆, 刘 光 斌, 赵 曦 晶, 等. 控 制 与 决 策, 2013, 28: 1121 [18] Brahim-Belhouari S, Bermak A. Computational Statistics & Data Analysis, 2004, 47: 705 [19] 刘 冬, 张 清 华. 测 绘 学 报, 2011, 40: 59 [20] 张 晓 红, 王 琪 洁, 朱 建 军, 等. 中 国 科 学 院 上 海 天 文 台 年 刊, 2011, 32: 147 [21] Akyilmaz O, Kutterer H, Shum C K, et al. Applied Soft Computing, 2011, 11: 837 The Prediction of Length-of-day Variations Based on Gaussian Processes LEI Yu 1,2,3 ZHAO Dan-ning 1,3 GAO Yu-ping 1,2 CAI Hong-bing 1,2 (1 National Time Service Center, Chinese Academy of Sciences, Xi an 710600) (2 Key Laboratory of Time and Frequency Primary Standards, Chinese Academy of Sciences, Xi an 710600) (3 University of Chinese Academy of Sciences, Beijing 100049) ABSTRACT Due to the complicated time-varying characteristics of the length-of-day (LOD) variations, the accuracies of traditional strategies for the prediction of the LOD variations such as the least squares extrapolation model, the time-series analysis model, and so on, have not met the requirements for real-time and high-precision applications. In this paper, a new machine learning algorithm the Gaussian process (GP) model is employed to forecast the LOD variations. Its prediction precisions are analyzed and compared with those of the back propagation neural networks (BPNN), general regression neural networks (GRNN) models, and the Earth Orientation Parameters Prediction Comparison Campaign (EOP PCC). The results demonstrate that the application of the GP model to the prediction of the LOD variations is efficient and feasible. Key words astrometry, time, methods: data analysis