一种中药网络药理深度分析方法与系统

技术领域

本发明涉及生物信息技术领域，尤其涉及中药网络药理深度分析方法与系统。

背景技术

中药及复方中药通过多成分多靶点协同作用发挥治疗疾病作用，虽然疗效显著，但中药活性成分多，与机体相互作用复杂，作用机制不清，成为中药现代化及中药被广泛推广应用的瓶颈问题。目前，研究中药及复方中药的作用机制方法很多，大体包括化学研究、靶点研究和系统生物学3种研究模式。这些方法虽在一定程度上从不同层次阐释了中药及复方中药的作用机制，但通过单成分-单通路-单靶点的模式研究中药及复方中药的作用机制仍存在局限性。

网络药理学是基于系统生物学的理论，对生物系统的网络分析，强调对信号通路的多途径调节，这与中医药的防治疾病的整体观念和系统性特点相契合，在一定程度上，能更准确地揭示中医药的整体观念和辨证论治理论，为中药及复方中药多组分多靶点协同作用机制的研究提供了全新思路。然而，在实际应用中，网络药理学分析方法尚存在数据量大、生成网络复杂问题，难以从中准确分析中药及复方中药作用机制。

目前网络药理学预测中药作用机理通常采用以下方法，从OMIM、DrugBank、Genecards 等数据库中搜索与疾病相关的靶点，利用TCMSP、Batman等数据库或者化合物成分预测的方法得到与中药成分相关的靶点。利用所搜集的疾病相关靶点以及中药成分相关作用靶点，分别结合String、HPRD、BIND、DIP、MINT、INTACT、BIOGRID蛋白质-蛋白质相互作用(PPI) 数据库，分别构建疾病相关PPI网络及中药成分调控PPI网络，进而利用Cytoscape软件对构建的疾病相关PPI网络和中药成分调控PPI网络进行intersection交集运算，得到的交集网络，即中药干预疾病网络。进而，利用Cytoscape软件对得到的中药干预疾病网络计算其拓扑特征，根据degree等值筛选得到网络的核心节点(核心靶点)，利用DAVID等数据库对筛选得到的核心靶点进行GO功能注释及KEGG通路富集分析。KEGG通路富集分析后根据得到的p值进行排序。但是利用p值进行排序，p＜0.05的通路往往较多，而与中药干预疾病特别相关的通路却可能排序靠后，此外，利用此方法富集得到的通路较多，实验验证准确率较低。

CN110880354 A公开了一种基于群体智能的药物-靶标相互作用预测方法，包括：1.从数据库获取药物与靶标的数据，分别计算药物相似性、靶标相似性，建立全局药物与靶标相互作用网络；2.采用随机行走的方法，从全局药物与靶标相互作用网络中提取用来分类的关键特征，将网络中的药物与靶标相互作用连接权重作为类编号；3.建立基于网络特征的分类模型；4.根据相似的药物与相似的靶标绑定的原则，筛选分类的负样本；采用极限学习机作为药物-靶标相互作用分类器；采用群体智能算法对极限学习机输入层权重进行优化，群体智能算法使用全局搜索与局部搜索相结合的方法，群体智能算法为优化的混合蛙跳算法。

CN104978474 A公开了一种基于分子网络的药效评价方法及系统，方法包括：1)通过多种途径获取数据；2)构建药物效应的分子网络和疾病的分子网络；3)将药物分子相应网络映射到疾病的分子网络上，产生药物治疗效应的分子网络；4)对疾病的分子网络进行中心性分析；5)以药物治疗效应的分子网络中分子在疾病的分子网络中的中心值分布定量评价药物对疾病影响程度及其各治疗机制对疾病影响程度。

CN111462814 A公开了一种基于网络药理学分析马钱子活性成分镇痛作用的方法，包括以下步骤：S1获取马钱子的化学成分与成分作用靶点：以马钱子为关键词，利用BATMAN-TCM 数据库检索获得符合条件的马钱子活性化合物及其靶向基因；S2识别疾病作用靶点：利用 DisGeNET数据库获得与疼痛相关的虚拟目标靶标基因，马钱子活性化合物及其靶向基因进行交集，得药物-疾病交集靶基因；S3构建马钱子与疾病网络，以及筛选核心节点；S4利用DAVID 对核心节点进行GO富集和KEGG通路注释分析，预测靶点功能分布，并根据富集因子分析核心通路富集程度，解析马钱子镇痛作用机制。

因此，对于中药尤其是复方中药来说，目前方法富集的通路多，准确性差，可操作性差，仍需更加有效的方法。

发明内容

本发明的目的是克服现有技术的不足，提供一种中药网络药理深度分析方法与系统，其能够有效分析中药包括中药复方的作用机制，本发明方法富集的通路更准确、集中，可操作性强，为中药作用机制的研究以及临床应用提供更加科学客观、快捷的分析方法及系统。

本发明的技术方案概述如下：

本发明提供一种中药网络药理深度分析方法，所述方法包括如下步骤：

(1)获取疾病和中药相关靶点的数据；

(2)根据步骤(1)获取的靶点数据，分别构建疾病的蛋白质-蛋白质相互作用网络(PPI 网络)和中药成分调控的蛋白质-蛋白质相互作用网络(PPI网络)，并进行交集运算得到交集网络，筛选交集网络的关键节点，得到核心靶点；

(3)根据步骤(2)得到的核心靶点进行生物学功能注释和KEGG通路分析；

(4)根据步骤(3)的结果，对网络核心节点的度(degree)值采用归一化法处理，分析 KEGG通路t的重要性值tsig，

将网络核心节点的度(degree)值的最大值作为1，设核心靶点i归一化后的度值为i(d)，所述y为核心靶点i所在通路t上富集到y个核心靶点，所述N为通路t上共有N个靶点，所述n为网络中共有核心靶点-通路一对一关系n条；

(5)根据步骤(4)通路的重要性值tsig对富集得到的通路进行排序，对通路上靶点的上下游及影响效应重新构建网络，最后选取富集到通路上影响与疾病相关效应的靶点作为中药干预疾病作用靶点。

所述步骤(1)中的中药相关靶点是指中药复方、中药材、中药有效部位、中药有效成分相关靶点。

所述步骤(2)中的交集运算为intersection交集运算。

其中步骤(2)所述关键节点的筛选为对交集网络拓扑特征进行分析，以Degree值大于2倍中位数，进一步优选大于3、4、5、6、7、8、9、10倍，Betweenness Centrality、Closeness Centrality 大于1倍中位数为参考标准，进一步优选大于2、3、4、5、6、7、8、9、10倍。

步骤(3)中所述生物学功能注释为GO注释，所述GO注释是根据基因功能将基因分类，其包括生物过程、分子功能、细胞组分分类。

步骤(5)对通路的重要性值tsig排序后，选取排序前20的通路进一步分析，优选排序前15、 10，最优选排序前15。

本发明还提供一种中药网络药理深度分析方法的系统，该系统包括

疾病和中药相关靶点获取模块；

核心靶点分析模块，根据疾病和中药相关靶点获取模块获取的靶点数据，分别构建疾病的蛋白质-蛋白质相互作用网络(PPI网络)和中药成分调控的蛋白质-蛋白质相互作用网络 (PPI网络)，并进行交集运算得到交集网络，筛选交集网络的关键节点，得到核心靶点并保存；

生物学功能注释和KEGG通路分析模块，用于核心靶点分析模块得到的核心靶点富集；

通路重要性值tsig生成模块，根据生物学功能注释和KEGG通路分析模块富集的核心靶点，对网络核心节点的度(degree)值采用归一化法处理，分析KEGG通路t的重要性值tsig，

将网络核心节点的度(degree)值的最大值作为1，设核心靶点i归一化后的度值为i(d)，所述y为核心靶点i所在通路t上富集到y个核心靶点，所述N为通路t上共有N个靶点，所述n为网络中共有核心靶点-通路一对一关系n条；

中药药效靶点生成模块，根据通路重要性值tsig生成模块，对通路重要性值tsig排序，选取该排序通路上的靶点上下游及影响效应重新构建网络，确定富集到通路上影响与疾病相关效应的靶点为中药药效靶点。

其中疾病和中药相关靶点获取模块，疾病相关靶点的获取，采用OMIM数据库和DrugBank数据库、Genecards数据库获取疾病靶点。

其中核心靶点分析模块，根据选自HPRD、BIND、DIP、MINT、INTACT、BIOGRID中多种蛋白质-蛋白质相互作用构建疾病PPI网络和中药成分调控PPI网络

其中所述关键节点的筛选为对交集网络拓扑特征进行分析，以Degree值大于2倍中位数，进一步优选大于3、4、5、6、7、8、9、10倍；Betweenness Centrality、ClosenessCentrality 大于1倍中位数为参考标准，进一步优选大于2、3、4、5、6、7、8、9、10倍。

所述生物学功能注释为GO注释。

本发明根据影响网络核心靶点通路富集的各环节，对每一环节中的影响因素进行综合考虑，开发出一种计算分析通路t的重要性值tsig方法，进而对tsig排序靠前的信号通路上靶点的上下游及影响效应重新构建网络，根据富集到通路上影响与疾病相关效应预测中药治疗疾病的作用靶点和影响通路。该方法既保留了网络药理学系统生物学理论强调对信号通路的多途径调节，又可从复杂网络中准确分析中药及复方中药的作用机制。

说明书附图：

图1网络核心节点GO功能注释分析结果

图2重要通路所影响效应图

图3口腔黏膜组织中关键基因的表达

有益效果：

(1)本发明的方法以及系统能快速、准确分析中药干预疾病的作用靶点，富集得到的通路排序更准确，实验验证准确率显著提高。

(2)中药包括中药复方成分复杂、干预疾病的机理复杂，作用机制不清，本发明的中药网络药理深度分析方法，通过快速、准确地预测其作用靶点及影响通路，阐释多组分多靶点中药作用机制，可加快中药现代化和中药的推广应用。

具体实施方式

下面结合具体实施例对本发明作进一步的说明。

实施例1

1.实验材料

1.1疾病相关靶点的收集

利用OMIM数据库(https://omim.org)和DrugBank数据库(https://go.drugbank.com/)、 Genecards数据库(https://www.genecards.org/)检索口腔溃疡相关靶点，筛选关键词设置为“Oral ulcer”，并筛除不同数据库检索到的重复靶点。

1.2牙痛停滴丸相关靶点的收集

利用TCMSP数据库(https://tcmspw.com/tcmsp.php)检索牙痛停滴丸的主要活性成分的作用靶点信息，以“丁香”、“荜茇”、“冰片”为关键词进行检索，以口服生物利用度(Oral Bioavailability，OB)大于等于30％，类药性(Drug Like index，DL)大于等于0.18为筛选参数进行筛选，并将相关文献报道的有生物活性和药理作用的成分也纳入，检索出相应的靶点。并利用Swiss Target Prediction数据库对成分进行检索，获得每个成分的预测靶点。

1.3牙痛停滴丸治疗口腔溃疡的网络构建以及核心节点的确定

根据1.1和1.2得到的疾病相关靶点及中药成分相关靶点，利用Cytoscape软件中的 Bisogenet插件中关联的HPRD、BIND、DIP、MINT、INTACT、BIOGRID等蛋白质-蛋白质相互作用(PPI)数据库，分别构建疾病PPI网络和中药成分调控PPI网络。利用Cytoscape 软件中merge模块对构建的疾病PPI网络和中药成分调控PPI网络进行intersection交集运算，得到交集网络，即牙痛停滴丸治疗口腔溃疡的网络。

利用Cytoscape软件中NetworkAnalyzer工具对上述得到的牙痛停滴丸治疗口腔溃疡网络进行拓扑特征分析，以Degree值大于2倍中位数，Betweenness Centrality、Closeness Centrality 大于1倍中位数为筛选标准，筛选得到网络的关键节点。

1.4 GO功能注释及KEGG通路富集分析

生物学信息注释数据库(DAVID，https://david.ncifcrf.gov/)是为大规模的基因或蛋白列表提供系统综合的生物功能注释信息，能够找出最显著富集的生物学注释。以人为物种限定条件，对1.3分析得到的网络的关键节点进行GO功能富集和KEGG通路分析，在GO富集分析和KEGG富集分析结果中，P值越小，显著性越明显。

1.5 KEGG通路重要性分析及中药可能作用机制的预测

由于degree在网络核心靶点分析中占有重要地位，我们根据得到的网络核心节点的度 (degree)值进行归一化处理，将最大值作为1，设核心靶点i归一化后的度值为i(d)，核心靶点i所在通路t上富集到y个核心靶点，即通路t上含核心靶点-通路一对一的关系共y条，通路t上共有N个靶点，网络中共有核心靶点-通路一对一关系n条，则通路t的重要性值(tsig) 为

根据通路的重要性值对富集得到的通路进行排序，选择排序靠前的重要通路进行更深层次地分析。查看富集到该通路上靶点的上下游，分析其影响的效应。对重要程度靠前的通路分别进行此分析，对靶点上下游及影响效应重新构建网络，有多个通路可影响的效应说明是该中药治疗疾病可能表现的主要药效形式，最后选取富集到通路上可影响与疾病相关效应的靶点作为中药干预疾病作用靶点。

1.6牙痛停滴丸治疗口腔溃疡作用机制的验证

利用PCR对1.5分析得到的中药干预疾病作用靶点进行验证。SD大鼠(雄性，240±20g， SPF级，购于北京维通利华实验动物技术有限公司，许可证号：SCXK(京)2016-0006)，水合氯醛麻醉，使用80％浓度冰醋酸，使用6mm滤纸片，每只30uL，贴敷30s，贴敷于唇侧黏膜/颊侧黏膜，24h后观察溃疡表面情况。造模成功后进行给药，于给药第4天时进行取材，取溃疡粘膜处组织，冻于液氮中。组织按说明书取约50mg，加500μL RNA裂解液。

1.6.1总RNA提取

(1)加入裂解液后，用高压灭菌过的剪刀剪碎组织后，将组织置于冰里用组织匀浆机破碎，功率70HZ，时间60s，破碎后组织加入500μL稀释液，用移液枪吹打3-4次，70℃水浴3min，12000rpm，4℃，离心5min，吸取上清液。

(2)加入0.5倍上清体积的无水乙醇，用移液枪吹打3-4次，将混合物转移到离心柱，4℃， 12000rpm，离心1min，弃滤液。

(3)离心柱中加600uLRNA洗液，4℃，12000rpm，离心45sec。弃滤液。

(4)离心柱中加50uLDNA酶Ⅰ孵育液至吸附膜中央，室温放置15min。

表1 DNA酶Ⅰ孵育液

(5)离心柱中加600uLRNA洗液，4℃，12000rpm，离心45sec。弃滤液。

(6)离心柱中加600uLRNA洗液，4℃，12000rpm，离心45sec。弃滤液。将离心柱重新安置于收集管中，4℃，12000rpm，离心2min。

(7)将离心柱转移至洗脱管上，在离心柱膜中央加入100μL无核酸酶水，室温静置2min， 4℃，12000rpm，离心2min。

(8)重复步骤7。将RNA保存在-80℃备用。滤液测RNA浓度后，进行下一步操作。

1.6.2反转录(以20uL体系为例)

(1)配置gDNA去除反应体系：20uL体系mRNA的n在50-2000ng之间，V

(2)配置反转录反应体系：

表2反转录反应体系混合液

(3)将该反应混合液加入gDNA去除反应体系中，混匀。42℃孵育15min，95℃孵育3min，后置于冰上得到cDNA用于后续操作。

1.6.3扩增(20uL反应体系)

(1)将cDNA和扩增试剂盒在冰上解冻，简短离心，混匀。

(2)在冰上配置RT-PCR反应体系：

表3 RealTime PCR反应体系

(3)八联排封盖后简短离心。将八联排置于Bio-Rad扩增仪中，进行Real time PCR反应。

建议采用两步法PCR反应程序。当出现模板浓度过低引起非特异性扩增，引物Tm值较低导致的扩增效率低下或扩增曲线重现性不佳等现象时，建议尝试进行三步法PCR扩增反应。

表4两步法反应程序

*先使用60℃32sec进行扩增。

数据分析：用2^(-ΔΔCt)法。

1.7统计方法

利用SPSS23.0统计软件使用单因素方差分析(one-way ANOVA)进行统计学分析，以均值±标准差(

2实验结果

2.1口腔溃疡相关靶点的搜集

利用OMIM数据库和DrugBank数据库、Genecards数据库检索口腔溃疡相关靶点，筛选关键词设置为“Oral ulcer”，筛除不同数据库检索到的重复靶点，得到3845个。

2.2牙痛停滴丸相关靶点的收集

通过TCMSP数据库以“丁香”、“荜茇”、“冰片”关键词检索牙痛停滴丸中的化合物，以OB≥30％，并且DL≥0.18为筛选标准，共得到丁香成分10个，荜茇成分18个，冰片成分 4个，共得到482个靶点。利用Swiss Target Prediction数据库对成分进行检索，预测得到成分靶点253个。

2.3牙痛停滴丸治疗口腔溃疡的网络构建以及核心节点的确定

根据2.1和2.2得到的疾病相关靶点及中药成分相关靶点，利用Cytoscape软件中的 Bisogenet插件中关联的HPRD、BIND、DIP、MINT、INTACT、BIOGRID等蛋白质-蛋白质相互作用(PPI)数据库，分别构建疾病PPI网络和中药成分调控PPI网络。利用Cytoscape 软件中merge模块对构建的疾病PPI网络和中药成分调控PPI网络进行intersection交集运算，得到交集网络，即牙痛停滴丸治疗口腔溃疡的网络。通过NetworkAnalyzer工具对其拓扑特征进行分析，以Degree值大于2倍中位数，Betweenness Centrality、Closeness Centrality大于 1倍中位数，为参考标准，筛选得到交集网络中的关键节点，得到261个核心靶点，如附图1。

2.4牙痛停治疗口腔溃疡可能作用机制的预测

对2.3得到的网络核心靶点，进行GO功能注释和KEGG通路分析。GO注释是根据基因功能将基因分类的一种分析方法，其包括生物过程(biological process，BP)、分子功能(molecular function，MF)、细胞组分(cell component，CC)。将上述获得的261个核心靶点进行GO功能注释分析，设置P<0.01进行筛选，共得到生物学过程245个，细胞组分64个，分子功能68个。生物学过程结果显示，关键靶点主要参与细胞周期蛋白酶活性，蛋白质分解代谢，信号传导相关，这些与溃疡创面愈合密切相关。细胞组分结果显示，关键靶点主要在细胞核、细胞质细胞膜发生反应。分子功能结果显示，这些靶点主要与蛋白质结合、酶结合等相关，结果见附图1。

2.4.1 p值排序法

富集得到KEGG通路后，根据p值进行排序，获得排名前15的通路，结果见表5

表5根据p值排序得到的排名前15的信号通路

2.4.2网络核心节点的度(degree)值归一化法

根据得到的网络核心节点的度(degree)值进行归一化处理，将最大值作为1，设核心靶点i归一化后的度值为i(d)，核心靶点i所在通路t上富集到y个核心靶点，即通路t上含核心靶点-通路一对一的关系共y条，通路t上共有N个靶点，网络中共有核心靶点-通路一对一关系n条，根据上述技术通路重要性值公式计算富集得到的通路的重要性值t

表6富集得到的通路重要性排序表

根据通路的重要性值对富集得到的通路进行排序，选择排序靠前的重要通路进行更深层次地分析。查看富集到该通路上靶点的上下游，分析其影响的效应。对重要程度靠前的通路分别进行此分析，对靶点上下游及影响效应重新构建网络，有多个通路可影响的效应说明是该中药治疗疾病可能表现的主要药效形式，最后选取富集到通路上可影响与疾病相关效应的靶点作为中药干预疾病作用靶点，结果见附图2。根据附图2，可知牙痛停滴丸可能主要通过 PI3K-Akt signaling pathway、cell cycle、Thyroid hormone signalingpathway、FoxO signaling pathway、p53 signaling pathway、HIF-1signaling pathway、ErbB signaling pathway、T cell receptor signaling pathway、B cell receptorsignaling pathway、Estrogen signaling pathway、TGF-beta signaling pathway等信号通路调控细胞周期、血管生成、免疫反应，发挥治疗口腔溃疡的作用。

2.5牙痛停滴丸治疗口腔溃疡作用机制的验证

根据2.4预测得到的牙痛停滴丸治疗口腔溃疡可能相关通路及影响效应，我们选取富集到该通路上的影响主要效应的关键靶点采用PCR方法进行验证。选取调控通路上关键靶点CDK2、CDK4、CDK6、FOS、HIF-1、MAPK8(JNK)、TGF-β、GRB2、ERK、TP53进行验证，结果见附图3。结果显示与正常对照组相比，模型组口腔黏膜组织中CDK2、CDK4、FOS、 HIF-1、MAPK8(JNK)、TGF-β表达增加，具有统计学意义(*P<0.05，**P<0.01)；与模型组比较，给予牙痛停滴丸治疗后上述基因表达均有不同程度的降低，具有统计学意义(#P<0.05， ##P<0.01)。与模型组比较，给予牙痛停滴丸治疗后GRB2、ERK表达有下降趋势，具有统计学意义(#P<0.05，##P<0.01)。与正常对照组相比，模型组口腔黏膜组织中CDK6、TP53表达增加，具有统计学意义。

富集得到KEGG通路后，按照常规方法根据p值进行排序，发现排名靠前的信号通路包含多个癌症和病毒感染相关通路，Ribosome、Viral carcinogenesis、Prostate cancer、Epstein-Barr virus infection、Proteasome、Hepatitis B、Chronic myeloid leukemia、Pancreatic cancer、HTLV-I infection、Proteoglycans in cancer、Endometrialcancer。根据此排序进行分析牙痛停滴丸可能通过抑菌、抗炎等发挥治疗口腔溃疡的作用。然而，根据后期实验验证结果，发现牙痛停滴丸主要通过影响PI3K-Akt signalingpathway、cell cycle、Thyroid hormone signaling pathway、 FoxO signaling pathway、p53 signaling pathway、HIF-1signaling pathway、ErbB signaling pathway、 T cellreceptor signaling pathway、B cell receptor signaling pathway、Estrogensignaling pathway、 TGF-beta signaling pathway等信号通路调控细胞周期、血管生成、免疫反应，发挥治疗口腔溃疡的作用，这与本发明所采用网络药理分析方法获得的通路相一致。