基于生物信息学分析和网络药理学探讨黄连解毒汤干预动脉粥样硬化斑块破裂的作用机制①

李晓辉 彭广操 李兴渊 王建茹（河南中医药大学第一附属医院，郑州 450000）

动脉粥样硬化（atherosclerosis，AS）作为一种慢性进行性血管疾病，其发生发展涉及内皮功能障碍、炎症损伤、脂质浸润和纤维帽破裂等多种病理因素。历经几十年的缓慢发展，动脉斑块可能会突然继发血栓，或斑块出血，常表现为急性冠状动脉综合征（acute coronary syndrome，ACS）。而导致ACS的罪魁祸首是不稳定斑块破裂。2004年共识声明中首先提出了导致ACS的血栓斑块为“斑块破裂”，并指出其是一种具有严重损伤的斑块，在纤维帽中具有真正的缺陷或间隙，将富含脂质的动脉粥样化核心与流动的血液分离，从而暴露出斑块的血栓形成核心［1］。研究发现，全球1 847例致命性动脉血栓中有73%存在斑块破裂［2］。目前导致斑块破裂的机制尚不明确，以斑块破裂为导向的治疗效果欠佳，所以识别不稳定斑块破裂的靶基因及机制，对于防治斑块破裂、形成血栓，预防急性心血管事件的发生具有重要意义。

中药复方具有多成分、多靶点、多途径和多效应的特点，在治疗复杂性疾病方面表现出疗效确切、毒副作用小的优势，但其作用机制研究是一直以来的重点和难点。黄连解毒汤（Huanglian Jiedu Decoction，HLJDD）出自《外台秘要》，由黄连、黄芩、黄柏、栀子4味药物组成。本方苦寒直折，具有清热解毒之功。研究证实，HLJDD联合常规西药能明显改善急性心肌梗死介入术后患者心绞痛症状，显著降低超敏肌钙蛋白Ⅰ（hypersensitivity troponinⅠ，hs-TnⅠ）、N末端B型脑钠钛前体（N-terminal Bbrain natriuretic peptide，NT-proBNP）、肌酸激酶MB同工酶（creatine kinase isoenzyme MB，CK-MB）、血小板聚集功能测定、血清磷脂酶A2及冠脉校正TIMI帧数水平，抑制巨噬细胞极化，有效减少心肌损伤改善心肌功能，且安全性较好［3-4］；还可显著改善ACS患者PCI术后内皮功能（VEGF、NO），降低炎症反应，促进侧支循环的建立［5］。提示HLJDD可延缓或抑制斑块破裂，但其作用机制不明显，限制了其在临床的应用与发展。

近年来，利用生物信息学方法对芯片基因表达数据和高通量测序数据进行整合分析已经广泛应用于各种疾病的研究。稳健排序整合（robust rank aggregation，RRA）方法可用于各种排序基因列表的比较，可有效避免在筛选差异表达基因（differentially expressed genes，DEGs）的过程中因样本数有限、技术平台差异等因素对研究结果产生偏倚的情况［6-7］。网络药理学是基于目标分子、生物功能和生物活性化合物产生复杂的相互作用网络，符合中药复方的天然特征，不仅能够从分子水平上系统地阐明中药复方的作用机制；还可以对中药复方做出整体性预判和前瞻性预测，以弥补分子生物学的局限性［8］。目前已有诸多学者成功地将其应用于阐释中药复方治疗疾病的物质基础和潜在作用机制的研究中。基于以上的认识，本研究利用网络药理学结合生物信息学分析的方法探索具有解毒功效的HLJDD延缓或抑制AS由稳定斑块到破裂过程的作用机制，以期为后续的研究提供一定参考；同时，以方测证，以证验因，在一定程度上可阐述AS中医“毒损脉络”“毒致易损斑块”的科学内涵。

1 资料与方法

1.1 芯片数据的来源从GEO数据库（http：//www.ncbi.nlm.nih.gov/geo/）和ArrayExpress数据库（https：//www.ebi.ac.uk/arrayexpress/）中下载含有人颈AS稳定和破裂斑块组织相关的数据集，详见表1。

1.2 数据预处理和DEGs的筛选 以E-MTAB-2055和GSE41571为基因芯片数据集，利用perl脚本对表达数据进行注释，将探针矩阵转化为基因矩阵；然后，利用avereps函数对基因对应多个探针取均值，利用normalizeBetweenArrays函数对数据进行矫正，若基因表达数值较大，则需取log2；最后，利用limma包筛选DEGs。GSE120521为转录组测序数据，利用avereps函数对重复基因名取均值后，利用wilcoxTest筛选DEGs。以|log（fold change，FC）|＞1和P＜0.05为阈值对3个数据集的表达数据进行差异分析。

表1 数据集的详细信息Tab.1 Details of data set

1.3 RRA法整合分析芯片数据利用稳健排序整合RRA包对3个数据集中的DEGs进行整合，以|log2FC|＞1和P＜0.05为条件筛选稳健DEGs，然后利用pheatmap包绘制P值Top 20的上调和下调基因的热图。

1.4 HLJDD相关活性成分及作用靶点的筛选在TCMSP数据库（http：//tcmspw.com/tcmsp.php）中，输入HLJDD中4味中药名称，即黄连、黄柏、黄芩、栀子，搜索得到各药化学成分，通过口服生物利用度（oral bioavailability，OB）≥30%及类药性（drug likeness，DL）≥0.18对4味中药所含的所有化合物进行筛选，获取其活性成分；同时，在该数据库中获取相应活性成分的作用靶点。在PubChem数据库（http：//pubchem.ncbi.nlm.nih.gov/）中获取HLJDD的活性成分的Canonical SMILES结构；然后，在Swiss Target Prediction数据库（http：//swisstargetprediction.ch/）中对HLJDD活性成分的作用靶点进行预测，并剔除Probability=0的靶点。整合所有活性成分的作用靶点，即为HLJDD的作用靶点。

1.5 HLJDD-AS交集DEGs的获取及HLJDD调控网络的构建将1.3中的稳健DEGs与1.4中的HLJDD作用靶点取交集，并绘制韦恩图，这些交集靶点即为HLJDD延缓或抑制AS斑块破裂的DEGs；然后，用Cytoscape3.7.2软件构建中药-成分-靶标-疾病网络。

1.6 功能富集分析为了探究稳健DEGs和HLJDD-AS交集DEGs可能的功能作用，利用cluster-Profiler和enrichplot包对其进行GO和KEGG富集分析，并利用ggplot2和GOplot包进行可视化。

1.7 交集DEGs的蛋白互作（protein-protein interaction，PPI）分析将HLJDD-AS交集DEGs上传至STRING数据库（https：//string-db.org/），互作得分设置为最高置信度≥0.400，进行PPI分析，并用Cytoscape3.7.2软件进行可视化。

1.8 Hub基因的筛选和Module分析利用MCODE插件 以 默认参数（degree cut-off=2，node score cut-off=0.2，Max depth=100，k-score=2）为标准对PPI网络进行Module分析，筛选最重要的子模块。Hub基因由CytoHubba插件中的MCC、DMNC、EPC、Degree、MNC算法排序前10的基因的重叠基因组成。

1.9 分子对接 将Hub基因与其所对应HLJDD的活性成分进行分子对接。首先，从PDB数据库（http：//www.rcsb.org/pdb/）中检索人源核心靶点的PDB ID号，用PyMOL去除水分子、配体等，Auto dock Tools加氢，然后转换成pdbqt格式。其次，从PubChem数据库中下载活性成分的分子结构，用ChemBio3D软件进行优化。最后，利用Autodock Vina软件进行分子对接，获取结合能值，对接结果采用Discovery Studio Visualizer 4.5软件进行可视化。结合能数值为负则提示活性成分与Hub基因可自由结合，数值越小说明对接体系越稳定，对接结果越可靠，本研究设定小于-5 kcal/mol的结合能代表结合具有显著性，以此佐证活性成分与核心靶点之间作用关系的可信度。

2 结果

2.1 各数据集DEGs的筛选差异分析的结果显示，E-MTAB-2055筛选出1 152个DEGs，其中上调572个，下调580个；GSE41571筛选出558个DEGs，其 中 上调126个，下调432个；GSE120521筛 选 出1 974个DEGs，其中上调1 152个，下调822个。图1为3个数据集的火山图和热图。

2.2 稳健DEGs的筛选及其富集分析RRA方法对3数据集进行整合分析后，共筛选出864个稳健DEGs，其中上调453个，下调411个。图2A为前20个上调和下调稳健DEGs的热图。GO分析结果显示，依据P＜0.05确定了793个GO条目；生物学过程（biological process，BP）条目为620条，主要涉及中性粒细胞的脱颗粒、活化、迁移、趋化作用、炎症反应的调节、中性粒细胞介导的免疫反应、白细胞迁移等；分子功能（molecular function，MF）条目为64条，主要涉及胶原结合区、受体配体活性、趋化因子的活性等；细胞组分条目（cellular component，CC）为109条，主要涉及含细胞外基质的胶原、溶酶体腔、三级颗粒膜等（图2B）。KEGG分析结果显示，依据P＜0.05共映射出56条信号通路，包括吞噬体、溶酶体、细胞因子受体相互作用、趋化因子信号通路、PPAR信号通路、细胞黏附分子等（图2C）。

2.3 HLJDD活性成分及其作用靶点在TCMSP数据库中以OB≥30%、DL≥18%为条件进行筛选，黄芩获得36个活性成分，黄连为14个，黄柏为37个，栀子为15个，合并去重后为84个（表2）。将HLJDD的84个活性成分在TCMSP数据库和Swiss Target Prediction数据库进行靶点预测，剔除重复项后共获得928个。

2.4 HLJDD-AS交集DEGs的获取及其调控网络的构建将RRA方法鉴定出的864个稳健DEGs，与网络药理学方法预测的HLJDD的928个作用靶点取交集，结果共有42个交集DEGs，其中上调32个，下调10个（图3）；即HLJDD通过抑制32个上调DEGs表达，促进10个下调DEGs表达，延缓或抑制AS由稳定斑块向破裂斑块进展（图4、表3）。

图1 各数据集的火山图和热图Fig.1 Volcano maps and heat maps of each data set

图2 稳健DEGs及其富集分析Fig.2 Robust DEGs and its enrichment analysis

表2 HLJDD“中药-活性成分-预测靶点”基本信息Tab.2 Basic information of HLJDD'Traditional Chinese Medicine-Active Ingredients-Predictive Target'

图3 HLJDD作用靶点与稳健DEGs交集基因的韦恩图Fig.3 Venn diagram of intersection genes between targets of HLJDD and robust DEGs

2.5 交集DEGs的功能富集分析GO分析结果显示，依据P＜0.05确定了392个GO条目，主要涉及胶原蛋白的分解代谢、组织重构、中性粒细胞的脱颗粒、活化、中性粒细胞介导的免疫反应、细胞外基质分解等（图5A）。KEGG分析结果显示，依据P＜0.05共映射出25条信号通路，主要包括溶酶体、钙信号通路、PPAR信号通路、补体与凝血级联反应等（图5B）。

2.6 交集DEGs的PPI网络构建、Hub基因的筛选和Module分析将交集DEGs导入STRING数据库中进行PPI分析后，再用Cytoscape软件绘制PPI网络。如图6A所示，该网络包括35个节点、1 058条边。通过MCC、DMNC、EPC、Degree、MNC算法，共筛选出7个Hub基因，模块分析共筛选出2个子模块，其中大部分Hub基因在模块1中（图6B～D）。

2.7 分子对接结果为了验证HLJDD调控交集DEGs的可能性，本研究将Hub基因与其所对应的活性成分进行了分子对接。结果表明，提取出来的36对“活性成分-作用靶点”关系对中，HLJDD活性成分与7个Hub基因均存在不同程度的结合，结合自由能均小于-5 kcal/mol，提示两者之间具有良好的亲和力，结果具有一定的可信度（表4）。现将结合能小于-10对接结果进行可视化展示，见图7。

图4 HLJDD干预AS的“中药-成分-靶标-疾病”网络图Fig.4 Network diagram of'Traditional Chinese Medicine-Ingredients-Target-Disease'of HLJDD Intervening AS

图5 交集DEGs的功能富集分析Fig.5 Functional enrichment analysis of intersection DEGs

图6 交集DEGs的PPI网络构建、Hub基因的筛选和Module分析Fig.6 PPI network construction，Hub gene screening and Module analysis of intersection DEGs

表3 HLJDD调控交集DEGs在AS进展中的作用Tab.3 Role of HLJDD in regulating and intersecting DEGs in progress of AS

表4 HLJDD活性成分与作用靶点分子对接的结果Tab.4 Results of docking of active components of HLJDD with target molecules

图7 对接结果模拟图Fig.7 Simulation diagram of docking result

3 讨论

AS是多种心脑血管疾病的共同病理基础。所以AS的发生发展机制一直是研究热点，尤其斑块破裂的发生机制。随着基因芯片和高通量测序技术的不断发展，对于斑块破裂的作用及机制的相关研究日渐深入。鉴于动脉斑块由稳定到破裂的发生机制的复杂性，筛选出与斑块破裂相关的有价值的差异基因对于识别诊断、治疗及预后尤为必要。

本研究通过RRA法整合3个数据集后，共筛选出864个稳健DEGs，其中上调453个，下调411个。这些稳健DEGs主要通过补体和凝血级联信号通路、铁死亡、IL-17信号通路、白细胞跨内皮迁移、溶酶体途径、PPAR信号通路、抗原处理和呈递、TGF-β信号通路等参与斑块破裂的病程变化。与斑块破裂相关的重要稳健DEGs中，目前有研究证实补体C1q子成分A（C1QA）、补体C1q子成分B（C1QB）等作为补体系统的一部分，当补体C1q水平的显著降低时可能是导致AS斑块不稳定或破裂的一个促成因素［9］，但其通过与ApoE形成复合物可以减轻无法解决的炎症［10］；血红素加氧酶-1（heme oxygenase-1，HMOX1）与纤维帽较薄和血栓形成破裂斑块有关［11］；基质金属蛋白酶表达水平与斑块稳定程度及斑块纤维帽厚度呈负相关［12］，其中基质金属蛋白酶9（matrix metalloproteinase 9，MMP9）通过炎症反应在介导斑块破裂中发挥重要作用［13］；纤维帽中G蛋白信号通路调节蛋白5（G protein signaling pathway 5，RGS5）的缺失促进了斑块破裂的进程［14］；脂肪细胞增强结合蛋白1（adipocyte enhancer binding protein 1，AEBP1）是一种新型巨噬细胞促炎介质，可以通过介导巨噬细胞浸润加速AS病程变化［15］。与斑块破裂相关的信号通路中，研究发现铁死亡是一种主要由铁依赖性脂质过氧化介导的调节性细胞死亡，抑制铁死亡能够减轻内皮细胞的脂质过氧化和内皮功能障碍，减轻AS程度［16］；脂质过氧化和铁沉积是晚期AS斑块的共同特征，故推测铁死亡参与了斑块不稳定的发展进程。IL-17是由T细胞或一个独特的辅助性T细胞亚群产生的细胞因子，可通过增加血管平滑肌细胞的凋亡水平，促进易损斑块的破裂［17］；过氧化物酶体增殖物激活受体（peroxisome proliferator activated receptor，PPAR）是核受体超家族的成员，作为配体诱导转录因子，研究发现通过AMPK-PPARα通路介导巨噬细胞极化参与了AS斑块稳定性调节［18］；TANG等［19］研究发现通过加剧TGF-β/ROCK1调节的内皮间质转化进展，亦可以导致颈动脉斑块的脆弱性。此外，有报道指出在类风湿关节炎和Ⅰ型糖尿病患者中，AS负担和斑块破裂的脆弱性增加［20-21］。总之，以上相关文献的报道从侧面证实了本研究结果的可靠性和科学性。

现代医学从微观层面可以阐释出斑块破裂的发生机制与重要相关基因，而传统医学有别于现代医学的精准分子致病理论。AS属于发生在脉管上的病变，因于“毒邪最易腐筋伤脉”，与AS斑块破溃、出血及大量炎症因子浸润等病理变化有关，故现代医家认为斑块破裂的致病因素乃是中医“毒”邪［22］。既有“毒者，邪气蕴结不解之谓”；又有“其久心痛者，是心之支别络脉，为风邪冷热所乘痛也，故成疹不死，发作有时，经久不瘥也。”可知毒邪有内毒与外毒之分，而现代医家偏向于AS的发生与内毒紧密相关。内毒主要因脏腑功能失衡，气血紊乱，使得机体内代谢产物不能及时重新利用或排出，累积日久，败坏形体而化生。如陈可冀院士根据中医关于瘀、毒致病特点，具有广泛性、顽固性、趋内性、趋本性、兼夹性等，类似于AS易损斑块炎性反应及血栓的形成，提出“瘀毒致易损斑块”的观点［23-24］。王永炎院士根据“毒损络脉”的理论学说及病程进展，将病络分为络弛、络结、络破三种状态，其中络破即为斑块破裂，指络壁受损，脉中血液溢出的病理状态［25-26］。此外，苏文全等［27］已提出了相似学说“脉郁毒损”，即脉郁生毒、毒损脉络，指出毒损是AS病程中发展变化的关键环节。不同医家根据AS的病机演变提出了瘀毒、痰毒、热毒学说，而热毒是促进斑块破裂的病机基础［28］。故有学者倡导在动脉斑块破裂之际应加重清热解毒通脉之品［29］。

将HLJDD与稳健DEGs取交集靶点，发现HLJDD通过抑制32个上调DEGs的表达，促进10个下调DEGs的表达，延缓或抑制AS由稳定斑块向破裂斑块的进展。经GO富集及KEGG富集分析，发现HLJDD可通过调节胶原代谢、组织重构、免疫反应、血液循环调节等，参与钙信号通路、细胞凋亡、流体剪应力与AS、PPAR信号通路及Apelin信号通路等信号途径。目前研究证实抑制KCa3.1通道可以有效延缓动脉斑块破裂的进展［30］；Apelin-APJ信号参与调节类风湿关节炎斑块稳定性介质水平的改变［20］。此外，研究发现与无斑块破裂动脉相比，有斑块破裂的动脉组织中存在高壁剪切应力（wall shear stress，WSS），其导致斑块破裂的机制可能与诱导内皮细胞行为的特异性改变及加剧炎症和刺激AS脂质核心有关［31-32］。这些不仅补充了HLJDD对斑块破裂的治疗靶点，还展示出其以多环节、多途径方式参与斑块破裂的调节作用，同时以方测证，以证验因，上述靶点和功能富集结果可能是AS中医“毒损脉络”“毒致易损斑块”的科学内涵，为进一步中医学从毒论治AS提供了客观依据。

经Hub基因的筛选分析，推测HLJDD主要通过下调基因MMP8、CTSD、SPP1、PLAU、CTSB、CTSS、CTSL，发挥抑制或缓解斑块破裂病理过程。其中基质金属蛋白酶8（matrix metalloproteinase 8，MMP8）不仅参与胶原降解，提高斑块易损指数，还被鉴定为防止斑块破裂的良好靶标［33-34］。CTSD、CTSB、CTSS、CTSL属于组织蛋白酶，通过凝血、免疫反应、补体激活、脂肪生成、凋亡等一系列途径，对改善AS进展或斑块破裂发挥至关重要的作用［35-36］。重组人分泌型磷蛋白1（secreted phosphoprotein 1，SPP1）又称骨桥蛋白，是一种分泌型多功能糖蛋白，主要由内皮细胞、血管平滑肌细胞和巨噬细胞分泌，参与矿物质沉积和斑块钙化生物过程，研究发现SPP1在破裂斑块中的表达比稳定斑块中更为显著［37］。尿激酶-纤溶酶原激活物（urokinase-plasminogen activator，PLAU），亦称为uPA，是一种丝氨酸蛋白水解酶，通过调控巨噬细胞浸润参与AS斑块破裂［38］。目前研究证实HLJDD通过降解胶原纤维、改善血液流变学紊乱、促进胆固醇流出及调节巨噬细胞炎症浸润等延 缓AS进 程［39-42］。基 于Hub基 因 的 筛 选 可 知，HLJDD不仅可以通过下调这些关键靶点改善斑块破裂病理变化，还进一步深化了HLJDD稳定斑块及延缓斑块破裂的分子机制。

通过分子对接技术进一步预测Hub基因与其所对应的活性成分的结合关系，结果显示36对“活性成分-靶点蛋白”关系对的结合自由能均小于-5.0 kcal/mol，这提示HLJDD活性成分与斑块破裂的潜在靶点之间具有良好的结合力；其中kihadanin A、芸香碱（rutaecarpine）、（R）-氢化小檗碱［（R）-Canadine］、（S）-氢化小檗碱［（S）-Canadine］、棕榈酸酯A（Palmidin A）与其对应靶点蛋白CTSB、CTSS、CTSL、MMP8的结合自由能均小于-10 kcal/mol，表现出最强的结合能力，具有进一步深入开发研究的价值。

本研究运用生物信息学的方法发现了在AS由稳定斑块发展至破裂斑块过程中的一些DEGs以及HLJDD发挥延缓或抑制功效的作用机制，但仍有一些局限性。本文虽运用RRA法分析了多个数据集，但由于数据库中可获得的数据集有限，同时每个数据集中满足研究条件的样本较少，最终导致本研究纳入的样本量有限；生信预测的结果虽进行了分子对接实验的验证，但其结果仍需要进一步在后续的研 究 中 进 行 佐 证。总 之，MMP8、CTSD、SPP1、PLAU、CTSB、CTSS、CTSL在AS破裂斑块中发挥着重要作用，HLJDD可通过抑制其表达来延缓或抑制AS由稳定斑块发展至破裂斑块，研究结果可为进一步探讨AS破裂斑块的分子机制以及HLJDD干预其进展提供新的思路和切入点，并为中医从毒论治AS易损斑块提供一定的科学依据。