还在用烂大街的单疾病模型？看别人如何用机器学习构建双疾病模型，寻找关联疾病的关键节点

过去几十年，人们普遍认为肾病(CKD)启动并加速了主动脉瓣钙化(CAVD)，这反过来又增加了CKD患者的死亡风险，但相关机制仍不清楚。

2023年6月发表在“Journal of Translational Medicine(IF7.4)“上的文章，利用多种生物信息学工具，从GEO数据库收集了3个CAVD数据集和2个CKD数据集，揭示了CKD相关CAVD的核心基因和潜在机制。

Screening of immune-related secretory proteins linking chronic kidney disease with calcific aortic valve disease based on comprehensive bioinformatics analysis and machine learning

基于综合生物信息学分析和机器学习筛选慢性肾脏疾病与钙化主动脉瓣疾病相关的免疫相关分泌蛋白

https://doi.org/10.1186/s12967-023-04171-x

一 背景

• 慢性肾脏疾病(CKD)影响全球超过8%的人口，不仅表现为肾功能下降，还表现为矿物质过度沉积、炎症级联和氧化应激[1-3]，这些都与动脉粥样硬化、主动脉瓣性心脏病等心血管疾病密切相关[4]。主动脉瓣钙化病(CAVD)是最常见的瓣膜疾病之一，可能导致心力衰竭和心源性猝死等[5]。近期研究表明，与一般人群相比，CKD中CAVD更为常见，CKD是CAVD预后的独立危险因素[1-2]，提示CKD患者发生CAVD的风险可能更高。然而，导致CKD相关CAVD的潜在分子机制仍不清楚。

• 研究表明，过量的内源性和外源性介质可诱导CKD炎症，释放多种促炎细胞因子，这与CKD的进展和心血管疾病的发生有关[6]。此外，CKD的特征是细胞过早衰老，并表现出衰老表型，炎症介质、Wnt/β-catenin信号相关配体等的分泌，导致肾脏和其他靶器官或组织的级联衰老[7]。需要指出的是，衰老显著参与了各种疾病的病理过程，尤其是血管钙化过程[8]。这些研究表明CKD可能导致包括CAVD在内的后续并发症，至少部分是通过分泌蛋白引起的。

二 方法

01 数据处理与差异表达分析

从GEO数据库下载CAVD和对照组的3个表达谱数据集GSE12644、GSE51472和GSE83453。通过GEO获得CKD患者外周血单个核细胞(PBMC)芯片数据集(GSE37171)和肾脏组织芯片数据集(GSE66494)。基于R(v4.2.1)包SVA 的combat功能，对3个CAVD数据集进行批量校正整合。

对CAVD整合数据集和CKD数据集(GSE37171和GSE66494)进行背景校正、归一化和基因ID转换。随后用R包Limma识别CAVD和CKD数据集中的DEGs(CAVD p≤0.05，|log2FC|≥1；CKD p≤0.05，|log2FC|≥0.585)。用R包ggplot2和pheatmap可视化。

02 分泌蛋白鉴定与蛋白相互作用(PPI)网络构建

分泌蛋白从人类蛋白质图谱数据库下载。从"SPOCTOPUS预测分泌蛋白"的蛋白类中下载编码分泌蛋白的基因共3970个。

基于STRING数据库建立连接CKD和CAVD的PPI网络(中等置信评分>0.4)，挖掘CKD相关分泌蛋白与CAVD关键基因间的相互作用。PPI网络通过Cytoscape软件(v3.8.2)可视化，用分子复合物检测(MCODE)插件检测有意义的模块，取得分前2的模块进一步分析。

03 功能富集与cMAP分析

导入DAVID数据库进行GO和KEGG通路富集分析，探索CKD相关致病基因的生物学功能和具体机制(p<0.05)。通过气泡图和circos图可视化。将得分最高的前2个显著模块中的上调基因纳入cMAP在线数据库，寻找潜在的治疗CAVD的小分子药物，鉴定出富集分数最高的前10个化合物。

04 机器学习与诊断模型构建

为了识别候选生物标志物并建立CAVD的诊断模型，本研究首先采用逻辑回归筛选变量的最小绝对收缩和选择算子(LASSO)算法，用glmnet包筛选候选生物标志物。然后用随机森林(RF)算法整合多棵树以获得更好的准确性，用randomForest包筛选候选生物标志物。将LASSO模型与RF模型中MeanDecreaseGini>2重叠的基因定义为hub基因，构建CKD相关CAVD诊断模型。

用rms包构建基于两个hub基因的列线图。绘制受试者工作特征(ROC)曲线下面积，评估各基因及列线图诊断CAVD的效能。校准曲线和决策曲线分析(DCA)评估列线图对CKD相关CAVD的预测效能。

05 免疫浸润分析

通过CIBERSORT包评估CAVD基因表达谱中免疫细胞浸润的数量。用Wilcoxon检验比较CAVD与对照样本间22种免疫细胞比例的差异(p<0.05)。用corrplot包显示侵袭免疫细胞间的关联。最后用Spearman秩相关系数分析生物标志物表达与浸润免疫细胞含量的相关性(p<0.05)。

06 外部验证数据采集

采集对照组、CAVD患者、CKD伴或不伴CAVD患者的血清样本，测定血清SLPI和MMP9水平。提取总RNA进行逆转录。实时荧光定量PCR(RT-qPCR)，据2−ΔΔCt法计算mRNA的相对表达量。

引物序列:

SLPI-F 5ʹ-GAGATGTTGTCCTGACACTTGTG-3ʹ

SLPI-R 5ʹ-AGGCTTCCTCCT TGTTGGGT-3ʹ

MMP9-F 5ʹ-ACGCAGACATCGTCATCCAGT-3ʹ

MMP9-R 5ʹ-GACCACAACTCGTCATCGTC-3ʹ

GAPDH-F 5ʹ-GAGTCAACGGATTGGTCGT-3ʹ

GAPDH-R 5ʹ-GACAAGCTTCCCGTTCTCAG-3ʹ

三 结果

01 数据处理与差异表达基因鉴定

对CAVD的3个原始数据集去除批次效应然后合并+归一化，包括34个钙化样本和23个对照样本。去除批次效应后，三个数据集间的差异显著减小(图2A-B)。整合的CAVD数据集中发现173个差异表达基因(DEGs)，其中119个上调，54个下调(图2C-D)。

图2 CAVD数据集的整合与差异表达基因

(A)批效应校正前三个原始CAVD数据集的主成分分析；(B)校正后CAVD数据集主成分分析；(C)综合CAVD数据集DEGs火山图；(D)CAVD前30个和下调DEGs热图。

02 WGCNA鉴定CAVD中的关键模块

进行加权基因共表达网络分析，识别钙化主动脉瓣样本中最相关的基因模块。选择软阈值为5(图3A)，生成14个模块并聚类(图3B-C)。CAVD与基因模块的相关性(图3D)表明，粉色模块与CAVD正相关(358个基因，r= 0.84，p=5e−16)，黄色模块与CAVD最负相关(769个基因，r=−0.72，p=2e−10)。将粉色和黄色模块作为关键模块。分别在粉色和黄色模块中发现模块成员和基因显著性之间的强关联(图3E, F)，共确定1127个与CAVD显著相关的关键基因。将CAVD样本中的DEGs与WGCNA中的关键基因取交集，共获得124个基因进行后续分析(图3G)。

图3 WGCNA筛选CAVD数据集的关键模块基因

(A)无尺度拓扑模型，平均连通性和尺度无关性图。(B)基因树状图和模块特征基因网络热图。(C)模块特征基因的聚类树图。(D)模块特征基因与CAVD状态关系。(E)粉色模块隶属度与基因显著性的相关图。(F)黄色模块隶属度与基因显著性的相关图。(G)关键模块基因与DEGs的交集。

03 肾病分泌蛋白鉴定与PPI网络构建

CKD可能加速CAVD的发生和进展。分析CKD患者外周血单个核细胞(PBMC)和肾脏的表达谱，在PBMC中共发现2681个DEGs，在肾组织中发现4111个DEGs(图4A-D，p≤0.05，|log2FC|≥0.585)。PBMC(图4E)和肾组织(图4F)中分别有376和607种差异表达的分泌蛋白。

图4 鉴定PMBC和肾脏组织中CKD相关分泌蛋白

(A)PMBC的DEGs。(B)PBMC数据集中前30个上调/下调的DEGs。(C)Kidney的DEGs。(D)Kidney数据集中前30个上调/下调的DEGs。(E)PMBC的DEGs与分泌蛋白基因交集。(F)Kidney的DEGs与分泌蛋白基因交集。

MCODE软件确定显著性最高的前2个模块，纳入的76个基因被确定为CDK相关致病基因。将致病基因导入DAVID在线数据库进行功能富集和KEGG分析。GO分析显示，CDK相关CAVD的致病基因主要富集在“炎症反应”和“免疫反应”中(图5C，生物过程(BP))；致病基因多位于“膜的组成部分”和“细胞外区”(图5D，细胞成分(CC))；“蛋白结合”和“相同蛋白结合”是致病基因中最相关的项目(图5E，分子功能(MF))。KEGG通路分析显示CKD相关的CAVD致病基因与“细胞因子-细胞因子受体相互作用”、“PI3K-Akt信号通路”和“NF-Kappa B信号通路”(图5F)有很强的相关性。

图5 CDK相关分泌蛋白和CAVD关键基因的PPI及富集分析

(A)MCODE分析得分最高的模块1基因的PPI网络。(B)模块2基因的PPI网络。(C-F)模块1和模块2的GO分析结果，基因的生物过程(C)、细胞成分(D)和分子功能(E)。(F)模块1和模块2中基因的KEGG分析结果。

04 鉴定治疗CAVD的候选小分子化合物

将CAVD中上调的CDK相关致病基因导入connectivity map (cMAP)数据库，预测可以逆转CDK相关致病基因在CAVD中表达改变的小分子化合物。美替拉酮、吉非替尼、地拉西普、氨基五胺、甲氧沙林、福斯可林、CGP-37157、IKK2抑制剂、阿糖腺苷和TG-101348是阴性评分最高的前10个化合物，是治疗CDK相关CAVD的潜在药物(图6A)。这10个化合物的靶向途径和化学结构描述如图6B-C。

图6 cMAP分析筛选治疗CAVD的潜在小分子化合物

(A)基于cMAP分析的10个细胞系中富集程度最高的前10个化合物的热图。(B)前10名化合物信息。(C)10种化合物的化学结构。

05 机器学习筛选hub基因并构建诊断模型

在CKD分泌蛋白与CAVD的关键基因和DEGs取交集鉴定出17个基因，构建区分CKD伴或不伴CAVD患者的CAVD诊断模型(图7A)。用LASSO回归算法从17个基因中筛选出8个潜在基因(图7B, C)。同时，用随机森林(RF)机器学习算法，根据每个基因的变量重要性对17个常见基因排序，提取出6个MeanDecreaseGini>2的基因(图7D)。

图7 机器学习鉴定CKD相关CAVD的潜在诊断标志物

两种方法只有两个hub基因重叠，分别是分泌型白细胞蛋白酶抑制剂(SLPI)和基质金属蛋白酶9(MMP9)(图7E)。通过logistics回归构建了基于这两个hub基因的列线图(图8A)。应用受试者工作特征(ROC)曲线评估各核心基因的曲线下面积(AUC)值和列线图。两个hub基因的AUC值均>0.9，且列线图的AUC值高于每个hub基因，提示列线图可能具有较强的诊断价值(图8B-D)。校准曲线、DCA分析显示，列线图诊断模型的预测概率与理想模型几乎相同(图8E)，可能有利于CKD相关CAVD的诊断(图8F)。在GSE51472数据集中，列线图在CKD硬化主动脉瓣患者中显示出理想的预测价值(图8G)，这意味着该列线图模型对CKD并早期CAVD患者具有良好的诊断效能。

图8 nomogram模型开发与疗效评估

06 CAVD免疫细胞浸润及与hub基因相关性分析

研究发现CAVD中CKD相关致病基因与炎症和免疫过程密切相关。用CIBERSORT算法分析CAVD免疫细胞的特征，探索CAVD的免疫调节机制以及诊断标志物与免疫细胞浸润的相关性。10个免疫细胞亚群在CAVD和对照样本间存在显著差异(图9A)。与对照组相比，CAVD组巨噬细胞M0、CD8 T细胞和调节性T细胞(Tregs)的比例较高，而幼B细胞、活化树突状细胞、巨噬细胞M2、活化肥大细胞、活化NK细胞、浆细胞和CD4初始T细胞的比例较低(图9B)。此外，对22种免疫细胞的相关性分析表明，CD4初始T细胞与Tregs显著正相关(r=0.57,p<0.05)，活化的肥大细胞与活化的树突状细胞负相关(r=−0.68,p<0.05)(图9C)。hub基因SLPI和MMP9在CAVD中均与免疫细胞聚集显著相关(图9D)。

图9 CAVD免疫细胞浸润分析

07 hub基因表达与列线图模型验证

验证生信分析准确性，检测两个hub基因在外部队列中的表达模式。与对照组相比，RT-qPCR证实CAVD样本中两个hub基因上调(图10A)。SLPI和MMP9在血清中的水平在CKD和CAVD患者及CKD合并CAVD患者中显著升高(图10B)。然后，基于验证队列构建CAVD诊断列线图模型A，预测从对照组和CAVD组中发生CAVD的可能性(图10C)。根据ROC曲线，与各标志物相比，列线图A在对照组和CAVD患者间的AUC最高(图10D)。校准曲线和评估列线图A的DCA显示，基于列线图A的决策可能有利于CAVD的预测(图10E-F)。此外，构建列线图模型B区分CKD患者是否伴有CAVD(图10G)。同样，ROC曲线和校准曲线以及DCA显示列线图B对CKD患者CAVD的预测价值较好(图10H-J)。

图10 CAVD样本hub基因表达模式的验证及列线图模型区分CAVD的诊断性能评估

总结

01 CAVD数据集中发现173个差异表达基因(DEGs)，其中119个上调，54个下调。与WGCNA分析取交集，获得124个关键基因。

02 CKD PBMC和肾组织中分别有376和607种差异表达的分泌蛋白。

03 CDK相关CAVD的致病基因主要富集在“炎症反应”和“免疫反应”中，与细胞因子-细胞因子受体相互作用、PI3K-Akt和NF-Kappa B信号通路强相关。

04 CDK分泌蛋白与CAVD关键基因取交集，得到17个基因，分别用LASSO和RF算法鉴定，得到2个共有的hub基因SLPI和MMP9。

05 两个hub基因的AUC值均>0.9，在CAVD中均与免疫细胞聚集显著相关。

06 RT-qPCR证实两个hub基因在CAVD样本中上调，在血清检测中CKD和CAVD患者及CKD合并CAVD患者中显著升高。

REFERENCES

[1]. Brandenburg VM, Schuh A, Kramann R. Valvular calcification in chronic kidney disease. Adv Chronic Kidney Dis. 2019;26:464–71.

[2]. Rattazzi M, et al. Aortic valve calcification in chronic kidney disease. Nephrol Dialysis Transplantation. 2013;28:2968–76.

[3]. Benz K, Hilgers K-F, Daniel C, Amann K. Vascular calcification in chronic kidney disease: the role of inflammation. Int J Nephrol. 2018. https://doi.org/10.1155/2018/4310379.

[4]. Go AS, Chertow GM, Fan D, McCulloch CE, Hsu CY. Chronic kidney disease and the risks of death, cardiovascular events, and hospitalization. ACC Curr J Rev. 2004. https://doi.org/10.1056/NEJMoa041031.

[5]. Driscoll K, Cruz AD, Butcher JT. Inflammatory and biomechanical drivers of endothelial-interstitial interactions in calcific aortic valve disease. Circul Res. 2021;128:1344–70.

[6]. Speer T, Dimmeler S, Schunk SJ, Fliser D, Ridker PM. Targeting innate immunity-driven inflammation in CKD and cardiovascular disease. Nat Rev Nephrol. 2022. https://doi.org/10.1038/s41581-022-00621-9.

[7]. van Deursen JM. The role of senescent cells in ageing. Nature. 2014. https://doi.org/10.1038/nature13193.

[8]. Sutton NR, et al. Molecular mechanisms of vascular health: insights from vascular aging and calcification. Arterioscler Thromb Vasc Biol. 2022. https://doi.org/10.1161/ATVBAHA.122.317332.

往期推荐

1、SangerBox生信分析工具使用过程的各种问题解答

2、【生信文章复现】视频教程集合

3、【sangerbox生信分析】在线解答关于GSEA小工具的各种使用问题

4、定制化气泡图绘制工具操作指南