基于网络药理学和分子对接探讨千层纸素A治疗椎间盘退变的分子机制

刘明阳，赵翔，张宇，张锴广，高延征

(河南省人民医院脊柱脊髓科，郑州大学人民医院，河南大学人民医院，河南 郑州 450000)

椎间盘退变是肌肉骨骼系统的一种高发疾病，常导致椎管狭窄、椎间关节不稳定、疼痛甚至畸形[1]。其作为椎间盘源性腰痛的重要原因，常常给家庭和社会带来巨大的经济负担[2]。椎间盘退变的发病机制与髓核细胞凋亡、衰老、炎症、氧化应激、细胞外基质降解以及异常血管形成密切相关[2-4]。目前，椎间盘退变的药物治疗旨在缓解症状，而不是阻止或逆转椎间盘受损的病理进展。千层纸素A是从黄芩中提取的一种黄酮类化合物，其化学全名为5，7-二羟基-6-甲氧基-2-苯基-4H-1-苯并呋喃-4-酮[5]。据报道，千层纸素A具有广泛的药理作用，其中抗炎和抗癌作用最受关注，其他作用包括抗病毒、神经保护和抗血管生成作用等[6]。千层纸素A对于包括骨性关节炎、类风湿关节炎和骨质疏松在内的多种肌肉骨骼系统的疾病均显示出良好的治疗作用[7-9]。考虑到其显著的抗炎作用，千层纸素A可能是治疗椎间盘退变的潜在选择。网络药理学将分子化合物与药理学、系统生物学与生物信息学结合起来构建“化合物-靶点-疾病”的网络预测模型，能够系统地阐述化合物对于疾病的潜在作用靶点和相关通路[10]。现阶段对千层纸素A治疗椎间盘退变的作用靶点以及其发挥作用的分子机制的研究仍然有限。因此，本研究使用网络药理学方法，进一步预测千层纸素A治疗椎间盘退变的作用靶点和其作用机制，为后续开展千层纸素A相关的药理研究及椎间盘退变防治提供参考依据。

1 材料与方法

1.1 千层纸素A靶点基因的筛选在PubChem数据库(https://pubchem.ncbi.nlm.nih.gov/)搜索关键字Oroxylin A并下载结构式和SMILES式，然后将结构式导入PharmMapper数据库(http://www.lilab-ecust.cn/pharmmapper/)，SMILES式导入SwissTarget Prediction(STR)数据库(http://www.swisstargetprediction.ch/)和Similarity ensemble approach(SEA)数据库(https://sea.bkslab.org/)以获取药物作用靶点。STR数据库使用Probability>0为过滤条件筛选数据。PharmMapper数据库使用norm fit score≥0.7为过滤条件筛选数据，并使用UniProt数据库(http://www.Uniprot.org/)将从PharMapper数据库下载的数据转换为相应的基因。SEA数据库不设置过滤条件。3个数据库物种选择均为人类。随后将3个数据库获取的基因合并，去除重复的基因，即得到千层纸素A的潜在的药物作用靶点

1.2 千层纸素A在椎间盘退变中相关靶点的筛选以“Intervertebral disc degeneration”为关键词，在GeneCards数据库(https://www.genecards.org/)和在线人类孟德尔遗传(Online Mendelian Inheritance in Man，OMIM)数据库(https://omim.org/)中搜索椎间盘退变相关基因靶点。将2个数据库的靶点合并得到椎间盘退变的潜在致病基因。将椎间盘退变相关靶点和千层纸素A的药物靶点取交集，获得二者的共同基因即为千层纸素A抗椎间盘退变的潜在关键靶点。

1.3 蛋白相互作用(protein-protein interaction，PPI)网络构建PPI网络的构建对于研究蛋白质功能至关重要，它可以更加系统、全面地研究疾病分子机制、发现新药靶点等。将千层纸素A对应的椎间盘退变相关基因导入在线STRING数据库，将物种设置为“人类(homo sapiens)”，可信度设置为0.4，选择“hide disconnected nodes in the network”删除无相互作用的蛋白质。下载蛋白质互作的tsv文件并将其导入到Cytoscape 3.8.0软件中，运用Analysis network工具进行网络分析。度值显示与单个节点连接的节点数，其可作为反映靶点基因关键程度的指标。根据度值的高低调整节点的大小和颜色的深浅。节点度值越高，则节点越大，颜色越深。边的粗细根据combine_score进行调整。combine_score值越高，则边越粗。最后利用R语言的“count R”包处理tsv文件，按照度值排列得到蛋白质互作网络相关核心基因的条形图。

1.4 基因本体(Gene Ontology，GO)及京都基因与基因组百科全书(Kyoto Encyclopedia of Genes and Genomes，KEGG)通路富集分析在R语言中导入千层纸素A抗椎间盘退变相关的潜在靶点，用“clusterProfile”包进行GO功能富集分析及KEGG通路富集分析并进行可视化展示，筛选条件设定为P<0.05以及q<0.05。得到千层纸素A抗椎间盘退变的生物学过程(biological process，BP)、分子功能(molecular founction，MF)和细胞组成(cellular component，CC)以及KEGG信号通路信息。通过“ggplot”包绘制气泡图，用Cytoscape 3.8.0软件对靶点-通路进行网络分析，筛选出度值高的靶点。将靶点-通路网络和PPI网络中筛选出的度值较高的靶点取交集，得到关键基因。

1.5 分子对接分子对接是一种用于预测生成结合模型的两个分子之间的相互作用的方法，可以估算出小分子配体和生物大分子受体的结合能大小，为药物设计提供新的思路。将千层纸素A作为配体，关键基因编码的蛋白作为受体。在PubChem数据库中下载配体的二维结构文件，利用Chem3D软件进行化合物结构优化并保存为mol2格式。从蛋白质结构数据库中获得大分子受体的蛋白结构并用PyMol软件对受体蛋白进行前处理。使用AutoDock软件给蛋白质分子加氢和加电子并将配体和受体文件转换为pdbqt格式，对受体活性口袋进行定义并保存为gpf文件。利用AutoDock vina进行分子对接。一般认为自由能在-4～0 kcal·mol-1之间代表结合力较弱，-7～-4 kcal·mol-1之间代表结合力中等；自由能小于-7 kcal·mol-1代表结合力较强[11]。最后通过PyMol软件对对接结果进行可视化。

2 结果

2.1 靶点筛选将从STR、SEA、PharmMapper数据库中得到的药物作用靶点去重后共得到202个靶点基因。然后将GeneCards数据库和OMIM数据库搜索到的椎间盘退变相关基因靶点去重得到1 212个致病靶点。利用R语言将千层纸素A的药物作用靶点与椎间盘退变的致病靶点进行匹配取交集，最终得到35个共同靶点。见图1。

图1 千层纸素A与椎间盘退变靶点图

2.2 蛋白互作网络根据节点的度值将靶基因分布为2个同心圆，外层靶基因的度值为1～15，内层靶基因的度值为16～25，得到千层纸素A与椎间盘退变相关主要靶点PPI网络图。经R语言的“count R”包得到按核心基因度值排列的条形图。其中排名靠前的核心靶点有白细胞介素-6(interleukin 6，IL-6)、血管内皮生长因子A(vascular endothelial growth factor A，VEGFA)、表皮生长因子受体(epidermal growth factor receptor，EGFR)、雌激素受体1(estrogen receptor 1，ESR1)、前列腺素内过氧化物合酶2(prostaglandin-endoperoxide synthase 2，PTGS2)、基质金属蛋白酶9(matrix metallopeptidase 9，MMP9)等。见图2。

A为千层纸素A治疗椎间盘退变的蛋白互作网络；B为度值排名前15的基因。

2.3 GO与KEGG分析GO分析发现千层纸素A治疗椎间盘退变的靶点在生物过程主要包括对炎症反应的调节、对缺氧的反应和细胞对化学应激的反应等。在细胞组成主要富集于染色体区、线粒体基质、树突棘等。在分子功能主要富集于金属蛋白酶活性、生长因子受体结合等。见图3。

图3 GO功能富集分析

KEGG分析显示相关信号通路主要涉及雌激素信号通路、PI3K/AKT信号通路、IL-17信号通路、TNF信号通路以及糖尿病并发症中的AGE-RAGE信号通路等。见图4。随后将千层纸素A抗椎间盘退变的预测靶基因以及根据 KEGG通路的富集分析结果筛选出来的通路导入Cytoscape 3.8.0软件进行“靶标-通路”分析，展现二者之间的联系。见图5。经过与PPI网络分析得到的关键靶基因对照，筛选出IL-6、VEGFA、EGFR和MMP9共4个关键靶基因用于后续的分子对接。

图4 KEGG功能富集分析

图5 靶点-通路分析

2.4 分子对接为进一步探索千层纸素A与关键靶点之间的相互作用关系，IL-6(PBD ID：1IL6)、VEGFA(PBD ID：1VPF)、EGFR(PBD ID：1XKK)和MMP9(PBD ID：2OVX)被用来与千层纸素A进行分子对接。选取最低结合能的组合，结果表明千层纸素A与IL-6的结合能为-7.2 kcal·mol-1，与VEGFA的结合能为-7.8 kcal·mol-1，与EGFR的结合能为-8.7 kcal·mol-1，与MMP9的结合能为-10.0 kcal·mol-1。分子对接的结果显示千层纸素A与靶基因均具有良好的结合能力。见图6。

图6 分子对接示意图

3 讨论

本研究通过网络药理学的研究方法，探讨千层纸素A抗椎间盘退变的潜在分子机制。既往研究显示，千层纸素A有较强的抗炎和抗血管生成作用[6]。Zhang等[7]发现千层纸素A可以通过抑制ERK和PI3K/AKT信号通路的激活来减轻IL-1β诱导的关节软骨细胞炎症反应。而另一项研究则证实千层纸素A可以通过抑制NF-κB和Wnt/β-catenin信号通路减轻关节软骨细胞的肥大和炎症反应，从而延缓大鼠模型中骨性关节炎的进展[12]。Gao等[13]发现，千层纸素A通过阻断VEGF诱导的KDR/Flk-1(VEGFR-2)磷酸化和相关下游信号分子(包括p38丝裂原活化蛋白激酶、ERK和Akt)抑制血管生成。有关千层纸素A治疗椎间盘退变的研究虽未见报道，但考虑到炎症与异常的血管生成与椎间盘退变的发生发展密切相关[3]，千层纸素A可能是治疗椎间盘退变的潜在选择。

GO富集分析结果显示生物过程主要包括对炎症反应的调节、对缺氧的反应和对化学应力的反应等。炎症反应直接参与椎间盘退变的进展，还可导致继发性腰痛和神经根性症状。炎症反应参与椎间盘退变的各种病理过程，包括细胞外基质降解、髓核细胞的凋亡和衰老[2]。椎间盘患者椎间盘组织中促炎细胞因子如TNF-α、IL-1β和IL-6的表达显著增加[14-15]。Kwon等[16]发现缺氧可以增强炎症反应环境下椎间盘的血管生成能力，并通过抑制细胞外基质的合成促进椎间盘退变的发展。KEGG富集分析结果显示与千层纸素A治疗椎间盘退变相关的通路可能包括PI3K/AKT信号通路、雌激素信号通路、IL-17信号通路、TNF信号通路以及糖尿病并发症中的AGE-RAGE信号通路等。PI3K/AKT信号通路被证明与椎间盘退变的疾病进程密切相关。Jiang等[17]发现白藜芦醇可以通过激活PI3K/Akt通路抑制IL-1β诱导的髓核细胞凋亡，而雌激素可以通过激活PI3K/Akt通路来减轻炎症反应，降低氧化应激水平从而抑制髓核细胞的凋亡。TNF和IL-17途径主要通过促进炎症因子的释放、髓核细胞的凋亡和细胞外基质的降解在椎间盘退变的进展中发挥协同作用[18]。此外，AGE-MG-H1的积累与椎间盘退变过程中的软骨细胞肥大和成骨分化有关，并且这些与RAGE直接相关，表明AGE-RAGE可能是一个潜在的治疗靶点[19]。

PPI网络分析和KEGG富集分析发现IL-6、VEGFA、EGFR和MMP9是千层纸素A抗椎间盘退变作用关联度最高的靶点。分子对接结果表明，千层纸素A与核心靶基因IL-6、VEGFA、EGFR和MMP9均具有良好的亲和力。多项研究显示IL-6的过表达与椎间盘退变密切相关[20-22]。作为一种重要的炎症因子，既往研究显示IL-6可通过PI3K/Akt和Ras通路调控下游基因的表达从而参与椎间盘退变的进程[20]，与本研究中KEGG富集分析的结果相符。VEGFA属于VEGF家族，是促进微血管形成的重要细胞因子[23]。罗坚等[24]发现VEGF可以促进椎间盘组织新生血管的形成和毛细血管的浸润从而导致椎间盘组织的吸收退变。本研究结果表明千层纸素A可能通过靶向作用于VEGFA从而抑制椎间盘组织异常的血管生成从而减缓椎间盘退变的进展。EGFR属于ErbB受体酪氨酸激酶家族，有研究证实EGFR可以通过Src/ERK/STAT3通路诱导髓核细胞的炎症反应[25]。而千层纸素A则可能通过作用于EGFR从而抑制下游的炎症反应。MMP9属于MMP家族，其主要水解底物为变性胶原以及Ⅳ、Ⅴ型胶原[26]。Li等[27]发现退变的椎间盘组织中MMP-9的水平高于正常水平，且升高的水平与退变程度呈正相关。千层纸素A可通过靶向作用于MMP抑制髓核细胞外基质的降解。

本研究通过网络药理学的方法分析千层纸素A与椎间盘退变靶点的相互关系。IL-6、VEGFA、EGFR和MMP9是潜在的关键靶点。结合文献资料，千层纸素A抗椎间盘退变的药理机制与抑制髓核细胞凋亡、抑制细胞外基质降解、抗炎和抗血管生成等相关通路密切相关。由于当前研究是基于数据库的预测分析，研究结果有存在偏差的可能。后续研究需结合预测的靶点及通路，进一步通过实验验证千层纸素A抗椎间盘退变的具体机制。尽管如此，本研究仍为探讨千层纸素A对于椎间盘退变治疗作用的实验及临床研究提供了新方向。