发展“时空分解+机器学习”新方法预测太阳活动黑子数和综合磁图的变化

太阳活动是日地因果关系的源头，深刻影响着地球系统与人类社会。太阳大气的局部区域时不时会发生爆发活动，主要表现为太阳耀斑和日冕物质抛射，导致地球空间环境的强烈扰动。从本质上来说，太阳爆发活动是日冕中积累的磁能的爆发性释放过程，因此太阳磁场的拓扑性质和复杂程度是理解太阳爆发活动、预测太阳活动水平的关键。

    由于观测手段的限制，目前对于太阳磁场的认识在很大程度上局限于光球磁场。对于太阳光球磁场演化的研究由来已久，对于太阳磁活动的预测更是一直以来的热点话题。现有研究中不乏对于太阳黑子数和太阳耀斑等物理对象的预测，提供了预测未来太阳活动的多维视角。然而，作为太阳活动的决定性因素，太阳磁场的全球时空分布如何演化、又该如何预测，至今仍是个悬而未决的挑战性课题。

    针对这一问题，何建森研究团队与合作者基于威尔科克斯天文台（Wilcox Solar Observatory, Stanford/WSO）的观测数据，以太阳综合磁图为切入口，分析了综合磁图及其球谐展开系数中蕴含的潜在规律，发展了“空间分布球谐分解+时间序列经验模分解+长短时记忆神经网络”的综合方法，实现了对第25活动周（SC25）的太阳黑子数和太阳综合磁图的预测。研究成果作为封面文章发表在最近出版的《Science China: Earth Science》2024年第8期上[1]（图1）。

图1.《Science China: Earth Science》2024年第8期封面：我们正迎来第25个太阳活动周的峰年。封面展示了太阳黑子数在近五个活动周内的时间序列和在未来12年内的预测趋势, 同时展示了在每个太阳活动峰年和谷年期间的太阳综合磁图及其预测结果。封面背景为太阳极紫外辐射图像和由势场源表面模型得到的磁场分布, 以及与太阳活动相联系的地球极光。

【发现光球磁场七阶球谐系数演化的新规律】

光球径向磁场Br对应的球谐展开形式可表为：

其中 是l次m阶球谐函数，是相应的球谐系数，是球谐展开的截断阶数。

首先，作者提取了1阶到9阶球谐系数的小波频谱，并将其分为连续型（C型）和间歇型（I型）两种模式。最明显的C型模式存在于系数  、、、和中，它们具有与太阳磁周期相符的22年周期性。而I型模式则间歇型地出现在活动峰年附近，表明光球磁场在峰年具有复杂的短周期扰动。

在此基础上，作者进一步分析了太阳黑子数与m=0的各个球谐系数（对应于轴向磁多极矩）之间的相位关系（图2）。最为有趣的规律是系数几乎总是与太阳黑子数同时达到极值，作者猜测这一结果可能与太阳的子午环流有关。

图2. ISN与m=0的各个球谐系数的李萨如图。用于绘图的球谐系数时间序列已预先通过13个月的移动窗口进行平滑处理。不同的颜色代表不同的太阳活动周, 在一个周期内颜色逐渐从深到浅。

【利用长短时记忆机器学习方法预测太阳黑子数变化】

用于太阳黑子数预测的原始数据是13个月平滑后的月度国际太阳黑子数（International Sunspot Number, ISN），时间跨度从1749年7月到2022年11月。

    预测模型的主体采用长短时记忆（Long-Short Term Memory, LSTM）的架构，可进行选择性的记忆和遗忘，在训练中逐轮更新各个单元的内部参数，优化“输入-输出”对的映射关系。LSTM网络结构为3层，输入维度为356，隐藏层维度为356，输出维度为1。原始数据被划分为长度为356的历史值和下一个值的“输入-输出”对，即根据前30年左右的月度太阳黑子数预测下一个月的黑子数。

    “输入-输出”对按照0.6:0.2:0.2的比例划分为训练集、验证集和测试集。在模型的训练过程中，选取均方误差作为损失函数，选取Adam优化器进行优化，并设计早停策略以防止过拟合，提高模型的泛化能力。训练完毕的模型将执行递归预测，以获得未来的太阳黑子数序列（图3）。

图3. 基于LSTM网络进行ISN时间序列的建模、训练和预测过程的示意图。

    在第24活动周（SC24）上评估模型的适用性和预测误差。将原始数据截去最后11年，执行对SC24的递归预测（图4a和4b）。对比实测数据得出，模型对于ISN峰值强度预测的相对误差为13.3%，对于ISN峰值时间预测的绝对误差为8个月左右。

    对SC25的太阳黑子数进行预测（图4c和4d）。结果表明，SC25的ISN峰值将出现在2024年6月，误差估计为8个月，峰值强度可能达到166.9±22.6，活动周长度约为11年零2个月。因此，SC25将强于SC24，而略弱于SC23。

图4. 基于LSTM模型的ISN预测。(a)和(b)为SC 24的预测结果，(c)和(d)分别为SC 25的预测结果。黑线表示13个月平滑的月度ISN序列, 红线表示未来12年的预测值, 蓝线表示训练集上的模型预测, 橙线表示验证集上的模型预测，绿线表示测试集上的模型预测。(b)和(d)底部的数字表示太阳活动周的编号。

【提出“时空分解+长短时记忆”方法预测光球综合磁图变化】

由于球谐系数中包含跨越多尺度的复杂周期性，首先采用经验模态分解方法将1阶到5阶的球谐系数分解为不同频段的分量序列，再对每个分量序列执行与上述相似的机器学习预测方法，之后相加得到预测的球谐系数时间序列。最后，通过球谐重构方法得到截断到5阶的预测磁图。

    图5定性展示了不同活动周相位的磁图预测结果。可以看到，预测低阶磁图与太阳磁场的极性规律相吻合，且极性反转的时间与预测的黑子数峰值时间基本一致。

图5. 截断到5阶的综合磁图预测结果。(a) CR2271, 当前活动周的上升相; (b) CR2278, 临近下一次极性反转; (c) CR2286, 完成下一次极性反转; (d) CR2361, 下一个太阳活动谷年。卡林顿周的选择基于对SC25的ISN峰值的预测, 即峰值将于2024年6月前后到达(CR2284)。

    进一步定量比较低阶的预测磁图和观测磁图之间的差异（图6）。对于预测的第一个卡林顿周（CR2259），磁场的分布和强度都高度匹配，绝大部分区域的预测误差低于1.5G，预测磁场和观测磁场之间的相关系数可达0.92。对于预测的第三个卡林顿周（CR2261），虽然大尺度磁场区域的中心位置与观测结果大致相近, 但它们的覆盖范围有较大差异，相关系数下降到0.77。预测准确性的降低可能来源于递归预测中误差的累积效应。

图6. CR2259(左列)和CR2261(右列)的预测磁图与观测磁图的比较(a)和(e)是预测磁图; (b)和(f)是观测的WSO磁图; (c)和(g)是用预测磁图减去观测磁图得到残差磁图, 其中实线标记了残差为±1.5G的等值线, 虚线标记了残差为±0.5G的等值线; (d)和(h)是观测磁场与预测磁场的散点分布。

    综上所述，作者基于光球磁场的球谐展开，将机器学习、模态分解和级数重构的方法相结合，发展出一套预测太阳黑子数和综合磁图的综合方法，填补了太阳磁场全球分布经验预测的空白, 对未来的太阳观测计划具有潜在的参考价值。当然，目前的预测准确度仍然有很大的提升空间，有待更长期的观测支持和更先进的模型辅助。

作者信息：

卓睿（第一作者，北京大学博士生），何建森（通讯作者，北京大学），段叠（国防科技大学）、林荣（北京大学）、吴子祺（北京大学）、闫丽梅（中国科学院地质与地球物理研究所）、魏勇（中国科学院地质与地球物理研究所）。

致谢：

相关工作得到了国家自然科学基金项目、国家重点研发计划项目、民用航天技术预研项目和中国科学院青年创新促进会的支持。感谢汪景琇院士对本研究工作的撰文评价与宝贵建议 [2]。

文章链接：

http://engine.scichina.com/doi/10.1007/s11430-023-1354-4

参考文献：

[1] Zhuo, R., He, J., Duan, D. et al. Prediction of solar activities: Sunspot numbers and solar magnetic synoptic maps. Sci. China Earth Sci. (2024).

[2] Wang, J. Predict the spatiotemporal distribution of solar magnetic field and understand the evolutionary feature of solar activity. Sci. China Earth Sci. (2024).