引用本文: 耿祺焜, 李吉利, 胡雲迪, 李润泓, 向涛, 张双弋, 李镭. 基于可解释机器学习的重症慢性阻塞性肺疾病的预后模型. 中国呼吸与危重监护杂志, 2024, 23(3): 153-159. doi: 10.7507/1671-6205.202311003 复制
慢性阻塞性肺疾病(简称慢阻肺)是世界发病率、病死率最高的疾病之一[1],全球慢阻肺患者已增加至3.84亿人[2]。根据世界卫生组织的数据,慢阻肺的病死率在全世界排名第三,每年在全球造成的经济负担超过1 000亿美元[3]。随着治疗方法的发展,尽管慢阻肺住院率有所下降,但呼吸衰竭、疾病恶化仍然会导致患者入住重症监护病房(intensive care unit,ICU),造成严重危害[4]。研究表明慢阻肺是可预防、可治疗的[5]。早期预测重症患者的死亡风险对制定合适的干预措施、合理分配有限的医疗资源至关重要[6-7]。BODE指数是一种广泛应用于慢阻肺患者死亡风险预测的工具,但其中的一些指标(如6分钟步行试验距离)对于重症患者难以测量,限制了其在临床工作中的应用[8]。一些传统的危重症评分系统[如急性生理学和慢性健康状况评价(acute physiology and chronic health evaluation,APACHE)、简化急性生理学评分Ⅱ(simplified acute physiology score Ⅱ,SAPS-Ⅱ)、死亡概率模型(mortality probability model,MPM)]已广泛用于预测重症患者的病死率,但其并非针对慢阻肺而设计[9]。此外,这些工具往往将每项评分结果直接相加,计算死亡风险,忽略了各项指标间的非线性关系[10]。
机器学习是计算机科学和数学相结合的一门学科,它可以捕捉变量之间的非线性关系。近年来,基于机器学习研发了多种应用于临床的预测模型[11-13],为ICU设计的模型中有27.1%用于死亡风险预测[14],其性能优于传统的危重症评分系统[15-17]。然而,许多研究只报道了模型的预测能力,而忽视了可解释性,难以获取患者和医护人员的信任,限制了其在临床实践中的应用[16]。随机森林(random forest,RF)是一种机器学习分类器,由决策树和单个树的输出组成[18]。在医学领域,RF在预测病死率、住院费用等方面性能优异[19-20],但其“黑箱”问题使医生难以信任,即只给出模型的预测结果,例如死亡风险的预测概率值,而无法提供具体的推理过程,导致难以解释每个指标对预测结果的具体影响。此外,机器学习模型的数据维度往往较高,导致预测的分界线是人类难以直观想象的超平面或超曲面,这种多维性使得人们难以理解机器学习模型的内核。近年来,SHAP(SHapley Additive exPlanations)方法被用于解释模型的预测结果,并已应用于预测术后并发症、疾病发展等[21-22]。在这项回顾性研究中,我们旨在基于美国多中心急诊重症监护病(emergency intensive care unit,eICU)数据库开发一种可解释的机器学习模型,用于预测重症慢阻肺患者的全因住院病死率。
1 资料与方法
1.1 临床资料
eICU数据库V2.0是一个多中心的开源数据库,包含2014—2015年间美国208家医院的335个ICU收治的139 367名患者的200 859份入院病例[23]。数据包含患者的性别、年龄、实验室指标、并发症、出院状态等。经过在线培训课程,获得数据库的访问权限(认证号:44589323),负责提取数据。所有程序均遵守1964年《赫尔辛基宣言》及其修正案的道德标准。患者的健康信息全部是匿名的,因此不需要进行伦理审查和患者个人同意。纳入标准:① 根据国际疾病分类第9版(ICD-9)编码,所有患有慢阻肺的患者[24];② ICU入住时间超过24 h且年龄超过18岁;③ 多次入住ICU选取首次入住ICU的数据。排除标准:数据缺失>30%的患者[25]。筛选流程见图1。

1.2 方法
1.2.1 预测变量
主要结局是全因住院病死率,定义为患者出院时的生存状况。通过MySQL 8.0.11软件(甲骨文公司),使用结构化查询语言(SQL)从eICU数据库中提取人口统计学资料、生命体征、实验室指标、并发症中的83个变量作为预测指标[26-28]。除了人口统计学资料外,其他变量都是从每次入住ICU的前24 h内提取。对于连续变量,分别计算其最大值、最小值和平均值,以找出所有可能的关联[25]。在提取变量后,采用LASSO回归筛选变量,避免过拟合[29-30]。
1.2.2 缺失值处理
eICU数据库包含缺失值,这可能会影响模型性能,甚至产生偏倚[31]。在数据提取过程中,我们排除了缺失值超过30%的变量。对于缺失值少于30%的变量,通过R 4.1.0软件,采用missForest非参数缺失值插补方法补全缺失值[32]。该方法可以迭代地使用RF处理缺失的数据,直到满足停止标准[33],已广泛应用于医学研究[34-35]。
1.2.3 可解释性
SHAP方法是一种基于博弈论和加特征归因方法的工具,为预测模型中的危险因素提供测试数据集上的解释,是一种统一、可视化的方法,使机器学习模型更透明、更可行,可用于临床实际的应用[36]。通过Python 3.8.0软件中的SHAP包,计算SHAP值并完成可视化分析,实现各变量重要性的排序以及对个体的死亡风险预测结果进行解释。在我们的研究中,阳性或阴性的SHAP值分别表明相应的特征导致更高或更低的死亡风险[22]。
1.3 统计学方法
比较生存患者和死亡患者的基线特征。通过R 4.1.0软件,取截断点为1%和99%,采用Winsorization方法替换生命体征的离群值。在使用Kolmogorov-Smirnov完成正态性检验后,将连续变量表示为中位数(四分位数)[M(P25,P75)],并使用Wilcoxon秩和检验进行比较。分类变量用数字和百分比表示,采用χ2检验进行比较。P<0.05为差异具有统计学意义。
通过Python 3.8.0软件,建立了RF、XGBoost、logistic回归(LR)、朴素贝叶斯(NB)和支持向量机(SVM)5种机器学习模型。设置特异性为85%,通过受试者操作特征(receiver operating characteristic,ROC)曲线下面积(area under curve,AUC)、敏感性、准确率、F1值、阳性预测值(positive predicative value,PPV)和阴性预测值(negative predicative value,NPV)等多种指标比较机器学习模型和APACHE Ⅳa评分系统的预测性能[37]。其中,AUC的比较通过DeLong检验进行。
2 结果
2.1 患者特征
如图1所示,共有8 088名患者符合我们的筛选标准。将数据集随机分为两组,70%(n=5661)作为训练集,30%(n=2427)作为验证集。训练集中的住院病死率为10.56%(598/5661),测试集中的住院病死率为10.59%(257/2427)。使用LASSO回归在训练集的83种变量中选择17种作为预测指标,每种变量在训练集与测试集中的差异均有统计学意义,见表1。

2.2 机器学习模型开发与评估
利用LASSO回归筛选的预测指标,在训练集中建立5种机器学习模型,分别是RF、XGBoost、LR、NB和SVM。如表2所示,测试集的AUC分别为0.830、0.815、0.814、0.811和0.698。5种机器学习模型中,RF模型的AUC最大[AUC=0.830,95%置信区间(confidential interval,CI)0.806~0.855],SVM模型的AUC最小(AUC=0.698,95%CI 0.659~0.738)。如图2所示,APACHE Ⅳa的AUC为0.782(95% CI 0.752~0.811),经过DeLong检验显著低于RF(P=0.002)。


2.3 用SHAP方法解释RF模型
通过SHAP方法确定每个预测指标与RF模型预测结果之间的相关性,为预测结果提供可视化的解释。变量重要性图(图3)显示了RF模型中排名前10位的危险因素。无创收缩压的最小值是预测全因住院病死率的最重要特征,其次是血尿素氮的平均值、体温的最低值、心率的平均值、是否进行机械通气、白细胞的最大值。图4显示了每个预测指标的SHAP值及其对模型预测结果的影响之间的关联。点的颜色表示该指标的值是高(红色)是低(蓝色)。例如,无创收缩压最小值的升高与病死率下降相关,而血尿素氮平均值的升高与病死率增加相关。

min:最小值;avg:平均值;max:最大值。

min:最小值;avg:平均值;max:最大值。
此外,我们将SHAP框架应用于RF模型,解释测试集中具体患者的预测结果。图5显示了两名患者的预测结果,其中患者A是死亡患者,患者B是生存患者。f(x)的值为概率预测值,来自所有指标的综合结果。每个指标都有不同长度的箭头,箭头越长,该指标对病死率预测的贡献就越大。红色箭头代表该指标使病死率增加,蓝色箭头代表该指标使病死率降低。对于死亡患者,其概率预测值较高(0.60);对于生存患者,其概率预测值较低(0.07)。

BUN:血尿素氮;nibp:无创血压;heart rate:心率;WBC:白细胞计数;SpO2:血氧饱和度;nibp_systolic:无创收缩压;respiratory rate:呼吸频率;mechanical ventilation:机械通气。
3 讨论
随着人工智能的飞速发展,机器学习已广泛应用于呼吸医学领域,尤其是慢阻肺患者[38-40]。与传统的预测方法相比,机器学习能够捕预测指标与预后之间的非线性关系[41- 42],输出最有价值的结果。重症慢阻肺患者的死亡风险预测因其对临床医生的指导意义而备受关注。在这项回顾性研究中,我们筛选出17个易于获取的与预后相关的指标,基于多中心eICU数据库建立了5个预测重症慢阻肺患者病死率的机器学习模型。
APACHE Ⅳa是基于142个变量开发的传统的危重症评分系统,但其中一些变量是多余的,并没有提高预测的准确率[43]。此外,在本项研究中,APACHE Ⅳa的AUC只有0.782,远远低于机器学习模型。近年来研究表明机器学习模型比传统的评分系统具有更好的预测能力[44- 45]。与这些研究一致,我们也发现,相比于传统的风险评分,机器学习模型对于重症慢阻肺患者具有更多的临床获益并且可以更好地预测住院病死率。
在我们的研究中,RF模型在5种算法中性能最好,对这一特定人群具有最佳的预测性能。这是因为LR模型、NB模型和SVM模型本质上都属于线性分类器,即决策边界为超平面,模型结构较为简单,而我们构建的原始数据集包含的变量数目多,数据结构更为复杂,变量之间存在交互作用,并且影响慢阻肺预后的相关因素也较为广泛,因此这三种结构简单的线性分类器很难取得较为满意的效果。XGBoost模型是boosting算法的一种,而RF模型采用的则是bagging思想,前者更善于降低偏差,后者则更善于降低方差;并且RF模型通过随机选取特征产生子树,因此子树彼此之间不同,提高了模型的多样性,更容易学习到特征之间的相互作用。以上两点使得与XGBoost相比,RF模型的泛化性强,过拟合的可能性低,因此取得了更好的分类性能。此外,我们使用SHAP方法为RF模型的预测结果进行解释,有助于临床医生更好地理解其决策过程。利用SHAP方法,我们还根据变量的重要性对危险因素进行了排序。在我们的研究中,无创收缩压最小值的下降具有最强的预测能力,这与以往的临床研究结论一致。Huang等[46]认为收缩压降低是重症慢阻肺患者死亡的独立的危险因素。Franco Palacios等[47]报道,高疾病负担的老年患者收缩压低于120 mm Hg(1 mm Hg=0.133 kPa)会导致不良的结局。一种被广泛接受的观点是,慢阻肺患者持续缺氧甚至高碳酸血症可能导致血管舒张作用,从而抵消交感神经兴奋的作用,减轻神经源性高血压的倾向[48]。对于老年人来说,血压调节机制发生了一定的变化,因此低血压可能增加电解质紊乱、器官灌注不足、晕厥和肾脏疾病恶化的风险[47]。
然而,本研究中所有患者的病死率仅为10.57%,RF模型的阳性预测值仅为0.322,这意味着当RF模型预测患者的病死率远远高于其他患者时,实际死亡的比例不到1/3。因此,该模型的临床适用性有待于在实际临床环境中进一步检验。尽管如此,由于其良好的性能和可解释性,我们认为RF模型可以为临床医生提供合理的指导,并在临床决策中发挥辅助作用。本研究也存在一定的限制性,包括未对慢阻肺的分期进行区分,只能进行早期的死亡风险预测而无法对整个住院期间进行实时的监测。此外,本项回顾性研究在未来的工作中还可以通过大规模的前瞻性研究来完善研究结果。
利益冲突:本研究不涉及任何利益冲突。
慢性阻塞性肺疾病(简称慢阻肺)是世界发病率、病死率最高的疾病之一[1],全球慢阻肺患者已增加至3.84亿人[2]。根据世界卫生组织的数据,慢阻肺的病死率在全世界排名第三,每年在全球造成的经济负担超过1 000亿美元[3]。随着治疗方法的发展,尽管慢阻肺住院率有所下降,但呼吸衰竭、疾病恶化仍然会导致患者入住重症监护病房(intensive care unit,ICU),造成严重危害[4]。研究表明慢阻肺是可预防、可治疗的[5]。早期预测重症患者的死亡风险对制定合适的干预措施、合理分配有限的医疗资源至关重要[6-7]。BODE指数是一种广泛应用于慢阻肺患者死亡风险预测的工具,但其中的一些指标(如6分钟步行试验距离)对于重症患者难以测量,限制了其在临床工作中的应用[8]。一些传统的危重症评分系统[如急性生理学和慢性健康状况评价(acute physiology and chronic health evaluation,APACHE)、简化急性生理学评分Ⅱ(simplified acute physiology score Ⅱ,SAPS-Ⅱ)、死亡概率模型(mortality probability model,MPM)]已广泛用于预测重症患者的病死率,但其并非针对慢阻肺而设计[9]。此外,这些工具往往将每项评分结果直接相加,计算死亡风险,忽略了各项指标间的非线性关系[10]。
机器学习是计算机科学和数学相结合的一门学科,它可以捕捉变量之间的非线性关系。近年来,基于机器学习研发了多种应用于临床的预测模型[11-13],为ICU设计的模型中有27.1%用于死亡风险预测[14],其性能优于传统的危重症评分系统[15-17]。然而,许多研究只报道了模型的预测能力,而忽视了可解释性,难以获取患者和医护人员的信任,限制了其在临床实践中的应用[16]。随机森林(random forest,RF)是一种机器学习分类器,由决策树和单个树的输出组成[18]。在医学领域,RF在预测病死率、住院费用等方面性能优异[19-20],但其“黑箱”问题使医生难以信任,即只给出模型的预测结果,例如死亡风险的预测概率值,而无法提供具体的推理过程,导致难以解释每个指标对预测结果的具体影响。此外,机器学习模型的数据维度往往较高,导致预测的分界线是人类难以直观想象的超平面或超曲面,这种多维性使得人们难以理解机器学习模型的内核。近年来,SHAP(SHapley Additive exPlanations)方法被用于解释模型的预测结果,并已应用于预测术后并发症、疾病发展等[21-22]。在这项回顾性研究中,我们旨在基于美国多中心急诊重症监护病(emergency intensive care unit,eICU)数据库开发一种可解释的机器学习模型,用于预测重症慢阻肺患者的全因住院病死率。
1 资料与方法
1.1 临床资料
eICU数据库V2.0是一个多中心的开源数据库,包含2014—2015年间美国208家医院的335个ICU收治的139 367名患者的200 859份入院病例[23]。数据包含患者的性别、年龄、实验室指标、并发症、出院状态等。经过在线培训课程,获得数据库的访问权限(认证号:44589323),负责提取数据。所有程序均遵守1964年《赫尔辛基宣言》及其修正案的道德标准。患者的健康信息全部是匿名的,因此不需要进行伦理审查和患者个人同意。纳入标准:① 根据国际疾病分类第9版(ICD-9)编码,所有患有慢阻肺的患者[24];② ICU入住时间超过24 h且年龄超过18岁;③ 多次入住ICU选取首次入住ICU的数据。排除标准:数据缺失>30%的患者[25]。筛选流程见图1。

1.2 方法
1.2.1 预测变量
主要结局是全因住院病死率,定义为患者出院时的生存状况。通过MySQL 8.0.11软件(甲骨文公司),使用结构化查询语言(SQL)从eICU数据库中提取人口统计学资料、生命体征、实验室指标、并发症中的83个变量作为预测指标[26-28]。除了人口统计学资料外,其他变量都是从每次入住ICU的前24 h内提取。对于连续变量,分别计算其最大值、最小值和平均值,以找出所有可能的关联[25]。在提取变量后,采用LASSO回归筛选变量,避免过拟合[29-30]。
1.2.2 缺失值处理
eICU数据库包含缺失值,这可能会影响模型性能,甚至产生偏倚[31]。在数据提取过程中,我们排除了缺失值超过30%的变量。对于缺失值少于30%的变量,通过R 4.1.0软件,采用missForest非参数缺失值插补方法补全缺失值[32]。该方法可以迭代地使用RF处理缺失的数据,直到满足停止标准[33],已广泛应用于医学研究[34-35]。
1.2.3 可解释性
SHAP方法是一种基于博弈论和加特征归因方法的工具,为预测模型中的危险因素提供测试数据集上的解释,是一种统一、可视化的方法,使机器学习模型更透明、更可行,可用于临床实际的应用[36]。通过Python 3.8.0软件中的SHAP包,计算SHAP值并完成可视化分析,实现各变量重要性的排序以及对个体的死亡风险预测结果进行解释。在我们的研究中,阳性或阴性的SHAP值分别表明相应的特征导致更高或更低的死亡风险[22]。
1.3 统计学方法
比较生存患者和死亡患者的基线特征。通过R 4.1.0软件,取截断点为1%和99%,采用Winsorization方法替换生命体征的离群值。在使用Kolmogorov-Smirnov完成正态性检验后,将连续变量表示为中位数(四分位数)[M(P25,P75)],并使用Wilcoxon秩和检验进行比较。分类变量用数字和百分比表示,采用χ2检验进行比较。P<0.05为差异具有统计学意义。
通过Python 3.8.0软件,建立了RF、XGBoost、logistic回归(LR)、朴素贝叶斯(NB)和支持向量机(SVM)5种机器学习模型。设置特异性为85%,通过受试者操作特征(receiver operating characteristic,ROC)曲线下面积(area under curve,AUC)、敏感性、准确率、F1值、阳性预测值(positive predicative value,PPV)和阴性预测值(negative predicative value,NPV)等多种指标比较机器学习模型和APACHE Ⅳa评分系统的预测性能[37]。其中,AUC的比较通过DeLong检验进行。
2 结果
2.1 患者特征
如图1所示,共有8 088名患者符合我们的筛选标准。将数据集随机分为两组,70%(n=5661)作为训练集,30%(n=2427)作为验证集。训练集中的住院病死率为10.56%(598/5661),测试集中的住院病死率为10.59%(257/2427)。使用LASSO回归在训练集的83种变量中选择17种作为预测指标,每种变量在训练集与测试集中的差异均有统计学意义,见表1。

2.2 机器学习模型开发与评估
利用LASSO回归筛选的预测指标,在训练集中建立5种机器学习模型,分别是RF、XGBoost、LR、NB和SVM。如表2所示,测试集的AUC分别为0.830、0.815、0.814、0.811和0.698。5种机器学习模型中,RF模型的AUC最大[AUC=0.830,95%置信区间(confidential interval,CI)0.806~0.855],SVM模型的AUC最小(AUC=0.698,95%CI 0.659~0.738)。如图2所示,APACHE Ⅳa的AUC为0.782(95% CI 0.752~0.811),经过DeLong检验显著低于RF(P=0.002)。


2.3 用SHAP方法解释RF模型
通过SHAP方法确定每个预测指标与RF模型预测结果之间的相关性,为预测结果提供可视化的解释。变量重要性图(图3)显示了RF模型中排名前10位的危险因素。无创收缩压的最小值是预测全因住院病死率的最重要特征,其次是血尿素氮的平均值、体温的最低值、心率的平均值、是否进行机械通气、白细胞的最大值。图4显示了每个预测指标的SHAP值及其对模型预测结果的影响之间的关联。点的颜色表示该指标的值是高(红色)是低(蓝色)。例如,无创收缩压最小值的升高与病死率下降相关,而血尿素氮平均值的升高与病死率增加相关。

min:最小值;avg:平均值;max:最大值。

min:最小值;avg:平均值;max:最大值。
此外,我们将SHAP框架应用于RF模型,解释测试集中具体患者的预测结果。图5显示了两名患者的预测结果,其中患者A是死亡患者,患者B是生存患者。f(x)的值为概率预测值,来自所有指标的综合结果。每个指标都有不同长度的箭头,箭头越长,该指标对病死率预测的贡献就越大。红色箭头代表该指标使病死率增加,蓝色箭头代表该指标使病死率降低。对于死亡患者,其概率预测值较高(0.60);对于生存患者,其概率预测值较低(0.07)。

BUN:血尿素氮;nibp:无创血压;heart rate:心率;WBC:白细胞计数;SpO2:血氧饱和度;nibp_systolic:无创收缩压;respiratory rate:呼吸频率;mechanical ventilation:机械通气。
3 讨论
随着人工智能的飞速发展,机器学习已广泛应用于呼吸医学领域,尤其是慢阻肺患者[38-40]。与传统的预测方法相比,机器学习能够捕预测指标与预后之间的非线性关系[41- 42],输出最有价值的结果。重症慢阻肺患者的死亡风险预测因其对临床医生的指导意义而备受关注。在这项回顾性研究中,我们筛选出17个易于获取的与预后相关的指标,基于多中心eICU数据库建立了5个预测重症慢阻肺患者病死率的机器学习模型。
APACHE Ⅳa是基于142个变量开发的传统的危重症评分系统,但其中一些变量是多余的,并没有提高预测的准确率[43]。此外,在本项研究中,APACHE Ⅳa的AUC只有0.782,远远低于机器学习模型。近年来研究表明机器学习模型比传统的评分系统具有更好的预测能力[44- 45]。与这些研究一致,我们也发现,相比于传统的风险评分,机器学习模型对于重症慢阻肺患者具有更多的临床获益并且可以更好地预测住院病死率。
在我们的研究中,RF模型在5种算法中性能最好,对这一特定人群具有最佳的预测性能。这是因为LR模型、NB模型和SVM模型本质上都属于线性分类器,即决策边界为超平面,模型结构较为简单,而我们构建的原始数据集包含的变量数目多,数据结构更为复杂,变量之间存在交互作用,并且影响慢阻肺预后的相关因素也较为广泛,因此这三种结构简单的线性分类器很难取得较为满意的效果。XGBoost模型是boosting算法的一种,而RF模型采用的则是bagging思想,前者更善于降低偏差,后者则更善于降低方差;并且RF模型通过随机选取特征产生子树,因此子树彼此之间不同,提高了模型的多样性,更容易学习到特征之间的相互作用。以上两点使得与XGBoost相比,RF模型的泛化性强,过拟合的可能性低,因此取得了更好的分类性能。此外,我们使用SHAP方法为RF模型的预测结果进行解释,有助于临床医生更好地理解其决策过程。利用SHAP方法,我们还根据变量的重要性对危险因素进行了排序。在我们的研究中,无创收缩压最小值的下降具有最强的预测能力,这与以往的临床研究结论一致。Huang等[46]认为收缩压降低是重症慢阻肺患者死亡的独立的危险因素。Franco Palacios等[47]报道,高疾病负担的老年患者收缩压低于120 mm Hg(1 mm Hg=0.133 kPa)会导致不良的结局。一种被广泛接受的观点是,慢阻肺患者持续缺氧甚至高碳酸血症可能导致血管舒张作用,从而抵消交感神经兴奋的作用,减轻神经源性高血压的倾向[48]。对于老年人来说,血压调节机制发生了一定的变化,因此低血压可能增加电解质紊乱、器官灌注不足、晕厥和肾脏疾病恶化的风险[47]。
然而,本研究中所有患者的病死率仅为10.57%,RF模型的阳性预测值仅为0.322,这意味着当RF模型预测患者的病死率远远高于其他患者时,实际死亡的比例不到1/3。因此,该模型的临床适用性有待于在实际临床环境中进一步检验。尽管如此,由于其良好的性能和可解释性,我们认为RF模型可以为临床医生提供合理的指导,并在临床决策中发挥辅助作用。本研究也存在一定的限制性,包括未对慢阻肺的分期进行区分,只能进行早期的死亡风险预测而无法对整个住院期间进行实时的监测。此外,本项回顾性研究在未来的工作中还可以通过大规模的前瞻性研究来完善研究结果。
利益冲突:本研究不涉及任何利益冲突。