Marine Hydrography

Using artificial intelligence for identifying the depth of upper-ocean mixed layer

  • ZHANG Kang 1, 2 ,
  • GUO Shuangxi 1, 3 ,
  • HUANG Pengqi 1, 2 ,
  • QU Ling 1, 3 ,
  • LU Yuanzheng 1, 3 ,
  • CEN Xianrong 1, 3 ,
  • YU Lusha 1, 2 ,
  • ZHOU Weidong 1, 3 ,
  • ZHOU Shengqi , 1, 3
Expand
  • 1. State Key Laboratory of Tropical Oceanography (South China Sea Institute of Oceanology, Chinese Academy of Sciences), Guangzhou 510301, China
  • 2. University of Chinese Academy of Sciences, Beijing 100049, China
  • 3. Institution of South China Sea Ecology and Environmental Engineering, Chinese Academy of Sciences, Guangzhou 510301, China
ZHOU Shengqi. E-mail: .

Copy editor: YIN Bo

Received date: 2018-12-14

  Request revised date: 2019-04-10

  Online published: 2019-10-09

Supported by

National Natural Science Foundation of China(91752108)

National Natural Science Foundation of China(41476167)

National Natural Science Foundation of China(41706029)

National Natural Science Foundation of China(41606010)

Natural Science Foundation of Guangdong Province(2016A030311042)

Natural Science Foundation of Guangdong Province(2016A030310114)

Guangzhou Science and Technology Program Key Project(201804020056)

Strategic Priority Research Program of Chinese Academy of Sciences(XDA11030302)

Institution of South China Sea Ecology and Environmental Engineering, Chinese Academy of Sciences(ISEE2018PY05)

Copyright

Copyright reserved © 2019. Office of Acta Agronomica Sinica All articles published represent the opinions of the authors, and do not reflect the official policy of the Chinese Medical Association or the Editorial Board, unless this is clearly specified.

Abstract

An artificial intelligence (AI) method for identifying upper-ocean mixed layer depth (MLD) is proposed. In this method, a linear model, whose coefficient and variance are made into a set of statistics to characterize the profile, is established between temperature (density) and pressure (or depth). A subjective priori distribution is set in an initial window. The maximum posterior probability estimate of the mean coefficient value is obtained when the window is moved down from the sea surface, by learning the new data through the Bayesian chain rule and the minimum description length principle. The existence and depth of the mixed layer are determined when the jump of the coefficient is found by using the F-distribution. Using the Argo buoy data measured in the Pacific Ocean in February 2017, and taking the value of the quality index (QI) to estimate the accuracy of the MLD results, we find that this AI method is superior to the gradient method, the threshold method, the Hybrid method, and the relative-variant method.

Cite this article

ZHANG Kang , GUO Shuangxi , HUANG Pengqi , QU Ling , LU Yuanzheng , CEN Xianrong , YU Lusha , ZHOU Weidong , ZHOU Shengqi . Using artificial intelligence for identifying the depth of upper-ocean mixed layer[J]. Journal of Tropical Oceanography, 2019 , 38(5) : 32 -41 . DOI: 10.11978/2018137

海洋上混合层(ocean upper mixed layer)是位于大气底边界层和深层海洋之间的温度、盐度、密度等在垂向上分布近似均匀的水体。它通常形成于海洋上层的湍流混合,在海洋和大气之间的淡水、动量以及热交换过程中起着重要作用(Kantha et al, 1999; Wunsch et al, 2004), 关于它的研究对建立平均气候态以及探索全球海洋变化极其重要(Deser et al, 1996; Hanawa et al, 2001; Yu et al, 2007; Oka et al, 2015)。混合层深度(Mixed Layer Depth, MLD)主要受到海表风应力、输入海洋的净热通量、蒸发降水等因素的影响, 其变化调节着海表温度(Sea Surface Temperature, SST)以及混合层热含量的演变等与热带气旋的形成、浮游植物的丰富程度以及气候变化等现象密切相关的过程(Carton et al, 2008; Dong et al, 2008; Vanin et al, 2009; Qiu et al, 2014)。风应力的剪切作用搅动海水, 加剧混合层内的湍流混合过程, 使MLD变深; 输入海洋的热通量与蒸发降水会改变海表的浮力通量, 进而改变MLD, 正的净热通量输入到上层海洋中会削弱湍流效应, 使MLD变浅, 反之造成一个更深的MLD; 蒸发(降水)导致海洋表层盐度变高(低), 使海表浮性频率降低(增加), 加强(减弱)湍流混合过程, 进而导致MLD的深度加深(变浅) (范聪慧, 2007; Talley et al, 2011; Qiu et al, 2016)。因而准确确定MLD对于建立海气热交换机制、研究海洋表层物理、生物过程以及气候变化等有着重要意义。
为了确定混合层深度, 前人做了大量的研究工作, 提出了一系列方法(Lukas et al, 1991; Brainerd et al, 1995; Obata et al, 1996; Kara et al, 2000; Lavender et al, 2002; Thomson et al, 2003; De Boyer Montégut et al, 2004; Lorbacher et al, 2006; Holte et al, 2009; Chu et al, 2011; Holte et al, 2017; Huang et al, 2018)。在这些方法中, 阈值类方法由于具有简单易于操作的特点成为最常用的一类。它们主要包括基于温度(密度)变化的阈值法、梯度变化的梯度法, 以及集合两者的混合法(Holte et al, 2009, 2017)。关于阈值的确定, Kara和Montegut等人在比较多种阈值的取值后, 总结出全球的温度和密度的最优阈值分别为0.2℃和0.03kg·m-3 (Kara et al, 2000; De Boyer Montégut et al, 2004; Holte et al, 2009)。在南大洋混合层的研究中, 密度和温度梯度的最优阈值分别被确定为0.0005~0.05kg·m-4和0.025℃·m-1 (Dong, 2008)。但是选择一个固定值作为最优阈值是比较主观的, 从大量温度和密度数据分析结果中发现, 不同区域(季节)需要不同的阈值, 其在热带与副热带区域为0.06kg·m-2, 在副极带区域为0.09kg·m-2 (Ohno et al, 2009)。另外, 选择一个合适的最优阈值也是比较困难的, 对于具有渐变密度跃层的廓线, 过小的阈值会使识别结果较实际情况偏浅, 而较大的阈值又会造成偏深的识别结果(Lorbacher et al, 2006)。
为了改进阈值法的这些不足, 学者提出了一些相对客观的方法。Thomson等(2003)对可变长度的密度廓线片段做线性最优拟合, 提出一种split-merge (分割-合并)方法, 但是他们发现这个方法在确定混合层深度的结果上和阈值法相近。Lorbacher等(2006)基于廓线的二阶导数提出了曲率法, 用全球海洋环流实验(World Ocean Circulation Experiment, WOCE)数据测试发现, 相比较于梯度阈值法, 该方法不仅在结果上更加精确, 而且可以更好地应用于不同分辨率的数据。然而该方法被验证难以应用于存在噪声的廓线数据中(Chu et al, 2011)。对水下滑翔机(seaglider)数据的研究中, Chu 等通过对温度(密度)廓线进行线性拟合提出了最优线性拟合法(Chu et al, 2010), 并且在此基础上进行改进提出了最大角度法(Chu et al, 2011), 分别对收集自西北大西洋和墨西哥湾流海域的滑翔机数据分析发现, 最优线性拟合法和最大角度法比阈值类方法具有更高的准确性。黄鹏起等人发现, 最大角度法和最优线性拟合法在分析WOCE数据中准确性不高, 他们通过构造廓线标准差和极差的比值—相对变量, 提出了相对变化法(Huang et al, 2018), 并且对WOCE数据进行了分析, 结果表明该方法比阈值类方法、曲率法和最大角度法等有着更好的表现。然而该方法仍存在些许不足, 在真实的海洋温度(位势密度)廓线中, 相对变量取得极大值的位置常常位于混合层下边界以下的某一深度上, 从该位置向上搜索混合层深度的过程中需要设定一些与廓线相关的参数。针对不同数据, 这些限定参数的最优取值是不同的。例如对于地转海洋学实时观测阵(Array for Real-time Geostrophic Oceanography, ARGO)数据, 限定参数的最优取值不同于其在WOCE数据中的取值。这些不足为使用该方法带来了不便, 因而提出一种既可以避开阈值类方法中最优阈值选取难题, 又可以弥补非阈值类方法的不足的方法是有意义的。
贝叶斯方法是Thomas Bayes在《An essay towards solving a problem in the doctrine of chances》中首次提出(Thomas Bayes, 1763), 并且经过De Fineti、Jefferys、Savage以及Lindley等人的努力, 逐渐发展出的一套完整的统计推断方法。该方法被广泛用于机器学习、人工智能等领域。本文基于贝叶斯方法, 提出一种识别海洋混合层深度的人工智能(Artificial Intelligence, AI)方法(下称AI法)。对温度(密度)和压强(或深度)建立线性模型, 并用模型系数和误差的方差因子描述温度(密度)及其梯度和涨落信息。在垂向梯度非常小的情况下, 把模型系数首次突变发生的位置识别为混合层下边界。使用2017年2月太平洋海域的ARGO浮标数据的测试结果表明新方法能够准确识别混合层深度, 与其他方法相比准确度更高, 并且在避开阈值类方法中最优阈值的选取难题的同时, 弥补了非阈值类方法的不足。

2 方法

2.1 模型的建立

设一条廓线具有n个有效数据点, 用n维列向量$\psi ={{\left( {{\psi }_{1}},{{\psi }_{2}},...,{{\psi }_{n}} \right)}^{\prime }}$和$P={{\left( {{p}_{1}},{{p}_{2}},...,{{p}_{n}} \right)}^{\prime }}$分别表示温度(位势密度)和对应的压强(或深度)。初始时, 为廓线设置一个长度为l的位于海表的压强(或深度)窗口。它将在数据学习过程中按照步长$\text{step}=1$向下移动。将下移k-1次后的窗口命名为wk, $k=1,2,...,n+1-l$, 并且把它内部的数据点记为${{S}_{k}}=\left( {{\psi }_{k}},{{P}_{k}} \right)$, 其中, ${{\psi }_{k}}={{\left( {{\psi }_{k}},{{\psi }_{k+1}},...,{{\psi }_{k+l-1}} \right)}^{\prime }}$ ${{P}_{k}}=$ ${{\left( {{p}_{k}},{{p}_{k+1}},\cdots ,{{p}_{k+l-1}} \right)}^{\prime }}$, 分别为窗口内的温度(位势密度)和压强(或深度)数据。
在窗口wk内建立线性模型(Koch, 2007):
$\begin{align}& {{\psi }_{k}}={{\beta }_{k,0}}\text{1}+{{\beta }_{k,1}}{{P}_{k}}+{{\varepsilon }_{k}} \\& {{\varepsilon }_{k}}\tilde{\ }N\left( 0,{{\tau }_{k}}^{-1}{{I}_{l\times l}} \right) \\\end{align}$
式中: ${{\mathbf{\varepsilon }}_{k}}={{\left( {{\varepsilon }_{k,1}},{{\varepsilon }_{k,2}},...,{{\varepsilon }_{k,l}} \right)}^{\prime }}$为窗口wk内的模型误差, 在此处认为它是正态的$N\left( 0,{{\tau }_{k}}^{-1}{{I}_{l\times l}} \right)$。Il×l是l×l的单位矩阵, 表示各点模型误差相互独立; 系数βk,0和βk,1蕴含温度(密度)及其垂向梯度信息, 方差因子τ-1描述数据的涨落信息, 他们可以被用于刻画数据的特征。记${{M}_{k}}=\left( 1,{{P}_{k}} \right)$, ${{\beta }_{k}}={{\left( {{\beta }_{k,0}},{{\beta }_{k,1}} \right)}^{\prime }}$。${{\tau }_{k}}^{-1}$和${{\beta }_{k}}$作为一组描述数据特征的统计量, 其先验概率密度函数$f(\tau,\beta)$是机器根据以往数据的信息, 给温度(密度)数据特征的预先认知。它的后验概率密度函数$f\left( {{\tau }_{k}},{{\beta }_{k}}|{{\psi }_{k}} \right)$ 是机器通过学习${{\psi }_{k}}$给出的对数据特征的新认知。它的似然函数
$f\left( {{\psi }_{k}}|{{\beta }_{k}},{{\tau }_{k}} \right)={{\left( 2\pi \right)}^{\frac{1}{2}}}{{\tau }^{\frac{1}{2}}}{{e}^{-\frac{\tau }{2}{{\left( {{\psi }_{k}}-{{M}_{k}}{{\beta }_{k}} \right)}^{\prime }}\left( {{\psi }_{k}}-{{M}_{k}}{{\beta }_{k}} \right)}}$
.表示在数据特征${{\tau }_{k}},{{\beta }_{k}}$取定的条件下, 温度(密度)取实测值的概率(Koch, 2007)。

2.2 数据学习过程

数据学习是更新机器对数数据学习是更新机器对数据特征的认知的过程。用正态-伽马分布${{N}_{\gamma }}\left( \mu ,V,b,\lambda \right)$描述对数据特征的认知, 其中$\mu ,V,b,\lambda $分别为该分布的均值、协方差、尺度参数和形状参数。这种分布具有共轭性, 即先验与后验概率密度函数都是正态-伽马分布(Degroot, 1970; Bernardo et al, 1994)。以窗口wk为例, τk和βk的先验和后验分布分别被记作${{\tau }_{k}},{{\beta }_{k}}\tilde{\ }{{N}_{\gamma }}\left( {{\mu }_{k}},{{V}_{k}},{{b}_{k}},{{\lambda }_{k}} \right)$和${{\tau }_{k}},{{\beta }_{k}}|{{\psi }_{k}}\tilde{\ }{{N}_{\gamma }}\left( \mu _{k}^{*},V_{k}^{*},b_{k}^{*},\lambda _{k}^{*} \right)$。学习过程分为更新和校正两个步骤。
更新步骤从前一窗口的后验分布${{\tau }_{k-1}},{{\beta }_{k-1}}$ $|{{\psi }_{k-1}}\tilde{\ }{{N}_{\gamma }}\left( \mu _{k-1}^{*},V_{k-1}^{*},b_{k-1}^{*},\lambda _{k-1}^{*} \right)$, k>1出发, 此时机器对数据的预先认知是前一窗口中的后验信息。增加数据点$\left( {{\psi }_{k+l-1}},{{p}_{k+l-1}} \right)$后, 根据贝叶斯链式法则更新认知(见附录公式A1)据特征的认知的过程。
$\begin{align} & f\left( {{\tau }_{k}},{{\beta }_{k}}|{{\psi }_{k-1}},{{\psi }_{k+l-1}} \right)= \\ & \frac{f\left( {{\tau }_{k-1}},{{\beta }_{k-1}}|{{\psi }_{k-1}} \right)f\left( {{\psi }_{k+l-1}}|{{\tau }_{k-1}},{{\beta }_{k-1}} \right)}{\iint{f\left( {{\tau }_{k-1}},{{\beta }_{k-1}}|{{\psi }_{k-1}} \right)f\left( {{\psi }_{k+l-1}}|{{\tau }_{k-1}},{{\beta }_{k-1}} \right)\text{d}{{\tau }_{k-1}}\text{d}{{\beta }_{k-1}}}} \\ \end{align}$
式中: $f\left( {{\tau }_{k}},{{\beta }_{k}}|{{\psi }_{k-1}},{{\psi }_{k+l-1}} \right)$是机器对数据信息τk和βk新的认知, 把它记为${{\tau }_{k}},{{\beta }_{k}}|{{\psi }_{k-1}},{{\psi }_{k+l-1}}\tilde{\ }$ ${{N}_{\gamma }}\left( \mu _{k}^{0},V_{k}^{0},b_{k}^{0},\lambda _{k}^{0} \right)$。
然而, 更新步骤得到的分布不能使模型(1)对窗口wk内数据的拟合达到最优化。校正过程则是从新的认知${{\tau }_{k}},{{\beta }_{k}}|{{\psi }_{k-1}},{{\psi }_{k+l-1}}$出发得到使拟合最优化的后验分布的过程。用g-先验(Zellner, 1986; Liang et al, 2008)对先验分布${{\tau }_{k}},{{\beta }_{k}}\tilde{\ }{{N}_{\gamma }}\left( {{\mu }_{k}},{{V}_{k}},{{b}_{k}},{{\lambda }_{k}} \right)$调整如下(附录公式A2):
$\begin{align} & {{\mu }_{k}}=\mu _{k}^{0} \\ & {{V}_{k}}=gV_{k}^{0} \\ & {{\lambda }_{k}}=\frac{1}{2} \\ \end{align}$
根据贝叶斯公式, 后验概率${{\tau }_{k}},{{\beta }_{k}}|{{\psi }_{k}}\tilde{\ }$ ${{N}_{\gamma }}\left( \mu _{k}^{*},V_{k}^{*},b_{k}^{*},\lambda _{k}^{*} \right)$为:
$f\left( {{\tau }_{k}},{{\beta }_{k}}|{{\psi }_{k}} \right)=\frac{f\left( {{\psi }_{k}}|{{\tau }_{k}},{{\beta }_{k}} \right)f\left( {{\tau }_{k}},{{\beta }_{k}} \right)}{\iint{f\left( {{\psi }_{k}}|{{\tau }_{k}},{{\beta }_{k}} \right)f\left( {{\tau }_{k}},{{\beta }_{k}} \right)\text{d}{{\tau }_{k}}\text{d}{{\beta }_{k}}}}$
调整过程中g和bk的取值, 可以由最小描述长度原理(Rissanen, 1978, 2009; Hansen et al, 2003)给出(附录公式A3、A4、A5)。数据学习过程中拟合直线如图1所示。
图1 数据学习过程示意图

图中粉色方框表示窗口

Fig. 1 An example of the learning process.

The blue line is the temperature profile. The magenta box represents the current window. The green, gray, and red dotted lines represent the fitting curves in the priori, the updating, and the posterior process, respectively

2.3 识别依据和检验方法

将${{\tau }_{k}},{{\beta }_{k}}|{{\psi }_{k}}\tilde{\ }{{N}_{\gamma }}\left( \mu _{k}^{*},V_{k}^{*},b_{k}^{*},\lambda _{k}^{*} \right)$沿着τk积分, 得到窗口wk内模型系数的概率密度函数为:
${{\beta }_{k}}\tilde{\ }t\left( \mu _{k}^{*},\frac{b_{k}^{*}V_{k}^{*}}{\lambda _{k}^{*}},2\lambda _{k}^{*} \right)$
这个分布的均值为$\mu _{k}^{*}$。其分布如图2所示:
图2 窗口k内拟合温度及其梯度的概率密度函数

Fig. 2 Probability density function of temperature and its gradient

当k>1时, 模型系数发生突变的数学表达为:
${{\beta }_{k}}={{\beta }_{k-1}}$
式中, ${{\beta }_{k-1}}$和${{\beta }_{k}}$分.别为窗口k-1和k内数据拟合的模型系数。在直观上, 这种突变表现为相邻两个窗口内概率分布的不同。如图3所示, 颜色表示两个相邻窗口内所计算的温度及其梯度的概率密度。图3左侧颜色条带表示后一窗口内的概率分布, 右侧