二仙汤防治骨质疏松症的网络药理学研究

王明亮 闫小龙 田思胜

1. 山东中医药大学，山东 济南 250355 2. 日照市中医医院，山东 日照 276800 3. 山东中医药大学附属医院，山东 济南 250014

随着世界人口老龄化程度的增高，骨质疏松症(osteoporosis，OP)的发病率也在逐年上升[1]。中医理论认为，肾藏精、主骨生髓，骨的生长发育有赖于肾中精气的滋养与推动[2-3]，肾精不足则骨缺乏滋养。二仙汤是临床上治疗OP的经典名方[4]，其中仙茅、淫羊藿、巴戟天温补肾阳，滋补肾精；佐以知母、黄柏滋补肾阴、泻肾火；当归养血活血，调理冲任为使药[5]，全方补肾壮阳药与滋阴泻火药并用，适合用于阴阳俱虚的复杂症侯[6]。临床和动物实验均证实，二仙汤能有效防治OP，改善OP患者以及OP模型动物骨密度、骨小梁厚度[7-8]，提高体外培养的成骨细胞增殖能力和骨髓间充质干细胞的成骨转化[9]。但是，中药复方化学成分多、结构复杂以及有效成分不确定等因素的存在，使中药复方作用机理的研究受到了限制，二仙汤治疗OP的具体分子机制尚未清楚。网络药理学是在系统生物学与计算机技术基础上通过网络分析，系统综合观察药物对疾病的干预与影响，与中医学整体观念理论相契合，能够对中药复方及方中化学成分进行系统分析[10-12]。本研究应用网络药理学方法，对中药经典名方二仙汤治疗OP的分子机制进行预测，为中药复方治疗机制的研究提供科学可靠的方法。

1 材料和方法

1.1二仙汤中活性成分筛选和靶点预测

中药系统药理学数据库及分析平台(Traditional Chinese Medicine Systems Pharmacology Database and Analysis Platform, TCMSP, 版本2.3, http://ibts.hkbu.edu.hk/LSP/tcmsp.php)能分析所含化合物的药物动力学(ADME：absorption，distribution，metabolism，excretion)的相关参数，包括口服吸收利用度(oral bioavailability，OB)、药物相似性(drug-likeness，DL)、Caco-2细胞渗透性、半衰期、血脑屏障、脂水分布系数等[13]。本研究分析二仙汤有效活性成分时，设定筛选条件为OB≥30%，DL≥0.18；同时应用TCMSP预测活性成分的作用靶点[14]。

1.2二仙汤“活性成分-靶点网络”的构建

将获得的二仙汤活性成分、活性成分作用靶点导入Cytoscape(版本3.7.1，https://cytoscape.org/)软件中，构建“活性成分-靶点网络”，应用network analyzer计算网络中节点参数以分析二仙汤中重要的成分和靶点。

1.3OP发病相关基因筛选

本研究以“Osteoporosis”为关键词检索网络已知公共疾病靶点数据库，包括TTD (Therapeutic Target Database，更新至2017.09.15, http://bidd.nus.edu.sg/group/cjttd/); DrugBank 数据库 (版本5.1.4, https://www.drugbank.ca/); OMIM数据库(Online Mendelian Inheritance in Man, https://www.omim.org/); GAD databases (Genetic Association Database,更新至2014.08.18, https://geneticassociationdb.nih.gov/)；PharmGKB数据库(https://www.pharmgkb.org/); DisGeNET 数据库 (版本6.0, http://www.disgenet.org/home/)，删除各数据库中重复靶点，获得已知的OP发病靶点。

GEO数据库(Gene Expression Omnibus, https://www.ncbi.nlm.nih.gov/geo/)是一个由基因表达数据、芯片和微阵列组成的公共功能基因组数据库。以“Osteoporosis”为关键词检索GEO数据库，筛选获得GSE56116数据集，包含肾阴虚组样品4个(GSM1356155，GSM1356157，GSM1356158，GSM1356161)，肾阳虚组样本3个(GSM1356160，GSM1356163，GSM1356164)，非肾虚组样本3个(GSM1356156，GSM1356162，GSM1356165)，对照组样品3个(GSM1356159，GSM1356166，GSM1356167)。利用R语言 limma包分析芯片数据的差异基因，筛选条件设定为P值<0. 05，差异倍数(fold change，FC)>2。

1.4PPI网络构建

构建蛋白质与蛋白质相互作用(protein-protein interaction, PPI)网络，有助于挖掘核心调控基因、深入了解疾病的发生和发展机制。本研究利用STRING数据库(版本11.0，https://string-db.org/)[15]对PPI网络进行预测，综合得分大于0.4的交互作用被认为具有统计学意义。为进一步分析网络中核心靶点，采用Cytoscape进行网络拓扑分析，运用自由度(degree，DC)、介度中心性(betweenness centrality，BC)、接近中心性(closeness centrality，CC)等指标进行筛选，其中DC值、BC值、CC值大于所有节点对应平均值的节点为网络中的重要节点。

1.5关键作用靶点KEGG和GO富集分析

表 1 预测出靶点活性成分的基本信息Table 1 Basic informaiton for predicting target active components

注释、可视化和集成发现数据库(Database for Annotation, Visualization and Integrated Discovery, DAVID, 版本6.8, https://david.ncifcrf.gov/)[16]是一个在线生物信息数据库，用于从高通量实验技术生成的大规模分子数据集中了解高级功能[17]和注释基因的生物学过程[18]。为了分析关键靶点的功能，采用DAVID数据库进行功能富集分析。

2 结果

2.1二仙汤活性成分筛选与ADME分析

通过TCMSP数据库筛选二仙汤中6味中药入血的活性成分，其中仙茅7个，淫羊藿23个，当归2个，巴戟天20个，黄柏27个，知母15个，删除重复成分，得到84个活性成分。二仙汤各味中药中符合筛选标准的成分比例图，见图1。

图1 二仙汤活性成分ADME参数分布Fig.1 ADME parameter distribution of active components in Erxian decoction

2.2二仙汤“活性成分-靶点网络”构建

应用Cytoscape构建二仙汤活性性成分和作用靶点网络以寻找网络中的关键节点，确定二仙汤有效成分中发挥重要作用的成分和靶点。结果如图2和表1所示，构建的网络中，共有84个活性成分，21个活性成分采用TCMSP预测出145个作用靶点。

图2 中药-活性成分-靶点网络图Fig.2 Traditional Chinese Medicine-active component-target network diagram

红色菱形代表中药，绿色倒三角代表预测出靶点的活性成分，绿色正方形代表未预测出靶点的活性成分，蓝色圆形代表作用靶点，圆形大小代表degree值的大小，圆形越大表示degree值越大

2.3OP相关致病靶点

疾病相关数据库中GAD、TTD、OMIM、DrugBank、PharmGKB和DisGeNET分别检索到275、27、94、48、5和426个靶点，去除重复项，共纳入691个与OP发生发展相关的已知靶点。

图3 基于GSE56116芯片OP致病靶点的聚类图和火山图Fig.3 Clustering map and volcano map of OP pathogenic targets based on GSE56116 chip

利用limma包分析GEO芯片数据库中GSE56116数据集，共找出84个变化最为明显的基因，其中上调表达的有70个，下调表达的有14个，这些基因可能与OP致病密切相关。差异表达基因见聚类图(图3 A)和火山图(图3B)所示。

2.4PPI网络的构建

通过二仙汤活性成分作用靶点与OP致病靶点直接映射，获得关键作用靶点40个，如图4所示。采用String构建关键靶点的PPI网络，结果如图5所示，共有40个节点，278个相互作用关系。通过DC值、BC值和CC值进行网络拓扑筛选，获得核心节点9个，见表2。

表2 二仙汤治疗OP核心靶点的基本信息Table 2 Basic information of Erxian decoction in treating OP core targets

图4 OP关键靶点的筛选图Fig.4 Screening of OP key targets

图5 关键靶点蛋白质相互作用网络Fig.5 Key target protein interaction network

图5中，节点的大小代表节点的自由度(degree)，节点越大表示degree值越大；线条表示各个靶点间的相互作用关系，条线粗细表示相互作用关系的紧密程度，线条越粗表示相互作用越紧密。

2.5关键靶点的KEGG和GO富集分析

应用DAVID对上述关键基因进行功能富集分析。GO分析结果表明，关键基因富集的生物学过程变化在RNA聚合酶II启动子对转录的正调节作用(positive regulation of transcription from RNA polymerase II promoter)、一氧化氮生物合成过程的正调控(positive regulation of nitric oxide biosynthetic process)、序列特异性结合转录因子活性的正调控(positive regulation of sequence-specific DNA binding transcription factor activity)等，见图6 A。关键基因富集细胞组分主要为细胞外间隙，胞外区，核浆和质膜等，见图6B。分子功能主要富集在酶结合(enzyme binding)、类固醇激素受体活性(steroid hormone receptor activity)、甾体结合(steroid binding)，见图6C。KEGG通路分析结果所示，关键靶点涉及的非疾病命名信号通路主要富集在TNF信号通路(TNF signaling pathway)，低氧诱导因子1信号通路(HIF-1 signaling pathway)，雌激素信号通路(Estrogen signaling pathway)，血管内皮生长因子信号通路(VEGF signaling pathway)等通路过程，见图6D和表3。

图6 关键靶点的富集分析图注：A 生物学进程，B 细胞组分，C 分子功能，D KEGG信号通路。Fig.6 Enrichment analysis plot of key targets

KEGG分析结果如表3所示，TNF信号通路富集基因数目及P值均较为显著，分析其可能为二仙汤治疗OP的重要途径。其具体作用机制可能是通过调控细胞凋亡、炎症活动、细胞外基质重塑、细胞粘附和细胞存活等过程，发挥治疗OP的作用，见图7。

表3 OP相关基因的KEGG分析Table 3 KEGG analysis of OP-related genes

图7 TNF信号通路中富集基因具体作用机制图Fig.7 Specific mechanism of enriched genes in TNF signaling pathway

3 讨论

由于全球人口老龄化，估计OP的发病率将在未来20年翻一番，给公共卫生服务造成沉重的负担[19]。目前西药治疗OP疗效显著[20]，但是药物治疗的不良反应和长期依从性差的问题未得到有效解决[21]。中草药治疗骨病已有几千年的历史，具有不良反应少，疗效显著的特点，受到越来越多的应用[22]，但是具体的分子机制尚未清楚。

本研究基于网络药理学的方法，从TCMSP数据库中筛选出二仙汤的活性成分及作用靶点，检索OP的致病靶点，构建PPI网络并映射核心靶点，进行核心靶点的KEGG分析和GO富集分析，结果显示，TNF信号通路、低氧因子诱导因子1信号通路和雌激素受体信号通路可能为二仙汤治疗OP的关键途径。

骨髓间充质干细胞是一种成体干细胞，具有与其他组织相同的遗传基础[23]，BMSCs的增殖和分化是成骨的关键[24]。炎性细胞因子如肿瘤坏死因子(TNF-α)、白细胞介素1(IL-1)、白细胞介素-6(IL-6)和干扰素-γ(IFN-γ)也可影响骨髓间充质干细胞的增殖和成骨分化[25]。TNF-α是一种重要的炎症细胞因子，通过激活多种信号通路和分子，调节细胞的增殖、分化和凋亡[26]。Lu等[27]证明，由于转录因子核因子κ-α(NF-κB)的激活，TNF-α抑制了成骨细胞与前体细胞的分化，降低了RUNX 2和成骨细胞相关转录因子(Osterix)的表达，抑制了维生素D刺激的转录。Hess等[28]研究表明，TNF-α通过NF-κB信号通路促进BMSCs中BMP-2的表达，从而促进细胞外基质的矿物质化。

缺氧对多种组织和细胞类型的影响主要由HIF-1介导[29]，其中HIF-1α在骨发育和再生过程中起着重要作用。HIF-1α通路是适应低氧反应的中枢调节因子，这一途径调控血管生成基因如VEGF、血管生成素等[30]。骨髓腔和骨骺生长板的骨内膜区都是缺氧区，HIF-1α信号通路控制着这些区域的软骨细胞和成骨细胞的功能[31]。研究结果表明，在发育中的成骨细胞中激活HIF-1α通路，增强了VEGF介导的骨血管形成，增加了骨的生成[32]。

4 结论

基于网络药理学的对临床常用中药复方二仙汤治疗OP的分子机制进行了预测，揭示了二仙汤可能是通过TNF信号通路和HIF-1信号通路治疗OP，为中药复方的推广应用以及中医药现代化研究提供了理论依据。