Data-driven Bayesian Optimization Method for Intermittent hypoxic Training Strategy Decision
-
摘要: 青藏地区快速的经济发展使得进入高原的群体数量日益增加, 随之而来的高原健康问题也愈发突出. 间歇性低氧训练(Intermittent hypoxic training, IHT)是急进高原前常使用的预习服方法, 一般针对不同个体均设置固定的开环策略, 存在方案制定无标准、系统化的理论指导缺乏、效果不明显等问题. 针对以上情况, 设计了一种小样本数据驱动的IHT策略贝叶斯闭环学习优化框架, 建立自回归结构的高斯过程血氧饱和度(Peripheral oxygen saturation, SpO2)预测模型, 并考虑高低风险事件对训练的影响, 设计与氧浓度变化方向和速率相关的风险不对称代价函数, 提出具有安全约束的贝叶斯优化方法, 实现IHT最优供氧浓度的优化决策. 考虑到现有仿真器无法反映个体动态变化过程, 依据“最优速率理论”设计了合理的模型自适应变化律. 所提出闭环干预方法通过该仿真器进行了可行性和有效性验证. 说明该学习框架能够指导个体提升高原适应能力, 减轻其在预习服阶段的非适应性不良反应, 为个性化IHT提供精准调控手段.Abstract: The rapid economic development of Qinghai-Tibet region has led to an increasing number of groups entering the plateau, and the consequent problem of high-altitude health has become increasingly prominent. Intermittent hypoxic training (IHT) is a commonly-used preacclimatization approach before rapidly going to the plateau. It is usually designed as fixed open-loop strategies for different individuals, which has several disadvantages such as no standard formulation, lack of systematic theoretical guidance and poor efficacy. In this paper, a data-driven Bayesian closed-loop learning optimization framework of IHT strategy is designed by using small samples, and a Gaussian process model with autoregressive structure of peripheral oxygen saturation (SpO2) is built for prediction. Based on the predictive model, a risk asymmetric cost function related to the oxygen concentration rate and its direction is developed. Finally, a Bayesian optimization method with safety constraints is proposed to enable the optimal decision of IHT oxygen concentration. Given that the existing simulator cannot reflect the process dynamics of individuals, a reasonable model adaptation law is designed according to the “optimal rate theory”. The feasibility and effectiveness of the proposed closed-loop intervention method are verified by the simulator. These results indicate that the proposed learning framework can help individuals to improve their adaptability to high-altitudes, reduce their non-adaptive adverse reactions in the pretraining stage, and provide precise control solution to personalized IHT.
-
随着高原地区经济的快速发展, 因军事救援、旅游观光、徒步登山等需求进藏或援藏的群体正不断增加, 高原健康问题也日益凸显[1-2]. 高原病是一种长期居住在平原的个体进入高海拔地区时因对低氧环境适应能力不全或失调而导致的综合症, 包括急性高原病、高原肺水肿等, 其引发的头晕、恶心、呕吐等不适症状严重影响人们的正常生活, 极度不适应者甚至会存在呼吸困难和死亡的危险[3-5].
为提升高原环境下的机体适应能力, 进行高原预习服至关重要. 目前, 阶梯习服和预训练被公认为是最有效的降低高原病易感性的两种预习服方法[6]. 特别的, 间歇性低氧训练 (Intermittent hypoxic training, IHT)通过模拟高原环境实现低氧适应能力的提升, 是一种广泛应用的预习服方法[7-9]. 在训练期间, 个体被要求吸入周期性变化的低氧和常氧(高氧), 保证在模拟高海拔环境以提高个体抗缺氧能力的同时, 降低机体长期暴露于低氧环境中的有害影响. 在合适的训练协议设置下, 该方法有助于机体建立适应性增强的生理机制, 加速习服过程, 降低患病风险. 在现有研究中, 策略一般在训练前预先固定设置, 受训者一般被要求暴露在氧浓度为10% ~ 15% 的常压或低压环境中, 持续3 ~ 8分钟, 周期性穿插浓度为21%的常氧或30% ~ 40%的高氧, 为期2 ~ 5分钟. 每1天或2天进行一次训练, 总时长从几十分钟到几小时不等. 具体的, 连续7天的1小时低氧通气被证明能有效地降低再次暴露于高原环境时急性高原病的发生率[10]; 文献[11]尝试将个体暴露在12%的氧浓度中, 每天持续4小时, 通过减少缺氧诱导的炎症和血脂紊乱来促进适应; 间歇性递增式的常压低氧暴露方案在文献[12]中被应用. 可以从上述研究分析发现, 虽然有关IHT对高原低氧适应能力改善效果的研究从训练时长、供氧浓度选择等方面做了许多的探索, 但在整个训练过程中, 现有IHT策略是单一且固定不变的, 无法根据当前适应状态完成自主调节与优化, 并在短时间内达到显著效果, 急进高原需求难以保障. 因此, 本工作旨在设计具有反馈式结构的闭环控制算法, 提供个体化最优训练策略, 实现精准高效的适应性能力增强.
合理准确的预测算法是设计个体化动态训练策略的关键. 在疾病诊断与治疗中, 常采用基于知识、基于机理模型和基于数据的三种方法实现预测和随后的决策. 基于知识的预测一般结合临床经验, 通过创建医疗数据库, 采用统计分析方法设计疾病预测和诊断的定量规则[13]. 然而, 考虑到经验方法依赖专家临床诊断的长期记录结果, 且可迁移应用于不同患者的能力差, 许多研究转而采用基于模型的预测方法, 通过第一性原理对系统动力学建模以揭示疾病作用机制, 并根据预测结果实现诊断与治疗控制[14]. 虽然机理建模提升了预测决策效果, 但过程复杂耗时, 参数众多, 调试存在困难. 随着物联网以及先进测量技术的快速发展, 过程数据能够被连续监控、存储和处理, 因而基于数据的预测算法受到了各领域学者的青睐[15-17]. 回归预测算法可大致分为传统时间序列预测以及机器学习预测两大类. 前者利用统计分析方法推测时序变量的变化方向与趋势[18], 动力学描述能力存在局限, 后者所包括的长短期记忆网络等方法对数据量要求高, 小样本问题难以适用[19-20]. 特别的, 高斯过程是一种具备坚实理论支撑的非参数学习方法, 能够对未知的非线性结构进行建模, 不仅对于小样本具有优异的回归预测性能, 还能够衡量模型的不确定度, 保证模型的可靠性. 许多研究在其基础之上考虑了一系列的变体形式, 并给出了相应的理论证明, 例如: 文献[21]针对计算复杂度大的问题, 推导了一种变分近似的稀疏高斯过程; 文献[22]提出了一种近似的高斯过程动态模型, 并分析了其控制性能. 此外, 高斯过程在各领域也以不同形式得到了广泛的应用. 文献[23]针对术后功能预测问题考虑目标变量取值范围限制, 设计了基于截断正态分布的高斯过程模型, 预测结果显著提升; 文献[24]提出了一种稀疏非平稳的高斯过程回归, 用于预测柴油机氮氧化物.
良好的预测效果可极大提升数据驱动闭环控制的性能. 文献[25]针对血糖控制问题利用餐后预测序列设计了数据驱动的餐前胰岛素优化决策方法,并通过FDA认证的UVA/Podava仿真器进行了验证; 文献[26]基于高斯过程回归与深度强化学习设计了分层人机协控制方法; 文献[27]总结了自治系统的数据驱动预测控制方法与有效应用等. 虽然数据驱动的预测控制在多个问题中得到了有效应用, 但有关高原预习服闭环决策的学习优化方法尚未得到充分研究, 目前的IHT策略固定单一, 缺乏个体化动态调整的有效机制. 针对这一问题, 本文设计了一种数据驱动的IHT策略贝叶斯学习优化框架, 主要贡献如下:
1) 构建了适用于小样本数据的自回归结构高斯过程模型, 实现IHT低氧段血氧饱和度(Peripheral oxygen saturation, $ {\rm{ SpO}}_2) $的滑动时窗式预测. 在此基础上, 考虑到高训练强度下机体非适应性恶化的高风险事件以及低训练强度下见效周期慢的低代价事件, 依据供氧浓度速率及其方向设计了风险不对称代价函数, 提出了具有安全约束的贝叶斯优化算法以实现IHT供低氧浓度的优化决策.
2) 由于现有仿真器无法模拟个体IHT前后的低氧应激系统动态变化过程, 结合“最优速率理论”对其进行了改进[28], 设计模型自适应变化律以建立仿真器模型[29]动力学特征与供氧浓度变化之间的关联关系, 使得虚拟训练者的$ {\rm{ SpO}}_2 $水平能在合理范围内波动并体现随不同训练速率的变化效果.
3) 采用改进的仿真器对优化决策算法进行仿真验证. 实验结果表明, 在闭环调控算法的作用下, 受试者的动态血氧指标(Dynamic $ {\rm{ SpO}}_2 $ index, DSI)[30]以及$ {\rm{ SpO}}_2 $均值较开环调控而言均有显著增加, 且标准差有所下降, 指示数据向高值集中, 说明所设计的闭环贝叶斯优化决策算法能有效提升高原低氧适应性能力, 可为个体化IHT提供精准有效的调控手段.
1. 基于高斯过程的IHT策略优化设计
本文探究具有不同低氧适应能力的机体进行个性化IHT训练策略的有效闭环调控方法. 首先采取高斯过程依次建立低氧段$ {\rm{ SpO}}_2 $的回归模型, 时窗滚动式地迭代获取一组$ {\rm{ SpO}}_2 $预测序列, 并基于该序列设计与预测值相关的风险不对称代价函数, 同时引入供氧浓度的硬安全约束, 最后基于贝叶斯优化进行最优策略求解, 整体算法流程图如图1所示.
图 1 高原适应性能力提升的IHT策略优化决策算法流程图Fig. 1 Flow chart of IHT optimization decision algorithm for high-altitude adaptability improvement1.1 高斯过程
高斯过程定义为在连续域上由无穷个高斯随机变量所组成的随机过程, 是一种概率性的机器学习方法, 因其能从少量历史数据中生成准确的预测值并衡量对应的不确定性而被广泛应用[31]. 令$ x\in $${\bf{R}} ^{n_x} $为输入向量, 对应输出为$ h(x)\in {\bf{R}}$, 则有高斯分布表示为:
$$ \begin{align} h(x) \sim {\cal{G}}{\cal{P}}(m(x),k(x,x')) \end{align} $$ (1) 其中, $ m(x)={\rm E} [h(x)] $ 和$k(x,x')={\rm E}[(h(x)-m(x))\times (h(x')-m(x'))]$分别表示均值函数和协方差函数. 考虑到实际测量时会存在干扰, 引入满足零均值正态分布的噪声项$ \epsilon \sim {\cal{N}}(0,\sigma_n^2) $, 则有观测输出$y= h(x)+\epsilon$. 式(1)可改写为
$$ \begin{align} y \sim {\cal{G}}{\cal{P}}(m(x),k(x,x')+\sigma_n^2\delta(x,x')) \end{align} $$ (2) 其中, $ \delta(x,x') $为克罗内克函数. 当采用具有$ n $组样本的训练集$ {\cal{D}}=\{(x_i,y_i)\vert i=1,\cdots,n\} $进行高斯过程学习时, 可以根据系统先验知识或者经验选择均值函数$ m(x) $和协方差函数$ k(x,x') $的形式. 零均值函数下, 通过最大化如下边际似然函数
$$ \begin{split} \lg p({\boldsymbol{y}} \vert X)=\;&-\frac{1}{2}{\boldsymbol{y}}^\mathrm{T}(K(X,X)+\sigma_n^2I_n)^{-1}{\boldsymbol{y}}\;-\\ &\frac{1}{2}\lg \vert K(X,X)+\sigma_n^2I_n\vert -\frac{n}{2}\lg 2\pi \\[-5pt]\end{split} $$ (3) 确定对应的超参数, 其中$ X=[x_1,\cdots,x_n]^\mathrm{T}\in $${\bf{R}} ^{n \times n_x} $, $ {\boldsymbol{y}}=[y_1,\cdots,y_n]^\mathrm{T}\in $${\bf{R}} ^{n} $, 自身协方差矩阵为$ K(X,X) $. 在给定$ n_* $个新的输入向量后, 构成矩阵$ X_* $, 基于所训练的模型可对相应的输出值$ h_* $进行预测. 其联合先验高斯分布为
$$ \begin{align} \left[ \begin{array}{lr} {\boldsymbol{y}}\\h_* \end{array} \right] ={\cal{N}} \left({{\bf{0}}}, \left[\begin{array}{lr} K(X,X)+\sigma_n^2I_n, & K(X,X_*) \\ K(X_*,X) ,& K(X_*,X_*) \end{array} \right] \right) \end{align} $$ (4) 随后可求得$ h_* $的后验分布
$$ \begin{split} &h_* \vert X,{\boldsymbol{y}},X_* \sim {\cal{N}}(\mu_*,\sigma_*)\\ &\mu_*=K(X_*,X)[K(X,X)+\sigma_n^2I_n]^{-1}{\boldsymbol{y}}\\ &\sigma_*^2=K(X_*,X_*)-K(X_*,X)[K(X,X)+\sigma_n^2I_n]^{-1}\cdot\\ &\;\;\quad\quad K(X,X_*)\\[-10pt] \end{split} $$ (5) 根据式(5), 可利用历史数据求得新时刻的预测值, 并得到估计的不确定性.
1.2 低氧段的$ {\bf{ SpO}}_{\boldsymbol{2}} $预测
考虑到$ {\rm{ SpO}}_2 $浓度变化主要受外界氧浓度变化影响, 氧气进入人体内后, 会经过一系列复杂的生理过程, 对各部分进行精准机理建模存在显著困难, 因而考虑数据驱动的建模方法. 另外, 在已知历史血氧序列段和可承受范围内供氧浓度的条件下, 当前血氧值必定在一定合理范围内波动. 因此, 可采用高斯过程模型, 以上述信息作为输入向量, 近似刻画当前血氧值的变化.
在IHT中, 高氧段的设置是为了缓解受试者长期处于低氧状态的不适症状, 尽可能避免机体损伤, 其对于适应性提升不会产生显著效果. 相反的, 氧浓度由高降低过程以及持续低氧过程的生理信号具有高原适应性能力的高度表征作用, 故筛选低氧段$ {\rm{ SpO}}_2 $用于建模与预测. 在一日IHT中, 一个供氧周期内的高氧浓度$ c_h $和低氧浓度$ c_l $的持续时长分别为$ t_h $和$ t_l $, 共包含$ n_s $个低氧段样本. 在这里, 我们将前$ n_s-1 $个样本作为训练集, 最后一个样本作为测试集. 为减轻个体深呼吸、摘下面罩等行为所带来的非高斯分布干扰, 在建模训练前首先采用滑动平均滤波对$ {\rm{ SpO}}_2 $样本数据进行预处理. 考虑到高氧相关变量对适应性提升无改善作用, 且用于评估的DSI指标是在固定低氧时长下提出的, 设定$ c_h=35 $%、$ t_h=180 $s和$ t_l=300 $s, 其数值遵从文献[28]中的临床试验设置, 均在常规IHT的合理范围内. 此外,将最为关键的低氧浓度$ c_l $作为输入调控变量.
记$ t $为一个低氧段的采样时刻, 并认为各低氧段中对应$ t $时刻的$ {\rm{ SpO}}_2 $数据属于同一自回归结构$ \Gamma_{t}(\cdot) $, 则有如下非线性自回归模型结构形式:
$$ \left\{\begin{aligned} &y_{t+1}=\Gamma_t(z_t)+\epsilon_{t}\\ &z_t=[c_l,y_{t-l},\cdots,y_t] \end{aligned} \right.$$ (6) 其中, $ t\geq l+1 $, $ l $为自回归模型阶数. 当采样周期为$ t_s $时, 低氧段数据点个数为$ n_t=t_l/t_s $, 则待预测的点个数为$ n_p=n_t-l-1 $. 值得注意的是, 各低氧段同一时刻的数据点均用同一高斯过程来进行描述, 即该步骤涉及到多个独立的高斯过程模型参与, 每一模型都对应实现某时刻$ {\rm{ SpO}}_2 $值的预测. $ t $时刻的输入$ z_t $由低氧浓度$ c_l $和前$ l+1 $个 $ {\rm{ SpO}}_2 $数据组成, 输出为$ y_{t+1} $. 在实际应用中, $ t_s $通常为1 s. 由于低氧段长度是固定的, 因此自回归结构中的历史数据长度决定着待预测序列长度, 选取$ l=99 $(即$ z_t $包含100个$ {\rm{ SpO}}_2 $数据点)使得历史序列包含足够丰富的特征信息的同时, 确保后段部分渐进稳态的$ {\rm{ SpO}}_2 $序列(更能反映出受试者的上高原后的真实适应情况)得到完全的预测. 在这里, 选取常用的零均值函数以及平方指数函数进行高斯过程训练
$$ \begin{align} m(z_t)&=0 \end{align} $$ (7) $$ \begin{align} k(z_t,z_t')&=\sigma^{2}\exp\left(-\frac{r^2}{2\ell^2}\right) \end{align} $$ (8) 其中$ r=\Vert z_t-z_t'\Vert $表示两输入变量的欧几里得距离, $ \sigma^{2} $和$ \ell $是待优化确定的超参数. 针对$ t+1,\cdots,t+n_p $时刻重复上述训练, 可获得一组描述$ {\rm{ SpO}}_2 $变化轨迹概率分布的高斯过程回归模型$ {\cal{G}}{\cal{P}}_1,\cdots,{\cal{G}}{\cal{P}}_{n_p} $. 接着, 利用上述模型对向量$ {\boldsymbol{y}}_p=[y_{t+1},\cdots,y_{t+n_p}]^\mathrm{T} $中的输出值进行逐一预测. 将测试集的数据代入$ {\cal{G}}{\cal{P}}_1 $获得预测值$ \hat{y}_{t+1} $, 随后, 对于时刻$ t+1 $, 将$ z_t $中的$ z_{t-l} $剔除, 并按时序加入预测值$ \hat{y}_{t+1} $组成新的输入向量$ z_{t+1} $. 将其代入到模型$ {\cal{G}}{\cal{P}}_2 $得到新的预测值$ \hat{y}_{t+2} $, 依次迭代后即可获得各元素均满足$ \hat{y}_i \sim {\cal{N}}(m_i,\sigma_i) $分布的预测向量 $ {\boldsymbol{\hat{y}}}_p=[\hat{y}_{t+1},\cdots,\hat{y}_{t+n_p}]^\mathrm{T} $(见图2). 该向量与待确定变量$ c_l $相关, 将设计于代价函数的首项中参与优化求解.
1.3 代价函数与优化问题设计
在高原健康研究中发现, 不同机体存在适应高原环境变化机制的个性化最优适应速率$ \Delta c_{op} $[22]. 若所处环境的低氧浓度变化快于个体最优的适应速率, 则会导致显著性高原反应, 因此, IHT所模拟环境的变化应当处于受试者可接受的适应性增强范围内, 使得产生有利于机体稳态水平提升的生理学低氧应激反应. 考虑上述适应性变化规律, 设计了如下的不对称代价函数:
$$ \min_{c_l} {\cal{L}}(c_l)\;\;\;\;\qquad\qquad\qquad\qquad\qquad\qquad $$ $$ \begin{align} {\rm{s.t.}}\; \; &\hat{y}_i \sim {\cal{N}}(m_i,\sigma_i), i=1,2, \cdots, n_p \end{align} $$ $$ c_{\rm min}\leq c_l\leq c_{\rm max} \qquad\qquad\qquad\qquad\qquad$$ (9) 其中, $ c_l $是待优化的低氧浓度, $ {\cal{L}}(c_l) $设计为:
$$ \begin{align} &{\cal{L}}(c_l):={\cal{L}}_s+w{\cal{L}}_{v} \end{align} $$ (10) $$ \begin{align} &{\cal{L}}_{s}:={\rm E}[({\boldsymbol{\hat{y}}}_{p}-{\boldsymbol{y}}_{r})^\mathrm{T} Q({\boldsymbol{\hat{y}}}_{p}-{\boldsymbol{y}}_{r})] \end{align} $$ (11) 其中, w为权重系数, $ {\cal{L}}_s $用于衡量预测向量与目标$ {\rm{ SpO}}_2 $向量$ {\boldsymbol{y}}_{r} $的距离. 在这里, 考虑概率分布中包含的预测值不确定性信息将其设计为随机变量的形式. $Q={\rm diag}\{q_1,\cdots, q_{n_p}\}$为惩罚矩阵, $ q_i $取值越大表明$ {\boldsymbol{\hat{y}}}_{p} $向量中第$ i $个$ {\rm{ SpO}}_2 $预测值相对于目标值的距离惩罚越大. 期望值是否能准确求出取决于矩阵$ Q $的具体形式, 当矩阵$ Q $的设计形式较为复杂, 则无法求得期望的准确值. 考虑到这种情况, 本文设计了一种通用的求解方式(具体详见第2节). 定义$ \Delta c_i=c_{l,i}-c_{l,i-1} $为第$ i $日和第$ i-1 $日IHT的供低氧浓度之差, 将其用于衡量训练中供氧环境的变化速率. 由于下述代价函数是基于一日IHT进行优化计算的, 故省略下标$ i $以便于表述, 则所设计的$ {\cal{L}}_{v} $具体表达式如下:
$$ {\cal{L}}_{v}= \left\{ \begin{array}{lll} &\ln(1+\Delta c),&\Delta c\geq 0 \\& 0,&-2<\Delta c<0\\ &-\dfrac{10}{1+0.1\exp(-3\Delta c-9)}\;+\\ &\quad\dfrac{10}{1+0.1\exp(-3)},&\text{其他} \end{array} \right. $$ (12) 对应的变化曲线如图3所示. 考虑到IHT氧浓度在不同速率调控下对机体的影响, 将该项设计成了风险不对称的形式. 具体的, 低氧浓度大幅度降低被判定为高安全风险的过强训练事件, 浓度较高则属于使得适应效果减弱的低代价事件. 尽管不同个体的$ \Delta c_{op} $存在差异, 但其总体属于一个分布密集的区间内, 且无法测量, 因此$ {\cal{L}}_v $给出了一个较大的适应性范围. 在这里, $ -2< \Delta c<0 $被认为是使得适应性增强的合适低氧下降速率范围, 因而不进行惩罚, $ {\cal{L}}_{v}=0 $; $ \Delta c\leq-2 $表明环境氧浓度处于快速变化状态, 机体可能产生高原不适反应, 且$ \Delta c $越小表明严重程度更高, 导致不可逆损伤的发生, 故在该情况下$ {\cal{L}}_{v} $项的惩罚力度随着$ \Delta c $的减小而显著增加; 当$ \Delta c\geq0 $时, 表明本日IHT氧浓度较上一日有所增加, 训练强度有所降低, 适应性改善力度减小, 然而, 其变化并不会引入更多的训练风险. 因此, 该情况下为$ {\cal{L}}_v $设置比$ \Delta c\leq -2 $低的惩罚力度. 除此之外, 设置了保证生理安全需求的硬约束条件$ c_{\rm min}\leq c_l\leq c_{\rm max} $以限制$ c_l $的可调控范围, 在受试者安全的前提下提升IHT策略的有效性.
图 3 所设计代价函数$ {\cal{L}}_{v} $部分的惩罚强度在不同$ \Delta c $下的变化Fig. 3 The penalty changes of designed $ {\cal{L}}_v $ term under different $ \Delta c $ values1.4 迭代可行性与稳定性
如上文所述, 人体呼吸系统运作时需要机体各子系统同步参与, 非线性程度强, 且受多源的内外部干扰, 因此对高原预习服闭环优化决策系统进行全面的理论分析存在很强的挑战性. 另一方面, 本文主要从实际应用角度出发, 针对特殊的被控对象(人体)设计控制器, 旨在提升高原预习服实际训练效果, 因此, 控制器结构在设计时具有明显的针对性和特殊性. 然而, 考虑到算法可行性和稳定性对系统安全和有效性的影响, 本文也对其进行了充分的考虑. 从IHT的实际应用情景不难得到, 该优化问题在初始时刻是可行的. 由于其内含关于优化变量的硬约束, 因此对于任意优化求解时刻$ T\in $N$ ^{+} $以及对应的输出值$ \hat{y}_i $的最优解均满足$ c_{l,T}^{*}\in [c_{\rm min},c_{\rm max}] $, 故该问题具有迭代可行性. 此外, 由于本工作考虑的系统为复杂人体, 其输出$ {\rm{ SpO}}_2 $值$ {\boldsymbol{y}}_{T} $必定有界, 故对于$ \forall T\geq0 $和常值目标向量 $ {\boldsymbol{y}}_{r} $, $ \exists \delta>0 $, 满足$ \Vert {\boldsymbol{y}}_{T}-{\boldsymbol{y}}_{r}\Vert < \delta $. 令$ c_{l,T} $和$ c_{l,T}^{*} $分别表示 $ T $时刻优化问题的可行解和最优解, 由上述迭代可行性分析可得, $ T $时刻的最优解必定为$ T+1 $时刻的可行解, 故若初始解$ c_{l,0}^{*} $对应的预测值满足$ \exists \delta_{1}>0 $使得$\Vert {\boldsymbol{y}}_{p,0}^{*}- {\boldsymbol{y}}_{0} \Vert < \delta_{1}$, 则有$ \exists \delta,\delta_{1}>0 $, $\Vert {\boldsymbol{y}}_{p,0}^{*}-{\boldsymbol{y}}_{r}\Vert=\Vert {\boldsymbol{y}}_{p,0}^{*}-{\boldsymbol{y}}_{0}\;+ {\boldsymbol{y}}_{0}-{\boldsymbol{y}}_{r} \Vert < \delta+\delta_{1}$. 针对$ T=1 $时刻, 令$ c_{l,1}=c_{l,0}^{*} $, 则有$ \exists\delta,\delta_1>0 $, $\Vert {\boldsymbol{\hat{y}}}_{p,1}-{\boldsymbol{y}}_{r} \Vert=\Vert {\boldsymbol{y}}_{p,0}^{*}-{\boldsymbol{y}}_{r} \Vert < \delta+ \delta_1$, 故在取最优解$ c_{l,1}^{*} $的情况下满足$\Vert {\boldsymbol{y}}_{p,1}^{*}-{\boldsymbol{y}}_{r}\Vert\leq\Vert {\boldsymbol{\hat{y}}}_{p,1}\,- {\boldsymbol{y}}_{r}\Vert$. 类似的, 对于$ \forall T\geq0 $, 有$\Vert {\boldsymbol{y}}_{p,T+1}^{*}-{\boldsymbol{y}}_{r} \Vert \leq \Vert {\boldsymbol{y}}_{p,T}^{*}\,- {\boldsymbol{y}}_{r}\Vert$. 此外, 随着训练增多, 样本数增加, 预测精度提升, 在$ \exists \delta_1>0 $使得$ \Vert {\boldsymbol{y}}_{p,0}^{*}-{\boldsymbol{y}}_{0} \Vert <\delta_1 $的初始条件下, 可得$ \forall T $, $ \Vert {\boldsymbol{y}}_{p,T}^{*}-{\boldsymbol{y}}_{T} \Vert <\delta_1 $, 则有$\lim_{T\rightarrow \infty}\Vert {\boldsymbol{y}}_{T}-{\boldsymbol{y}}_{r} \Vert= $ $\lim_{T\rightarrow \infty}\Vert {\boldsymbol{y}}_{T}-{\boldsymbol{y}}_{p,T}^{*}\;+\;{\boldsymbol{y}}_{p,T}^{*}\;-\;{\boldsymbol{y}}_{r}\Vert \; < \; \delta_1+\lim_{T\rightarrow \infty}$ $\Vert{\boldsymbol{y}}_{p,T}^{*}-{\boldsymbol{y}}_{r}\Vert < \delta+2\delta_1$, 故系统是有界稳定的. 特别地, 当满足$ \Vert {\boldsymbol{y}}_{p,T+1}^{*}-{\boldsymbol{y}}_{r}\Vert<\Vert {\boldsymbol{y}}_{p,T}^{*}-{\boldsymbol{y}}_{r}\Vert $时, 有$\lim_{T\rightarrow \infty} \Vert {\boldsymbol{y}}_{T}- {\boldsymbol{y}}_{r}\Vert < \delta_1$, 预测误差上界$ \delta_1 $随着预测精度的提升而减小, 稳定误差越小.
2. IHT供氧浓度贝叶斯优化决策及高原适应性评估
考虑到式(9)中代价函数存在不对称风险惩罚项, 且当权重矩阵设计复杂时无法求得数学期望, 难以获取函数梯度信息, 通过传统方法求得待优化低氧浓度$ c_l $的闭式解, 故采用贝叶斯优化方法, 并结合蒙特卡洛算法计算求解上述优化问题. 另外, 为衡量每日IHT后的高原适应能力改善效果, 采用指标DSI进行算法的有效性评估.
2.1 贝叶斯优化
贝叶斯优化是一种目标函数优化算法, 可通过机器学习方法对黑箱函数进行猜测, 并基于此求解该函数的最小值或最大值. 该方法无需提供优化函数的具体形式以及梯度信息, 仅通过采样数据即可实现迭代求解. 贝叶斯优化具备两大核心过程, 分别为: 构建贝叶斯统计模型用于目标函数建模以及利用采集函数决策采样点. 在本方法中, 采取高斯过程对所设计的不对称代价函数进行建模, 利用初始化样本得到高斯过程模型后, 选择改善期望(Expected improvement, EI)采集函数确定每一步迭代优化的采样点, 求得采样点对应的观测值后, 将该样本加入到数据集用于更新GP后验分布, 重复上述步骤以获取令代价函数最小的供低氧浓度取值, 下面进行详细说明.
2.1.1 代价函数建模
由于代价函数(9)形式较为复杂, 难以获取$ c_l $和$ {\cal{L}}(c_l) $的映射关系, 采用高斯过程回归对其进行拟合. 设定均值函数为常数形式, 协方差函数为平方指数形式. 超参数是基于$ n_o $个代价函数观测样本$ {\cal{D}}_{1:n_o} = \{c_{l,1:no} ,{\cal{L}}(c_{l,1:no})\} $, 通过最大化对数似然函数(3)确定的. 在预测步骤中, 新低氧浓度$ c_{l}^{*} $对应的代价函数值$ {\cal{L}}(c_{l}^{*}) $通过后验分布(5)获得, 随后将被用于构造采集函数以确定下一次待评估的供氧浓度. 在得到待评估氧浓度后, 采用蒙特卡洛算法对$ {\cal{L}}_s $随机变量项进行估计. 根据概率分布$ y_p^{′} \sim {\cal{N}}(m_i, \sigma^2_i) $生成1000条低氧段$ {\rm{ SpO}}_2 $轨迹样本并计算对应期望值, 在加上风险不对称项的对应取值后, 得到最终的$ {\cal{L}}(c_l) $观测值.
2.1.2 采集函数
采集函数能够反映新样本给优化带来的收益, 被用于确定贝叶斯优化搜索空间中的采样点, 对优化方法具有关键的作用. 该函数需要权衡开发和探索之间的关系, 兼顾当前可行域和高预测不确定性区域的影响, 保证在新增较少采样点的同时实现函数最优化. 考虑到在常见的采集函数中, 置信区间上界采集函数和改善概率采集函数分别存在对权重参数敏感, 无法考虑改善的幅度, 易陷入局部最优的问题, 在本工作中选择EI采集函数, 其不仅能衡量改进的概率, 还能通过数学期望反映改进量大小, 具体表达式如下:
$$\left\{ \begin{aligned} &I(c_{l}^{*}):=({\cal{L}}_m-\hat{{\cal{L}}}(c_{l}^{*})){\boldsymbol{1}}({\cal{L}}_m-\hat{{\cal{L}}}(c_{l}^{*}))\\ &\alpha_{EI}(c_{l}^{*}):={\rm E}[I(c_{l}^{*})] \end{aligned}\right. $$ (13) 其中, 截止到当前时刻的代价函数最小观测值用$ {\cal{L}}_m $表示, 待评估低氧浓度$ c_{l}^{*} $的代价函数预测值$ \hat{{\cal{L}}}(c_{l}^{*}) $服从正态分布$ \hat{{\cal{L}}}(c_{l}^{*}) \sim {\cal{N}}(m_{*}(c_{l}^{*}), \sigma_*^2(c_{l}^{*})) $. 相应的, $ {\cal{L}}_m-\hat{{\cal{L}}}(c_{l}^{*}) $ 和 $ {\boldsymbol{1}}({\cal{L}}_m-\hat{{\cal{L}}}(c_{l}^{*})) $分别表示最小值的预期改善量及其发生的概率. 在利用高斯过程模型进行拟合的情况下, 式(13)可解析地表示为如下形式:
$$ \begin{split} &\alpha_{EI}= \left\{ \begin{array}{lll} &({\cal{L}}_m-m_{*}(c_{l,m}))\Psi(C)\;+&{}\\ &\qquad\sigma_*^2(c_{l}^{*})\psi(C),&m_{*}(c_{l,m})\geq 0 \\ &0,&\text{其他} \end{array} \right.\\ &C=\frac{{\cal{L}}_m-m_{*}(c_{l,m})}{\sigma_*^2(c_{l}^{*})}\\[-15pt] \end{split} $$ (14) 在这里, $\Psi(\cdot)$和$ \psi $分别表示标准正态分布的累积分布函数和概率密度函数. 在最大化$ \alpha_{EI} $后即可确定下一次待评估的低氧浓度值$ c_{l}^{*} $.
$$ \begin{split} &\max_{c_{l}^{*}} \; \alpha_{EI} (c_{l}^{*})\\&\; {\rm{s.t.}} \; \; c_{\rm min}\leq c_l\leq c_{\rm max},\\ &\;{\cal{L}}(c_l):={\cal{L}}_s+w{\cal{L}}_{v} \end{split} $$ (15) 具体的计算步骤见算法1.
算法1. IHT供低氧浓度的贝叶斯优化决策算法
输入. 初始化样本集$ {\cal{D}} $
输出. 最优低氧浓度$ c_{l}^{*} $
1) while 迭代次数$ \leq n_o $ do
2) 利用数据集$ {\cal{D}} $训练高斯过程
3) 通过式(4), (5)计算待评估低氧浓度$ c_{l}^{*} $对应的代价函数预测值
4) 最大化采集函数, 确定待评估低氧浓度$ c_{l}^{*} $
5) 结合蒙特卡洛方法获取$ c_{l}^{*} $对应的代价函数观测值$ {\cal{L}}(c_{l}^{*}) $
6) 将样本$ (c_{l}^{*},{\cal{L}}(c_{l}^{*})) $加入数据集
7) end while
8) 选择令代价函数最小的观测值作为最终求解的低氧浓度值$ c_{l}^{*} $
2.2 低氧适应能力评估
为评价所设计的算法性能, 本文采用指标DSI对每日训练成效进行评估[28]. 令$ \{\eta_{\kappa}^{i}: \kappa=1,\cdots,n_t\} $为一日IHT内第$ i $个低氧段序列, 则单个低氧段的DSI指标$ \tau_i $被定义为从供氧浓度降低起, $ {\rm{ SpO}}_2 $持续保持在$ 85\% $以上的时长. 收集一日IHT内所有低氧段的DSI并将其按从小到大排列, 即$\{\tau_i:i= 1, \cdots,n_s\}$, $ \tau_i<\tau_j $, $ i<j $. 基于直方图选取较低20%百分位数作为一日IHT的DSI指标. 在IHT结束之后, 利用新获得的数据计算指标DSI, 用于评估经前期训练后机体适应能力是否得到改善.
3. 结果
3.1 高原适应能力评估的仿真器改进
本文的算法通过文献[27]中的仿真器进行验证. 考虑到现有仿真器是基于线性时不变的三阶自回归滑动平均模型(Auto-regressive moving average with exogenous inputs, ARMAX)结构辨识所得, 仅能在不同参数选择下, 仿真不同虚拟受试者单日IHT中的$ {\rm{ SpO}}_2 $波动情况, 不能描述经过变化供氧浓度作用后的机体低氧应激动态变化过程, 因此对仿真器做了合理的修正. 原始模型的表达式如下:
$$ \begin{align} A(q)y_{k}=B(q)u_{k}+C(q)e_{k} \end{align} $$ (16) 其中, $ u_k $为供氧浓度, $ y_k $为$ {\rm{ SpO}}_2 $输出. 在这里, 采样变量$ k $描述的是受试者在整个IHT过程中的采样时刻, 不同于上述高斯过程模型中的$ t $变量仅用于描述低氧段, 故选取不同采样时刻变量用以区分. 由于仿真器中表征的机体低氧适应特性与模型中的三个极点存在关联$(A(q)=1+a_1q^{-1}+a_2q^{-2}\;+ a_{3}q^{-3}= 0$的解), 将通过调节各极点数值模拟虚拟受试者不同情况下的动力学过程. 鉴于独立地改变$ a_1 $、$ a_2 $和$ a_3 $较为复杂, 采取仅通过系数$ a_1 $进行调节的方法改变三个极点的取值. 由“最优速率理论”可知, 当机体所处环境的氧浓度降低的变化速率快于个体最优适应速率$ \Delta c_{op} $, 则产生显著性高原反应; 若慢于$ \Delta c_{op} $, 则个体会在该环境下产生有利于系统稳态水平提升的生理学低氧应激反应. 此外, 当氧浓度增加时, 适应提升强度会逐渐降低. 遵循以上规律, 将$ a_1 $与训练期间的氧浓度变化信息进行关联, 设计了如下的模型变化律:
$$ a_{1,i}= \left\{ \begin{array}{lll} &a_{1,i-1}+a_T,& a_S \geq a_T \\ &a_{1,i-1}-a_T,& a_S \leq -a_T \\ &a_{1,i-1}+\displaystyle \sum\limits_{i=1}^{n_d}\Delta a_{1,i},& \text{其他} \end{array} \right. $$ (17) 其中, $ a_{1,i} $为第$ i $天虚拟受试者对应的$ a_1 $取值, $i= 1,\cdots,n_d$, 且有$ a_S=\sum_{j=1}^{i}\Delta a_{1,i} $, $ a_T=4\Delta a_{1,\max} $. $ \Delta a_{1,i} $通过下式确定:
$$ \Delta a_{1,i}= \left\{ \begin{array}{lll} &0.3(6-\Delta c_i) \times 10^{-8},& \Delta c_i \geq 0 \\ &-3\Delta c_i \times 10^{-8},& \Delta c_{op}\leq \Delta c_i < 0 \\ &0.5\Delta c_i\times 10^{-8},& \text{其他} \end{array} \right. $$ (18) 当$ \Delta c_i >0 $时, 以较慢的增长速度提升$ \Delta a_{1,i} $值, 进而增加$ a_{1,i} $, 使$ {\rm{ SpO}}_2 $水平得以提升; 当$\Delta c_{op}\leq \Delta c_i \leq 0$时, 高原适应能力提升成效较$ \Delta c_i >0 $时更好, 且随着$ \Delta c_i $越靠近$ \Delta c_{op} $, $ {\rm{ SpO}}_2 $的恢复水平越好. 然而, 当$ \Delta c_i < \Delta c_{op} $时, 机体会由于氧浓度下降过快而产生非预期反应. 考虑到在实际上高原时, 氧浓度随海拔的增加而降低, 设置$ \Delta c_{op}<0 $以模拟个体缓慢上高原的过程. 此外, 该式对$ a_1 $参数的调试设置了界限, 并考虑了每次训练后对个体影响的累加效应. 虽然不同受试者存在一定的个体间差异, 但均属于正常的生理变化范围, 故对参数$ a_1 $的选取设定了上下阈值, 阈值大小由$ \Delta a_{1,\max} $决定, 取值为$ 10^{-7} $. 该值是在10名虚拟受试者低氧应激系统的基础上调试得出的, 调试准则为血氧饱和度响应不得超过合理上下限范围, 即SpO2上限为100%, 下限由于个体化差异有所不同, 但大多在70% ~ 85%之间波动, 设置阈值能防止仿真器输出异常生理值, 保证仿真的正确性. 在结束每一次IHT后, 产生的效果均有所不同, 可能会使得机体性能改善或者恶化. 例如: 在训练前期, 当决策算法尚未充分学习到适应性变化规律时, 可能会产生不利的操作, 从而导致习服能力下降. 基于上述认识, 设计了$ a_S= \sum_{j=1}^{i}\Delta a_{1,i} $项, 该项反映了截止到当前训练时的总成效, 并综合体现在了$ a_1 $的变化里, 其中$ \Delta a_{1,i} $表示每日IHT的效果, 当$ \Delta a_{1,i} >0 $时, 说明当天的训练增强了机体的低氧适应性, 反之则减弱. 另外, 对虚拟受试者的初始适应性做了如下的设定: 若尚未进行IHT前, 虚拟受试者的初始适应性浓度$ c_{o} $高于训练给定的初始氧浓度, 那么机体会出现适应能力恶化的情况, 直至IHT氧浓度提升至$ c_{o} $以上的水平后, 将根据上述设计规律进行不同程度的能力提升或减弱. 值得注意的是, 虽然在代价函数不对称风险项以及仿真器中模型变化规律均遵循“最优速率理论”进行设计和改进, 但两者具有本质的不同. 在优化决策中, 由于$ \Delta c_{op} $不明确, 因此$ {\cal{L}}_{s} $中不同区域的惩罚设计是模糊的, 在已知大区间范围变化趋势的前提下通过调试后确定. 然而, 在仿真器使用时, 需要对虚拟受试者设置不同的$ \Delta c_{op} $参数以及初始适应情况以模拟个体间的高原适应性差异, 并将设计的模型变化规律作为真值, 对算法学习该规律的能力进行评价验证.
3.2 仿真验证
在经过上述修正之后, 通过仿真器对所设计的优化决策算法进行了有效性评估. 在第一日IHT的初始适应浓度设置为$ c_o=11 $, 后续训练的浓度待确定, 总时长共64分钟. 在此期间, 以$ t_s=1 $s的采样周期连续测量$ {\rm{ SpO}}_2 $数据, 用于算法的训练与测试. 由于DSI指标是基于氧浓度为11%的IHT数据设计的, 因此在一日IHT结束后, 会将虚拟受试者所处环境重新返回至11%氧浓度进行DSI指标的计算, 待确定DSI值后, 再继续进行优化求解步骤. 这在实际应用中也具有一定的可操作性, 受训后个体可进行11%氧浓度的短暂测试以计算DSI, 随后再继续进行IHT. 考虑到2500 m是目前急性高原病发病的普遍最低点[6, 32], 故选取2500 m海拔高度所对应的氧浓度(约15%)作为约束的上阈值. 此外, 下阈值遵循文献[33]中的高强度实验设置, 以保证满足安全要求的同时, 在最大可操作性范围内进行探索与决策. 其他具体参数取值见表1, $ H_i $表示含全1元素的行向量, $I $表示单位矩阵.
表 1 相关参数取值Table 1 Related parameters参数 含义 取值 $ n_s $ 一日IHT的低氧段总数 8 $ n_t $ 一段低氧段预测点总数 200 $ n_p $ 一段低氧段采样点总数 300 $ \sigma^{2} $ 平方指数核函数超参数 10 $ \ell $ 平方指数核函数超参数 10 $ \boldsymbol{y}_{r} $ 目标$ {\rm{ SpO}}_2 $向量 $ [95H_{50},90H_{50}, 85H_{100}]^\mathrm{T} $ $ Q $ 代价函数惩罚矩阵 136I $ w $ 代价函数权重系数 1000 $ c_{l,\min} $ 供氧浓度下阈值 (%) 10 $ c_{l,\max} $ 供氧浓度上阈值 (%) 15 为了较为详细地分析所设计算法的性能, 首先对虚拟受试者个例进行了对比验证. 图4和图5展现了虚拟受试者在不同氧浓度变化速率$ \Delta c_{op} $和初始适应浓度$ c_o $设定下的低氧适应变化. 可以从中看出, 两图中在采用开环策略和闭环策略下进行IHT的$ {\rm{ SpO}}_2 $曲线具有明显的差异. 当进行开环策略时, $ {\rm{ SpO}}_2 $持续降低, 这是因为在训练初期产生了大幅度的环境变化, 引起了个体的非适应性反应. 相反的, 在闭环情况下, 尽管IHT首日的氧浓度设定值低于个体的最初适应浓度, 但经过基于新样本的高斯过程学习后, 预测模型不断更新, 优化决策性能也逐步得到改善, 从而个体的适应性效果也有所提升. 下面将针对两种设置情况进行具体数值对比分析.
图 4 虚拟受试者1和2采取开环和闭环策略进行IHT的$ {\rm{ SpO}}_2 $曲线Fig. 4 The $ {\rm{ SpO}}_2 $ curves of simulated subject 1 and 2 that perform IHT based on traditional open-loop strategy and proposed closed-loop strategy
图 5 虚拟受试者3和4采取开环和闭环策略进行IHT的$ {\rm{ SpO}}_2 $曲线Fig. 5 The $ {\rm{ SpO}}_2 $ curves of simulated subject 3 and 4 that perform IHT based on traditional open-loop strategy and proposed closed-loop strategy在图4中, 虚拟受试者1和2被设置为具有相同的最优适应速率$ \Delta c_{op}=-1.5 $, 而未训练前的初始低氧适应浓度有所不同, 分别为$ c_o=13 $和$ c_o=12 $(表2). 在这两种情况下, 持续性的高强度开环策略均恶化了个体的高原习服能力, 虚拟受试者1和2的DSI由最初的137 s分别降至了84 s和90 s, $ {\rm{ SpO}}_2 $均值从90.9%降至了84.9%和87.8%. 相反的, 将闭环优化策略作用于受试者后, DSI性能均有大幅度提升, 分别从初始值137 s提升至了256 s和300 s, 且$ {\rm{ SpO}}_2 $均值从90.9%提升至了93.6%和96.2%. 此外, 采用开环和闭环策略下的标准差较初始状态均有所降低, 结合$ {\rm{ SpO}}_2 $变化趋势来看可知, 其数据分别向低值和高值集中. 虽然不同个体设定下的适应性能均有所增强, 但$ c_o = 12 $设定下的加强效果整体优于$ c_o = 13 $设定下的情况, 出现这种结果的原因可能是前者的初始适应性更好, 这意味着在进行供氧浓度优化决策时可给定更高训练强度, 加速个体习服, 且在氧浓度限制范围内, 当变化速率超出最优适应速率$ \Delta c_{op} $时, 非适应性症状更轻.
表 2 虚拟受试者1和2采取开环策略进行IHT和采取闭环策略进行IHT的效果对比Table 2 Comparison results of simulated subject 1 and 2 trained by using traditional open-loop strategy and proposed closed-loop strategy个体设定 指标 初始状态 开环策略 闭环策略 $ c_{o}=13 $,
$ \Delta c_{op}=-1.5 $DSI (s) 137 84 256 SpO2平均值(%) 90.9 84.9 93.6 SpO2标准差(%) 5.5 5.4 4.8 $ c_{o}=12 $,
$ \Delta c_{op}=-1.5 $DSI (s) 137 90 300 SpO2平均值(%) 90.9 87.8 96.2 SpO2标准差(%) 5.5 5.4 3.6 在图5中, 虚拟受试者3和4具有相同的适应性变化趋势, 其个体设定分别为$c_o = 12$, $ \Delta c_{op} = -1 $和$ c_o = 12 $, $ \Delta c_{op} = -1.5 $(表3). 可以看出, 在闭环策略下, 机体的低氧适应性能力均有所增强, 其DSI从初始的154 s提升为203 s和300 s, 分别增长了49 s和146 s, $ {\rm{ SpO}}_2 $均值从91.4%分别增加到94.7%和96.4%. 然而, 在开环情况下, DSI分别降低了18 s和73 s, 从最初的154 s减至136 s和81 s, 降幅达11.6%和49.3%. $ {\rm{ SpO}}_2 $均值均由91.4%降至89.5%, 两种设定下闭环策略标准差明显降低, 说明$ {\rm{ SpO}}_2 $向高水平上升显著. 此外, 在$ c_o= 12 $均相同的前提下, $ \Delta c_{op}=-1.5 $时的适应性效果优于$ \Delta c_{op}=-1 $的情况, 因为更快的最优适应速率给予了闭环策略更宽泛的决策范围, 使得控制器能够采取氧浓度更低的训练方式, 保证个体在固定时间范围内获得更好的低氧适应性改善效果.
表 3 虚拟受试者3和4采取开环策略进行IHT和采取闭环策略进行IHT的效果对比Table 3 Comparison results of simulated subject 3 and 4 trained by using traditional open-loop strategy and proposed closed-loop strategy个体设定 指标 初始状态 开环策略 闭环策略 $ c_{o}=12 $,
$ \Delta c_{op}=-1 $DSI (s) 154 136 203 SpO2平均值(%) 91.4 89.5 94.7 SpO2标准差(%) 6.1 6.5 5.0 $ c_{o}=12 $,
$ \Delta c_{op}=-1.5 $DSI (s) 154 81 300 SpO2平均值 (%) 91.4 89.5 96.4 SpO2标准差(%) 6.1 6.6 4.1 进一步地, 为了对控制算法进行综合评价, 通过10名虚拟受试者做了整体的性能分析. 图6展示了10个不同虚拟受试者的训练效果对比结果, 对应的量化指标分别在表4和表5列出. 可以看出, 闭环和开环训练后的效果具有显著差异. 在初始情况下, DSI、$ {\rm{ SpO}}_2 $平均值和标准差分别为137.1 s、91.7%和4.97%, 开环情况下各指标数值分别为56.3 s、86.3%和4.95%, 较未训练前水平有明显降低. 然而, 在闭环策略作用下, DSI提升了114.6 s、$ {\rm{ SpO}}_2 $平均值增加了3.6%, 且标准差降低1.09%. 上述针对单个受试者和整体的分析结果均表明, 所提出的IHT低氧浓度优化决策算法对个体高原低氧适应性能力具有明显的改善作用, 说明了该算法用于预习服性能提升的可行性和有效性.
图 6 10名虚拟受试者采取开环和闭环策略进行IHT的$ {\rm{ SpO}}_2 $曲线Fig. 6 The $ {\rm{ SpO}}_2 $ curves of 10 simulated subjects that perform IHT based on traditional open-loop strategy and proposed closed-loop strategy表 4 10名虚拟受试者采取开环策略进行IHT和采取闭环策略进行IHT的效果对比Table 4 Comparison results of 10 simulated subjects trained by using traditional open-loop strategy and proposed closed-loop strategy受试者 指标 初始状态 开环策略 闭环策略 1 DSI (s) 137 32 300 SpO2平均值(%) 90.9 84.9 96.4 SpO2标准差(%) 5.46 5.37 3.53 2 DSI (s) 140 53 300 SpO2平均值(%) 92.7 86.1 95.6 SpO2标准差(%) 4.72 4.85 3.73 3 DSI (s) 138 31 300 SpO2平均值(%) 92.3 84.8 97.4 SpO2标准差(%) 4.81 4.70 2.88 4 DSI (s) 271 153 300 SpO2平均值(%) 94.6 89.9 98.5 SpO2标准差(%) 4.04 4.10 2.29 5 DSI (s) 138 35 184 SpO2平均值(%) 92.3 87.1 96.3 SpO2标准差(%) 4.32 4.33 3.51 6 DSI (s) 138 20 279 SpO2平均值(%) 92.3 84.8 94.1 SpO2标准差(%) 4.81 4.79 4.34 7 DSI (s) 42 35 78 SpO2平均值(%) 89.0 85.1 91.5 SpO2标准差(%) 5.44 5.52 5.19 8 DSI (s) 91 84 176 SpO2平均值(%) 89.9 87.6 92.1 SpO2标准差(%) 6.17 6.16 5.63 9 DSI (s) 138 66 300 SpO2平均值(%) 91.3 88.2 95.1 SpO2标准差(%) 5.01 4.86 4.02 10 DSI (s) 138 44 300 SpO2平均值(%) 91.3 84.8 96.0 SpO2标准差(%) 4.94 4.78 3.64 表 5 10名虚拟受试者的对比结果Table 5 Comparison results of 10 simulated subjects指标 初始状态 开环策略 闭环策略 DSI (s) 137.1 56.3 251.7 SpO2平均值(%) 91.7 86.3 95.3 SpO2标准差(%) 4.97 4.95 3.88 4. 结论
本文设计了一种面向高原适应力增强的IHT供低氧浓度闭环优化决策算法. 利用各低氧段少量样本数据建立了高斯过程预测模型, 并将其应用到具有安全约束设计的高低风险不对称代价函数中, 保证在不影响个体生命健康的同时实现快速预习服. 考虑到优化问题中具有随机变量, 采用贝叶斯优化算法迭代求解每一日IHT中的最优供氧浓度. 此外, 针对现有仿真器无法描述个体不同训练前后动态变化过程的问题, 对仿真器进行了修正. 引入了调试界限设置以及个体训练累加效应, 并结合“最优速率理论”合理设计了模型变化律. 仿真验证结果表明, 所设计的闭环干预算法相较于开环策略而言具有较明显的改善效果, 说明了其用于加速高原适应能力的可行性和有效性. 未来的工作将进一步探讨针对供低氧浓度和低氧持续时长的双变量优化决策问题, 并对算法进行更广泛的验证.
-
图 1 高原适应性能力提升的IHT策略优化决策算法流程图
Fig. 1 Flow chart of IHT optimization decision algorithm for high-altitude adaptability improvement
图 3 所设计代价函数$ {\cal{L}}_{v} $部分的惩罚强度在不同$ \Delta c $下的变化
Fig. 3 The penalty changes of designed $ {\cal{L}}_v $ term under different $ \Delta c $ values
图 4 虚拟受试者1和2采取开环和闭环策略进行IHT的$ {\rm{ SpO}}_2 $曲线
Fig. 4 The $ {\rm{ SpO}}_2 $ curves of simulated subject 1 and 2 that perform IHT based on traditional open-loop strategy and proposed closed-loop strategy
图 5 虚拟受试者3和4采取开环和闭环策略进行IHT的$ {\rm{ SpO}}_2 $曲线
Fig. 5 The $ {\rm{ SpO}}_2 $ curves of simulated subject 3 and 4 that perform IHT based on traditional open-loop strategy and proposed closed-loop strategy
图 6 10名虚拟受试者采取开环和闭环策略进行IHT的$ {\rm{ SpO}}_2 $曲线
Fig. 6 The $ {\rm{ SpO}}_2 $ curves of 10 simulated subjects that perform IHT based on traditional open-loop strategy and proposed closed-loop strategy
表 1 相关参数取值
Table 1 Related parameters
参数 含义 取值 $ n_s $ 一日IHT的低氧段总数 8 $ n_t $ 一段低氧段预测点总数 200 $ n_p $ 一段低氧段采样点总数 300 $ \sigma^{2} $ 平方指数核函数超参数 10 $ \ell $ 平方指数核函数超参数 10 $ \boldsymbol{y}_{r} $ 目标$ {\rm{ SpO}}_2 $向量 $ [95H_{50},90H_{50}, 85H_{100}]^\mathrm{T} $ $ Q $ 代价函数惩罚矩阵 136I $ w $ 代价函数权重系数 1000 $ c_{l,\min} $ 供氧浓度下阈值 (%) 10 $ c_{l,\max} $ 供氧浓度上阈值 (%) 15 
下载: 导出CSV
表 2 虚拟受试者1和2采取开环策略进行IHT和采取闭环策略进行IHT的效果对比
Table 2 Comparison results of simulated subject 1 and 2 trained by using traditional open-loop strategy and proposed closed-loop strategy
个体设定 指标 初始状态 开环策略 闭环策略 $ c_{o}=13 $,
$ \Delta c_{op}=-1.5 $DSI (s) 137 84 256 SpO2平均值(%) 90.9 84.9 93.6 SpO2标准差(%) 5.5 5.4 4.8 $ c_{o}=12 $,
$ \Delta c_{op}=-1.5 $DSI (s) 137 90 300 SpO2平均值(%) 90.9 87.8 96.2 SpO2标准差(%) 5.5 5.4 3.6
下载: 导出CSV
表 3 虚拟受试者3和4采取开环策略进行IHT和采取闭环策略进行IHT的效果对比
Table 3 Comparison results of simulated subject 3 and 4 trained by using traditional open-loop strategy and proposed closed-loop strategy
个体设定 指标 初始状态 开环策略 闭环策略 $ c_{o}=12 $,
$ \Delta c_{op}=-1 $DSI (s) 154 136 203 SpO2平均值(%) 91.4 89.5 94.7 SpO2标准差(%) 6.1 6.5 5.0 $ c_{o}=12 $,
$ \Delta c_{op}=-1.5 $DSI (s) 154 81 300 SpO2平均值 (%) 91.4 89.5 96.4 SpO2标准差(%) 6.1 6.6 4.1
下载: 导出CSV
表 4 10名虚拟受试者采取开环策略进行IHT和采取闭环策略进行IHT的效果对比
Table 4 Comparison results of 10 simulated subjects trained by using traditional open-loop strategy and proposed closed-loop strategy
受试者 指标 初始状态 开环策略 闭环策略 1 DSI (s) 137 32 300 SpO2平均值(%) 90.9 84.9 96.4 SpO2标准差(%) 5.46 5.37 3.53 2 DSI (s) 140 53 300 SpO2平均值(%) 92.7 86.1 95.6 SpO2标准差(%) 4.72 4.85 3.73 3 DSI (s) 138 31 300 SpO2平均值(%) 92.3 84.8 97.4 SpO2标准差(%) 4.81 4.70 2.88 4 DSI (s) 271 153 300 SpO2平均值(%) 94.6 89.9 98.5 SpO2标准差(%) 4.04 4.10 2.29 5 DSI (s) 138 35 184 SpO2平均值(%) 92.3 87.1 96.3 SpO2标准差(%) 4.32 4.33 3.51 6 DSI (s) 138 20 279 SpO2平均值(%) 92.3 84.8 94.1 SpO2标准差(%) 4.81 4.79 4.34 7 DSI (s) 42 35 78 SpO2平均值(%) 89.0 85.1 91.5 SpO2标准差(%) 5.44 5.52 5.19 8 DSI (s) 91 84 176 SpO2平均值(%) 89.9 87.6 92.1 SpO2标准差(%) 6.17 6.16 5.63 9 DSI (s) 138 66 300 SpO2平均值(%) 91.3 88.2 95.1 SpO2标准差(%) 5.01 4.86 4.02 10 DSI (s) 138 44 300 SpO2平均值(%) 91.3 84.8 96.0 SpO2标准差(%) 4.94 4.78 3.64
下载: 导出CSV
表 5 10名虚拟受试者的对比结果
Table 5 Comparison results of 10 simulated subjects
指标 初始状态 开环策略 闭环策略 DSI (s) 137.1 56.3 251.7 SpO2平均值(%) 91.7 86.3 95.3 SpO2标准差(%) 4.97 4.95 3.88
下载: 导出CSV
-
[1] Andrew M L, Peter H H. Medical conditions and high-altitude travel. New England Journal of Medicine, 2022, 386(4): 364-373 doi: 10.1056/NEJMra2104829 [2] Joshua C T, Philip N A. Global and country-level estimates of human population at high altitude. Proceedings of the National Sciences, 2021, 118(18): e2102463118 doi: 10.1073/pnas.2102463118 [3] Cobb A B, Levett D Z H, Mitchell K, Aveling W, Hurlbut D, Gilbert-Kawai E, et.al. Physiological responses during ascent to high altitude and the incidence of acute mountain sickness. Physiological reports, 2021, 9(7): e14809 [4] Gudbjartsson T, Sigurdsson E, Gottfredsson M, Bjornsson O M, Gudmundsson G. High altitude illness and related diseases - A review. Laeknabladid, 2019, 105(11): 499-507 [5] Victor S, Jan C P, and Katarína K. Manifestation of intracranial lesions at high altitude: Case report and review of the literature. High Altitude Medicine & Biology, 2021, 22(1): 87-89 [6] Fulco C S, Beidleman B A, Muza S R. Effectiveness of preacclimatization strategies for highaltitude exposure. Exercise and Sport Sciences Reviews, 2013, 41(1): 55-63 doi: 10.1097/JES.0b013e31825eaa33 [7] Ambroży T, Maciejczyk M, Klimek A T, Wiecha S, Stanula A, Snopkowski P, et.al. The effects of intermittent hypoxic training on anaerobic and aerobic power in boxers. International Journal of Environmental Research and Public Health, 2020, 17(24): 9361 doi: 10.3390/ijerph17249361 [8] Wille M, Gatterer H, Mairer K, Philippe M, Schwarzenbacher H, Faulhaber M, et.al. Short-term intermittent hypoxia reduces the severity of acute mountain sickness. Medicine & Science in Sports, 2012, 22(5): e79-e85 [9] 刘园园. 高原健康理论框架下的渐进型间歇性低氧预习服训练研究 [博士学位论文], 山东大学, 中国, 2014Liu Yuan-Yuan. Short-Term Intermittent Hypoxia Reduces the Severity of Acute Mountain Sickness [Ph.D. dissertation], Shandong University, China, 2014 [10] Treml B, Kleinsasser A, Hell T, Knotzer H, Wille M, Burtscher M. Carry-over quality of pre-acclimatization to altitude elicited by intermittent hypoxia: A participant-blinded, randomized controlled trial on antedated acclimatization to altitude. Frontiers in Physiology, DOI: 10.3389/fphys.2020.00531 [11] Gangwar A, Pooja, Sharma M, Singh K, Patyal A, Bhaumik G, et.al. Intermittent normobaric hypoxia facilitates high altitude acclimatization by curtailing hypoxia-induced infammation and dyslipidemia. Pflugers Archiv, 2019, 471(7):949-959 doi: 10.1007/s00424-019-02273-4 [12] 杨军, 俞梦孙, 曹征涛, 吴峰, 张宏金, 王海涛, 等.间歇性递增式常压低氧暴露训练对高原习服效果的研究. 中华航空航天医学杂志, 2012, 3: 161-164Yang Jun, Yu Meng-Sun, Cao Zheng-Tao, Wu Feng, Zhang Hong-Jin, Wang Hai-Tao, et.al. Study on the effect of increasing intermittent hypoxia exposure on altitude acclimatization. Chinese Journal of Aerospace Medicine, 2012, 3: 161-164 [13] Kwiatkowska M, Atkins M S, Ayas N T, Ryan C F. Knowledge-based data analysis: First step toward the creation of clinical prediction rules using a new typicality measure. IEEE Transactions on Information Technology in Biomedicine, 2007, 11(6):651-660 doi: 10.1109/TITB.2006.889693 [14] Sakellarios A I, Räber L, Bourantas C V, Exarchos T P, Athanasiou L S, Pelosi G, et.al. Prediction of atherosclerotic plaque development in an In Vivo coronary arterial segment based on a multilevel modeling approach. IEEE Transactions on Biomedical Engineering, 2017, 64(8):1721-1730 doi: 10.1109/TBME.2016.2619489 [15] 喻勇, 司小胜, 胡昌华, 崔忠马, 李洪鹏. 数据驱动的可靠性评估与寿命预测研究进展:基于协变量的方法. 自动化学报, 2018, 44(2): 216-227Yu Yong, Si Xiao-Sheng, Hu Chang-Hua, Cui Zhong-Ma, Li Hong-Peng. Data driven reliability assessment and life-time prognostics: A review on covariate models. Acta Automatica Sinica, 2018, 44(2): 216-227 [16] 李天梅, 司小胜, 刘翔, 裴洪. 大数据下数模联动的随机退化设备剩余寿命预测技术. 自动化学报, 2022, 48(9): 2119-2141 doi: 10.16383/j.aas.c201068Li Tian-Mei, Si Xiao-Sheng, Liu Xiang, Pei Hong. Data-model interactive remaining useful life prediction technologies for stochastic degrading devices with big data. Acta Automatica Sinica, 2022, 48(9): 2119-2141 doi: 10.16383/j.aas.c201068 [17] 蒋珂, 蒋朝辉, 谢永芳, 潘冬, 桂卫华. 基于动态注意力深度迁移网络的高炉铁水硅含量在线预测方法. 自动化学报, DOI: 10.16383/j.aas.c210524"> 10.16383/j.aas.c210524 Jiang Ke, Jiang Zhao-Hui, Xie Yong-Fang, Pan Dong, Gui Wei-Hua. Online prediction method for silicon content of molten iron in blast furnace based on dynamic attention deep transfer network. Acta Automatica Sinica, DOI:10.16383/j.aas.c210524"> 10.16383/j.aas.c210524 [18] Box G E, Jenkins G M, Reinsel G C, Ljung G M. Time Series Analysis: Forecasting and Control. Hoboken: John Wiley & Sons, 2015. [19] Xie J, Wang Q. Benchmarking machine learning algorithms on blood glucose prediction for type I diabetes in comparison with classical time-series models. IEEE Transactions on Biomedical Engineering, 2020, 67(11): 3101-3124 doi: 10.1109/TBME.2020.2975959 [20] Moniri A, Terracina D, Rodriguez-Manzano J, Strutton P H, Georgiou P. Real-time forecasting of sEMG features for trunk muscle fatigue using machine learning. IEEE Transactions on Biomedical Engineering, 2021, 68(2): 718-727 doi: 10.1109/TBME.2020.3012783 [21] Michalis K T. Variational learning of inducing variables in sparse Gaussian processes. In: Proceedings of the 12th International Conference on Artificial Intelligence and Statistics. Florida, USA: PMLR, 2009. 567−574 [22] Beckers T, Hirche S. Prediction with approximated gaussian process dynamical models. IEEE Transactions on Automatic Control, 2022, 68: 6460-6473 [23] Lee S I, Mortazavi B, Hoffman H A, Lu D S, Li C, Paak B H, et.al. A prediction model for functional outcomes in spinal cord disorder patients using gaussian process regression. IEEE Transactions on Biomedical Engineering, 2016, 20(1): 91-99 [24] Huang H, Song Y, Peng X, Ding S X, Zhong W, Du W, et.al. A sparse nonstationary trigonometric gaussian process regression and its application on nitrogen oxide prediction of the diesel engine. IEEE Transactions on Industrial Informatics, 2021, 17(12): 8367-8377 doi: 10.1109/TII.2021.3068288 [25] 史大威, 蔡德恒, 刘蔚, 王军政, 纪立农. 面向智能血糖管理的餐前胰岛素剂量贝叶斯学习优化方法. 自动化学报, DOI: 10.16383/j.aas.c210067Shi Da-Wei, Cai De-Heng, Liu Wei, Wang Jun-Zheng, Ji Li-Nong. Bayesian learning based optimization of meal bolus dosage for intelligent glucose management. Acta Automatica Sinica, DOI: 10.16383/j.aas.c210067 [26] 金哲豪, 刘安东, 俞立. 基于GPR和深度强化学习的分层人机协作控制. 自动化学报, 2022, 48(9): 1-11Jin Zhe-Hao, Liu An-Dong, Yu Li. Hierarchical human-robot cooperative control based on GPR and DRL. Acta Automatica Sinica, 2022, 48(9): 1-11 [27] Rosolia U, Zhang X, Borrelli F. Data-driven predictive control for autonomous systems. Annual Review of Control, Robotics, and Autonomous Systems, 2018, 1(1): 259-286 doi: 10.1146/annurev-control-060117-105215 [28] Yu M. Human-performance engineering at high altitude. Science Supp, 2014: 7−8 [29] Chen J, Xiao R, Wang L, Zhu L, Shi D. Unveiling interpretable key performance indicators in hypoxic response: a system identification approach. IEEE Transactions on Industrial Electronics, 2022, 69(12): 13676-13685 doi: 10.1109/TIE.2021.3137618 [30] Chen J, Tian Y, Zhang G, Cao Z, Zhu L, Shi D. IoT-enabled intelligent dynamic risk assessment of acute mountain sickness: The role of event-triggered signal processing. IEEE Transactions on Industrial Informatics, 2023, 19(1): 730−738 [31] Williams C K I, Rasmussen C E. Gaussian Processes for Machine Learning. Cambridge: The MIT Press, 2006. [32] Hackett, Peter H. and Roach, Robert C. High-altitude llness. New England Journal of Medicine, 2001, 345(2): 107-114 doi: 10.1056/NEJM200107123450206 [33] Levine B D, Stray-Gundersen J. Dose-response of altitude training: how much altitude is enough? Advances in Experimental Medicine and Biology, 2006, 69: 233-247 期刊类型引用(2)
1. 孙鹏翔，蔡勇元，蔡承希. 基于MWBO-XGBoost-ERT的森林火灾火焰辐射功率预测. 计算机时代. 2025(02): 24-30 . 
百度学术2. 安泉旭，邵林海. 间歇性低氧训练在运动训练中的应用研究. 拳击与格斗. 2024(20): 31-33 . 百度学术其他类型引用(0)
-
下载:
计量
- 文章访问数: 791
- HTML全文浏览量: 162
- PDF下载量: 186
- 被引次数: 2
