一种基于相似度算法的地磁Kp指数预报方法及系统

摘要

本发明涉及空间物理和人工智能交叉学科，尤其涉及一种基于相似度算法的地磁Kp指数预报方法及系统。所述方法包括：接收当前太阳风事例，读取预先建立的历史太阳风事例库，通过相似度算法计算当前太阳风事例与历史太阳风事例的相似度，根据相似度排序得到相似度最高的若干个历史太阳风事例，由历史太阳风事例对应的Kp指数，得到当前太阳风事例对应的地磁Kp指数预报。本发明将推荐思想运用在地磁指数预报上，是一次全新的尝试；采用本发明的方法使得预报结果生动、友好、可以加入人工预报经验；并且采用本发明的方法使得预报输入参数简单、预报事件提前量领先，随着历史数据量不断积累，预报结果会更加精准。

说明书

技术领域

本发明涉及空间物理和人工智能交叉学科，尤其涉及一种基于相似度算法的地磁Kp指数预报方法及系统，并应用于空间天气预报的现实业务中。

背景技术

太阳风是太阳高层大气向外流动所形成的超声速等离子体流。当由日冕物质抛射事件或冕洞产生的太阳风到达地球时，会向磁层注入更多能量和粒子，可引起磁层扰动，引发地磁暴。地磁暴会严重影响运行在其中的卫星等飞行器的性能和安全，从而影响着整个空间和地面的技术系统。因此，对太阳风造成的地磁扰动进行及时的预报具有重要的实用价值。

当前，机器学习的技术不断突破，其自适应、自学习的强大并行信息处理能力，在很多方面收获了突破性进展。同时，空间环境监测也有了长期的积累，作为衡量近地空间全球磁扰强度的重要指标之一的地磁Kp指数，已有长达80年无间断的数据；太阳风方面也积累了几十年的卫星监测数据。结合这两方面的优势，可以让我们更有效的进行地磁扰动预报。

相似度计算是机器学习无监督聚类问题和推荐算法的核心思路，目前相似度算法种类繁多，且对于各类问题都有与之相应的、比较适合的算法。相似推荐的实际应用中，很多场景都取得了理想的结果。

现有的地磁Kp指数预报模型主要基于Kp指数与太阳风参数的统计关系、Kp指数与磁层内部的统计关系、Kp指数与能量耦合函数的关系和人工神经网络，还没有基于推荐的思想进行Kp指数预报的先例。且现有模型中，因为物理原理和物理过程的复杂，所以普遍所需参数种类多，准确率参差不齐，预报的时间提前量较短，又因为人工神经网络的“黑盒”过程，预报员只能得出单一的、生硬的预报的数值结果，无法在预报中加入人工预报经验。

发明内容

本发明的目的在于克服现有技术缺陷，提出了一种基于相似度算法的地磁Kp指数预报方法及系统。

为了实现上述目的，本发明的实施例1提出了一种基于相似度算法的地磁Kp指数预报方法，所述方法包括：

接收当前太阳风事例，读取预先建立的历史太阳风事例库，通过相似度算法，计算当前太阳风事例与历史太阳风事例的相似度，根据相似度排序，得到相似度最高的若干个历史太阳风事例，由历史太阳风事例对应的Kp指数，得到当前太阳风事例对应的地磁Kp指数预报。

作为上述方法的一种改进，所述方法还包括历史太阳风事例库的建立步骤；具体包括：

获取历史太阳风参数数据和地磁Kp指数；

对太阳风参数数据进行数据清洗及归一化处理，对地磁Kp指数进行临近数值转换；

采用XGBoost算法进行太阳风特征参数选取，筛选出对地磁Kp指数影响最大的两个特征参数；

采用单层感知机计算两个特征参数对地磁Kp指数影响的权重；

以固定步长的滑动窗口将历史太阳风参数数据划分为若干个太阳风事例，将时间上紧随其后的Kp指数定义为与该太阳风事例对应的Kp指数。

作为上述方法的一种改进，所述采用XGBoost算法进行太阳风特征参数选取，筛选出对地磁Kp指数影响最大的两个特征参数；具体包括：

以太阳风速度V、太阳风质子密度N、太阳风温度T、行星际磁场强度B、行星际磁场By、Bz分量各自1小时的平均值作为输入，以相应的Kp指数作为标签，通过XGBoost模型迭代建立决策树，对输入的六个参数进行了打分，取分数最高的两个特征参数。

作为上述方法的一种改进，所述采用单层感知机计算两个特征参数对地磁Kp指数影响的权重；具体包括：

以两个特征参数作为单层感知机的输入层参数，时间分辨率为1小时，输出层为每小时相应的Kp指数，以均方差函数为损失函数，以Adam算法为优化器，得到两个特征参数对地磁Kp指数影响的权重。

作为上述方法的一种改进，所述滑动窗口的大小按Kp指数的时间分辨率进行划分，滑动窗口步长在0.5至1小时范围内。

作为上述方法的一种改进，所述通过相似度算法，计算当前太阳风事例与历史太阳风事例的相似度；具体包括：使用动态时间规整算法或使用分段线性表示算法计算当前太阳风事例与历史太阳风事例的相似度。

作为上述方法的一种改进，所述使用动态时间规整算法计算当前太阳风事例与历史太阳风事例的相似度，具体包括：

通过对当前太阳风事例与每个历史太阳风事例的时间序列进行局部非线性缩放，使两个序列形态对齐后的欧氏距离累计值为两个时间序列的最小距离，从而得到当前太阳风事例与每个历史太阳风事例的相似度。

作为上述方法的一种改进，所述使用分段线性表示算法计算当前太阳风事例与历史太阳风事例的相似度，具体包括：

对当前太阳风事例与每个历史太阳风事例的时间序列采用分段线性表示算法，各得到若干个模式，然后进行等模式化，使得当前太阳风事例与历史太阳风事例的模式序列长度相等，再使用绝对值距离或欧式距离计算方法得到每个模式的距离值，对模式距离值进行加权计算得到当前太阳风事例与每个历史太阳风事例的相似度。

本发明的实施例2提出了一种基于相似度算法的地磁Kp指数预报系统，所述系统包括：历史太阳风事例库、太阳风事例接收模块、相似度计算模块和地磁Kp指数预报输出模块；其中，

所述太阳风事例接收模块，用于接收当前太阳风事例；

所述相似度计算模块，用于读取预先建立的历史太阳风事例库，通过相似度算法计算当前太阳风事例与历史太阳风事例的相似度；

所述地磁Kp指数预报输出模块，用于根据相似度排序得到相似度最高的若干个历史太阳风事例，由历史太阳风事例对应的Kp指数，得到当前太阳风事例对应的地磁Kp指数预报。

与现有技术相比，本发明的优势在于：

1、本发明将推荐思想运用在地磁指数预报上，是一次全新的尝试；

2、采用本发明的方法使得预报结果生动：不再是单一的地磁指数作为输出，而是几个相似的历史太阳风事例作为输出(包括历史太阳风各特征的时序变化图像以及其造成的地磁影响(后续地磁Kp指数))；

3、采用本发明的方法使得预报结果友好：对于新预报员，即使不了解以前典型的太阳风事件，也可以像“百度搜索”一样，给出推荐的太阳风结果；

4、采用本发明的方法使得预报结果可以加入人工预报经验：以往的模型受限于“黑盒”过程，只返回单一的结果，而现在的模型不光提供预测值，还可以给出历史的太阳风特征参数时序变化情况，与当前输入太阳风特征参数时序变化对比，预报员可以从图像中加入自己的预报经验对预测值做选择或调整；

5、采用本发明的方法使得预报输入参数简单：以往的模型输入参数种类多，而现在可以根据XGboost特征选择结果降低输入特征种类数；

6、采用本发明的方法使得预报事件提前量领先：历史事例已经发生，我们可以看见其后续足够长的影响，只要我们能够找到最相似的事例，其历史的后续影响就可以作为当前预报的参考；

7、采用本发明的方法使得随着历史数据量不断积累，预报结果会更加精准。

附图说明

图1是本发明的方法流程图；

图2(a)是两个时序构建的矩阵网格；

图2(b)是动态时间规整算法原理；

图3是分段线性表示中两个时序的等模式数化。

具体实施方式

本方法首次从推荐的角度出发，进行地磁Kp指数的预报。通过推荐获得的预报结果，不光有(历史相似太阳风事例相应的地磁Kp指数作为)Kp指数的预测值，还会输出历史相似太阳风事例的演化过程。这让预报结果更加生动形象，且预报员可以根据太阳风历史情况和现有情况的对比，在地磁指数预报中加入自身的预报经验。

下面结合附图和实施例对本发明的技术方案进行详细的说明。

实施例1

如图1所示，本发明的实施例1提供了一种基于相似度算法的地磁Kp指数预报方法。包括以下步骤：

接收当前太阳风事例，读取预先建立的历史太阳风事例库，通过相似度算法，计算当前太阳风事例与历史太阳风事例的相似度，根据相似度排序，得到相似度最高的若干个历史太阳风事例，由历史太阳风事例对应的Kp指数，得到当前太阳风事例对应的地磁Kp指数预报。

1数据

太阳风及行星际数据来源于OMNI网站。OMNI网站整合了来自IMP 8、Geotail、Wind和Ace卫星的太阳风数据。本模型选取的太阳风数据的时间范围是从1998年至2019年，为ACE卫星的监测数据，数据时频为1分钟。由于太阳风数据存在缺失，对于缺失数据(缺省值为-9999.9)，做了数据清洗，将个别的缺省数据取值为其前1分钟对应的数据，保证了数据的完整。OMNI网站提供的太阳风参数包括太阳风速度V、太阳风质子密度N、太阳风温度T、GSM坐标下行星际磁场By、Bz分量的值、行星际磁场强度B。为了便于训练和学习，对太阳风参数数据进行了归一化处理，消除不同量纲对计算的影响。

地磁Kp指数来源于OMNI网站整合自德国地球科学研究中心的结果。为了计算分析方便，将Kp指数转化为临近的数值，即Kp指数3

本模型将世界时整点时刻后的3个小时内的太阳风数据作为一个太阳风事例，将时间上紧随其后的Kp指数定义为与这一太阳风事例对应的Kp指数。太阳风事例时间上领先其对应的Kp指数1-3小时。

2 模型方法

2.1 基本原理

太阳风是造成地磁暴的直接原因，因此特征一致的太阳风引起的地磁扰动也应该相近。如果历史数据足够多，总可以在历史数据中

寻找到与当前太阳风相似的历史太阳风事例。历史太阳风事例引起的地磁变化，就可以作为当前Kp指数预报的参考。

本模型主要利用机器学习中的相似度算法来构建寻找相似太阳风的推荐。具体的步骤分为3步。第1步，特征参数筛选。太阳风特征参数主要有速度、密度、温度、行星际磁场等6个；在相似计算中，过多的特征参数会导致推荐结果混乱，需要选取对地磁扰动起主要作用的参数作为相似计算的特征参数。第2步，参数权重计算。在地磁扰动中，太阳风各种特征参数的影响作用不同，需要计算不同特征参数的权重值。第3步，推荐模型构建。以选定的太阳风特征参数作为输入，基于机器学中的相似度算法，寻找太阳风特征时间序列变化的相似，构建筛选相似太阳风事例的推荐模型。

2.2特征参数选取

本模型采用XGBoost算法进行太阳风特征参数选取。XGBoost是一个决策树训练模型，该模型依据输入的特征参数在迭代建立决策树中的作用和价值，给每个特征参数打分。

以归一化后的1998年至2019年太阳风速度V、太阳风质子密度N、太阳风温度T、行星际磁场强度B、行星际磁场By、Bz分量1小时的平均值作为输入，以相应的Kp指数作为标签，通过XGBoost模型迭代建立决策树，对输入的六个参数进行了打分。评分高低即为XGBoost模型认为不同参数对地磁扰动(Kp指数)影响的重要程度。为了避免过多特征参数造成推荐结果的混乱，以及增加模型运行的时间成本，选择对地磁扰动影响最重要的两个参数作为推荐模型的输入参数。

2.3参数权重计算

本模型采用单层感知机评估筛选出的特征参数对地磁扰动影响的权重。以特征参数作为单层感知机的输入层参数，时间分辨率是1h，输出层为每小时相应的Kp指数，损失函数为均方差函数，优化器是Adam(adaptive moment estimation)算法

2.4推荐模型构建

本模型以机器学习中的相似度算法构建太阳风相似推荐模型，主要采用了动态时间规整算法。

欧式距离，如公式1，是指m维空间中两个点之间的真实距离，或者向量的自然长度，在二维和三维空间中的欧氏距离就是两点之间的实际距离。当其输入为两组特征参量时，可以反映两组特征参量之间的综合差异。

动态时间规整(dynamic time warping，简称DTW算法)算法是流行的寻找时间序列相似的算法之一，其原理是给定两个序列Q＝{Q1,Q2,…,Qi,…,Qn}和C＝{C1,C2,…,Cj,…,Cm}，其长度分别是n和m，构造一个n*m的矩阵网格，如图2(a),以欧氏距离D为标准，矩阵元素(i,j)为Qi和Cj两个点的距离D(Qi,Cj)。从序列起始段所在的矩阵角为边界条件，在满足连续性和单调性约束的同时，通过动态规划求得距离累积值最小的路径即为最佳路径。

如图2(b)，动态时间规整算法认为，通过对时间序列进行局部非线性缩放，使两个序列形态尽可能对齐后的欧氏距离累计值才是两个时间序列的最小距离。

将推荐思想运用在地磁指数预报上，是一次全新的尝试，主要优点：(1)预报结果生动：不再是单一的地磁指数作为输出，而是几个相似的历史太阳风事例作为输出(包括历史太阳风各特征的时序变化图像以及其造成的地磁影响(后续地磁Kp指数))(2)预报结果友好：对于新预报员，即使不了解以前典型的太阳风事件，也可以像“百度搜索”一样，给出推荐的结果(3)预报结果可以加入人工预报经验：以往的模型受限于“黑盒”过程，只返回单一的结果，而现在的模型不光提供预测值，还可以给出历史的太阳风特征参数时序变化情况，与当前输入太阳风特征参数时序变化对比，预报员可以从图像中加入自己的预报经验对预测值做选择或调整(4)预报输入参数简单：以往的模型输入参数种类多，而现在可以根据XGboost特征选择结果降低输入特征种类数(5)预报事件提前量领先：历史事例已经发生，我们可以看见其后续足够长的影响，只要我们能够找到最相似的事例，其历史的后续影响就可以作为当前预报的参考(6)当历史数据量不断积累，一定可以找到与输入事例相似的历史事例，所以预报结果未来会越来越精准。

可替换点：除了使用动态时间规整算法，还可以使用PLR(分段线性表示)算法，寻找时间序列变化相似。也可以使用两种方法进行混合推荐。PLR：(piecewise linearrepresentation)分段线性表示。如图3所示。一个时间序列，无非有三种状态：上升、下降和不变，我们将这种状态对应表示为[1，-1，0]。假设有某个长度为S的序列，我们将其等分为K段。对于每一段计算一个斜率，斜率为正表示上升，为负表示下降，为0表示不变。那么我们就可以将这个序列表示为[1,1,0,-1...]这样的序列，将相邻的相同模式进行合并，我们得到[1,0,-1...]的序列。这就是PLR算法。因为我们将相邻的相同模式进行了合并，所以我们得到的模式序列一定是1,-1,0间隔排列的，每一个模式可能跨越了不同的时间长度，在两个时间序列合并模式后，序列S1可能有N个模式，S2可能有M个模式。现在我们要将他们等模式数化。通俗的说就是让他们使用共同的分割点，以获得最后长度相等的模式序列。接下来我们计算距离。这里我们使用绝对值距离，当然使用欧式距离也是可以的。模式距离公式为：D＝|m

太阳风是造成地磁暴的直接原因，因此特征一致的太阳风引起的地磁扰动也应该相近。如果历史数据足够多，总可以在历史数据中寻找到与当前太阳风相似的历史太阳风事例。历史太阳风事例引起的地磁变化，就可以作为当前Kp指数预报的参考。

通过计算得到历史上与当前输入最相似的几个太阳风，并以相似的历史结果作为推荐结果。预报结果不光包含地磁指数预测值，还包含推荐的历史太阳风与当前输入的太阳风的演化对比，让预报员可以更好的结合自身预报经验。

推荐的思想避免了目前尚未完全清楚的太阳风参数与地磁指数复杂的非线性关系式，而是直接从相似出发，结合预报员的人工经验，所以预报结果更友好、直观、可解释，且预报结果能够看到更多的信息，预报时间提前量和准确率也会比较领先。

本方法主要利用机器学习中的相似度算法来构建寻找相似太阳风的推荐模型。

具体的步骤为：

准备：

1)对历史数据进行清洗，保证数据完整；

2)将数据归一化，消除不同参数的量纲影响；

3)在多个可能造成地磁影响的太阳风参数中进行特征选择，筛选出最重要的两个；

4)对分钟精度的历史数据按滑动窗口划分成无数个太阳风事例(滑动窗口大小可以按Kp指数的时间分辨率——三小时划分，固定0.5或1小时的步长)；

推荐：

1)输入事例通过相似度算法与历史太阳风事例遍历计算(特征参数赋权重后)相似度，并按相似度排序(可用分布式计算，提升效率)

2)推荐(输出)

实施例2

本发明的实施例2提出了一种基于相似度算法的地磁Kp指数预报系统。所述系统包括：历史太阳风事例库、太阳风事例接收模块、相似度计算模块和地磁Kp指数预报输出模块；其中，

所述太阳风事例接收模块，用于接收当前太阳风事例；

所述相似度计算模块，用于读取预先建立的历史太阳风事例库，通过相似度算法计算当前太阳风事例与历史太阳风事例的相似度；

所述地磁Kp指数预报输出模块，用于根据相似度排序得到相似度最高的若干个历史太阳风事例，由历史太阳风事例对应的Kp指数，得到当前太阳风事例对应的地磁Kp指数预报。

最后所应说明的是，以上实施例仅用以说明本发明的技术方案而非限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行修改或者等同替换，都不脱离本发明技术方案的精神和范围，其均应涵盖在本发明的权利要求范围当中。