2026, Vol. 56
2. 安徽建筑大学电子与信息工程学院，安徽 合肥 230601
海表温度(Sea surface temperature, SST)是研究海洋气候变化、海洋环流及其对全球气候影响的关键数据之一，其变化直接影响到全球气候系统和生态环境[1-2]。准确的海表温度数据有助于研究海洋-大气相互作用、气候变暖、厄尔尼诺和拉尼娜现象[3-5]等。因而，获取准确的海表数据对全球气候研究至关重要。SST数据的获取方法包括现场原位测量。近年来，随着全球Argo浮标观测网络的不断完善和广泛应用，海洋表层温度数据的采集在时空覆盖范围和观测精度方面都取得了显著提升。然而，基于浮标、科考船等传统观测手段获取的SST数据存在明显的局限性：一方面，这些离散点状数据在空间分布上呈现不均匀性，难以实现全球海洋的连续覆盖；另一方面，现场观测通常需要耗费大量的人力物力资源，观测成本居高不下。因此，如何将离散观测数据重构为连续的网格数据，且确保重构的网格数据具有高分辨率和高精确度成为当前研究中的一个重要问题。
为了实现从离散观测数据到网格数据的重构，许多研究者尝试使用多种技术。目前，海洋温度场的重构方法主要有动力学方法[6]、统计学法[7-8]、人工智能法[9-10]、变分法[11-13]。高志斌等[14]利用数据插值经验正交函数(Data interpolating empirical orthogonal functions, DINEOF)方法对MODIS卫星2003—2018年赤道东太平洋SST数据进行重构，能够填补一定的观测空缺；Yan等[15]采用DINEOF方法对东海缺失数据进行重构。然而，在数据密度特别低的区域，该方法往往难以满足精度要求，尤其是在非均匀分布的数据区域。范冬林等[16]使用多种机器学习算法如：多层感知机(Multi-layer perceptron, MLP)、随机森林回归(Random forest regression, RFR)和支持向量回归(Support vector regression, SVR)对海表温度进行反演。尽管机器学习算法在海表温度反演中展现了巨大潜力，但它们通常需要大量的训练数据，对参数选择敏感，且泛化能力(对新的或未见数据的适应性)可能不佳。Chao等[17]利用二维变分资料同化(2-Dimensional variational，2D-Var)重构海表温度，但未能充分考虑海洋物理约束的作用，如平流约束和海表面高度约束。在实际海洋环境中，洋流的输运作用导致温度沿流向及周边区域扩散[18]，而海表面高度(Sea surface height, SSH)反映了海洋温度的热胀冷缩效应[19]。因此，在温度重构过程中合理引入海洋的物理约束至关重要。
针对上述方法的不足，本文提出了一种利用变分同化分析温度距平的Argo海表温度重构方法，并在此基础上引入了多种约束条件以优化重构精度。该方法的核心思想如下：首先，在Argo浮标观测点提取温度数据，并进行距平处理，得到离散的温度距平数据。随后，采用基于多约束条件的二维变分(Multi-constraint 2-dimensional variational, MC2D-Var)同化方法生成连续距平场，并在同化代价函数中分别引入平流约束、SSH约束和平流加SSH联合约束三种方案(见表 1)，以增强温度场的物理一致性和细节表现。最终，通过将气候态SST叠加到距平场中，实现海表温度的重构。
|
表 1 三种重构方案及主要约束条件 Table 1 Three reconstruction schemes and the main constraints |
采用变分法，通过融合东向海水流速和北向海水流速(Eastward sea water velocity and northward sea water velocity, UV)、海表面高度等多源观测数据，实现了全球海表温度的高精度重构。重构方法的整体步骤(见图 1)如下：
|
图 1 海表温度重构及评估的总体流程图 Fig. 1 Overall flowchart of sea surface temperature reconstruction and evaluation |
(1) 对某月全球Argo数据集进行筛选、插值等一系列预处理，利用Akima插值将其统一到WOA18标准层深度，并提取深度层为0的海表温度数据。
(2) 随机筛选该月70%的Argo数据作为同化实验数据集xo(离散观测)，利用WOA18月均海表温度数据线性插值计算出实验数据集空间位置处的温度值作为历史温度数据 x，将实验数据集xo减去历史温度数据x得到离散距平值xanomaly。
(3) 利用海表流速UV数据构建平流约束，海表面高度数据构建SSH约束，通过融合多约束条件变分法分析得到连续的距平场。
(4) 将不同重构方案得到的连续距平场与WOA18数据进行叠加，实现海表温度重构。
(5) 将步骤(2)中剩余30%的Argo数据作为验证数据集，对不同方案重构的SST和最优插值海表温度(Optimum interpolation sea surface temperature, OISST)的业务化数据分别进行误差分析，并对比误差结果，以评估重构方案的有效性。
1.2 观测数据预处理距平(Anomaly)通常用于描述某一变量(如温度、降水等)相对于其长期平均状态的偏离程度。它反映了该变量相较于基准状态(如多年平均值)的变化情况。距平的计算方法是将原始观测数据与某个基准值(通常为多年平均观测值或气候均值)进行比较，得到偏差量，可表示为
| $ x_{\text {anomaly }}=x_{\mathrm{o}}-\bar{x}。$ | (1) |
式中：xanomaly为距平即偏离值；xo为Argo某一时间点(或空间位置)的实际观测值; x为该观测位置处的历史平均值，通常是多年的平均观测值或气候均值。WOA18数据是基于长期海洋观测数据统计计算的气候态数据集，本文采用WOA18温度数据，通过线性插值计算出Argo表层空间位置处的温度值，以此作为历史温度数据x。
1.3 基于多约束条件的二维变分同化传统二维变分数据同化是一种将观测数据与数值模型背景场结合的技术，广泛应用于气象学和海洋学等领域[20-22]，以改进数值预报精度，其代价函数为
| $ \begin{gathered} J(\boldsymbol{x})=\underbrace{\frac{1}{2}\left(\boldsymbol{x}-\boldsymbol{x}_b\right)^{\mathrm{T}} \boldsymbol{B}^{-1}\left(\boldsymbol{x}-\boldsymbol{x}_b\right)}_{J_1}+ \\ \underbrace{\frac{1}{2}\left(\boldsymbol{H} \boldsymbol{x}-\boldsymbol{y}_{\mathrm{o}}\right)^{\mathrm{T}} \boldsymbol{R}^{-1}\left(\boldsymbol{H} \boldsymbol{x}-\boldsymbol{y}_{\mathrm{o}}\right)。}_{J_2} \end{gathered} $ | (2) |
式中：x为待分析场数据向量；将所有观测数据重新分组为向量yo；将背景场数据重新分组为向量xb；H为观测算子，主要作用是将模式状态向量映射或插值到观测向量；R为观测误差协方差矩阵，为对角矩阵且对角线值为观测误差。J1项和J2项分别为背景约束项和观测约束项，其作用确保结果尽量符合观测数据同时保证同化结果不会偏离背景场太远。B为背景误差协方差矩阵，在高维情况下，直接构建背景误差协方差矩阵B并求其逆通常计算成本较高。采用Barth等[23]提出的范数离散化方法直接构建背景误差协方差逆矩阵B－1：
| $ \boldsymbol{B}^{-1}=\frac{1}{c}\left(\alpha_0 \boldsymbol{I}+\alpha_1 \boldsymbol{D}_1^{\mathrm{T}} \boldsymbol{D}_1+\alpha_2 \boldsymbol{D}_2^{\mathrm{T}} \boldsymbol{D}_2+\cdots+\alpha_m \boldsymbol{D}_m^{\mathrm{T}} \boldsymbol{D}_m\right) 。$ | (3) |
式中：c为归一化常数；I为单位矩阵；α为权重系数，取二项式系数；Dn是第n阶差分算子矩阵, 1≤n≤m。
| $ \begin{gathered} c=\frac{(4 \pi)^{\frac{n}{2}} \Gamma(m)|\boldsymbol{L}|}{\Gamma\left(m-\frac{n}{2}\right)} ; \\ \alpha_k=\binom{m}{k}=\frac{m!}{k!(m-k)!}, (k=0, 1, \cdots, m) ; \\ \boldsymbol{L}=\left(\begin{array}{ccc} \boldsymbol{L}_x & 0 & \\ 0 & \boldsymbol{L}_y & \\ & & \ddots \end{array}\right) 。\end{gathered} $ | (4) |
式中：n为空间维度；m为高阶导数的阶数；Γ(·)为Gamma函数；L为缩放相关长度。本文研究海表温度数据，只考虑二维空间，取n为2，m最高阶数取2即可，由于空间分辨率为0.25°×0.25°，Lx和Ly通过缩放均取1即可，则此时B－1为
| $ \begin{gathered} \boldsymbol{B}^{-1}=\frac{1}{4 \pi \boldsymbol{L}_x \boldsymbol{L}_y}\left(\boldsymbol{I}+2\left(\boldsymbol{L}_x^2 \boldsymbol{D}_{x x}+\boldsymbol{L}_y^2 \boldsymbol{D}_{y y}\right)+\right. \\ \left.\left(\boldsymbol{L}_x^2 \boldsymbol{D}_{x x}+\boldsymbol{L}_y^2 \boldsymbol{D}_{y y}\right)^2\right)，\end{gathered} $ | (5) |
式中：Dxx、Dyy为x轴(经向)和y轴(纬向)二阶导数通过中心差分法近似得到的离散形式。
| $ \begin{aligned} & \left(\boldsymbol{D}_{x x} \phi\right)_{i, j} \approx \frac{\phi_{i+1, j}-2 \phi_{i, j}+\phi_{i-1, j}}{\Delta x^2}, \\ & \left(\boldsymbol{D}_{y y} \phi\right)_{i, j} \approx \frac{\phi_{i, j+1}-2 \phi_{i, j}+\phi_{i, j-1}}{\Delta y^2} ; \end{aligned} $ | (6) |
| $ \boldsymbol{T}_N=\frac{1}{\Delta x^2}\left[\begin{array}{ccccc} -2 & 1 & 0 & \cdots & 0 \\ 1 & -2 & 1 & \cdots & 0 \\ \vdots & \ddots & \ddots & \ddots & \vdots \\ 0 & \cdots & 1 & -2 & 1 \\ 0 & \cdots & 0 & 1 & -2 \end{array}\right]；$ | (7) |
| $ \boldsymbol{D}_{x x}=\frac{1}{\Delta x^2}\left(\boldsymbol{I}_M \otimes \boldsymbol{T}_N\right), \boldsymbol{D}_{y y}=\frac{1}{\Delta y^2}\left(\boldsymbol{T}_M \otimes \boldsymbol{I}_N\right) 。$ | (8) |
多约束条件的二维变分同化方法(MC2D-Var)在式(2)的基础上可以引入任意数量的附加约束。这些约束由对称的正定矩阵Qi、矩阵Ci和向量zi定义，
| $ J_c(\boldsymbol{x})=\sum\limits_i\left(\boldsymbol{C}_i \boldsymbol{x}-\boldsymbol{z}_i\right)^{\mathrm{T}} \boldsymbol{Q}_i^{-1}\left(\boldsymbol{C}_i \boldsymbol{x}-\boldsymbol{z}_i\right) \text { 。} $ | (9) |
这些额外约束可被观测算子H、方差矩阵R以及观测数据组成的向量yo重新吸收组合：
| $ \left(\begin{array}{c} \boldsymbol{y}_0 \\ \boldsymbol{z}_1 \\ \boldsymbol{z}_2 \\ \vdots \end{array}\right) \rightarrow \boldsymbol{y}_{\circ},\left(\begin{array}{c} \boldsymbol{H} \\ \boldsymbol{C}_1 \\ \boldsymbol{C}_2 \\ \vdots \end{array}\right) \rightarrow \boldsymbol{H},\left(\begin{array}{cccc} \boldsymbol{R} & & & 0 \\ & \boldsymbol{Q}_1 & & \\ & & \boldsymbol{Q}_2 & \\ 0 & & & \ddots \end{array}\right) \rightarrow \boldsymbol{R}。$ | (10) |
然后，将重组后的观测算子H、方差矩阵R以及观测数据组成的向量yo重新计算，以使得最终得到的分析场更加符合约束条件，从而提升其物理一致性与精确度。
由于平流和SSH对SST有重要影响，在同化代价函数中分别引入平流约束、SSH约束及平流+SSH联合约束三种方案，分别如下：
方案1：式(2)引入SSH约束。首先构建SSH约束JSSH：
| $ J_{\mathrm{SSH}}(\boldsymbol{x})=\frac{1}{2}\left(\boldsymbol{y}_{\mathrm{SSH}}-\boldsymbol{H}_{\mathrm{SSH}} \boldsymbol{x}\right)^{\mathrm{T}} \boldsymbol{R}_{\mathrm{SSH}}^{-1}\left(\boldsymbol{y}_{\mathrm{SSH}}-\boldsymbol{H}_{\mathrm{SSH}} \boldsymbol{x}\right) \text { 。} $ | (11) |
式中：ySSH为海表面高度数据；HSSH为提取算子，用于将模型中的分析场状态映射到SSH观测空间；RSSH为SSH约束的对角协方差矩阵，对角元素表示SSH观测数据的误差。
将JSSH(x)加入式(2)中，可得代价函数表达式：
| $ \begin{array}{r} J(\boldsymbol{x})=\underbrace{\frac{1}{2}\left(\boldsymbol{x}-\boldsymbol{x}_b\right)^{\mathrm{T}} \boldsymbol{B}^{-1}\left(\boldsymbol{x}-\boldsymbol{x}_b\right)}_{J_1}+ \\ \underbrace{\frac{1}{2}\left(\boldsymbol{H} \boldsymbol{x}-\boldsymbol{y}_{\mathrm{o}}\right)^{\mathrm{T}} \boldsymbol{R}^{-1}\left(\boldsymbol{H} \boldsymbol{x}-\boldsymbol{y}_{\mathrm{o}}\right)}_{J_2}+ \\ \underbrace{\frac{1}{2}\left(\boldsymbol{y}_{\mathrm{SSH}}-\boldsymbol{H}_{\mathrm{SSH}} \boldsymbol{x}\right)^{\mathrm{T}} \boldsymbol{R}_{\mathrm{SSH}}^{-1}\left(\boldsymbol{y}_{\mathrm{SSH}}-\boldsymbol{H}_{\mathrm{SSH}} \boldsymbol{x}\right)}_{J_{\mathrm{SSH}}}。\end{array} $ | (12) |
方案2：式(2)引入平流约束项
| $ J_c(\boldsymbol{\varphi})=\int_D\left[\left(\boldsymbol{u} \frac{\partial \boldsymbol{\varphi}}{\partial \boldsymbol{x}}\right)^2+\left(\boldsymbol{v} \frac{\partial \boldsymbol{\varphi}}{\partial \boldsymbol{y}}\right)^2\right] \mathrm{d} D \text { 。} $ | (13) |
式中：u和v分别是速度场在东西方向和南北方向的分量；
| $ \begin{gathered} J_c(\boldsymbol{x})=\frac{1}{2}\left(\boldsymbol{y}_{\mathrm{U}}-\boldsymbol{H}_{\mathrm{U}} \boldsymbol{x}\right)^{\mathrm{T}} \boldsymbol{R}_{\mathrm{U}}^{-1}\left(\boldsymbol{y}_{\mathrm{U}}-\boldsymbol{H}_{\mathrm{U}} \boldsymbol{x}\right)+ \\ \frac{1}{2}\left(\boldsymbol{y}_{\mathrm{V}}-\boldsymbol{H}_{\mathrm{V}} \boldsymbol{x}\right)^{\mathrm{T}} \boldsymbol{R}_{\mathrm{V}}^{-1}\left(\boldsymbol{y}_{\mathrm{V}}-\boldsymbol{H}_{\mathrm{V}} \boldsymbol{x}\right) 。\end{gathered} $ | (14) |
式中：yU、yV分别为海面经向流速和纬向流速；HU和HV为提取算子，用于将模型中的分析场状态映射到流速观测空间；RU和RV为平流约束的协方差矩阵，对角元素表示流速观测数据的误差。
将式(14)代入到代价函数中，此时式(2)可以表示为
| $ \begin{gathered} J(\boldsymbol{x})=\underbrace{\frac{1}{2}\left(\boldsymbol{x}-\boldsymbol{x}_b\right)^{\mathrm{T}} \boldsymbol{B}^{-1}\left(\boldsymbol{x}-\boldsymbol{x}_b\right)}_{J_1}+ \\ \underbrace{\frac{1}{2}\left(\boldsymbol{H} \boldsymbol{x}-\boldsymbol{y}_{\mathrm{o}}\right)^{\mathrm{T}} \boldsymbol{R}^{-1}\left(\boldsymbol{H} \boldsymbol{x}-\boldsymbol{y}_{\mathrm{o}}\right)}_{J_2}+ \\ \underbrace{\frac{1}{2}\left(\boldsymbol{y}_{\mathrm{U}}-\boldsymbol{H}_{\mathrm{U}} \boldsymbol{x}\right)^{\mathrm{T}} \boldsymbol{R}_{\mathrm{U}}^{-1}\left(\boldsymbol{y}_{\mathrm{U}}-\boldsymbol{H}_{\mathrm{U}} \boldsymbol{x}\right)+\frac{1}{2}\left(\boldsymbol{y}_{\mathrm{V}}-\boldsymbol{H}_{\mathrm{V}} \boldsymbol{x}\right)^{\mathrm{T}} \boldsymbol{R}_{\mathrm{V}}^{-1}\left(\boldsymbol{y}_{\mathrm{V}}-\boldsymbol{H}_{\mathrm{V}} \boldsymbol{x}\right)。}_{J_{\mathrm{uv}}} \end{gathered} $ | (15) |
方案3：在代价函数式(2)中综合添加平流约束和SSH约束，更新后的代价函数如下：
| $ \begin{gathered} J(\boldsymbol{x})=\underbrace{\frac{1}{2}\left(\boldsymbol{x}-\boldsymbol{x}_b\right)^{\mathrm{T}} \boldsymbol{B}^{-1}\left(\boldsymbol{x}-\boldsymbol{x}_b\right)}_{J_1}+ \\ \underbrace{\frac{1}{2}\left(\boldsymbol{H} \boldsymbol{x}-\boldsymbol{y}_{\mathrm{o}}\right)^{\mathrm{T}} \boldsymbol{R}^{-1}\left(\boldsymbol{H} \boldsymbol{x}-\boldsymbol{y}_{\mathrm{o}}\right)}_{J_2}+ \\ \underbrace{\frac{1}{2}\left(\boldsymbol{y}_{\mathrm{U}}-\boldsymbol{H}_{\mathrm{U}} \boldsymbol{x}\right)^{\mathrm{T}} \boldsymbol{R}_{\mathrm{U}}^{-1}\left(\boldsymbol{y}_{\mathrm{U}}-\boldsymbol{H}_{\mathrm{U}} \boldsymbol{x}\right)+\frac{1}{2}\left(\boldsymbol{y}_{\mathrm{V}}-\boldsymbol{H}_{\mathrm{V}} \boldsymbol{x}\right)^{\mathrm{T}} \boldsymbol{R}_{\mathrm{V}}^{-1}\left(\boldsymbol{y}_{\mathrm{V}}-\boldsymbol{H}_{\mathrm{V}} \boldsymbol{x}\right)}_{J_{\mathrm{UV}}}+ \\ \underbrace{\frac{1}{2}\left(\boldsymbol{y}_{\mathrm{SSH}}-\boldsymbol{H}_{\mathrm{SSH}} \boldsymbol{x}\right)^{\mathrm{T}} \boldsymbol{R}_{\mathrm{SSH}}^{-1}\left(\boldsymbol{y}_{\mathrm{SSH}}-\boldsymbol{H}_{\mathrm{SSH}} \boldsymbol{x}\right)。}_{J_{\mathrm{SSH}}} \end{gathered} $ | (16) |
平流约束确保了分析场与背景速度场之间的一致性，特别是确保海洋流场与温度、盐度等海洋状态变量的变化一致，避免了模型结果中的非物理性行为。而SSH约束则提供了海面变化的直接观测，反映了海洋表面和垂直结构的动力学特性，能够有效地引导模型的海表温度和其他变量向实际观测数据靠拢。
2 资料来源实验所采用的数据如表 2所示，具体数据说明如下。
|
|
表 2 实验数据资料来源 Table 2 Sources of experimental data |
(1) 观测数据。观测数据采用中国科学院海洋数据中心全球海洋科学数据库(CAS-Ocean Data Center-Global Ocean Science Database，CODC-GOSD)中的Argo浮标数据。Argo数据是CODC-GOSD的重要组成部分，提供了全球海洋的温度和盐度垂直廓线观测。数据库中的Argo数据经过严格的质量控制和偏差订正[24]，确保了数据的高精度和可靠性。本文分析的Argo数据来自2023年1—6月，涵盖了全球海洋多个区域的观测记录。每月进行重构时，随机筛选70%的当月观测数据进行同化，剩余30%作为验证集用于评估重构精度。
(2) 先验背景场数据。先验背景场数据采用WOA18气候态数据。该数据集是基于长期海洋观测数据统计计算得到的气候态数据集。它整合了全球海洋历史观测数据(如Argo浮标、CTD、XBT等)，通过客观分析方法计算月、季节和年尺度的温、盐等海洋环境要素的气候态平均值(Climatology)，代表了多年(通常为30年或以上)的平均状态。
(3) 约束条件数据。本文提出的三种约束方案的约束条件数据，包括SSH、经向流速(V)和纬向流速(U)均采用哥白尼数据中心提供的卫星观测数据，该数据空间分辨率为0.125°×0.125° (https://data.marine.copernicus.eu/product/SEALEVEL_GLO_PHY_L4_MY_008_047/description)，时间分辨率为逐日，将每月逐日数据进行平均统计获得逐月数据。本文利用线性插值将约束条件数据全部统一到0.25°×0.25°的分辨率网格中。
(4) 对比数据。为了评估不同重构方案得到的SST数据的精度，本文采用NOAA 0.25°逐日最优插值海表温度(OISST)数据[25]作为对比基准。具体方法是分别利用上述两类数据与当月30%的Argo观测数据计算误差，并进行对比分析，以评估各重构方案相对于OISST数据的精度表现。OISST数据为逐日数据，因此需要将每月逐日数据进行平均统计获得逐月数据。
3 实验及结果分析 3.1 2023年6月全球SST重构实验及其结果分析利用不同的重构方案对2023年6月全球海表温度进行重构，并分析不同方法的优劣。实验过程如下：首先，选取6月的Argo数据作为观测数据并随机筛选70%作为实验数据集，30%作为验证数据集。WOA18月均温度数据作为背景场数据，对实验数据进行距平分析，计算得到每个观测点的距平值。随后，采用本文提出的三种重构方法，分别对这些离散距平值进行重构，得到具有0.25°×0.25°分辨率的连续距平场。最后，将重构得到的距平场与WOA18数据进行叠加，最终生成全球0.25°×0.25°的海表温度数据。
2023年6月的Argo浮标观测数据共计约15 850个，随机筛选70%(约11 095个)作为实验数据集(见图 2)。剩余30%的现场观测数据(约4 755个)作为真值对重构结果进行精度验证，数据集如图 3所示。经观测数据预处理后得到的离散距平值如图 4所示。
|
图 2 2023年6月全球70%的Argo观测浮标站点分布图 Fig. 2 The distribution map of the 70% of Argo observation buoy stations globally in June, 2023 |
|
图 3 2023年6月全球剩余30%的Argo观测浮标站点分布图 Fig. 3 The distribution map of the remaining 30% of Argo observation buoy stations globally in June, 2023 |
|
图 4 2023年6月全球70%的Argo温度距平值 Fig. 4 The global Argo temperature anomaly values for 70% of the Argo stations in June, 2023 |
将不加任何约束条件得到的距平场作为对照组(见图 5(a))。将离散距平值分别采用三种约束方案进行融合多约束条件的变分分析，得到连续的距平场，如图 5(b)—(d)所示。
|
( 图(a)为传统2D-Var得到的温度距平场，图(b)为引入SSH约束得到的温度距平场，图(c)为引入平流约束项得到的温度距平场，图(d)为综合SSH和平流约束得到的温度距平场。Figure (a) shows the temperature anomaly field obtained by the traditional 2D-Var method, figure (b) shows the temperature anomaly field obtained with the inclusion of SSH constraints, figure (c) shows the temperature anomaly field obtained with the inclusion of advection constraints, and figure (d) shows the temperature anomaly field obtained with the combined SSH and advection constraints. ) 图 5 2023年6月不同方案重构出的温度距平场 Fig. 5 Temperature anomaly fields reconstructed by different schemes in June, 2023 |
由图 5可知，传统2D-Var方法未引入额外物理约束，重构结果主要依赖观测数据的分布，这导致如南海等观测稀疏区域缺少异常值，且异常分布较为集中，空间扩散能力有限。加入SSH同化后，异常值可在更广泛区域中出现，尤其在南海和西北太平洋。引入平流约束后，异常值虽不如SSH方案广泛，但更贴近实际海洋动力特征，如赤道区域的热量输运趋势更明显且在70°W—50°W、30°N—50°N海域内异常值减小。在同时加入SSH和平流约束的方案中，异常分布兼具空间广度与动力合理性，相较方案1，在180°—70°W、20°S—20°N赤道区域，温度分布更符合实际洋流对热量输运的影响，且相比于方法2，温度异常分布更广泛。
图 6(a)表示OISST数据，图 6(b)—(d)分别表示基于不同约束方案的变分重构方法所得的SST，其中SST由海表温度距平场(Sea surface temperature anomaly, SSTA)与WOA18海表数据叠加得到。
|
( 图(a)为OISST的温度场，图(b)为方案1重构的温度场，图(c)为方案2重构的海表温度，图(d)为方案3重构的海表温度数据。Figure (a) shows the OISST temperature field, figure (b) shows the SST reconstructed by method 1, figure (c) shows the SST reconstructed by method 2, and figure (d) shows the SST reconstructed by method 3. ) 图 6 不同方案重构后的全球海表温度 Fig. 6 Global sea surface temperature reconstructed by different schemes |
在OISST数据和重构数据中添加温度为0、5、25和30 ℃的等温线进行直观对比，与OISST数据(见图 6(a))相比，方案1和方案3重构的SST数据在60°E—120°E、15°S—15°N海域内呈现出更复杂的等温线结构，能够刻画更多中小尺度温度变化特征；而方案2的整体等温线结构整体较为平滑，复杂程度低于OISST数据，但是在小尺度海域内(如150°E—180°，15°S—15°N)，对温度细节的刻画优于OISST数据。
3.2 2023年1—6月实验及其结果分析进一步采用三种重构方案对2023年1—6月期间不同月份的全球海表温度进行重构实验。每月随机筛选剩余30%的Argo数据作为真值，并与不同重构结果以及当月OISST产品数据进行均方根误差(Root mean square error, RMSE)、平均绝对误差(Mean absolute error, MAE)和相关系数(r)计算，以进一步验证本文提出的不同重构方案的精度。
| $ R M S E=\sqrt{\frac{1}{N} \sum\limits_{i=1}^N\left(T_{\mathrm{Argo}, i}-T_{\mathrm{rec}, i}\right)^2} ; $ | (17) |
| $ M A E=\frac{1}{N} \sum\limits_{i=1}^N\left|T_{\mathrm{rec}, i}-T_{\mathrm{Argo}, i}\right| ; $ | (18) |
| $ r=\frac{\sum\limits_{i=1}^n\left(T_{\mathrm{rec}, i}-\overline{T_{\mathrm{rec}, i}}\right)\left(T_{\mathrm{Argo}, i}-\overline{T_{\mathrm{Argo}, i}}\right)}{\sqrt{\sum\limits_{i=1}^n\left(T_{\mathrm{rec}, i}-\overline{T_{\mathrm{rec}, i}}\right)^2 \sum\limits_{i=1}^n\left(T_{\mathrm{Argo}, i}-\overline{T_{\mathrm{Argo}, i}}\right)^2}} 。$ | (19) |
式中：TArgo, i为Argo数据，用来作为真值数据；Trec, i分别用重构SST数据、OISST数据以及WOA18温度数据来与真值数据计算相对应的误差；N是温度数据点的总数。三种方案重构的SST、OISST和WOA18的精度对比见表 3和图 7。
|
|
表 3 不同方案重构结果与WOA18海温数据以及OISST数据误差分析对比 Table 3 Comparison of error analysis between different reconstruction schemes, WOA18 sea surface temperature data, and OISST data |
|
图 7 三种方案重构的RMSE和MAE对比 Fig. 7 Comparison of RMSE and MAE for the reconstruction of three methods |
综合表 3与图 7可知，WOA18数据在各月份均表现出较大的误差，6月的RMSE高达1.386 ℃，MAE亦达1.012 ℃。相比之下，OISST数据的RMSE主要稳定在0.6~0.7 ℃之间，其中5月和6月的RMSE分别达到0.865和0.902 ℃，MAE则基本保持在0.4~0.5 ℃的范围内。本文提出的三种约束方案在RMSE和MAE方面均较OISST数据有所优化，其中方案3的重构结果最精确，1—6月的RMSE维持在0.5~0.7 ℃之间，6月的最大误差为0.864 ℃，仍低于OISST数据。同时，方案3的MAE也稳定在0.3~0.5 ℃之间，其中3月的重构精度最高，RMSE和MAE分别为0.566和0.369 ℃，均优于同月的OISST数据。此外，方案3在所有月份的相关系数最高，表明其重构的海表温度数据与真实值具有更强的一致性。
4 结语本文提出了一种基于二维变分融合多约束条件的SST重构方法，并设计了三种不同约束方案以进行对比分析。研究结果表明，方案3在各项误差指标上均表现出较高的重构精度。具体而言，方案3的均方根误差在1—6月间维持在0.5~0.7 ℃范围内，6月的最大误差(0.864 ℃)低于OISST数据(0.902 ℃)；其平均绝对误差亦稳定在0.3~0.5 ℃之间，并在3月达到最佳重构精度(RMSE为0.566 ℃，MAE为0.369 ℃)，优于OISST数据。此外，方案3在所有月份的相关系数均高于其他方案，表明其重构的SST数据与真实值具有更高的一致性。综合各项指标，方案3在重构精度和一致性方面优于其他两种方案，能够较好地融合Argo资料与卫星观测数据。本文为提高Argo网格化数据的精度提供了有效方法。
| [1] |
O'Carroll A G, Armstrong E M, Beggs H M, et al. Observatio-nal needs of sea surface temperature[J]. Frontiers in Marine Science, 2019(6): 00420. (  0) |
| [2] |
Melet A, Teatini P, Le Cozannet G, et al. Earth observations for monitoring marine coastal hazards and their drivers[J]. Surveys in Geophysics, 2020, 41(6): 1-46. ( 0) |
| [3] |
Wang X, Lu Q, Zhong S, et al. Differential impacts of the 2015-2020 El Niño/El Niño Modoki on seasonal ozone levels across China[J]. Atmospheric Pollution Research, 2025, 16(5): 102449. DOI:10.1016/j.apr.2025.102449 ( 0) |
| [4] |
Zhang M, Cheng Y, Zhang H, et al. Spatiotemporal variability of air-sea CO2 fluxes in response to El Niño-related marine heatwaves in the tropical Pacific Ocean[J]. Marine Environmental Research, 2025, 204: 106949. DOI:10.1016/j.marenvres.2025.106949 ( 0) |
| [5] |
Ruby S, Dibakar G. Analysis of pre-El Niño and La Niña events using climate network approach[J]. Chaos, Solitons and Fractals, 2025, 191: 115781. DOI:10.1016/j.chaos.2024.115781 ( 0) |
| [6] |
Mike C, Keith H. Altimetric assimilation with water property conservation[J]. Journal of Geophysical Research: Oceans, 1996, 101(1): 1059-1077. ( 0) |
| [7] |
Paul C. Fiedler, surface manifestations of subsurface thermal structure in the California Current[J]. Journal of Geophysical Research: Oceans, 1988, 93(5): 4975-4983. ( 0) |
| [8] |
王喜冬, 韩桂军, 李威, 等. 利用卫星观测海面信息反演三维温度场[J]. 热带海洋学报, 2011, 30(6): 10-17. Wang X D, Han G J, Li W, et al. Inversion of three-dimensional temperature fields using satellite-observed sea surface information[J]. Journal of Tropical Oceanography, 2011, 30(6): 10-17. ( 0) |
| [9] |
Ali M M, Swain D, Weller A R. Estimation of ocean subsurface thermal structure from surface parameters: A neural network approach[J]. Geophysical Research Letters, 2004, 31(20): L20308. ( 0) |
| [10] |
Wu X B, Yan X H, Young-Heon J, et al. Estimation of subsurface temperature anomaly in the North Atlantic using a self-organizing map neural network[J]. Journal of Atmospheric and Oceanic Technology, 2012, 29(11): 1675-1688. DOI:10.1175/JTECH-D-12-00013.1 ( 0) |
| [11] |
Fujii Y, Kamachi M. A reconstruction of observed profiles in the sea east of Japan using vertical coupled temperature-salinity EOF modes[J]. Journal of Oceanography, 2003, 59(2): 173-186. DOI:10.1023/A:1025539104750 ( 0) |
| [12] |
Helber R, Townsend L, Barron C, et al. Validation test report for the improved synthetic ocean profile (ISOP) system, part Ⅰ: Synthetic profile methods and algorithm[J]. Naval Research Laboratory, 2013, 7320-13-9364. ( 0) |
| [13] |
He Z K, Wang X D, Wu X R, et al. Projecting three-dimensional ocean thermohaline structure in the North Indian Ocean from the satellite sea surface data based on a variational method[J]. Journal of Geophysical Research: Oceans, 2021, 126(1): e2020JC016759. DOI:10.1029/2020JC016759 ( 0) |
| [14] |
高志斌, 王云涛, 柴扉. 基于海表温度遥感资料重构的热带不稳定波特征分析[J]. 厦门大学学报(自然科学版), 2024, 63(3): 512-525. Gao Z B, Wang Y T, Chai F. Analysis of tropical instability wave characteristics reconstructed from sea surface temperature remote sensing data[J]. Journal of Xiamen University (Natural Science), 2024, 63(3): 512-525. ( 0) |
| [15] |
Yan Z, Yingying G, Junxiao L. Research on Sea Surface temperature Reconstruction from long-term MODIS data[J]. IOP Conference Series: Earth and Environmental Science, 2021, 631(1): 012032. DOI:10.1088/1755-1315/631/1/012032 ( 0) |
| [16] |
范冬林, 杨鑫, 曾优, 等. Himawari-8卫星云下海表温度反演的机器学习方法比较研究[J]. 遥感信息, 2023, 38(5): 39-48. Fan D L, Yang X, Zeng Y, et al. Comparative study of machine learning methods for sea surface temperature retrieval under cloudy conditions using Himawari-8 satellite data[J]. Remote Sensing Information, 2023, 38(5): 39-48. ( 0) |
| [17] |
Chao Y, Li Z, Farrara J D, et al. Blending sea surface temperatures from multiple satellites and in situ observations for coastal oceans[J]. Journal of Atmospheric and Oceanic Technology, 2009, 26(7): 1415-1426. DOI:10.1175/2009JTECHO592.1 ( 0) |
| [18] |
Chen L, Zhang R H, Gao C. Effects of temperature and salinity on surface currents in the equatorial pacific[J]. Journal of Geophysical Research: Oceans, 2022, 127(4): e2021JC018175. DOI:10.1029/2021JC018175 ( 0) |
| [19] |
Marčelja S. The timescale and extent of thermal expansion of the global ocean due to climate change[J]. Ocean Science, 2010, 6(23): 179-184. ( 0) |
| [20] |
Dimet L F, 马怀存. 气象观测资料分析和同化的变分方法: 理论部分[J]. 气象科技, 1987(5): 14-22. Dimet L F, Ma H C. Variational methods for the analysis and assimilation of meteorological observations: Theoretical aspects[J]. Meteorological Science and Technology, 1987(5): 14-22. ( 0) |
| [21] |
Philippe L, Angela B, Peter B, et al. Experimental 2D-Var assimilation of ARM cloud and precipitation observations[J]. Quarterly Journal of the Royal Meteorological Society, 2006, 132(617): 1325-1347. DOI:10.1256/qj.05.24 ( 0) |
| [22] |
Lv S R, Lin W M, Wang Z X, et al. Blending sea surface winds from the HY-2 satellite scatterometers based on a 2D-Var method[J]. Remote Sensing, 2022, 15(1): 193. DOI:10.3390/rs15010193 ( 0) |
| [23] |
Barth A, Beckers J M, Troupin C, et al. DIVAnd-1.0:n-dimensional variational data analysis for ocean observations[J]. Geoscientific Model Development, 2014, 7(1): 225-241. DOI:10.5194/gmd-7-225-2014 ( 0) |
| [24] |
Tan Z T, Cheng L J, Gouretski V, et al. A new automatic quality control system for ocean profile observations and impact on ocean warming estimate[J]. Deep-Sea Research Part I, 2023, 194: 103961. DOI:10.1016/j.dsr.2022.103961 ( 0) |
| [25] |
Huang B, Liu C, Banzon V E, et al. Improvements of the daily optimum interpolation sea surface temperature (DOISST) version 2.1[J]. Journal of Climate, 2020, 34: 2923-2939. ( 0) |
2. School of Electronics & Information Engineering, Anhui Jianzhu University, Hefei 230601, China