整合多维组学数据鉴定登革热感染关键生物分子标记及潜在治疗药物

徐琪，谢树仁，王嘉琪，才让杰，黄佳雨，陈强，张娅，李永生

海南医学院生物医学信息与工程学院，海南 海口 571199

绝大部分传染性疾病源于热带地区，近年来SARS、禽流感、埃博拉、寨卡热的流行和新型冠状病毒感染对人类健康都造成了重大影响[1]。登革热(dengue fever)是由登革热病毒感染、伊蚊传播的一种急性传染病，无法直接经过呼吸道、消化道或接触在宿主之间传播。理解宿主-病原体间复杂的关系是开发新的诊断、预防和治疗模型及建立动态防控体系的关键前提。随着高通量测序技术的发展，基于多组学生物大数据已经发现了一些与传染病相关的关键基因及药物靶点。例如，基于外周血的转录组测序，Nikolayeva等[2]识别了18个基因生物标记能够用于识别严重的登革热患者。Robinson等[3]通过整合多套基因表达数据识别了20个基因标记物用于预测登革热患者的疾病进展。然而对于登革热感染关键分子标记的识别仍然具有挑战。

此外，登革热感染经常会伴随着人类其他复杂疾病的产生，患者治愈后期会有多种并发症。解析登革热感染与人类其他复杂疾病的关系对于登革热感染后期并发症的早期预测具有十分重要的意义。同时，目前已经开发了一些登革热治疗的方法，但是仍然需要新的治疗方案的提出。如何基于多维组学数据整合分析识别新的候选药物也是亟待解决的关键问题之一。本研究通过整合登革热感染相关的多层次信息、挖掘潜在的分子标记，并基于网络分析技术识别登革热感染与其他疾病间的关联，优选潜在的治疗药物，对于登革热感染等疾病的预防和治疗具有重要的理论意义。

1 资料与方法

1.1 微阵列芯片数据的获取与预处理首先从GEO和ArrayExpress数据库中共获得了4套登革热感染相关基因表达数据集，如表1所示。通过对芯片数据进行对数化、探针过滤、补缺失值、标准化、探针注释等预处理获得基因表达矩阵。

表1 本研究使用数据集情况Table 1 Datasets used in this study

1.2 差异表达基因的筛选与整合根据样本分组信息，基于Limma包对数据集进行差异表达分析。首先，输入基因表达矩阵，分组矩阵并构建差异比较矩阵；然后将数据拟合到具有lmFit函数的线性模型；使用contrasts.fit函数根据对比度模型计算差异；进一步使用eBayes函数进行贝叶斯测试；最后针对具有topTable功能的所有基因生成测试结果，使用Benjamini和Hochberg矫正方法校正P值。基于阈值FDR＜0.05，得到各个登革热感染相关数据集的差异表达基因，并绘制对应的火山图与热图。

由于实验样本的异质性，使用不同的检测平台和数据处理方法均会导致结果不一致。因此，进一步利用“RobustRankAggreg”[4]方法将来自不同数据集的分析结果整合，设阈值FDR＜0.05，fold change＜1.2，得到整合的差异表达基因。

1.3 基于随机游走模型优选生物分子标记物随机游走是给定源节点开始从当前节点到随机选择的邻居迭代walker转换[5]。以RRA整合得到的生物分子标记物中FDR＜0.01的基因作为种子节点，映射到HPRD数据库的人类蛋白质相互作用网络进行随机游走分析。提取得分前5%的基因作为优选的生物分子标记物，利用Cytoscape从人类蛋白质互作网络中提取了439个优选生物分子标记物的PPI子网，并删除了网络中没有互作的节点。最后将随机游走得到的节点得分进行排序后作为输入，进行基因集合功能富集分析(GSEA)。

1.4 登革热病毒与人类蛋白之间的相互作用蛋白质相互作用网络的结构和性质是系统生物学中重要的研究内容，通过网络分析来评估相关疾病基因的特性是研究疾病最简单方法之一，这有利于找到病毒潜在靶向的蛋白，并揭示有关其病因的重要线索。通过已有的登革热公共数据库(包括DenHunt[6]和DenvInt[7])和文献[8-9]，收集到了登革热病毒蛋白和人类蛋白质之间的互作关系，得到8 526个PPIs，涉及11种登革热病毒蛋白和2 863种人类蛋白。将随机游走优选后的节点映射到网络中并使用“clusterprofiler”包对病毒蛋白靶向的人类蛋白进行功能富集分析。

1.5 基于PPI的登革热感染与复杂疾病的关联为了揭示登革热感染与人类复杂疾病之间的关系，通过查找文献收集了多种复杂疾病相关的基因[10]，共涉及299种人类复杂疾病的3 173个基因，将这些疾病进一步划分为10个大类。随后基于网络之间的重叠评估了疾病相关基因和登革热感染相关基因之间的相似度Svb值[11]，公式如下：

dvv和dbb分别代表登革热感染相关蛋白网络和人类复杂疾病蛋白网络的平均最短路径，dvb表示疾病b和登革热v之间的成对平均最短距离，Svb＜0则表明登革热v感染相关蛋白质和与疾病b相关蛋白质之间存在基于网络的重叠模块。

1.6 构建药物-靶点网络优选潜在抗病毒药物利用DrugBank数据库(v4.3)[12]、治疗靶点数据库(TTD)[13]、PharmGKB数据库[14]、ChEMBL(v20)[15]、BindingDB[16]和IUPHAR/BPS药理学指南[17]收集药物-靶点相互作用信息。并从DrugBank[12]中提取每种药物的化学结构，采用SMILES格式。在这里，使用了符合以下三个标准的药物-靶相互作用：(i)结合亲和力，包括Ki、Kd、IC50、和EC50均小于等于10μmol/L；(ii)该目标在UniProt数据库[18]中标记为“已审查”；(iii)人类靶标由唯一的UniProt登录号表示。本课题组最近的研究[19-21]提供了经实验验证的药物-靶点网络的详细信息。

进一步基于互作网络对候选药物进行了优先级排序。基于与病毒互作的宿主基因编码蛋白质V和药物靶标T，计算了V与每个候选药物的靶标蛋白质T的网络邻近度。计算公式如下：

d(i,j)是人类蛋白质相互作用网络中蛋白质i(病毒互作蛋白)和j(药物靶标蛋白)之间的最短距离。基于随机检验将网络邻近度结果进一步转换为Z-score：

2 结果

2.1 转录组分析识别登革热感染的关键基因使用limma R包对四套数据集进行差异分析，卡阈值FDR＜0.05筛选差异表达的基因，各个数据集差异表达基因数目如表2所示。并对四套数据集分别绘制了火山图和热图(图1)。

表2 不同数据集识别的登革热感染上下调基因情况Table 2 Numbers of differentially expressed genes

图1 登革热感染四套数据集差异分析结果Figure 1 Differential expression of genes in dengue infections

通过热图观察到基因表达模式具有异质性，在GSE84331数据集中，转录组分析可以很明显地将登革热疾病样本与正常样本区分开，而在GSE18090和GSE51808的热图中，正常样本和登革热早期(DF)阶段的疾病样本之间并未检测到明显的差异，却在DF到DHF发展的阶段疾病转录组发生了急剧的变化(图1)。在疾病恢复期，转录组水平与严重登革热相比并未发生明显改变。因此，转录组分析可区分早期登革热患者和登革出血热患者。

2.2 RRA整合分析识别稳健的生物分子标记随后基于RRA方法整合了四个数据集的DEGs，并基于logFC的大小进行排序，取FDR＜0.01的691个上调和1 071个下调DEGs绘制热图(图2A)。CD1C阳性树突状细胞是真皮树突状细胞所有亚群中最丰富的，研究表明，人类皮肤树突状细胞可能是登革热病毒感染的重要目标，它们会将传染性病毒从皮肤转移到引流淋巴结，从而有效的激活全身免疫反应，CD1C可能是引发早期适应性抗登革病毒T细胞反应的重要因素，这进一步确定了CD1C细胞是登革热的潜在靶标[22]。随后通过观察在四个数据集里均差异表达的基因CASP71(图2B)，以及在三个数据集中显著差异的基因FEZ1、PITPNC1、SEMA4C、TXLNA(图2C～E)的表达情况，发现差异基因在不同数据集中的差异程度都有所不同。

图2 稳健差异表达基因的鉴定Figure 2 Robust differentially expressed genes

2.3 基于随机游走优选登革热分子标记蛋白质-蛋白质相互作用(PPI)网络由无向图表示(图3A)。进一步基于专门存储经过实验验证人类蛋白质互作信息的HPRD数据库中[23]的网络进行随机游走分析，并取经过RRA整合后FDR＜0.01的1 762个差异基因作为随机游走模型的种子节点，通过迭代得到每个节点的在随机游走模型中的得分，排序后提取前5%共439个节点，使用Cytoscape3.8.0将其映射到HPRD网络来构建优选后的蛋白互作网络(图3A)。该网络删除了没有互作的节点之后包括304个节点(node)，726条边(edge)，同时对网络的拓扑属性进行了分析如表3所示。

图3 优选登革热分子标记的功能富集Figure 3 Functional enrichment of prioritized biomarkers for dengue infection

表3 PPI网络拓扑属性Table 3 Topological features of PPI network

TP53作为常见的与细胞周期控制密切相关的基因在该网络中具有最高的连通度为36，已有研究证明，细胞周期停滞是p53激活后观察到的主要生物学结果，可防止受损DNA的积累以及基因组不稳定性[24]。其次CDK1基因编码的蛋白质作为Ser/Thr蛋白激酶家族的成员是高度保守的蛋白激酶复合物催化亚基，被称为M期促进因子(MPF)，它对真核细胞周期的G1/S和G2/M相变必不可少，这表明在登革热疾病感染过程出现了细胞周期的改变[25]。

为了分析登革热疾病中的分子生物学功能，将随机游走的得分排序后作为GSEA的输入，并对其进行KEGG和GO功能富集分析得到了一些免疫反应激活信号转导通路，例如，体液和细胞适应性免疫是人体抵御病毒感染的重要途经，它在DENV感染期间被激活，众多研究表明人体细胞被病毒感染后，免疫系统功能异常[26-28]，有研究通过对严重登革热过度免疫激活进行研究充分证实了免疫反应在登革热发病机制中的关键作用[29]。同时还发现差异基因主要富集在DNA复制、细胞周期蛋白体酶等通路中(图3B～3C)。且已有研究表明，病毒感染细胞时，病毒编码的蛋白或DNA会扰乱细胞周期通路，促进细胞向G1/S期转化或者使细胞静息于G2/M期，这些现象表明在登革热疾病感染过程出现了细胞周期的改变[30]。

2.4 登革热相关基因与病毒蛋白高频率互作接下来利用已有的登革热公共数据库(包括DenHunt[6]和DenvInt[7])和文献[8-9]收集了病毒与人类蛋白的互作信息，将随机游走优选后的304个节点映射到网络中，生成病毒与人类蛋白互作网络(图4A)。由于节点的拓扑性质可以通过其在网络中的位置来显示，因此分别比较受病毒蛋白靶向的人类蛋白和HPRD网络中其他蛋白的介数、紧密度、连通度三个拓扑属性(图4B～4D)。发现病毒靶向的人类蛋白介数、紧密度、连通度均明显高于其他蛋白，表明病毒与人类蛋白互作网络具有很强的模块性，且位于网络的中心位置，进一步证实病毒可能会针对位于网络中心或在信息传播中起重要作用的蛋白质进行复制。随后利用“clusterprofiler”包对病毒靶向的人类蛋白进行功能富集分析，发现大部分显著富集的功能都与病毒感染相关，而且差异蛋白主要富集在蛋白酶体蛋白质分解代谢过程、肽链内切酶复合体、单链DNA结合等通路。已有研究表明，在登革热感染过程中，自噬有助于病毒复制，而病毒可以通过蛋白酶体降解逐渐降低自噬受体的水平，说明蛋白酶体分解代谢通路与登革热病毒的发展密切相关[31]。

图4 登革热相关基因与病毒蛋白互作Figure 4 PPI interactions among dengue proteins and human proteins

2.5 登革热感染与多种复杂疾病关联随后从文献中收集了涉及299种人类复杂疾病的3 173个基因，来探索登革热感染与人类复杂疾病之间的关系。通过计算疾病相关基因和登革热感染相关基因之间基于网络的的相似度Svb值，发现了登革热感染可能与癌症、消化系统肿瘤、哮喘、类风湿性关节炎，以及一些心血管疾病相关(图5A)。特别是已有研究证实类风湿性关节炎[32]、消化系统肿瘤[33-34]、脑血管疾病[35]与登革热感染相关。通过对这三类疾病进行分子网络可视化发现，MIF、TP5、HNRNPD等基因不仅靶向登革热病毒还与复杂疾病有关。其中巨噬细胞迁移抑制因子MIF(图5B)是一种多效性炎症细胞因子，在先天性和适应性免疫反应的调节中非常重要，它与自身免疫性疾病的发病机制有关[36]，并且该基因的中和会导致DENV感染巨噬细胞的迁移能力提高，与DENV感染的严重程度有关[37]。有研究表明，p53(图5C)介导的干扰素刺激基因在体外和体内都能抑制病毒的复制[38]，并且TP53的突变可作为患者的预后预测指标[39]。HNRNPD(图5D)作为一种通用的黄病毒宿主因子与AUF1蛋白结合可以促使登革热病毒复制[40]。

图5 登革热与人类复杂疾病之间的关系Figure 5 Relationships of dengue infection and human diseases

2.6 登革热感染潜在药物的挖掘根据ZdVT对优选药物进行排序，共获得39个与登革热疾病相关的药物。经筛选后获得了23种至少靶向两个登革热蛋白的药物(图6A)。已有研究表明，DENV可能导致人体多种疾病，包括癌症、哮喘、溃疡性结肠炎等[41]，其中哮喘可以作为预测登革热发展的标志[42]，而药物关联性分析获得的oxtriphylline已被FDA批准用于治疗哮喘[43](图6B)。此外，登革热病毒感染常伴有关节炎发生[44]，Tenoxicam、Acemetacin、Mefenamic acid、Etoricoxib、Loxoprofen作为常见的非甾体抗炎药(图6C)，常用于缓解骨关节炎的症状[45-49]。同时，Omega-3 fatty acids(图6D)可通过降低甘油三酯血症来降低动脉粥样硬化性心血管疾病的发病风险[50]。Romidepsin(图6E)作为一种组蛋白去乙酰化酶抑制剂，对多种实体瘤具有抗肿瘤作用[51]。此研究可以为DENV的治疗及潜在药物的发掘提供参考。

图6 潜在药物的挖掘Figure 6 Potential drugs for dengue infection

3 讨论

登革热感染作为一种传染性极强、发病率极高的蚊媒传染性疾病，目前还没有发现有效的抗病毒疗法可以预防和治疗，因此进一步研究登革热的分子生物学机制，探索更加有效的预防和治疗方法是热带疾病研究领域亟待解决的关键。本研究旨在探讨登革热感染相关的生物分子标记、潜在的药物组合疗法的挖掘，为登革热感染发病机制的探索与诊疗提供了新的思路。

本文通过RRA方法将从四个数据集中获取的差异基因信息整合到一起，再基于随机游走模型优选出关键生物分子标记物，并建立全面的登革热感染相关功能调控网络，经过功能富集分析发现大部分差异基因都富集在DNA复制、细胞周期相关的通路上。这表明在登革热感染能影响细胞周期相关通路，感染登革热病毒时会影响细胞活动，此外这些基因还富集在免疫相关通路上，而患者的免疫反应与后期疾病发展的严重程度相关。

另外，对于登革热(DF)和登革出血热(DHF)两组也做了差异分析，发现与DF和DHF相关基因主要富集在T细胞受体信号通路、适应性免疫应答、T细胞激活、调节淋巴细胞活化等免疫相关通路，进一步表明人体细胞被病毒感染后，会出现免疫功能失调，导致免疫系统功能异常。最后通过分析登革热与299种复杂疾病的关系进一步寻找与登革热相关的潜在药物，在探究登革热的预防和治疗道路上取得新的进展。

综上所述，本研究通过建立全面的登革热感染相关功能调控网络，优选登革热感染相关的基因分子标记，解析了与病毒感染相关的功能调控通路。通过探究登革热与其他复杂疾病的关联，和潜在药物的挖掘，为登革热感染预防和诊疗提供了新的思路。