基于网络药理学和分子对接筛选黄芪治疗高血压心室重构的生物标记物

代华磊，胡成成，张桂敏，陶四明，陈建昆

（云南大学附属医院心血管病中心，云南 昆明 650021）

高血压是心血管疾病共同的危险因素，高血压的长期压力超负荷必将引起心脏结构和功能的改变，即心室重构。心室重构是心肌细胞凋亡、心脏肥大和心肌纤维化的复杂病理过程，主要表现为左室肥厚和慢性心力衰竭，其发生率会随着高血压的进程逐渐上升，是人类健康的一大威胁[1]。心室重构会进一步损害心功能，最终引发心律失常、心衰、心肌梗死、猝死等严重并发症[2]。

黄芪（Astragalus membranaceus）是一种传统的中草药，其主要成分为黄芪多糖、皂苷、黄酮类化合物等。研究表明，黄芪对心肌具有良好的保护作用，可在一定程度内控制高血压[3-4]，但其作用的分子机制尚不明确。因此，借助网络药理学的方法挖掘黄芪治疗高血压心室重构作用的靶点与通路，从微观的角度科学论证黄芪治疗高血压心室重构的作用机制具有重要意义。

网络药理学以生物学与药理学的交叉学科理论为基础，从整体的角度探索药物与疾病间的关联性，运用各组学、高通量筛选、网络分析等多种前沿技术，揭示 “疾病-基因-靶点-药物”之间复杂的网络关系，具备整体性和系统性的特点。与中药及其方剂多组分、多靶点干预及系统调控的原理观念基本一致，能够多维度了解疾病的分子基础，预测药物潜在的药理作用机制[5-6]。本研究基于网络药理学原理，结合分子对接，全面分析黄芪治疗高血压心室重构的分子机制，筛选出关键生物标记物，为黄芪治疗高血压心室重构提供更多潜在的生物标记物以供后续的研究。

1 数据来源与方法

1.1 数据来源

利用TCMSP（http://tcmspw.com/tcmsp.php）数据库搜索关键词黄芪的有效成分，共获得87 个活性药物。使用ETCM（http://www.tcmip.cn/ETCM/index.php/Home/）数据库获得黄芪的有效成分，共获得20 个活性药物。在TCMSP 数据库及ETCM 数据库中筛选每种药物活性成分对应的靶点。其中，TCMSP 数据库共获得药物靶点181 个，ETCM 数据库共获得药物靶点254 个。对两个数据库获得的药物靶点取并集，共获得390 个药物靶点。利用Genecards（https://www.genecards.org/）数据库获取高血压心室重构的疾病靶点，仅保留Category为Protein Coding 的基因，共获得3 281 个相关靶点。利用GEO 数据库（https://www.ncbi.nlm.nih.gov/geo/）下载数据集GSE74144 的高血压心室重构及正常样本全血的mRNA 测序数据及样本信息。

1.2 高血压心室重构相关模板及其基因筛选

将高血压心室重构作为WGCNA 的性状数据，使用R 包“WGCNA”筛选高血压心室重构相关的关键模块及其基因。对样本进行聚类，软阈值β=12，构建共表达网络，并且使用动态树切割算法获得模块并分析各模块与高血压心室重构的相关性，筛选出关键模块。利用limma 包（version3.46.0）比较基因表达水平的差异性，使用R 语言ggplot2（version3.3.3）绘制火山图来展示基因差异表达的情况。利用R 语言VennDiagram 包（version1.6.20）对差异表达基因及相关关键模块基因取交集，获得关键模块基因[7-8]。

1.3 关键靶点的筛选及其调控网络构建

使用R 语言VennDiagram 包对药物靶点、疾病靶点及差异的关键模块基因取交集，并使用R语言ClusterProfiler 包进行GO 功能/KEGG 通路富集分析[9]。利用Cytoscape（version3.8.2）软件对关键靶点-功能、关键靶点-通路的调控网络图进行可视化[10]。

1.4 蛋白互作（protein-protein interaction，PPI）网络构建

针对关键靶点利用STRING（https://string-db.org）网站置信度为0.4（Confidence=0.4）绘制PPI网络，并采用cytoscape 进行可视化。

1.5 生物标记物的筛选及其表达分析

利用LASSO[11]及RF 构建高血压心室重构的诊断模型并进行评估，并使用R 语言Venn-Diagram 包对模型基因取交集获得生物标记物。对生物标记物进行Pearson 相关性分析，并利用R 语言ggpubr 包绘制散点图进行可视化[12]。

1.6 中药-活性成分-生物标记物网络

对生物标记物、药物活性成分绘制中药药理调控网络，并使用Cytoscape 进行可视化。

1.7 分子对接

利用 PubChem（https://pubchem.ncbi.nlm.nih.gov/）数据库下载药物活性成分的3D 结构。并从PDB（https://www1.rcsb.org/）数据库得到关键靶点的蛋白结构，通过AutoDockTools（version1.5.6）完成蛋白质加氢、计算电荷。从PubChem（https://pubchem.ncbi.nlm.nih.gov/）数据库下载活性成分结构，使用AutoDockvina 进行分子对接，最后利用PyMol（version2.5）软件进行可视化以及美化。

2 结果

2.1 高血压心室重构关键模块及其关键基因的筛选

使用R 包“WGCNA”构建共表达网络（图1），筛选高血压心室重构相关的关键模块及基因。通过聚类分析筛选表达谱相似的基因模块（图1A，1B）。为保证基因间相互作用最大限度符合无尺度分布，选择软阈值β 为12（图1C），以合并动态剪切树算法分析获得11 个模块（图1D，1E）。通过模块与性状相关性分析（图1F，1G），可以知道MEgrey60 模块与高血压心室重构具有最高的基因显著性（Cor=0.78，P=3e-4），且相关系数的绝对值最大，MEgrey60 基因与高血压心室重构性状的相关系数为0.64（P＜ 0.05），MEgrey60 模块共4 278 个基因，故我们选择模块MEgrey60 作为高血压心室重构的关键模块进行后续分析。

图1 高血压心室重构关键模块及其关键基因的筛选Fig.1 Screening of key modules of hypertensive ventricular remodeling and their key genes

通过R 语言ggplot2（version3.3.3）绘制火山图来展示差异基因及其表达情况（图1H）。共存在5 073 个差异表达的基因，其中3 132 个基因表达上调，1941 个基因表达下调。

R 语言VennDiagram 包（version1.6.20）对高血压心室重构vs 对照组之间差异表达基因及WGCNA 获得的高血压心室重构相关关键模块基因取交集（图1I），结果共获得2 103 个差异表达的关键模块基因。

2.2 关键靶点的筛选及其调控网络

使用R 语言VennDiagram 包对上述获得的差异表达的关键模块基因、药物靶点及高血压心室重构疾病靶点取交集，共获得24 个关键靶点（图2A）。

图2 关键靶点筛选及其调控网络Fig.2 Screening of key targets and their regulatory networks

使用R 语言“enrichplot”（version1.10.2）对关键靶点基因进行功能富集分析，绘制结果图。GO 功能富集分析结果条形图（图2B）。在生物过程方面，共获得288 个Terms，关键靶点基因与细胞对非生物刺激的反应、细胞对化学应激的反应及对氧化应激的反应显著相关。在分子功能方面，共获得15 个Terms，关键靶点基因与异构酶活性、转氨酶活性、有机酸结合、激酶调节活性功能等显著相关。在细胞组成方面，共获得29个Terms，关键靶点基因与膜区、膜微区、质膜结合细胞投射细胞质等显著相关。KEGG 通路富集分析结果气泡图（图2C）。结果显示，关键靶点基因与长寿调节途径、精氨酸生物合成、酪氨酸代谢、甲状腺激素信号通路等显著相关。

将GO_BP、GO_CC、GO_MF 的Terms 其 对应的的关键靶点提出来，构建了关键靶点-功能的调控网络图，并进行可视化（图2D），图中包含30 个Term，22 个关键靶点，共131 个关系对。将KEGG TOP50 Pathways 及其对应的的关键靶点提出来，构建关键靶点-通路的调控网络图，并进行可视化（图2E），网络包含20 条KEGG Pathways，16 个关键靶点，共73 个关系对。

利用STRING（https://string-db.org）网站对24个关键靶点进行PPI 网络构建（图2F）。置信度为0.4（Confidence=0.4），去除离散的蛋白，得到21个蛋白的互作网关系，包括21 个节点，88 条边。

2.3 高血压心室重构诊断模型的构建

合并后数据集的所有样本（n=22）以5∶5 的比例随机分为训练集（n=11）及验证集（n=11）用来进行诊断模型的构建及验证（图3）。LASSO 分析工具为R 语言glmnet 包（version4.1-1），得到了LASSO 回归常见的有两个图形（图3A），一个是基因系数的图形，一个是交叉验证的误差图。使用R 语言pROC（version1.17.0.1）绘制训练集及验证集的ROC（receiver operating characteristic curve）曲线，并计算AUC（Area Under Curve）值，一般AUC 值越大，说明预测的越准确（图3B）。图中训练集及验证集的AUC 均大于0.9，说明该诊断模型对高血压心室重构具有较高的预测能力。

图3 高血压心室重构最佳诊断模型Fig.3 Best diagnostic model for hypertensive ventricular remodeling

基于得到的24 个关键靶点在GSE74144 数据库（高血压心室重构=14，Normal=8）中各个样本的表达值，结合样本的分组信息，其中以样本分组作为响应变量，24 个关键靶点作为解释变量，我们采用R 语言“caret”包（version6.0-86）构建RF 模型，并使用R 语言DALEX 包（version2.3.0）包的explain 函数对RF 模型进行解释性分析，并用plot 函数对模型表现分布进行可视化，分别绘制累积残差分布图和箱线图分布图（图3C，D），接着分析在RF 模型中，不同变量对于模型预测的相对重要性程度（图3E），由表1 可知，在RF模型中，CSNK2B，MAPK1，GHR，MBL2，SELE，IL2 这6 个变量对响应变量（group，score=0.136 6）的预测值有较大的影响（响应变量以上的基因被选为特征基因），因此本研究将这6 个基因作为候选诊断标志物用作下一步分析。本研究对RF 模型的诊断价值进行了评估，结果RF 模型的AUC 为0.955，说明RF 模型对高血压心室重构具有较高的预测能力（图3F）。

表1 基因在RF 模型中的重要性排序Tab.1 Ranking of importance of genes in the RF model

2.4 关键基因的筛选及其表达分析

使用R 语言VennDiagram 包（version1.6.20）对LASSO 及RF 算法获得的基因取交集，结果共获得4 个候选基因（MAPK1，IL2，CSNK2B，SELE）（图4A），并且本研究对4 个基因绘制单基因ROC 曲线（图4B），并计算AUC 值，所有交集基因的AUC 值均在0.8 以上。因此，本研究将MAPK1，IL2，CSNK2B，SELE 作为生物标记物进行后续分析。

图4 关键基因的筛选及其表达分析Fig.4 Screening and expression analysis of key genes

对4 个生物标记物在高血压心室重构样本及正常样本之间的表达情况进行了研究，并使用R 语言ggpubr 包（version0.4.0）及ggplot2 包（version3.3.3）绘制散点图进行可视化（图4C）。由图可知，MAPK1 和CSNK2B 在高血压心室重构样本中上调表达；而IL2 和SELE 在高血压心室重构样本中下调表达。同时，本研究对生物标记物之间进行了Pearson 相关性分析，并绘制相关性热图（图4D）。可知，生物标记物之间存在较强的相关性。

2.5 生物标记物的中药药理调控网络

将上述获得4 个生物标记物对应的活性成分提出来，构建了药物-活性成分-关键靶点基因网络关键，并用Cytoscape 进行可视化（图5），网络包含1 个药物（黄芪），3 个药物活性成分，4 个生物标记物，共12 个关系对。

图5 生物标记物的中药药理调控网络Fig.5 TCM pharmacological regulatory network of biomarkers

2.6 分子对接

提取了生物标记物的靶向活性分子，分别进行了分子对接（图6）。SELE 与Kaempferol 进行分子对接（图6A），其中ASN-83、LYS-55、GLU-71 及LYS-74 残基与Kaempferol 分子有氢键相互作用。活性分子与蛋白之间的对接亲和力为-5.9 kcal/mol。MAPK1 与quercetin 进行分子对接（图6B），其中LYS-114、MET-108、ASP-167、LYS-54 及ILE-31 残基与quercetin 分子有氢键相互作用。活性分子与蛋白之间的对接亲和力为-8.3 kcal/mol。IL2 与quercetin 进行分子对接（图6C），其中GLU-52，SER-127 及GLN-126残基与quercetin 分子有氢键相互作用。活性分子与蛋白之间的对接亲和力为-6.8 kcal/mol。CSNK2B 与Kumatakenin 进行分子对接（图6D），其中 GLN-96，ARG-92 及 LYS-100 残基与Kumatakenin 分子有疏水键相互作用。活性分子与蛋白之间的对接亲和力为-6.2 kcal/mol。

图6 分子对接结果图Fig.6 Molecular docking result

3 讨论

高血压是心血管疾病发生的基础，其与心血管疾病的发病率呈线性相关，长期的高血压可诱导心室重构发生，进而增加患者死亡率[13]。特别在顽固性高血压人群中，血压往往难以控制，需使用的降压药物种类多，增加了药物不良反应[14]。黄芪能够通过抑制通路的活化改善心脏的功能、减轻心室重构[15]，但大剂量使用会引起胃肠不适、腹泻、恶心和呕吐等症状。因此，探索黄芪治疗高血压心室重构新的治疗靶点和方法显得尤为重要。

利用网络药理学原理筛选出黄芪治疗高血压心室重构的关键生物标记物主要为MAPK1、IL2、CSNK2B、SELE。丝裂原活化蛋白激酶（mitogenactivated protein kinase，MAPK）是一组丝氨酸-苏氨酸蛋白激酶，其介导调节多种细胞过程的关键信号通路，与细胞的增值、分化、凋亡有关。MAPK 通路能够促进平滑肌增殖和胶原合成，改变细胞外基质（extracellular matrix，ECM）的构成，从而使血管重塑加重[16]。MAPK1（也称为细胞外信号调节蛋白激酶2，ERK2）是MAPK 家族的成员，研究表明，ERK1/2 信号通路可调节心脏偏心性和向心性生长之间的平衡[17]。MEK-ERK1/2信号传导有促心室重构的作用，可对大多数已知诱导肥厚生长的压力刺激作出反应，而抑制MEKERK1/2 信号通路对压力超负荷所致的心脏肥大具有保护作用[18]。

炎症与高血压等心血管疾病的发展有关[19]。先天性和适应性的免疫反应都是通过与肾脏、血管和大脑相关的炎症变化，从而促进了高血压的病理发展[20]。炎症介质通过增加血管通透性，释放细胞因子、ROS、NO 和MMP 促进高血压的发生发展。细胞因子的释放会导致阻力血管管腔直径缩小，血管阻力和硬度进而增加[21]。IL2 是一种多效细胞因子，对T 细胞增殖、效应细胞和记忆细胞的产生有至关重要的作用[22]。有研究表明，由IL2 和抗IL2 mAb（JES2-6）组成的IL2 免疫复合物可抑制炎症反应并减轻心肌梗死后的心室重塑[23]。

CSNK2B 编码络蛋白激酶Ⅱ（casein kinase2，CK2）的β 亚基，可调节信号转导、代谢途径、复制、转录和翻译[24]。CSNK2B 参与炎性细胞因子的表达，促使内质网应激[25]。内质网应激在高血压及其并发疾病中发挥重要作用，其通过损伤血管内皮介导的血管舒张功能参与高血压的形成[26-27]，内质网应激还参与了心室重构过程，导致心肌细胞肥大和心脏纤维化[28-29]。

高血压的发生还与内皮细胞功能障碍有关[30]。SELE 编码的E-选择素，存在于细胞因子刺激的内皮细胞中，是一种内皮细胞粘附分子，对内皮细胞活化具有特异性，被认为是介导细胞与血管壁粘附从而导致白细胞在炎症部位积聚的原因[31]。SELE 在mRNA 和蛋白质水平表达的增加与高血压有显著的相关性[32]。而在中国人群中SELE 基因多态性（A561C）的C 等位基因携带者易患高血压[33]。

综上所述，结合我们的研究结果，所筛选出的4 个生物标记物MAPK1、IL2、CSNK2B、SELE具有较高生物学活性，阐明了黄芪对高血压心室重构治疗可能的分子机制，可以作为黄芪治疗高血压心室重构的潜在生物标志物和关键治疗靶点，为高血压心室重构治疗提供了新的治疗靶点和思路方法。