三氮脒处理后犬吉氏巴贝斯虫转录组测序分析

冯笑笑，张佩豪，张 悦，田思瑾，刘亚娟，姚大伟

(南京农业大学 动物医学院，江苏 南京 210095)

三氮脒(diminazene aceturate，DA)作为一种广谱抗血液原虫的芳香双脒类衍生物药物，因其特效以及价格低廉成为国内治疗犬巴贝斯虫病的首选药物[1]。但是，近年来随着耐药虫株的不断出现，临床上DA很难彻底清除体内的巴贝斯虫，导致治疗周期较长，而且复发率极高。目前，尚不清楚犬吉氏巴贝斯虫(Babesiagibsoni)对DA产生耐药性的机理，对于DA药物作用机理的研究也相对较少。转录组学测序(RNA sequencing，RNA-Seq)技术[2]可从整体水平上对细胞中基因转录情况以及转录调控规律进行研究。目前转录组学测序已在B.bigemina[3]、B.microti[4]、B.bovis[5]以及B.divergens[6]的研究中得到应用。但是,关于犬吉氏巴贝斯虫转录组学的研究报道较少，并没有公开可用的转录组信息。为了更好的了解DA药物与犬吉氏巴贝斯虫之间的相互作用，本试验对犬吉氏巴贝斯虫进行体外培养，DA药物处理24 h后，提取样品中总RNA进行转录组学测序分析，探索、筛选与DA药物作用相关的犬吉氏巴贝斯虫基因，以提高对DA与巴贝斯虫相互作用机制的了解，寻求潜在的药物作用靶点。

1 材料与方法

1.1 虫株来源犬吉氏巴贝斯虫由本实验室保存，接种于试验犬保虫。

1.2 主要试剂胎牛血清(FBS)、RPMI 1640培养液购自Gibco；青链霉素混合液(100×)、红细胞裂解液购自Solarbio；DA购自上海麦克林生化科技有限公司。

1.3 带虫血样的培养与处理采集保虫犬血液2 mL，EDTA抗凝，3 000 r/min离心10 min。弃掉上层血清和白细胞层，收集下层的红细胞。用PBS洗涤3次，红细胞泥用10倍体积的完全培养液悬浮(RPMI1640培养基含40% FBS、1%青链霉素混合液(100×))。无菌采取健康犬血液2 mL EDTA抗凝，处理方法同上。将保虫犬红细胞悬液与健康犬红细胞悬液等比例混合，放至24孔细胞培养板内，每孔2 mL。每24 h吹打混匀红细胞1次，取5 μL制备血涂片，经瑞氏吉姆萨染色后在油镜下观察虫体，计算红细胞染虫率。培养3 d后弃去上清800 μL。在试验组添加含DA药物的完全培养液800 μL，药物终浓度为800 nmol/L，对照组添加完全培养液800 μL，每组3个重复。

1.4 红细胞裂解处理DA药物处理24 h后，取1.5 mL培养的红细胞，加入3倍体积的红细胞裂解液，涡旋混匀后，置于冰上15 min。期间涡旋混匀2次。4℃、450×g离心10 min，弃上清液。向沉淀中加入2倍体积的红细胞裂解液，涡旋混匀。于4℃、450×g离心10 min，弃上清液。用500 μL DEPC水重悬，放至液氮迅速冷冻后，-80℃ 备用。对照组样品处理方法相同。

1.5 转录组测序流程将-80℃保存的样品送至深圳华大基因科技服务有限公司完成转录组测序。采用DNBSEQ测序平台，DNBSEQ测序流程为mRNA文库构建、测序、数据分析。由于数据库中没有关于犬吉氏巴贝斯虫的转录组信息，因此本试验采用无参比对。测序得到的原始测序序列用过滤软件SOAPnuke v1.4.0进行统计，并使用trimmomatic v0.36进行质量过滤，获得Clean reads，将样品的Clean reads与组装得到的Unigene库进行序列比对。

1.6 差异基因转录水平验证采集保虫犬血液和临床患犬吉氏巴贝斯虫病的病犬的血液样品2 mL，制备红细胞悬液，然后用DA处理，处理方法同1.3。取处理后的血液样品100 μL，加入1 mL RNAiso plus，按照说明书方法提取样品中总RNA，10 μL反转录体系中加入500 ng总RNA按照PrimeScript RT Master Mix(Perfect Real Time)说明书方法将上述RNA反转录合成cDNA。

筛选转录水平差异较大的基因CL118，根据测得的基因序列，采用Primer Primer (Version 5.0)软件设计引物，上游引物CL 118-F：5′-CCGGATTACGGGACAGAGTTA-3′，下游引物CL 118-R：5′-TAGTTCCAATGGAGCCGTCC-3′，扩增片段大小为72 bp。根据犬吉氏巴贝斯虫18S rRNA基因保守序列设计引物，上游引物B.com 151-F：5′-TGTTGCAGTTAAAAAGCTCGTAGTTG-3′，下游引物B.com 151-R：5′-TCCATGCTGAAGTATTCAAGACACA-3′，扩增片段大小为152 bp。以上引物序列由生工生物工程(上海)股份有限公司合成。在25 μL反应体系中加入2.5 μL cDNA，按照TB Green®Premix Ex Taq试剂说明书进行荧光定量PCR反应，反应在7300 Real-time PCR System中进行。以CL118作为靶基因，犬吉氏巴贝斯虫18S rRNA基因作为内参基因，采用2-ΔΔCt法进行相对定量分析，以未用药物处理组作为对照，检测DA药物作用后犬吉氏巴贝斯虫CL118基因转录水平。

2 结果

2.1 测序数据与组装结果的比对统计从Clean reads质量统计表(表1)可以看出：6个样品过滤后的序列数在3 643～4 256万条之间，测序质量值大于20(Q20)的碱基占所有碱基的92%以上，测序比对效率均值在86%左右，测序数据质量良好，满足后续分析要求。

表1 Clean reads质量统计

2.2 转录组功能注释结果

2.2.1GO数据库分析 通过转录组测序后的转录组基因数据进行筛选，根据Nt数据库注释共筛选出77个有关巴贝斯虫属的基因，将筛选后的基因进行GO数据库注释分类，基因注释信息如图1所示。属于生物过程的有74个，其中与细胞过程(39个)和代谢过程(29个)的unigenes较多，而生物之间的种间相互作用(2个)、生物调节(1个)、刺激反应(1个)、定位(1个)以及生物调节过程(1个)的unigenens较少。属于分子功能有77个，其中与结合活性(31个)、结构分子活性(22个)、催化活性(17个)、翻译调节活动(6个)的unigenes较多，而与蛋白质标签(1个)相关的unigenens较少。属于细胞组分有92个，二级分类分别为细胞结构实体(43个)、细胞内部组分(38个)和含蛋白质复合物(11个)。

图1 犬吉氏巴贝斯虫药物处理后差异表达基因GO分类

2.2.2KEGG数据库分析 将筛选后的基因进行KEGG代谢通路注释，主要聚集在细胞过程、环境信息处理、遗传信息处理、代谢和有机系统类别中。如图2所示，其中参与细胞过程的相关基因共有15个，其二级分类分别为真核细胞群落(2个)、细胞生长与死亡(5个)、运输与分解代谢(6个)以及细胞运动(2个)；参与环境信息处理代谢途径的相关基因有7个，均属于信号传导(7个)类；参与遗传信息处理的相关基因共有36个，其二级分类分别属于翻译(29个)、转录(2个)、蛋白质折叠、分类以及降解(5个)；与代谢相关的基因共有3个，其二级分类分别为概述代谢图(1个)、能量代谢(1个)以及碳水化合物代谢(1个)；与有机系统相关的基因共有22个，其二级分类的代谢途径分别为内分泌系统(8个)、免疫系统(5个)、成熟过程(3个)、环境适应(3个)、敏感度系统(2个)以及循环系统(1个)。

图2 吉氏巴贝斯虫药物处理后差异表达基因KEGG分类

2.3 不同样本差异基因表达分析本试验在筛选的过程中，将错误发现率FDR(flase discovery rate)<0.01、差异倍数(flod change,FC)≥2为筛选标准。其中，FC表示2个样本(组)间表达量的比值。测序结果显示，DA作用后犬吉氏巴贝斯虫CL118基因转录水平显著下调。该基因在目前有记载的犬吉氏巴贝斯虫基因数据库中尚无相关报道和注释，该基因与GenBank数据库中已知的线虫Parastrongyloidestrichosuri(LM523569)基因序列同源性约90%，该基因与GenBank数据库中已知的盘尾丝虫Onchocerca(VDN00566)氨基酸序列同源性约71%。

2.4 基因表达验证结果为验证转录组测序数据的可靠性，采用相对实时荧光定量PCR方法进一步验证，检测犬吉氏巴贝斯虫CL118基因相对18S rRNA内参基因的转录水平，结果如图3样品A组所示，试验室保存虫株经DA处理后的CL118基因转录水平显著性下调，该结果与RNA-Seq的分析结果相一致，说明CL118基因在DA处理后基因转录水平确实下调。

A.实验室保存虫株；B～D.临床收集待测虫株;与对照组比较，*示P<0.05;**示P<0.01

采集临床其他患犬吉氏巴贝斯虫病的病犬血液，检测DA药物处理后CL118基因的转录水平，所有临床样品CL118基因的转录水平都下调(降低为未经药物处理组的22.6%～61.2%)。结果表明，犬吉氏巴贝斯虫在DA处理后虫体CL118基因转录水平显著下调。

3 讨论

转录组学是从整体水平对基因在蛋白合成中所发挥的功能以及基因的结构进行研究，以揭示在生物应对特定的疾病或其他因素时在分子水平上的变化规律。这对于研究寄生虫的致病机制、疫苗的开发以及药物靶标的筛选都发挥重要的作用。有研究报道，由罗氏巴贝斯虫所引起的犬巴贝斯虫病与人的恶性疟原虫病的临床症状相似[7]，因此参考疟疾的研究模型来用作犬巴贝斯虫的研究的同时[8]，也有研究将巴贝斯虫与疟原虫转录组以及基因组进行比对[9-10]，以发现两者的异同。

目前,DA是普遍公认有效的抗虫药，还具有抗炎[11]和抗血压[12]等作用。该药作为抗原虫常用药已有60多年的历史，长时间的使用已有耐药虫株的出现。目前对其耐药性的作用机制研究较少。有研究表明[14]，通过体外培养可获得犬吉氏巴贝斯虫DA耐药虫株，而且耐药虫株与原虫株有能量和糖酵解等代谢的差异。具有DA抗性的犬吉氏巴贝斯虫分离株比自然感染的野生型虫株具有更高的ATP浓度,表明其糖酵解途径的活性并没有改变。为此，本试验通过对犬吉氏巴贝斯虫对体外培养扩增后用DA进行药物处理，以探究药物对于虫体基因转录水平的差异。

对于犬吉氏巴贝斯虫并没有该虫株可查的全基因组数据参考，已知的生物信息量又很匮乏，很多并未公开数据，这给转录组数据的分析带来了困难。目前对没有参考基因组的物种的研究，主要采取的方法是通过无参转录组测序进行组装，将获得的数据与已知的蛋白数据库(NR、Nt、Swiss-Prot、GO、KEGG)进行比对，依托于强大的生物信息学平台作为支撑，并根据“基因结构相似，功能同源”的原理，对基因的功能进行注释，从而从基因水平上预测未知基因的功能并进行有方向性的检测。本试验将犬属Unigene筛除后，将获得的9 699条Unigene数据进行7大功能数据库注释(KEGG、GO、NR、NT、SwissProt、Pfam和KOG)，共有9 332条Unigene获得了基因注释，占All-Unigene的96.22%。被所有基因库注释的基因共有131条Unigene(1.35%)，有3.78%的Unigene未被注释。对于没有得到注释的Unigene，有可能是吉氏巴贝斯虫经DA处理后有新基因，或由于数据库现有的基因资源有限，基因功能注释信息不丰富，从而造成部分序列暂时无法获得对应的功能注释信息。

对转录组测序后的Nt进行有关巴贝斯虫属注释的筛选，共筛选出77个表达基因，通过比较其表达基因GO数据库注释，发现注释富集基因最多的类别是细胞组分、分子功能和生物进程。而差异表达基因在KEGG通路注释中，显示参与翻译通路有关的基因较多，提示了DA处理后关于巴贝斯虫属的基因表达变化的方向和代谢途径。

通过对基因转录水平显著与否进行基因筛选，发现CL118基因的转录水平显著性下调。对该基因进行核酸数据库和蛋白数据库比对后发现其与Onchocercaochengi(VDN00566.1)氨基酸序列存在71%的相似度，与Parastrongyloidestrichosuri(LM523569.1)基因序列存在90%的同源性。对CL118基因转录水平进行qPCR验证，试验结果发现犬吉氏巴贝斯虫虫株在DA处理后，CL118基因的转录水平均出现明显下降，该结果与RNA-Seq的分析结果相一致。通过临床验证发现，经DA处理后的临床样本虫株的CL118基因的转录水平均出现下降，试验结果与RNA测序结果相一致。其下调的原因可能是与DA药物作用有关，有可能是DA药物作用的靶点。