Bayesian Learning Based Optimization of Meal Bolus Dosage for Intelligent Glucose Management
-
摘要: 餐前胰岛素剂量精准决策是改善糖尿病患者血糖管理的关键. 临床治疗中胰岛素剂量调整一般在较短时间内完成, 具有典型的小样本特征; 数据驱动建模在该情形下无法准确学习患者餐后血糖代谢规律, 难以确保胰岛素剂量的安全和有效决策. 针对这一问题, 设计一种临床经验辅助的餐前胰岛素剂量自适应优化决策框架, 构建高斯过程血糖预测模型和模型有效性在线评估机制, 提出基于历史剂量和临床经验决策约束的贝叶斯优化方法, 实现小样本下餐后血糖轨迹的安全预测和餐前胰岛素注射剂量的优化决策. 该方法的安全性和有效性通过美国食品药品监督管理局接受的UVA/Padova T1DM平台测试结果和1型糖尿病患者实际临床数据决策结果充分验证. 可为餐前胰岛素剂量智能决策及临床试验提供方法基础和技术支持, 也为中国糖尿病患者血糖管理水平的有效改善, 提供了精准医学治疗手段.Abstract: Precise decision of preprandial insulin bolus is of crucial importance to achieving enhanced glucose management for patients with diabetes. Insulin dosage adjustment in current clinical practice normally needs to be completed in a relatively short period of time and thus the data size is typically small. It would be consequently challenging to accurately learn the postprandial glucose metabolism through data-driven modelling techniques and to further ensure the safe and efficient insulin dosage decision. For this problem, a clinical-experience-assisted adaptive preprandial insulin bolus decision framework is proposed in this work. To achieve safe prediction of the postprandial glucose traces and optimal decision of the preprandial insulin dosage with limited training data, a Gaussian process based glucose prediction model and an online model efficiency assessment mechanism are constructed, and a Bayesian optimization method with historical data exploitation and clinical-experience guided decision constraints is proposed. The safety and effectiveness of the proposed method are extensively validated through in silico results of the US Food and Drug Administration accepted UVA/Padova T1DM simulator and advisory mode analysis based on clinical data from a subject with type 1 diabetes. The obtained results provide methodological and technical support for intelligent meal bolus decision and the forthcoming clinical studies, and introduce a precision medicine solution to effectively improved glucose management for Chinese patients with diabetes mellitus.
-
糖尿病是一种影响人类健康的严重慢性疾病, 主要包括1型(胰岛素完全不分泌)和2型(胰岛素分泌相对不足). 长期良好的血糖控制是避免出现糖尿病慢性并发症(如糖尿病视网膜病变、糖尿病、肾病)的关键, 同时也是降低心血管死亡风险的重要因素. 中国拥有全世界最多的糖尿病患者, 在预防和控制糖尿病方面面临巨大挑战[1-3]. 因而, 研究开发适合中国糖尿病患者的先进治疗技术具有重要意义.
胰岛素治疗是控制血糖达标的重要治疗手段, 其中每日多次胰岛素皮下注射和胰岛素泵治疗是实施胰岛素强化治疗的两种主要策略[4]. 随着血糖监测等硬件设备的发展, 人工胰腺系统成为一种新型的治疗方式, 它是一类闭环给药智能系统[5-7], 可全天候控制患者血糖. 然而, 在上述治疗方法中, 餐后血糖管理是关键步骤[8]. 现阶段主要采用前馈控制方法来控制餐后血糖水平, 即患者在就餐前提前输注大剂量胰岛素, 以控制餐后血糖的波动.
目前, 餐前胰岛素剂量主要根据标准计算方法确定, 其常用于患者的临床治疗中[9-10]. 具体地, 标准计算方法形式如下:
$$u_{s} = \frac{\text{CHO}}{\text{CR}}+\frac{G_{c}-G_{sp}}{\text{CF}}-u_{\text{IOB}}$$ (1) 式中, CHO为摄入膳食的碳水化合物含量; $ G_{c} $为餐前血糖水平; $ G_{sp} $为血糖控制目标; $ u_{\text{IOB}} $为体内残余活性胰岛素(Insulin on board, IOB)量; CR和CF为个体化生理参数, 分别为碳水化合物系数和胰岛素敏感系数. 现已有不少研究聚焦于提升餐前胰岛素剂量决策性能, 例如文献[11]利用血糖测量信息, 通过计算餐后血糖最小值与目标值之间的差值来调整CR和CF系数, 改善餐后血糖管理. 文献[12]提出一种基于胰岛素灵敏度指标优化CR系数的方法, 该指标利用患者的血糖监测和胰岛素输注数据估计得到. 文献[13]基于患者历史血糖管理数据, 利用案例推理(Case based reasoning, CBR)方法调整CR和CF系数. 另外, 也有一些研究致力于改善人工胰腺中的餐后血糖控制, 例如文献[14]结合伯格曼最小模型和无迹卡尔曼滤波自动确定餐前胰岛素剂量, 无需患者参与. 文献[15]结合CBR和R2R (Run-to-run)方法自适应确定餐前胰岛素剂量. 文献[16]提出一种贝叶斯优化辅助学习方法, 可利用患者历史血糖监测和胰岛素输注信息调整CR系数, 提升餐后血糖控制水平.
基于以上分析发现, 现有方法大多通过调整CR和CF系数来改善餐后血糖管理. 这种方式不能有效地将餐后血糖动力学考虑进餐前剂量决策. 随着血糖监测设备的普及和提升, 数据驱动控制方法有望为餐前剂量决策提供一种崭新的方式. 数据驱动控制无需系统模型知识, 仅利用系统的输入/输出数据, 即可完成建模和控制[17], 近年来被广泛应用于复杂的工程领域. 文献[18]采用神经网络, 提出数据驱动的赤铁矿磨矿过程运行优化控制方法. 文献[19]针对机理建模预测和控制性能差的问题, 建立基于多元时间序列的高炉铁水硅含量数据驱动预测模型. 文献[20]提出了一种基于迭代学习模型预测控制(Iterative learning model predictive control, ILMPC)的闭环血糖控制方法, 该方法利用患者历史数据迭代更新控制目标, 改善血糖控制, 并通过临床试验[21-22]验证了算法的有效性与安全性. 与该方法所研究的实时闭环血糖控制问题不同, 本文重点关注有限血糖和胰岛素数据训练样本情形下单次餐前胰岛素剂量决策问题, 旨在通过单次剂量决策改善餐后血糖控制; 在研究方法方面, 与文献[20-22]采用的ILMPC方法不同, 本文基于患者有限历史数据, 建立非参数化的餐后血糖高斯过程预测模型, 利用风险敏感(Risk sensitive, RS)优化框架, 完成考虑血糖预测不确定性的胰岛素剂量鲁棒优化决策, 并通过引入专家经验辅助决策, 确保学习优化过程的安全性, 以达到改善餐后血糖管理水平的目的. 文献[23]结合高斯过程和风险敏感控制, 提出一种数据驱动餐前胰岛素剂量决策方法. 然而该方法采用离线建模的方式, 所需样本数量大, 不适用于患者短期住院或小样本情况下餐前剂量的决策. 在小样本情况下, 数据驱动建模无法有效地学习患者血糖代谢规律, 难以确保餐前剂量的安全和有效决策.
针对以上问题, 本文提出临床经验辅助的餐前剂量贝叶斯学习优化方法, 通过评估数据驱动建立预测模型的有效性, 决定是否实施临床经验决策辅助, 避免小样本情况下贝叶斯优化决策失误. 同时, 实时收集样本, 更新预测模型, 形成餐前剂量自适应优化决策框架. 其中考虑到中国现有糖尿病患者大多采用固定的餐前胰岛素剂量[24], 缺乏系统的决策规则, 本文结合医生临床决策经验, 建立基于患者历史血糖管理数据的系统性临床经验决策规则. 另外, 为了在小样本情况下, 快速确定合适的餐前剂量, 利用迭代更新的思想, 设计包含上次剂量惩罚项的代价函数, 约束贝叶斯优化基于上次剂量进行求解. 本文使用美国食品药品监督管理局接受的UVA/Padova T1DM血糖代谢仿真器[25]对所设计的方法进行仿真验证, 同时利用1型糖尿病患者临床治疗数据, 基于顾问模式分析方法[26]评估所设计方法的决策性能. 实验结果表明, 本文设计的方法可为餐前胰岛素剂量智能决策及临床试验提供方法基础和技术支持, 也为中国糖尿病患者血糖管理水平的有效改善, 提供了精准医学治疗手段.
1. 基于高斯过程的餐前剂量优化设计
首先, 利用高斯过程建立餐后血糖预测模型; 然后, 利用预测的餐后血糖设计惩罚不对称风险敏感代价函数, 同时设计包含上次剂量的惩罚项; 最后, 形成决策餐前剂量的约束优化问题. 设计方法的总体框架如图1所示.
图 1 餐前胰岛素剂量贝叶斯学习优化方法框架图Fig. 1 Block diagram of the Bayesian learning based optimization method for meal bolus1.1 高斯过程
高斯过程是在机器学习领域发展形成的概率化建模方法[27], 与传统建立确定的参数化模型的方法不同, 其假设系统在任意的输入$ {\boldsymbol x}\in{\bf R}^n $下, 输出$ f({\boldsymbol x})\in{\bf R} $服从联合高斯分布, 记为:
$$ f({\boldsymbol x})\sim \mathcal{GP}(m({\boldsymbol x}), k({\boldsymbol x}, {\boldsymbol x'})) $$ (2) 式中, $ m({\boldsymbol x}) $和$ k({\boldsymbol x}, {\boldsymbol x'}) $分别表示均值函数和协方差函数, 可根据对系统的先验知识进行设计, 其中待定的参数称为超参数. 在实际观测中, 系统的输出可能会受到噪声污染, 假设噪声变量服从$ \omega\sim\mathcal{N}(0,\sigma_\omega^2) $且观测方程为$ y = f({\boldsymbol x})+\omega $, 则观测输出$ y $仍是一个高斯过程:
$$ y\sim \mathcal{GP}(m({\boldsymbol x}), k({\boldsymbol x}, {\boldsymbol x'})+\sigma_\omega^2\delta({\boldsymbol x},{\boldsymbol x'})) $$ (3) 式中, $ \delta({\boldsymbol x},{\boldsymbol x'}) $为克罗内克函数, 当且仅当 $ {\boldsymbol x} = {\boldsymbol x'} $时, 为1; 其他, 为0. 假定给定$ N $个训练输入样本${\boldsymbol X} = [{\boldsymbol x}_1, \;{\boldsymbol x}_2, \cdots , {\boldsymbol x}_N]^{\rm{T}}$对应的输出观测样本 ${ \boldsymbol Y} = [y_1, y_2, \cdots, y_N]^{\rm{T}}$, $ y $的后验超参数可通过最大似然法确定[27].
确定超参数后的高斯过程可预测新的输入$ \boldsymbol x_* $对应的观测值$ y_* $. 由联合高斯分布的性质可知, $ y_* $服从高斯分布$ y_*\sim\mathcal{N}(m_*,\sigma_*^2) $, 满足:
$$ m_*= m(\boldsymbol x_*)+\boldsymbol K_*^{\rm{T}}(\boldsymbol K+\sigma_\omega^2 \boldsymbol I)^{-1}(\boldsymbol Y-m(\boldsymbol X)) $$ (4) $$ \sigma_*^2 = k(\boldsymbol x_*,\boldsymbol x_*)-\boldsymbol K_*^{\rm{T}}(\boldsymbol K+\sigma_\omega^2 \boldsymbol I)^{-1}\boldsymbol K_* $$ (5) 式中, $\boldsymbol K_* = [k(\boldsymbol x_*,\boldsymbol x_1),\cdots,k(\boldsymbol x_*,\boldsymbol x_N)]^{\rm{T}}$和$ \boldsymbol K $是协方差矩阵, 其元素$ K_{ij} = k(\boldsymbol x_i,\boldsymbol x_j) $, $\boldsymbol I $为单位矩阵. 可见, 对于一个新的输入, 高斯过程可给出一个确定的预测值(式(4))和对该预测值不确定性的评估, 以方差进行衡量(式(5)).
1.2 餐后血糖预测
现有研究大多采用具有生理学意义的房室模型描述人体血糖代谢过程[25, 28]. 这些模型由多个耦合的微分方程组成, 变量繁多且包含大量个体化参数, 有时仅利用血糖测量数据难以完成参数辨识. 另外, 考虑到人体血糖变化的连续性, 在已知过去的血糖值、胰岛素输注剂量和碳水化合物摄入量的条件下, 当前血糖值必定在一定合理范围内变化, 因此在上述信息作为输入向量条件下, 可利用高斯分布去近似刻画当前血糖的变化. 基于以上分析, 本文利用高斯过程对餐后血糖动力学进行建模. 特别地, 通过自回归输入/输出信号作为模型的回归变量, 高斯过程被广泛应用于非线性控制系统的建模[29-30].
为此, 本文利用自回归模型描述餐后血糖动力学, 每一时刻$t+i,\; i = \{1,2,\cdots,n\}$都由带有加性噪声的非线性函数表示:
$$ \left\{\begin{aligned} &P_{t+i} = f_{t+i-1}({\boldsymbol z}_{t+i-1})+\omega_{t+i-1}\\ &{\boldsymbol z}_{t+i-1} = [P_{t+i-1-l},\cdots,P_{t+i-1},u,d]^{\rm{T}} \end{aligned}\right.$$ (6) 式中, $ t $表示摄入膳食的时刻, $ \omega $表示高斯白噪声, $ P $为血糖测量值, $ l $表示自回归模型的阶数, $ u $为餐前剂量, $ d $为摄入的碳水化合物CHO含量. 为了尽可能地包含血糖态势信息, 本文考虑$ l = 7 $和采样时间$ T = 15 $min, 时间长度为2 h. 另外, 考虑$n = 8$, 时间长度为2 h的餐后血糖轨迹.
高斯过程的应用分为建模和预测两个阶段. 在建模阶段, 基于模型(6), 为每种餐型构建预测模型, 餐型包括早餐、午餐和晚餐. 具体地, 针对某一餐型, 分别利用高斯过程对餐后每一时刻的血糖动力学进行建模. 各时刻的建模方式类似, 以$ t+1 $时刻为例, 收集${\boldsymbol z}_t = [P_{t-7},\cdots,P_t,u,d]^{\rm{T}}$作为训练输入, 其中血糖测量值利用最小/最大归一化法, 进行归一化, 同时为减小预测不确定性, 收集差值$\Delta P_t = P_{t+1}-P_t$作为对应的训练输出[29]. 均值函数考虑为线性形式, 协方差函数采用常用的平方指数形式:
$$ m({\boldsymbol z}_t) = \boldsymbol a^{\rm{T}} {\boldsymbol z}_t+b$$ (7) $$ \begin{split} k_{\Delta P_t}({\boldsymbol z}_t, {\boldsymbol z}_t') =\;& \sigma_f^2\exp \left( -\frac{1}{2}({\boldsymbol z}_t - {\boldsymbol z}_t')^{\rm{T}}{\boldsymbol \Omega}^{-1}({\boldsymbol z}_t - {\boldsymbol z}_t') \right)\;+\\ &\delta({\boldsymbol z}_t, {\boldsymbol z}_t')\sigma_\omega^2 \\[-10pt] \end{split} $$ (8) 式中, ${\boldsymbol \Omega} = \text{diag}\{l_1^2, l_2^2,\cdots,l_{10}^2\}$, $ \sigma_f^2 $为信号方差, $ l_j $为尺度参数. $\theta = \{\boldsymbol a, b, {\boldsymbol \Omega}, \sigma_f^2, \sigma_\omega^2\}$为超参数集合.
在预测阶段, 给定摄入的CHO含量和注射的胰岛素剂量, 根据餐前血糖水平, 利用上阶段训练得到的高斯过程, 逐时刻预测餐后血糖水平, 得到餐后2 h (8步)的血糖变化轨迹概率分布. 其中预测输入中的餐后血糖值采用血糖预测值. 由于建模阶段, 训练输出采用差值的形式, 因此在得到原始的预测值后, 需要将此值与上一时刻的血糖测量(预测)值相加, 确定为当前时刻的血糖预测值, 预测不确定性(方差)保持不变. 考虑到糖尿病患者一般具有规律的饮食习惯(如相似的用餐时间和膳食摄入量), 因此, 可将膳食摄入对餐后血糖变化的影响, 近似成为一个常值扰动, 利用高斯过程鲁棒性建模的特点, 建立无需膳食信息的餐后血糖预测模型, 进一步减轻患者估算CHO摄入量的负担, 实现无需膳食信息的剂量优化决策.
1.3 不对称风险敏感代价函数
血糖控制中存在高低血糖不对称风险; 与长时间高血糖所引发的并发症不同, 即时的严重低血糖会直接造成患者死亡. 因此在决策剂量纠正患者餐后高血糖的同时, 需要格外小心低血糖事件的发生. 现有研究中一个可行的解决方案是构建不对称代价函数[31-32]. 受文献[32]启发, 本文对高于和低于血糖目标值的偏离量实施不对称惩罚; 此外, 考虑到预测信息存在的不确定性, 本文以风险敏感的形式来设计代价函数. 记高斯过程提供的8步血糖预测量分别为$ G_i\sim\mathcal{N}(m_i,\sigma_i^2) $, $i\in \{1,2,\cdots,8\}$, 并组成一个向量状态$\boldsymbol G = [G_1, G_2, \cdots , G_8]^{\rm{T}}$, 该状态描述了餐后2 h血糖变化轨迹的概率分布. 进一步地, 根据风险敏感代价函数的定义[33], 设计不对称风险敏感代价函数如下:
$$ \mathcal{L}_{A} = -\frac{2}{\gamma}\ln {\rm{E}} \left[\exp\left(-\frac{\gamma}{2}(\boldsymbol G- \boldsymbol G_{ r})^{\rm{T}} \boldsymbol Q(\boldsymbol G- \boldsymbol G_{ r})\right)\right] $$ (9) 式中, ${\rm{E}}$为数学期望; $\boldsymbol G_{ r}$为餐后血糖控制目标; $ \boldsymbol Q $为惩罚矩阵, 其为对角形式, 用于对预测的高低血糖值施加不对称惩罚权值; $ \gamma<0 $表示风险敏感系数, 决定着优化者对预测不确定性的态度[34-35], $ \gamma< 0 $使得优化者倾向于选择减少预测不确定性的控制量. 基于上述不对称风险敏感代价函数, 优化者能够充分利用数据样本中的决策经验, 同时形成自我决策能力.
对于惩罚矩阵$ \boldsymbol Q $, 记$ \boldsymbol Q $的第$ i $个对角元素为$ Q_i $, 设计不对称惩罚如下:
$$ \left\{\begin{aligned} & Q_i = \left\{\begin{aligned} &Q_b^i\left(\frac{c_1}{1+\exp\{\alpha (\beta-|G_i-G_{ri}|)\}}+c_2\right),\\ &\qquad\qquad{}e_i > 0\\ &Q_b^i,\,\qquad {\text{其他}} \end{aligned}\right.\\ &e_i = \mathbf{1}(G_i-G_{ri}<0) \end{aligned}\right. $$ (10) 式中, $ \mathbf{1}(\cdot) $为示性函数, $ G_{ri} $为$ \boldsymbol G_{r} $的第$ i $个元素, $Q_b^i > 0$为基准惩罚系数, $ \Gamma : = \{\alpha,\beta,c_1,c_2\} $为决定惩罚强度的元素集合. 记$h(G_i,\ \Gamma) = {c_1}/(1 + \exp\{\alpha(\beta - (G_{ri}\;- G_i))\}+ c_2)$, 可得在$ G_i-G_{ri}<0 $时, $ h(G_i,\ \Gamma) $的变化曲线如图2所示. 由图2可知, 惩罚强度的上/下界分别由$ c_1+c_2 $和$ c_2 $决定, 中间值由$ \beta $决定, 增长速度由$ \alpha $控制. 总之, 在高于血糖目标值时, 设计常值的惩罚系数$ Q_b^i $可使得优化者对高血糖保持相对合理的响应; 在低于血糖目标值时, 设计随偏离量增大而指数式增大的惩罚系数可使得优化者在偏离血糖目标值不多时保持保守的响应, 但同时能及时快速地对大偏离量做出响应, 保证安全性.
根据上述设计原则, 本文利用美国食品药品监督管理局接受的血糖代谢仿真器确定参数集合$ \Gamma $的数值. 本文方法的参数确定选用仿真器中的“平均成人患者”完成, 性能验证利用10个虚拟成人患者完成. 在参数确定过程中, 首先选取平均成人患者进行8 h的仿真, 模拟仿真从9:00开始, 假定患者在12:00进食70 g碳水化合物, 观察不同餐前胰岛素剂量下患者餐后血糖的变化, 得到发生低血糖事件与餐后血糖偏离控制目标$\boldsymbol G_{ r}$程度之间的大致关系. 基于此关系, 首先根据图2选取确定合适的$ \beta $和$ c_2 $, 使得预测血糖在控制目标附近时, 惩罚程度大致不变; 然后, 设计$ \alpha $和$ c_1 $确定惩罚增长速度和最大惩罚强度, 使得预测血糖偏离控制目标的程度超过$ \beta $时, 惩罚急速增大, 对预测低血糖高敏感性; 最后, 利用所得的参数组进行不同含量碳水化合物进食测试, 确保参数组安全和有效.
1.4 优化问题
基于上述不对称风险敏感代价函数$ \mathcal{L}_{A} $, 本文将餐前剂量决策设计成带约束的随机优化问题, 具体形式如下:
$$ \begin{array}{l} \min\limits_{u}\; \mathcal{L}(u): = \mathcal{L}_{A}+R(u-u_p)^2 \end{array} $$ (11) $$ \begin{split} &\text{s.t.} \; \\ &G_i\sim\mathcal{N}(m_i,\sigma_i),\; i\in\{1,2,\cdots,8\} \end{split} $$ (12) $$ \begin{array}{l} \; -2 \leq( u - u_p) \leq 2 \end{array} $$ (13) $$ \begin{array}{l} \; 0 \leq u \leq u_{\max} \end{array} $$ (14) 式中, $ u $为待优化的餐前剂量, $ u_p $为上一次相同餐型餐前输注的剂量, $ R $为权重系数. 设计惩罚项$ R(u-u_p)^2 $引入迭代更新的方式, 使得在小样本情况下, 优化者能基于上次剂量, 快速确定合适的剂量. 其中$ R $可用于调整对上次剂量的信任程度, 可使得优化者偏向于选择上次剂量附近的剂量.
2. 临床经验辅助的贝叶斯优化决策
由式(11)可知, 随机变量的性质和不对称惩罚矩阵的设计使得难以计算该代价函数对待优化剂量$ u $的梯度. 因而, 本文结合贝叶斯优化[36]和蒙特卡洛法求解上述随机优化问题. 贝叶斯优化首先利用高斯过程拟合代价函数(11), 然后利用采集函数确定下次评估的剂量, 同时利用蒙特卡洛法计算对应的代价函数观测值, 收集观测样本更新高斯过程, 重复上述操作, 最后确定优化问题的解. 另外, 考虑到在小样本情况下, 数据驱动建模利用少量的输入/输出数据, 无法准确学习患者餐后血糖代谢规律, 本文设计模型有效性在线评估机制, 评估血糖预测模型的有效性, 进而实施临床经验决策辅助, 避免小样本情况下贝叶斯优化决策失误; 同时实时收集样本, 更新预测模型, 形成自适应优化决策.
2.1 贝叶斯优化
2.1.1 代价函数建模
首先, 利用高斯过程拟合式(11)中决策变量$ u $与代价函数$ \mathcal{L}(u) $之间的映射关系. 均值函数考虑为常数形式, 协方差函数采用平方指数形式. 给定$ N $个代价函数的观测样本, $ \mathcal{D}_{1:N} = \{u_{1:N},\mathcal{L}(u_{1:N})\} $, 通过最大化后验概率确定高斯过程的超参数. 接着, 利用预测式(4)、式(5), 预测新剂量$ u_* $对应的代价函数值, 该值将被用于构建采集函数(见第2.1.2节)来确定下一次评估的剂量. 在确定评估的剂量后, 利用蒙特卡洛法估算对应的代价函数观测值. 具体地, 首先, 给定一剂量, 根据第1.2节方法计算餐后血糖变化轨迹的概率分布; 然后, 独立生成1000条餐后血糖轨迹样本; 最后, 计算这些样本对应的代价函数的均值, 该均值即为该剂量下代价函数的观测值.
2.1.2 采集函数
采集函数是贝叶斯优化中最重要的一个环节, 其直接决定优化方向, 影响求解速度. 具体地, 基于建模阶段提供的预测信息, 采集函数通过权衡探索新的搜索空间和利用现有的可行区域确定下一次评估点. 目前, 相关研究已提出了多种采集函数[36], 本文采用改善期望(Expected improvement, EI)采集函数[37]来确定每次的评估剂量.
为最小化代价函数(11), 改善期望采集函数$ \alpha_{EI}(u_*) $, 设计如下:
$$ \begin{array}{l} I(u_*): = (\mathcal{L}^p_{m}-\hat{\mathcal{L}}(u_*))\mathbf{1}(\mathcal{L}^p_{m}>\hat{\mathcal{L}}(u_*)) \end{array} $$ (15) $$ \begin{array}{l} \alpha_{EI}(u_*): = {\rm{E}}(I(u_*)) \end{array} $$ (16) 式中, $ \mathcal{L}^p_{m} $表示截至目前代价函数观测到的最小值, $ \hat{\mathcal{L}}(u_*)\sim\mathcal{N}(m_*(u_*),\sigma^2_*(u_*)) $为待评估剂量$ u_* $对应的代价函数预测量. 式中$ \mathcal{L}^p_{m}-\hat{\mathcal{L}}(u_*) $表示最小值的预期改善量, $\mathbf{1}(\mathcal{L}^p_{m}>\hat{\mathcal{L}}(u_*)) $表示该改善量发生的概率. 计算上述期望, 可得:
$$\left\{\begin{aligned} &\alpha_{EI}(u_*)=\left\{\begin{aligned} &(\mathcal{L}^p_{m}-m_*(u_*))\Phi(U)+\sigma_*(u_*)\phi(U),\\ &\qquad\quad{}m_*(u_*)>0\\ &0,\qquad{\text{其他}} \end{aligned}\right.\\ &U=\frac{\mathcal{L}^p_{m}-m_*(u_*)}{\sigma_*(u_*)} \end{aligned}\right.$$ (17) 式中, $ \Phi(\cdot) $为标准正态分布的累积分布函数, $\phi(\cdot)$为标准正态分布的概率密度函数. 最大化$ \alpha_{EI} $, 即可确定下一次评估剂量$ u_* $的值.
$$ \begin{array}{l} \max\limits_{u_*}\; \alpha_{EI}(u_*) \end{array} $$ (18) $$ \begin{split} & \text{s.t.} \;\\ &\mathcal{L}(u_*) = \mathcal{L}_{A}+R(u_*-u_p)^2 \end{split} $$ (19) $$ \begin{array}{l} -2 \leq u_* - u_p \leq 2 \end{array} $$ (20) 式中, 约束条件由设计的不对称风险敏感代价函数引入, 使得采集函数能实时利用更新的血糖管理样本, 基于上次剂量$ u_p $快速确定下一次评估的剂量.
贝叶斯优化的求解步骤见算法1, 在进行$ M $次顺序优化求解后, 得到优化问题(11) ~ (14)的解, 记为$ u_b $. 为保证安全性, 利用体内残余活性胰岛素约束修正该剂量, IOB是指上次给药一段时间后在体内仍然残留的有降糖作用的胰岛素, 用来约束当前胰岛素输注, 以防止过量注射胰岛素导致低血糖, 其计算方法可参见文献[31]. 记IOB约束为$ u_{\text{IOB}} $, 则最终确定的餐前胰岛素剂量为:
$$ \begin{array}{l} \widetilde{u}_b = u_b-u_{\text{IOB}} \end{array} $$ (21) 算法1. 餐前剂量贝叶斯优化决策算法
输入. 初始化数据集${\cal{D}} $.
输出. 最优胰岛素剂量$u_* $.
1)$ \mathcal{D}\leftarrow $初始化. $ \{u_{1:5},\mathcal{L}(u_{1:5})\} $, 其中$ u_{1:5} $为在式(20)区间等距离选择的5个评估剂量;
2) while 迭代次数$ \leq M $ do;
3)利用数据集$ \mathcal{D} $训练高斯过程;
4)计算待评估剂量$ u_* $对应的代价函数预测值(式(4)、式(5));
5)最大化采集函数(18)确定待评估剂量$ u_* $的值;
6)计算确定剂量$ u_* $对应的$ \mathcal{L}(u_*) $观测值(见第2.1.1节);
7)将样本$ \{u_{*},\mathcal{L}(u_*)\} $加入数据集$\mathcal{D};$$\;{} $
8) end while;
9)选择代价函数所有观测值的最小值对应的剂量作为最后求解的结果.
2.2 预测模型有效性评估
餐后血糖预测模型影响着贝叶斯优化的决策方向. 考虑到只有符合生理代谢规律的预测模型才能保证贝叶斯优化能根据预测的高低血糖正确地决策剂量, 因而本文依据胰岛素具有降糖作用这一生理代谢规律设计模型有效性在线评估机制, 防止贝叶斯优化决策失误. 具体地, 该机制利用训练得到的高斯过程, 根据此次餐前血糖监测水平, 基于上次相同餐型的餐前剂量, 分别预测增加$ \Delta u $、保持不变和减少$ \Delta u $剂量情况下的餐后血糖轨迹, 并依次记为$ {\boldsymbol P}_e^+ $、$ {\boldsymbol P}_e $ 和$ {\boldsymbol P}_e^- $, 其都为由8步血糖预测值组成的向量, 进而设计如下的判断条件, 比较确定模型的有效性:
$$ \begin{array}{l} J = \left\{\begin{array}{ll} 1,&{} {\boldsymbol P}_e^+<{\boldsymbol P}_e<{\boldsymbol P}_e^-\\ 0,&{\text{其他}} \end{array}\right. \end{array} $$ (22) 式中, $ J $为判断结果, 若为1, 则认为预测模型有效, 采用贝叶斯优化决策推荐的剂量; 若为0, 则认为预测模型无效, 采用临床经验辅助决策推荐的剂量(见第2.3节), 避免贝叶斯优化决策失误. $“ < ”$表示一一比较两向量相同位置上元素的大小, 若都满足小于条件, 则返回逻辑1; 否则, 返回逻辑0.
2.3 临床经验辅助决策
临床经验辅助决策可约束在预测模型不可信时的贝叶斯优化, 推荐合适的餐前剂量, 保证安全性. 根据发布的中国1型糖尿病胰岛素治疗指南[24], 国内患者普遍规律饮食, 缺乏估算膳食CHO含量的经验, 通常采用固定的餐前剂量, 之后根据血糖监测结果进行调整. 一般地, 餐前剂量基于全天胰岛素总量的40% ~ 60%, 根据早、中、晚三餐按1/3、1/3、1/3或1/5、 2/5、2/5分配. 这种方式不能有效地利用餐后血糖管理历史数据, 同时没有形成系统的剂量实时调整方法, 容易受到医生或患者主观意向的影响. 鉴于此, 本文将医生临床决策经验量化, 形成系统的临床经验决策规则. 具体地, 利用历史数据, 分3档评估上次相同餐型餐后血糖的控制效果, 再结合考虑本次餐前血糖态势(包括血糖水平和变化趋势), 基于上次剂量调整确定本次剂量, 所设计的临床经验决策规则见图3.
图3中$ G_h $和$ G_l $分别表示评估餐后血糖控制水平的高低阈值, $G_{\min}^{I}$表示上一次相同餐型餐后到下一餐餐前血糖监测的最小值, $ \bar{G}_c $表示此次餐前半小时的血糖监测平均值, $ \bar{G}_p $表示上一次相同餐型餐前半小时的血糖监测平均值, $ \Delta G $表示血糖均值的增量, 上升、平稳和下降表示此次餐前血糖水平的变化趋势. 临床经验辅助决策, 首先通过查询上次相同餐型餐后血糖控制效果和此次餐前血糖态势符合的决策规则, 推荐图中对应的剂量$ u_b $; 然后, 同样利用IOB约束, 修正此次剂量.
3. 结果
本文方法参数利用美国食品药品监督管理局接受的UVA/Padova T1DM血糖代谢仿真器调试得到, 同时根据临床经验, 考虑其是否符合人体血糖代谢实际情况, 最后加以确认, 具体数值见表1. 进一步地, 设计不同仿真测试方案, 利用仿真器中的成年患者组(内含10个模拟成年人)验证本文方法的性能. 另外, 收集1型糖尿病患者的临床治疗数据, 利用顾问模式分析方法[26], 评估本文方法的决策性能.
表 1 方法参数Table 1 The parameters of the proposed method变量 含义 值 $ \gamma $ 风险敏感参数 −2 $ R $ 权重系数 10 $ u_{\max} $ 最大剂量约束值 20 $ M $ 最大优化次数 25 $ {Q\mathit{\boldsymbol{}}}_b $ 基准惩罚矩阵 $\rm{diag}\{0.0025,\cdots, 0.0025,$ $0.005,0.005\}$ $ \Gamma $ 不对称惩罚强度 $ \{1,10,5,1\} $ $ {G\mathit{\boldsymbol{}}}_r $ 餐后血糖控制目标 $[110,130,150,170, 160,$ $150,150,150]^{\rm{T} }$ $ \Delta u $ 评估增量 $ 2 $ $ G_h $ 餐后血糖高阈值 $ 120 $ $ G_l $ 餐后血糖低阈值 $ 80 $ $ \Delta G $ 餐前血糖均值增量 $ 30 $ 基于以下场景对本文方法进行仿真. 模拟仿真假定患者采用每日多次胰岛素皮下注射治疗方案, 即注射基础胰岛素和餐前大剂量胰岛素. 仿真从第1天的5:00开始, 持续5天, 共120 h. 考虑到糖尿病患者通常具有规律的日常活动(如类似的用餐时间和CHO摄入量), 但为了模拟生活方式存在的随机性, 假设患者早餐、午餐、晚餐的CHO摄入量服从正态分布(均值和标准差分别为[50, 75, 75] g和[3, 4, 4] g); 用餐时间分别在[7:00, 9:00]、[11:00, 13:00] 、[18:00, 20:00]之间, 服从均匀分布. 假定5天的胰岛素基础率保持正常值不变, 前3天的餐前胰岛素剂量由标准计算方法调整得到:
$$ u_{\text{bolus}} = \left(\frac{\text{CHO}}{\text{CR}}+\frac{G_{c}-G_{sp}}{\text{CF}}-u_{\text{IOB}}\right)\times r $$ (23) 式中, $ G_{sp} $通常选为140 mg/dL; $ r $为剂量调整系数, 可用于仿真前期不同的剂量决策情况; 后2天的餐前胰岛素剂量利用本文方法确定.
对于仿真器而言, 虚拟患者最优的CR和CF参数是已知的, 标准计算方法可在仿真器中决策近似最优的餐前剂量, 实现良好的餐后血糖控制. 因此, 为仿真现实生活中存在的剂量决策情况(偏多、偏少和正常), 本文将标准计算方法决策的剂量乘以不同调整系数$(r = 1.5、 r = 0.5$和$r = 1)$确定为前期决策的剂量, 生成本文方法前期所需的训练样本, 后2天的餐前剂量利用本文方法确定, 并与保持调整系数不变的前期决策方法对比, 以此说明本文方法的学习与改善血糖控制的能力. 特别地, 当本文方法能实现与标准计算方法$(r = 1)$类似的控制性能时, 可在一定程度上说明, 利用本文方法可实现无需CR和CF参数的餐前胰岛素剂量优化决策. 在下文图表中, 本文提出的方法用“本文方法”表示, 不同调整系数的标准计算方法用“标准方法”表示, 同时在表头和图注标明调整系数.
由于仿真实验中, 每天的用餐时间和CHO摄入量都类似, 本文方法分别收集前3天的早餐、午餐和晚餐样本, 各3个, 为每个餐型构建无需膳食信息的餐后血糖预测模型, 并在后2天依次收集样本进行更新. 仿真结果主要使用正常血糖70 ~ 180 mg/dL时间比率、低血糖低于70 mg/dL时间比率、严重低血糖低于54 mg/dL时间比率、高血糖高于250 mg/dL 时间比率和血糖浓度平均值大小等, 如表2、表3所示的7种评价指标, 对方法性能进行评估. 表中数据以中位数(四分位距)的形式展现, $ p $值由仿真的患者数量计算得到, $ p < 0.05 $被认为具有统计学差异性, 并在表中以粗体标出. 同时, 餐后血糖调节及剂量决策比较结果见图4 ~ 图6, 图中曲线包括5%、25%、75%和95%的四分位数曲线和中位数曲线, 相关讨论见第3.1.1节.
表 2 前期餐前剂量正常情形下$r = 1 $血糖控制评估结果 (%)Table 2 Evaluation results of the glucose control in the case of normal boluses in the early stage when$r=1 $ (%)评价指标 标准方法 本文方法 $ p $ 值 < 54 mg/dL (%) 0 (0) 0 (0) 1.000 < 70 mg/dL (%) 0 (0) 0 (0) 1.000 70 ~ 180 mg/dL (%) 93.6 (7.5) 92.3 (6.6) 0.188 > 250 mg/dL (%) 0 (0) 0 (0) 1.000 血糖平均值(mg/dL) 133.2 (9.4) 132.5 (14.0) 1.000 标准差(mg/dL) 30.2 (7.0) 30.3 (4.4) 0.037 7:00血糖平均值(mg/dL) 109.8 (11.0) 109.3 (8.0) 0.195 表 3 前期餐前剂量偏少/偏多情形下血糖控制评估结果(%)Table 3 Evaluation results of the glucose control in the case of under-estimated/over-estimated boluses in the early stage (%)评价指标 前期餐前剂量偏少($ r = 0.5 $) 前期餐前剂量偏多($ r = 1.5 $) 标准方法 本文方法 $ p $ 值 标准方法 本文方法 $ p $ 值 < 54 mg/dL (%) 0 (0) 0 (0) 1.000 1.6 (3.8) 0 (0) 0.031 < 70 mg/dL (%) 0 (0) 0 (0) 1.000 6.5 (7.8) 0 (1.4) 0.008 70 ~ 180 mg/dL (%) 64.8 (23.6) 77.1 (18.4) 0.004 91.6 (6.9) 95.6 (2.6) 0.232 > 250 mg/dL (%) 1.0 (4.5) 0 (0) 0.031 0 (0) 0 (0) 1.000 血糖平均值 (mg/dL) 162.1 (29.2) 148.7 (16.8) 0.004 114.9 (9.0) 125.6 (10.5) 0.002 标准差 (mg/dL) 37.9 (7.6) 35.8 (5.6) 0.010 29.4 (6.2) 30.9 (5.3) 0.049 7:00血糖平均值 (mg/dL) 118.0 (20.0) 113.0 (9.5) 0.008 104.8 (10.0) 108.0 (12.5) 0.039 
图 4 前期剂量正常$r=1$情形下, 餐后血糖调节和剂量决策比较结果Fig. 4 Comparison of postprandial glucose regulation and dosage decision in the case of normal boluses in the early stage when$r=1$
图 5 前期餐前剂量偏少$ r = 0.5 $情形下, 餐后血糖调节和剂量决策比较结果Fig. 5 Comparison of postprandial glucose regulation and dosage decision in the case of under-estimated boluses in the early stage when$ r = 0.5 $
图 6 前期餐前剂量偏多$ r = 1.5 $情形下, 餐后血糖调节和剂量决策比较结果Fig. 6 Comparison of postprandial glucose regulation and dosage decision in the case of over-estimated boluses in the early stage when$ r = 1.5 $为进一步说明本文方法的学习能力和鲁棒决策能力, 本文利用仿真测试实际生活场景下本文方法的血糖控制性能, 并与R2R方法[38]进行比较. 具体地, 模拟仿真从第1天的5:00开始, 持续5天, 共120 h. 考虑到患者在实际住院调整餐前胰岛素剂量过程中, 每天要求摄入固定碳水化合物含量的膳食, 餐前剂量由固定的碳水化合物含量确定, 但患者执行始终存在偏差(少食或过食). 因此为模拟实际生活场景, 假设患者早餐、午餐和晚餐的CHO摄入量分别在标准量[50, 75, 75] g附近随机变动, 变动范围为[−15, 15] g; 用餐时间分别在[7:00, 9:00]、[11:00, 13:00]、[18:00, 20:00], 服从均匀分布. 餐前剂量由标准计算方法按CHO摄入标准量确定, 但三餐分别有着不适当的CR参数, CF参数始终保持正常不变. 在后续4天, 利用R2R方法优化三餐的CR参数, 其中R2R方法形式如下:
$$ \begin{split} &{CR}_i^{-1}(k+1) = \\ &\quad\left\{ \begin{aligned} &{CR}_i^{-1}(k)-{CR}_i^{-1}(0)k_1T_b(k),\;\;\; {\text{}}T_b(k) > 0\\ &{CR}_i^{-1}(k)+{CR}_i^{-1}(0)\Bigg(k_2T_a(k)\;+\\ &\qquad k_3\frac{G_m-G_T}{G_T}\Bigg),\;\;\;\;\;\;\;\;\;\;\qquad {\text{其他}} \end{aligned}\right. \end{split} $$ (24) 式中, $ k $表示第$ k $次R2R优化, 下标$ i\in\{1,2,3\} $分别对应早餐、午餐和晚餐, $ {CR}_i^{-1}(0) $表示CR初始值的倒数, $k_1、 k_2$和$ k_3 $表示R2R增益, $T_b(k)、 T_a(k)$和$ G_m $为餐后血糖控制评价指标, 分别表示低血糖低于70 mg/dL时间比率, 高血糖高于180 mg/dL时间比率和血糖浓度平均值, 其统计时长最长为5 h, 或从本次餐后到下次餐前. 式中参数数值与文献[38]一致. 为确保比较的公平性, 本文方法收集R2R方法前3天的数据进行建模, 并在后续两天进行决策与R2R方法比较, 比较结果如表4和图7所示, 相关讨论见第3.1.2节.
表 4 膳食扰动下血糖控制评估结果(%)Table 4 Evaluation results of the glucose control in the case of meal disturbance (%)评价指标 R2R 本文方法 $ p $ 值 $< $54 mg/dL (%) 0 (0.7) 0 (0) 0.250 $< $70 mg/dL (%) 3.6 (8.3) 0 (0) 0.016 70 ~ 180 mg/dL (%) 87.7 (13.0) 89.1 (5.6) 0.625 $>$250 mg/dL (%) 0 (0) 0 (0) 1.000 血糖平均值(mg/dL) 123.8 (10.6) 131.3 (7.0) 0.006 标准差(mg/dL) 36.0 (9.4) 34.9 (7.5) 0.432 7:00血糖平均值(mg/dL) 102.8 (12.0) 103.8 (8.5) 0.004 
图 7 膳食扰动下, 本文方法与R2R方法餐后血糖调节和剂量决策比较结果Fig. 7 Comparison of postprandial glucose regulation and dosage decision in the case of meal disturbance between the proposed method and the R2R method顾问模式分析方法可将临床试验中获得的血糖监测数据重新提供给本文方法, 使得本文方法在对应的用餐时刻, 根据餐前血糖监测数据, 确定相应的胰岛素剂量, 并可与同时刻医生建议的胰岛素剂量进行比较, 评估方法的决策性能. 本文方法决策的胰岛素剂量只是用于比较, 并不会对历史血糖监测数据造成因果性影响. 评估结果见图8、图9, 相关讨论见第3.2节.
图 9 基于临床数据的餐前剂量决策评估结果Fig. 9 The evaluation results for the meal bolus decision based on the clinical data3.1 仿真验证
3.1.1 前期不同剂量决策情况
1)前期餐前剂量正常情形
取剂量调整系数$ r = 1 $, 仿真前期餐前剂量正常的情形. 将后2天的血糖控制效果进行比较, 统计结果如表2所示. 后3天的餐后血糖调节及剂量决策比较结果如图4所示. 由表2和图4可看出, 在前期训练样本良好时, 本文方法可在无需膳食信息下实现与标准方法类似的血糖控制效果, 表现在正常血糖70 ~ 180 mg/dL时间比率(92.3%对93.6%, $p = 0.188)$和血糖平均值(132.5对133.2, $p = 1.000)$等评价指标, 两方法均未出现低血糖情况(0对0, $p = 1.000)$.
2)前期餐前剂量偏低情形
取剂量调整系数$ r = 0.5 $, 仿真前期餐前剂量偏低的情形. 该情形用于测试本文方法在前期样本餐前剂量偏低时, 是否可快速地确定合适的餐前剂量, 安全有效地改善患者出现的餐后高血糖情况. 将仿真的后2天血糖控制效果进行比较, 比较结果如表3和图5所示. 可以看出, 在前期剂量决策偏少时, 与继续采用剂量偏少的决策方法相比, 在后2天的决策中, 本文方法可基于前1天的剂量, 快速确定合适的剂量, 有效缓解继续高血糖情况, 这可以从正常血糖70 ~ 180 mg/dL时间比率(77.1%对64.8%, $p = 0.004)$、高血糖高于250 mg/dL时间比率(0%对1%, $p = 0.031)$和血糖平均值(148.7对162.1, $p = 0.004)$等评价指标看出.
3)前期餐前剂量偏高情形
取剂量调整系数$ r = 1.5 $, 仿真前期餐前胰岛素剂量偏高的情形. 该情形用于测试本文方法在前期样本餐前剂量偏高时, 是否可以快速地确定合适的餐前剂量, 安全有效地改善患者出现的餐后低血糖情况. 将仿真的后2天血糖控制效果进行比较, 比较结果如表3和图6所示. 可以看出, 在前期剂量决策偏多时, 与上一种情况相似, 在后2天的决策中, 本文方法可基于前一天的剂量, 快速确定合适的剂量, 有效缓解继续采用剂量偏多的决策方式所带来的低血糖情况, 体现在低血糖低于54 mg/dL时间比率(0%对1.6%, $p = 0.031)$、低于70 mg/dL时间比率(0%对6.5%, $p = 0.008)$和血糖平均值(125.6对114.9, $p = 0.002 )$等评价指标.
3.1.2 膳食扰动场景
由表4和图7可看出, R2R方法容易受到膳食实际摄入量与标准量不匹配的影响, 决策过多或过少的餐前剂量, 导致餐后低血糖事件或持续高血糖事件的发生. 由于患者每天摄入类似CHO含量的膳食, 本文方法利用历史数据建立包含预测不确定性的高斯过程预测模型, 有效降低了膳食信息对决策性能的影响, 从而在改善血糖控制性能的同时, 对膳食扰动具有更好的鲁棒性. 这可以从低血糖情况的改善(0%对3.6%, $p = 0.016)$和血糖平均值的提升(131.3对123.8, $p = 0.006)$等评价指标的改善看出.
3.2 临床数据验证
本文收集在北京大学人民医院采用每日多次胰岛素皮下注射治疗方案的某位1型糖尿病患者临床数据, 并利用顾问模式评估方法性能. 在患者治疗过程中, 餐前胰岛素剂量仍由医生决策确定, 本文方法决策的剂量只用于后续的比较分析, 并未实际应用于患者治疗中. 患者住院期间采用瞬感血糖监测技术(Flash glucose monitoring, FGM), 每15 min监测一次皮下血糖浓度. 该数据收集通过北京大学人民医院伦理审查委员会的批准, 并获得所有参与者的书面知情同意[39]. 由于患者住院时间普遍在一周左右且在住院期间饮食规律, 因此利用前3天的血糖监测数据和医生餐前决策剂量训练得到无需膳食信息的餐后血糖预测模型, 并利用后4天临床数据验证系统性能. 患者前3天的历史血糖管理数据如图8所示, 后4天的顾问模式分析结果如图9所示, 图中医生决策的剂量以“医生建议”表示, 本文方法决策的剂量中由临床经验决策规则决策的以“决策规则”表示, 由贝叶斯优化决策的以“贝叶斯优化”表示.
由图8和图9可以看出, 本文方法可基于上次相同餐型的剂量, 快速确定合适的餐前剂量. 其中临床经验辅助决策考虑到上次相同餐型的餐后血糖控制水平高于高阈值$ G_h $且本次餐前血糖水平较高, 依据逻辑规则推荐高于医生建议的剂量, 有效地利用患者血糖历史管理数据进行剂量决策. 由图9可以看出, 在临床经验辅助决策下, 本文方法不断收集样本实现贝叶斯学习优化, 持续学习患者血糖代谢规律, 使得贝叶斯优化决策模块可基于预测的餐后高血糖, 优化确定高于医生建议的剂量. 这些高于医生建议的剂量, 有望缓解患者餐后的高血糖情况.
4. 结束语
本文设计一种临床经验辅助的餐前胰岛素剂量贝叶斯学习优化方法. 该方法通过评估高斯过程预测模型是否符合血糖代谢规律, 决定是否实施临床经验辅助决策, 避免小样本情况下贝叶斯优化决策的失误. 同时, 不断收集样本, 持续更新预测模型, 实现自适应优化决策. 仿真和临床数据验证结果表明, 本文所设计的临床经验决策规则, 可有效地利用患者血糖历史管理数据辅助决策, 使得贝叶斯优化决策可持续有效地学习患者血糖代谢规律, 根据预测的高/低血糖水平, 基于上次相同餐型的剂量, 快速确定合适的剂量. 本文方法辅助医生决策的临床试验申请已通过北京大学人民医院伦理审查(批件号: 2020PH338-01), 未来将通过随机对照临床试验, 进一步评估方法的安全性和有效性.
-
图 1 餐前胰岛素剂量贝叶斯学习优化方法框架图
Fig. 1 Block diagram of the Bayesian learning based optimization method for meal bolus
图 4 前期剂量正常$r=1$情形下, 餐后血糖调节和剂量决策比较结果
Fig. 4 Comparison of postprandial glucose regulation and dosage decision in the case of normal boluses in the early stage when$r=1$
图 5 前期餐前剂量偏少$ r = 0.5 $情形下, 餐后血糖调节和剂量决策比较结果
Fig. 5 Comparison of postprandial glucose regulation and dosage decision in the case of under-estimated boluses in the early stage when$ r = 0.5 $
图 6 前期餐前剂量偏多$ r = 1.5 $情形下, 餐后血糖调节和剂量决策比较结果
Fig. 6 Comparison of postprandial glucose regulation and dosage decision in the case of over-estimated boluses in the early stage when$ r = 1.5 $
图 7 膳食扰动下, 本文方法与R2R方法餐后血糖调节和剂量决策比较结果
Fig. 7 Comparison of postprandial glucose regulation and dosage decision in the case of meal disturbance between the proposed method and the R2R method
图 9 基于临床数据的餐前剂量决策评估结果
Fig. 9 The evaluation results for the meal bolus decision based on the clinical data
表 1 方法参数
Table 1 The parameters of the proposed method
变量 含义 值 $ \gamma $ 风险敏感参数 −2 $ R $ 权重系数 10 $ u_{\max} $ 最大剂量约束值 20 $ M $ 最大优化次数 25 $ {Q\mathit{\boldsymbol{}}}_b $ 基准惩罚矩阵 $\rm{diag}\{0.0025,\cdots, 0.0025,$ $0.005,0.005\}$ $ \Gamma $ 不对称惩罚强度 $ \{1,10,5,1\} $ $ {G\mathit{\boldsymbol{}}}_r $ 餐后血糖控制目标 $[110,130,150,170, 160,$ $150,150,150]^{\rm{T} }$ $ \Delta u $ 评估增量 $ 2 $ $ G_h $ 餐后血糖高阈值 $ 120 $ $ G_l $ 餐后血糖低阈值 $ 80 $ $ \Delta G $ 餐前血糖均值增量 $ 30 $ 
下载: 导出CSV
表 2 前期餐前剂量正常情形下$r = 1 $血糖控制评估结果 (%)
Table 2 Evaluation results of the glucose control in the case of normal boluses in the early stage when$r=1 $ (%)
评价指标 标准方法 本文方法 $ p $ 值 < 54 mg/dL (%) 0 (0) 0 (0) 1.000 < 70 mg/dL (%) 0 (0) 0 (0) 1.000 70 ~ 180 mg/dL (%) 93.6 (7.5) 92.3 (6.6) 0.188 > 250 mg/dL (%) 0 (0) 0 (0) 1.000 血糖平均值(mg/dL) 133.2 (9.4) 132.5 (14.0) 1.000 标准差(mg/dL) 30.2 (7.0) 30.3 (4.4) 0.037 7:00血糖平均值(mg/dL) 109.8 (11.0) 109.3 (8.0) 0.195
下载: 导出CSV
表 3 前期餐前剂量偏少/偏多情形下血糖控制评估结果(%)
Table 3 Evaluation results of the glucose control in the case of under-estimated/over-estimated boluses in the early stage (%)
评价指标 前期餐前剂量偏少($ r = 0.5 $) 前期餐前剂量偏多($ r = 1.5 $) 标准方法 本文方法 $ p $ 值 标准方法 本文方法 $ p $ 值 < 54 mg/dL (%) 0 (0) 0 (0) 1.000 1.6 (3.8) 0 (0) 0.031 < 70 mg/dL (%) 0 (0) 0 (0) 1.000 6.5 (7.8) 0 (1.4) 0.008 70 ~ 180 mg/dL (%) 64.8 (23.6) 77.1 (18.4) 0.004 91.6 (6.9) 95.6 (2.6) 0.232 > 250 mg/dL (%) 1.0 (4.5) 0 (0) 0.031 0 (0) 0 (0) 1.000 血糖平均值 (mg/dL) 162.1 (29.2) 148.7 (16.8) 0.004 114.9 (9.0) 125.6 (10.5) 0.002 标准差 (mg/dL) 37.9 (7.6) 35.8 (5.6) 0.010 29.4 (6.2) 30.9 (5.3) 0.049 7:00血糖平均值 (mg/dL) 118.0 (20.0) 113.0 (9.5) 0.008 104.8 (10.0) 108.0 (12.5) 0.039
下载: 导出CSV
表 4 膳食扰动下血糖控制评估结果(%)
Table 4 Evaluation results of the glucose control in the case of meal disturbance (%)
评价指标 R2R 本文方法 $ p $ 值 $< $54 mg/dL (%) 0 (0.7) 0 (0) 0.250 $< $70 mg/dL (%) 3.6 (8.3) 0 (0) 0.016 70 ~ 180 mg/dL (%) 87.7 (13.0) 89.1 (5.6) 0.625 $>$250 mg/dL (%) 0 (0) 0 (0) 1.000 血糖平均值(mg/dL) 123.8 (10.6) 131.3 (7.0) 0.006 标准差(mg/dL) 36.0 (9.4) 34.9 (7.5) 0.432 7:00血糖平均值(mg/dL) 102.8 (12.0) 103.8 (8.5) 0.004
下载: 导出CSV
-
[1] 纪立农. 活的更长、活得更好——糖尿病与“健康中国2030”. 中国糖尿病杂志, 2019, 27(01):1-2 doi: 10.3969/j.issn.1006-6187.2019.01.001Ji Li-Nong. Live longer and live better---diabetes and "Healthy China 2030". Chinese Journal of Diabetes, 2019, 27(01):1-2 doi: 10.3969/j.issn.1006-6187.2019.01.001 [2] 纪立农. 胰岛素无针注射器的中国研究:FREE研究. 中华糖尿病杂志, 2020, 12(9):674-676 doi: 10.3760/cma.j.cn115791-20200811-00497Ji Li-Nong. A Chinese study on needle-free injection in type 2 diabetes mellitus: FREE Study. Chinese Journal of Diabetes, 2020, 12(9):674-676 doi: 10.3760/cma.j.cn115791-20200811-00497 [3] Liu W, McGuire H C, Kissimova-Skarbek K, Zhou X, Han X, Wang Y, et al. Factors associated with resistance to complications in long-standing type 1 diabetes in China. Endocrine Research, 2020, 45(1):1-8 doi: 10.1080/07435800.2019.1624567 [4] Schiffrin A, Belmonte M M. Multiple daily self-glucose monitoring: its essential role in long-term glucose control in insulin-dependent diabetic patients treated with pump and multiple subcutaneous injections. Diabetes Care, 1982, 5(5):479-484 doi: 10.2337/diacare.5.5.479 [5] Quiroz G. The evolution of control algorithms in artificial pancreas: A historical perspective. Annual Reviews in Control, 2019, 48:222-232 doi: 10.1016/j.arcontrol.2019.07.004 [6] Kovatchev B. A century of diabetes technology: signals, models, and artificial pancreas control. Trends in Endocrinology and Metabolism, 2019, 30(7):432-444 doi: 10.1016/j.tem.2019.04.008 [7] 史大威, 杨肖, 蔡德恒, 牟治宇, 刘蔚, 纪立农. 基于胰岛素基础率估计的人工胰腺系统自抗扰控制. 自动化学报, 2021, 47(5):1043-1057Shi Da-Wei, Yang Xiao, Cai De-Heng, Mou Zhi-Yu, Liu Wei, Ji Li-Nong. Active disturbance rejection control for artificial pancreas system based on insulin basal rate estimation. Acta Automatica Sinica, 2021, 47(5):1043-1057 [8] El Fathi A, Raef Smaoui M, Gingras V, Boulet B, Haidar A. The artificial pancreas and meal control: An overview of postprandial glucose regulation in type 1 diabetes. IEEE Control Systems, 2018, 38(1):67-85 doi: 10.1109/MCS.2017.2766323 [9] Klonoff D. The current status of bolus calculator decision-support software. Journal of Diabetes Science and Technology, 2012, 6(5):990-994 doi: 10.1177/193229681200600501 [10] Shashaj B, Busetto E, Sulli N. Benefits of a bolus calculator in pre- and postprandial glycaemic control and meal flexibility of paediatric patients using continuous subcutaneous insulin infusion (CSII). Diabetic Medicine, 2008, 25(9):1036-1042 doi: 10.1111/j.1464-5491.2008.02549.x [11] Herrero P, Pesl P, Bondia J, Reddy M, Oliver N, Georgiou P, et al. Method for automatic adjustment of an insulin bolus calculator: In silico robustness evaluation under intra-day variability. Computer Methods and Programs in Biomedicine, 2015, 119(1):1-8 doi: 10.1016/j.cmpb.2015.02.003 [12] Schiavon M, Dalla Man C, Cobelli C. Insulin sensitivity indexbased optimization of insulin to carbohydrate ratio in silico study shows efficacious protection against hypoglycemic events caused by suboptimal therapy. Diabetes Technology and Therapeutics, 2018, 20(2):98-105 doi: 10.1089/dia.2017.0248 [13] Torrent-Fontbona F, Lopez B. Personalized adaptive CBR bolus recommender system for type 1 diabetes. IEEE Journal of Biomedical and Health Informatics, 2019, 23(1):387-394 doi: 10.1109/JBHI.2018.2813424 [14] Turksoy K, Hajizadeh I, Samadi S, Feng J, Sevil M, Park M, et al. Real-time insulin bolusing for unannounced meals with artificial pancreas. Control Engineering Practice, 2017, 59:159-164 doi: 10.1016/j.conengprac.2016.08.001 [15] Herrero P, Bondia J, Adewuyi O, Pesl P, El-Sharkawy M, Reddy M, et al. Enhancing automatic closedloop glucose control in type 1 diabetes with an adaptive meal bolus calculator - in silico evaluation under intra-day variability. Computer Methods and Programs in Biomedicine, 2017, 146:125-131 doi: 10.1016/j.cmpb.2017.05.010 [16] Shi D, Dassau E, Doyle III F J. Multivariate learning framework for long-term adaptation in the artificial pancreas. Bioengineering and Translational Medicine, 2019, 4(1):61-74 doi: 10.1002/btm2.10119 [17] 侯忠生, 许建新. 数据驱动控制理论及方法的回顾和展望. 自动化学报, 2009, 35(6):650-667 doi: 10.3724/SP.J.1004.2009.00650Hou Zhong-Sheng, Xu Jian-Xin. On data-driven control theory: the state of the art and perspective. Acta Automatica Sinica, 2009, 35(6):650-667 doi: 10.3724/SP.J.1004.2009.00650 [18] 代伟, 柴天佑. 数据驱动的复杂磨矿过程运行优化控制方法. 自动化学报, 2014, 40(9):2005-2014Dai Wei, Chai Tian-You. Data-driven optimal operational ccontrol of complex grinding processes. Acta Automatica Sinica, 2014, 40(9):2005-2014 [19] 郜传厚, 渐令, 陈积明, 孙优贤. 复杂高炉炼铁过程的数据驱动建模及预测算法. 自动化学报, 2009, 35(6):725-730 doi: 10.3724/SP.J.1004.2009.00725Gao Chuan-Hou, Jian Ling, Chen Ji-Ming, Sun You-Xian. Data-driven modeling and predictive algorithm for complex blast furnace ironmaking process. Acta Automatica Sinica, 2009, 35(6):725-730 doi: 10.3724/SP.J.1004.2009.00725 [20] Wang Y, Dassau E, Doyle III F J. Closed-loop control of artificial pancreatic beta-cell in type 1 diabetes mellitus using model predictive iterative learning control. IEEE Transactions on Biomedical Engineering, 2010, 57(2):211-219 doi: 10.1109/TBME.2009.2024409 [21] Wang Y, Zhang J, Zeng F, Wang N, Chen X, Zhang B, et al. "Learning" can improve the blood glucose control performance for type 1 diabetes mellitus. Diabetes Technology and Therapeutics, 2017, 19(1):41-48 doi: 10.1089/dia.2016.0328 [22] Song L, Liu C, Yang W, Zhang J, Kong X, Zhang B, et al. Glucose outcomes of a learning-type artificial pancreas with an unannounced meal in type 1 diabetes. Computer Methods and Programs in Biomedicine, 2020, 191(105416):1-9 [23] Cai D, Liu W, Ji L, Shi D. Bayesian optimization assisted meal bolus decision based on Gaussian processes learning and risk-sensitive control. Control Engineering Practice, 2021, 114(104881):1-11 [24] 中华医学会糖尿病学分会. 中国1型糖尿病胰岛素治疗指南. 中华糖尿病杂志, 2016, 8(10):591-597 doi: 10.3760/cma.j.issn.1674-5809.2016.10.005Chinese Diabetes Society. Insulin therapy guidelines for type 1 diabetes in China. Chinese Journal of Diabetes, 2016, 8(10):591-597 doi: 10.3760/cma.j.issn.1674-5809.2016.10.005 [25] Dalla Man C, Micheletto F, Lv D, Breton M, Kovatchev B, Cobelli C. The UVA/PADOVA type 1 diabetes simulator: New features. Journal of Diabetes Science and Technology, 2014, 8(1):26-34 doi: 10.1177/1932296813514502 [26] Gillis R, Palerm C C, Zisser H, Jovanovic L, Seborg D E, Doyle III F J. Glucose estimation and prediction through meal responses using ambulatory subject data for advisory mode model predictive control. Journal of Diabetes Science and Technology, 2007, 1(6):825-833 doi: 10.1177/193229680700100605 [27] Rasmussen C E, Williams C K I. Gaussian Processes for Machine Learning. Cambridge: The MIT Press, 2006. [28] Hovorka R, Canonico V, Chassin L J, Haueter U, Massi-Benedetti M, Federici M O, et al. Nonlinear model predictive control of glucose concentration in subjects with type 1 diabetes. Physiological Measurement, 2004, 25(4):905-920 doi: 10.1088/0967-3334/25/4/010 [29] Deisenroth M P, Fox D, Rasmussen C E. Gaussian processes for data-efficient learning in robotics and control. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2015, 37(2):408-423 doi: 10.1109/TPAMI.2013.218 [30] Jain A, Nghiem T X, Morari M, Mangharam R. Learning and control using Gaussian processes. In: Proceedings of the 9th ACM/IEEE International Conference on Cyber-physical Systems. Porto, Portugal: 2018. 140−149 [31] Gondhalekar R, Dassau E, Doyle III F J. Periodic zone-MPC with asymmetric costs for outpatient-ready safety of an artificial pancreas to treat type 1 diabetes. Automatica, 2016, 71:237-246 doi: 10.1016/j.automatica.2016.04.015 [32] Lee J B, Dassau E, Gondhalekar R, Seborg D E, Pinsker J E, Doyle III F J. Enhanced model predictive control (eMPC) strategy for automated glucose control. Industrial and Engineering Chemistry Research, 2016, 55(46):11857-11868 doi: 10.1021/acs.iecr.6b02718 [33] Whittle P, Kuhn J. A hamiltonian formulation of risk-sensitive linear/quadratic/Gaussian control. International Journal of Control, 1986, 43(1):1-12 doi: 10.1080/00207178608933445 [34] Yang X, Maciejowski J. Risk-sensitive model predictive control with Gaussian process models. IFAC-PapersOnLine, 2015, 48(28):374-379 doi: 10.1016/j.ifacol.2015.12.156 [35] Pan Y, Boutselis G, Theodorou E. Efficient reinforcement learning via probabilistic trajectory optimization. IEEE Transactions on Neural Networks and Learning Systems, 2018, 29(11):5459-5474 doi: 10.1109/TNNLS.2017.2764499 [36] Shahriari B, Swersky K, Wang Z, Adams R P, De Freitas N. Taking the human out of the loop: A review of bayesian optimization. Proceedings of the IEEE, 2015, 104(1):148-175 [37] Jones D R, Schonlau M, Welch W J. Efficient global optimization of expensive black-box functions. Journal of Global Optimization, 1998, 13(4):455-492 doi: 10.1023/A:1008306431147 [38] Toffanin C, Visentin R, Messori M, Palma F D, Magni L, Cobelli C. Towards a run-to-run adaptive artificial pancreas: In silico results. IEEE Transactions on Biomedical Engineering, 2018, 65(4):479-488 [39] Liu W, Chen J, He L, Cai X, Zhang R, Gong S, et al. Flash glucose monitoring data analyzed by detrended fluctuation function on β cell function and diabetes classification. Diabetes Obesity and Metabolism, 2021, 23(3):774-781 doi: 10.1111/dom.14282 期刊类型引用(2)
1. 林燕明，王慧璇，王少芳. 二甲双胍联合德谷门冬双胰岛素在2型糖尿病合并原发性肝癌患者中的应用价值分析. 糖尿病新世界. 2024(20): 56-59 . 
百度学术2. 陈向阳，李舍予. 非内分泌病房住院糖尿病患者的血糖管理. 中国全科医学. 2023(15): 1799-1803 . 百度学术其他类型引用(0)
-
下载:
计量
- 文章访问数: 1245
- HTML全文浏览量: 716
- PDF下载量: 201
- 被引次数: 2
