溴氰虫酰胺暴露对赤子爱胜蚓的转录组学的影响

乔治华，刘 翔，姚向峰，杜庆志，张建业，孙石昂，张风文，姜兴印*

(1.山东农业大学植物保护学院，泰安 271018；2.农药毒理与应用技术省级重点实验室，山东泰安 271018；3.中国农业科学院烟草研究所，山东青岛 266101)

在现代农业生产中，农药被广泛应用于防治病虫草害[1]，溴氰虫酰胺作为第二代双酰胺类杀虫剂，已得到广泛应用。溴氰虫酰胺具有触杀和内吸作用，当接触害虫并被虫体吸收后，激活靶标昆虫的鱼尼丁受体（RyRs），造成昆虫细胞内外钙离子浓度平衡破坏，最后肌肉过度收缩导致昆虫死亡[2-3]。目前关于溴氰虫酰胺的主要研究集中于对靶标害虫的防治。研究表明，溴氰虫酰胺可以很好地防治鳞翅目[4]、半翅目[5]和缨翅目[6-7]害虫。通常，溴氰虫酰胺在蔬菜和果树上可直接喷雾使用或者灌根防治害虫[8-10]，在玉米田中可以直接进行种子处理和土壤混施等方式来防治小地老虎和双委夜蛾[11-12]，先正达于2020年首次登记溴氰虫酰胺种子处理悬浮剂，用来防治玉米田小地老虎。这些防治方法都是利用了其优良的内吸活性，来达到防治害虫的效果。随着使用时间的延长，农药的弊端也逐渐暴露出来，农药对土壤环境、水环境的污染逐渐受到人们的关注[13-14]。溴氰虫酰胺的长期使用，在土壤中发生富集，也必然会对土壤环境带来问题。

蚯蚓广泛分布于土壤中，占土壤动物生物量的60%～80%，是控制陆地生态系统物质和能量循环的重要组成部分[15-16]。蚯蚓因为具有对环境污染物敏感，体积较大，常见，易于繁殖等特点，常作为监测土壤污染的指示生物，特别是用于确定土壤中污染物的生态毒理学[17-18]。

随着分子生物学的发展，一些现代化的手段应用于蚯蚓的生态毒理学研究，如转录组学、基因组学、蛋白组学、代谢组学等[19-21]。通过现代分子生物学手段研究对蚯蚓微观指标的影响成为目前生态毒理学研究的热点，通过转录组学可以研究污染物对生物的转录水平的影响，从而探究其分子机理。

目前溴氰虫酰胺对环境生物的影响研究较少，主要集中于急性毒性和生理生化水平的亚慢性毒性研究，而在转录组学水平上的研究尚未报道。本研究的主要目的是研究溴氰虫酰胺对于土壤生物蚯蚓的生态毒性，基于本研究前期溴氰虫酰胺对蚯蚓的基础生理生化指标的研究结果，进一步通过转录组学的方式研究溴氰虫酰胺对蚯蚓转录水平的影响，及其对蚯蚓使用的分子机理，同时评价溴氰虫酰胺对于土壤环境的影响，为溴氰虫酰胺的合理使用提供理论基础。

1 材料与方法

1.1 试剂与仪器

94%溴氰虫酰胺原药，上海杜邦农化有限公司；TRIzol试剂，上海生工生物技术有限公司；RNA Nano 6000检测试剂盒，安捷伦科技（中国）有限公司。QuantStudio 5 荧光定量PCR 仪，美国ABI 公司；BIO-BEST140E紫外凝胶成像仪，美国SIM公司；MyCyclerTMThermal Cycler梯度PCR基因扩增仪，美国Bio-rad公司；NanoDrop One分光光度计，赛默飞世尔科技中国有限公司。

1.2 赤子爱胜蚓

赤子爱胜蚓（Eisenia fetida）购于江苏省句容市王军蚯蚓养殖场。在天然的土壤中，20±1℃下驯养2 w后使其适应人工养殖环境。试验前，取环带明显且大小较一致，重量在300～500 mg的健康成蚓，用0.9%的生理盐水将蚯蚓体表洗净，放置在培养皿中，皿内铺有去离子水湿润的滤纸。20±1℃下放置一夜，目的是清除蚯蚓肠道内的杂质等[22]，保证试验过程中蚯蚓的死亡率不超过1%。

1.3 试验土壤

蚯蚓的试验土壤采用OECD土[23]。高岭土和苔藓泥炭土均购于普天试剂公司。

表1 OECD 土配方

1.4 试验方法

1.4.1 蚯蚓染毒

溴氰虫酰胺在土壤中的初始残留浓度为0.40～3.6 mg/kg，半衰期为9.2～20.8 d[14,24]，根据我们先前的研究结果[25-26]，转录组试验设计溴氰虫酰胺浓度为0.1、10 mg/kg干土和100 mg/kg干土，以丙酮处理为对照组，后面的转录组分析0.1 mg/kg干土处理组编号为Cya.L，10 mg/kg干土处理组编号为Cya.M，100 mg/kg干土处理组编号为Cya.H，丙酮对照组编号为CK。

称取0.212 8 g溴氰虫酰胺原药用丙酮溶解于100 mL的容量瓶中配制成2 000 mg/L的母液，然后将母液稀释为系列浓度。在50 g的人工土中加入适量的药液配制成所需要浓度的土壤[23]，对照组加入25 mL丙酮，搅拌均匀，放于通风处使丙酮挥发干净，再和450 g人工土混合均匀。加入一定量的无菌去离子水，使土壤含水量达到土壤质量的40%，放在1 L的烧杯中。

从低浓度到高浓度依次进行染毒，每个处理3个重复。每个烧杯放入15条经过清洗大小一致的蚯蚓，每周于土壤表面加入定量的牛粪和水，用皮筋将保鲜膜封住烧杯口并扎孔保持空气的流通[27]。放在生化培养箱，温度为20±1℃，L//D=12 h∶12 h[22]。于染毒后14 d取样，样品用液氮及时冷冻，并放入-80℃冰箱冻存待检测。

1.4.2 蚯蚓转录组学检测与分析

用TRIzol试剂按照说明书提取蚯蚓样品中的总RNA。RNA的完整性通过Agilent Bioanalyzer 2100系统（Agilent Technologies,CA,US）的RNA Nano 6000检测试剂盒检测。RNA的纯度和浓度由NanoDrop One 分光光度计测定（Thermo Fisher Scientific,Wilmington,DE）。样品经过检测合格之后，由南京集思慧远生物科技有限公司进行转录组测试。本试验参考基因组信息详见网址https://bigd.big.ac.cn/search/?dbId=gwh&q=GWHACBE00000000。对差异表达基因在7个数据库（Swiss-Prot，GO，KEGG，COG，KOG，Pfam，NR）进行比对。

采用DESeq2 软件进行差异筛选，将Fold Change≥2且FDR＜0.05作为筛选标准。使用R（v3.6.2）ggplot2程序包绘制差异基因火山图，使用R（v3.6.2）pheatmap绘制基因表达聚类热图。利用clusterProfiler R软件包实现差异表达基因（DEGs）的基因本体（GO）富集分析，使用Kobas软件检测[28]KEGG通路中差异表达基因的富集。

1.4.3 荧光定量PCR

通过RT-qPCR技术对转录组数据准确性进行验证，选择10个与差异显著相关的基因进行荧光定量表达分析，基因及其引物序列设计列表见表2。反应体系见表3、4。每处理3次生物学重复，每个重复3次技术重复，采用2-ΔΔCT法进行数据分析计算。

表2 用于实时荧光定量PCR 的蚯蚓基因引物序列设计

表3 实时荧光定量PCR 反应体系配制

表4 实时荧光定量PCR 反应条件

2 结果与分析

2.1 RNA质量检验

所有用于转录组的RNA质量检测如图1所示，通过电泳条带可以看出18SRNA条带清晰无拖尾，由于大多数的无脊椎生物中28SRNA 有一个“hidden break”的结构[29]，恰好位于28S的中央导致电泳时与18S重合。通过Nanodrop定性定量分析，各个处理的OD260/280值均在1.8～2.2之间，检测结果均为A，说明RNA的纯度和完整性较好，RNA浓度均大于400 ng/μL。则12个样品的RNA质量合格，符合建库的要求。

图1 RNA 提取电泳结果

2.2 差异表达基因筛选与分析

图2可以快速准确地比较出1、10 mg/kg干土溴氰虫酰胺处理的蚯蚓与CK的差异基因数目。图中一个点代表一个基因，其中点偏离横坐标的0点越远，说明该点对应的基因的表达量与CK差异越大；点的纵坐标越大，说明该点代表的基因在样本间差异表达越显著。因此试验中应尽量选择偏离横坐标原点且纵坐标值大的点所代表的基因作为功能验证基因。

图2 1 mg/kg 干土（A）、10 mg/kg 干土（B）和100 mg/kg 干土（C）溴氰虫酰胺处理14 d 后蚯蚓差异表达基因（DGEs）的数量火山图及差异基因韦恩图（D）

我们可以发现1 mg/kg干土的处理在第14 d时与CK相比，显著差异基因有960个，其中上调表达438个，下调表达522个（图2A）。10 mg/kg干土的处理在第14 d时与CK相比，筛选出显著差异基因1 212个，其中上调表达632个，下调表达580个（图2B）。100 mg/kg干土处理在第14 d时与CK相比，筛选出显著差异基因1 081个，其中上调表达608个，下调表达473个（图2C），同时从1 mg/kg干土到10 mg/kg干土差异基因显著增多，10 mg/kg干土到100 mg/kg干土差异基因有所下降。由图2D中的韦恩图可知，蚯蚓在3 个处理组中共有3 253 条差异表达基因（differentially expressed gene，DEG），不同差异表达基因集中的DEG有重叠，其中93个基因在不同处理组间均呈现出差异表达。

对差异表达的基因进行聚类分析（图3），其中横坐标表示不同样品对应的基因。颜色不同表示基因在样品中的表达量水平不同。2个处理间所有差异表达基因共聚类成6组，组内基因的表达水平基本相似。通过结果我们可以看出CKvsCya.L（图3A）中聚类结果将CK和Cya.L很好地分开，每个处理的3个重复聚类到一起，说明重复性良好。同样CKvsCya.M、CKvsCya.H也是如此。

图3 CK vs Cya.L（A）、CK vs Cya.M（B）、CK vs Cya.H（C）差异表达基因聚类图

图4 CK vs Cya.L（A）、CK vs Cya.M（B）、CK vs Cya.H（C）的差异表达基因GO 注释分类统计图

2.3 差异表达基因功能富集分析

2.3.1 差异基因注释结果

对筛选的差异基因与7个数据库进行比对分析，得出在各个基因库中注释的基因数量如表5所示。可以看出CKvsCya.L在7个数据库中比对注释的基因数为586个，其中在GO注释基因数为288个，KEGG数据库注释基因数为65个，COG数据库注释基因数86个；CKvsCya.M在7个数据库中比对的注释基因数为760个，其中在GO注释基因数为438个，KEGG数据库注释基因数为125个，COG数据库注释基因数为114个；CKvsCya.H在7个数据库中比对的注释基因数为669个，其中在GO注释基因数为413个，KEGG数据库注释基因数为111个，COG数据库注释基因数为113个。

表5 注释的差异表达基因数量统计表

2.3.2 差异基因GO富集分析

共有25 709条Unigene被成功注释到三大类的55个二级功能条目中。3个处理在细胞组分共富集了19个功能条目，在分子功能共富集了13个功能条目，在生物学过程共富集了23个功能条目。细胞组分上1、10、100 mg/kg干土处理与CK相比，胞外区、细胞、膜部分、突触部分、细胞组分等二级分类的功能条目显著增加。分子功能上1、10、100 mg/kg干土与CK相比，酶的催化活性、分子功能调节转录调节活性等二级分类的功能条目显著增加。在生物学过程上，1、10、100 mg/kg干土处理与CK处理相比，繁殖、细胞进程、繁殖过程、发育过程、生长、刺激反应、生物调节、解毒等二级分类的功能条目显著增加。

2.3.3 差异基因KEGG富集分析

在图5中，3个浓度的处理与CK之间的差异基因富集的KEGG一级代谢通路均为细胞进程、遗传信息处理、人类疾病、环境信息处理、代谢作用、有机系统六大通路。在细胞进程一级代谢通路中，各个处理的差异基因主要富集于细胞的生长与死亡、真核生物的细胞群落、运输和分解代谢3个二级代谢通路中，同时高浓度（10、100 mg/kg干土）的处理在细胞进程代谢通路富集的差异基因显著高于低浓度（1 mg/kg干土）的处理。在环境信息处理一级代谢通路中，低浓度（1 mg/kg干土）的处理的差异基因主要集中于信号转导、信号分子的相互作用二级代谢通路。高浓度（10、100 mg/kg干土）处理相比于低浓度又增加了膜信号二级代谢通路，同时在其他2个二级代谢通路富集的差异基因也显著增加。在遗传信息处理的一级代谢通路中，各个处理的差异基因主要富集于折叠、分类和降解，高浓度（10、100 mg/kg干土）处理均有差异基因在翻译二级代谢通路富集。在新陈代谢一级代谢通路中，低浓度（1 mg/kg干土）处理的差异基因在糖代谢二级通路富集最多，在其他6条二级通路富集差异基因较少。在高浓度（10、100 mg/kg干土）处理中，在脂质代谢二级通路富集的差异基因最多，100 mg/kg干土的处理组在糖代谢二级通路中与低浓度处理相比，富集的差异基因显著增多。同时高浓度处理组中差异的二级通路也比低浓度处理增多。在有机系统一级代谢通路中，低浓度的处理差异基因主要富集在内分泌系统和免疫系统。在高浓度的处理中差异基因显著增多，同时富集的二级通路也显著增加。主要表现在原有的二级通路富集基因数增加，同时增加了在循环系统、环境适应、消化系统等二级通路的富集基因。

图5 CK vs Cya.L（A）、CK vs Cya.M（B）、CK vs Cya.H（C）的差异表达基因KEGG 分类图

2.3.4 RT-qPCR对关键基因的验证

分别选取了10 个差异显著的基因进行了RT-qPCR验证，通过计算出相对表达量与CK的相对表达量倍数得出图6。可以看出经过RT-qPCR后的结果与转录组学的结果基本一致。结果显示控制谷胱甘肽S-转移酶基因evm.TU.Chr02.2134表达下调，控制热激蛋白70基因evm.TU.Chr11.1027上调，控制细胞色素P450 18 a1基因evm.TU.Chr01.306和evm.TU.Chr04.1250均下调。控制钙调蛋白-腺苷酸环化酶反应的基因evm.TU.Chr11.1424下调。控制调节幼虫生长发育的基因evm.TU.Chr09.1717下调，证明本次转录组测序结果可靠性较高。

图6 RT-qPCR 对关键基因的转录水平验证

3 讨论与结论

基于Illumina HiSeq测序技术探究溴氰虫酰胺对赤子爱胜蚓转录水平的影响。通过差异基因的韦恩图和火山图可以看出，溴氰虫酰胺1、10 mg/kg干土和100 mg/kg干土处理与CK相比，差异基因显著增加，这表明溴氰虫酰胺对蚯蚓转录水平产生显著影响。Zhu等[20]研究土壤中的TiO2纳米粒子对蚯蚓的影响也是如此，从5 mg/kg到500 mg/kg与CK相比，差异基因也显著增加。

GO数据库是2000年构建的一个结构化标准生物学注释系统，目的是建立一个基因及其产物适用于各个物种的知识标准词汇体系[30]。GO注释系统是一个有向无环图，其中包括3个主要分支，分别为生物学过程、细胞组分和分子功能[31]。通过GO富集分析可以看出，暴露于溴氰虫酰胺后，蚯蚓在细胞组分上与神经信号传导的相关二级分类条目增加，这可能与溴氰虫酰胺作用于鱼尼丁受体，导致细胞内外钙离子的变化有关。在分子功能的二级条目的基因富集说明蚯蚓暴露于溴氰虫酰胺后，蚯蚓体内的各种抗氧化酶及其他酶的活性会发生变化，同时蚯蚓体内开始对于溴氰虫酰胺的毒性做出反应调节，因此差异基因富集于酶的催化活性、分子功能调节和转录调节活性上。同时我们之前的研究也表明溴氰虫酰胺能够导致蚯蚓体内抗氧化酶活性的改变[25]。但是蚯蚓暴露于溴氰虫酰胺后差异基因最多的还是富集于生物学过程上，可以看出随着浓度的增加，溴氰虫酰胺对蚯蚓的繁殖、生长的生物学过程的影响也在加强。刺激反应、生物调节、解毒二级功能条目主要是对外界的刺激条件做出反应，同时进行自身的调节作用和对污染物的解毒作用。随着浓度的增加，这些二级功能条目的富集基因增加，说明蚯蚓对于溴氰虫酰胺的毒性应激反应增加。Liu等[32]研究氯虫苯甲酰胺对于蚯蚓生长繁殖的影响，试验结果表明，随着浓度的增加，蚯蚓的体重降低，繁殖量减少。同时也有研究表明蚯蚓暴露于污染物后体内的抗氧化酶活性、解毒酶活性、抗逆基因以及生殖基因等相关基因的表达量发生改变[33-35]。

作为有关Pathway的主要公共数据库[36]，KEGG主要提供代谢途径的查询，其中包括碳水化合物、核苷、氨基酸等代谢途径及有机物的生物降解途径，而且对催化各步反应的酶进行了全面的注解。KEGG是对生物体内代谢分析和代谢网络研究的有效工具。通过KEGG差异基因富集分析可以看出，低浓度的溴氰虫酰胺处理下主要影响环境信息处理和代谢途径一级代谢通路，差异基因主要富集于这2个通路。而高浓度处理下，则影响环境信息处理、代谢途径和有机系统一级代谢通路，同时富集的差异基因显著增多，差异基因富集到的二级通路也显著增加。这可能是溴氰虫酰胺浓度增加对蚯蚓的影响增大引起的。在代谢途径一级通路中，低浓度的溴氰虫酰胺影响蚯蚓的差异基因主要富集于糖代谢途径，而在高浓度处理下差异基因还更多的富集于脂质代谢途径。这可能是溴氰虫酰胺浓度的增高不仅对于糖代谢途径影响加大，同时也开始影响脂质代谢途径。内分泌系统和免疫系统是机体内起主导作用的调节系统，它们密切联系，互相配合，维持内环境相对稳定[37]。当随着浓度的升高，溴氰虫酰胺开始影响蚯蚓的消化系统和循环系统。我们认为溴氰虫酰胺可能首先导致蚯蚓的内分泌系统和免疫系统相关基因发生改变，然后导致内环境的紊乱，影响了蚯蚓的各种代谢途径，最终导致循环系统和消化系统的相关基因表达发生改变。陈婷[38]研究重金属铅对蚯蚓的生态毒性表明，铅能够导致蚯蚓能量代谢、氨基酸代谢等代谢途径失衡。

通过转录组学结果我们可以看出，低浓度的溴氰虫酰胺对于蚯蚓的转录水平影响较小，影响的代谢通路较少。随着浓度的升高，溴氰虫酰胺对蚯蚓影响的差异基因增多，差异基因富集的代谢通路在增加，说明溴氰虫酰胺对蚯蚓转录水平影响增大，可能对蚯蚓的毒性增大。