基于转录组测序的卤虫卵生和卵胎生相关基因表达研究

欧阳雪梅 郑毓毓 韩学凯 许如意 隋丽英

(天津科技大学海洋与环境学院, 亚洲区域卤虫参考中心, 天津 300457)

小型甲壳动物卤虫(Artemia)广泛分布于沿海盐池和内陆盐湖高盐水体中, 是卤水生态系统食物链的重要组成部分和生物调节者[1]。卤虫具有两种生殖方式, 分别为孤雌生殖和两性生殖, 同时雌性卤虫以卵生(Oviparous, OVI)、卵胎生(Ovoviviparous, OVOVI)两种方式进行繁育, 即雌性卤虫以产休眠卵或以生产无节幼体的方式来繁殖后代[2]。

卤虫的繁殖模式可能受内在遗传因子与外部环境因素的影响[3], 对于其繁殖模式形成的分子机制研究, 主要集中在休眠卵形成相关功能基因克隆和表达方面。目前已发现多种卤虫繁殖相关基因,如BAP1是一种参与调节G1/S期和休眠卵形成的分子, 有助于维持休眠卵的周期停滞[4]; 还鉴定了参与卤虫休眠卵形成过程的胁迫应答转录辅因子p8[5]。有研究报道, 在卤虫的卵生途径中Afr-AMPKα基因的mRNA表达水平显著降低, 这可能与卤虫的不同繁殖方式有关[6]。此外, 在卤虫中, 染色质组蛋白的积累可能导致细胞周期阻滞, 影响胚胎发育, 并参与卤虫休眠胚胎形成过程[7]。然而目前哪些基因在卤虫产卵或产幼过程中起关键作用还不清楚, 其不同繁殖方式产生的内在分子机制仍需进一步探究。

近年来转录组测序技术已被广泛的应用于甲壳动物发育、繁殖和免疫等方面的研究, 并用于发掘参与其各种生理过程的关键基因及信号通路[8]。在卤虫转录组学研究方面, 有报道采用转录组测序技术进行卤虫耐盐相关基因挖掘, 为研究卤虫独特的耐盐机制提供了基础资料[9]; Diego等[10]构建了雌雄差异转录组, 并开发了雌雄相关SNP标记, 为卤虫和其他甲壳类动物性别偏向基因的表达提供了重要信息。但迄今为止, 有关卤虫卵生和卵胎生不同繁殖模式的转录组学研究还未见报道, 相关的繁殖模式差异基因有待进一步发掘。因此本实验利用Illumina HiSeq 2500测序平台构建了孤雌生殖卤虫在卵生和卵胎生两种繁殖模式下的转录组文库,并对差异表达基因进行功能注释, 筛选与卤虫生殖相关的候选基因并开展表达研究, 为揭示其不同繁殖模式的发生机制奠定理论基础。

1 材料与方法

1.1 实验材料

本实验选取新疆艾比湖孤雌生殖卤虫, 在28℃和光照条件下孵化培养。在6个1 L锥形瓶中每个饲养200只卤虫并投喂新鲜培养的盐生杜氏藻Dunaliella salina。根据王智超[11]在对光诱导卤虫卵生卵胎生不同繁殖模式的研究, 得到相应繁殖模式的卤虫, 以备后用。收集晚期胚胎的卵巢组织,解剖后迅速置于液氮中, –80℃冰箱中保存, 用于分析候选生殖相关基因在不同繁殖模式的表达情况。

1.2 文库构建与测序

从卵巢组织样品中提取总RNA, 检验其浓度、纯度及完整性。在样品检测合格后, 将mRNA反转录为双链cDNA, 并将其纯化。经末端修复和接头连接后通过PCR富集得到最终的cDNA文库, 并对其质量进行检测, 然后在Illumina HiSeq 2500平台进行转录组测序[12]。

1.3 从头组装和功能注释

为获取后续分析所需的参考序列, 将测序得到的原始序列进行过滤, 首先对原始数据进行测序数据质量评估, 将带接头的、低质量的数据过滤得到测序读段, 后将获得的测序读段用Trinity[13]软件进行序列拼接, 将Trinity拼接得到的转录本序列, Corset层次聚类并去冗余后得到最长Cluster序列进行后续的分析。为获得全面的基因功能信息, 将拼接所得的序列与NR、NT、PFAM、KOG、KO、SwissProt、GO和KEGG等数据库进行比对注释[14]。

1.4 差异表达基因

将测序读段运用Bowtie 2进行比对, 利用DEG-seq得到卤虫卵生和卵胎生两种繁殖模式间的差异表达基因[14], 在筛选过程中, 以qvalue<0.005且log2|fold_change|>1为筛选阈值。为了获得关于差异表达基因更多信息, 通过KEGG通路富集和GO功能对其进行分析[15]。

1.5 qRT-PCR分析

从转录组筛选的差异基因中, 参照转录组注释并结合文献查找筛选得到Serine/threonine-like kinase(SLK)、Double sex and mab-3 related transcription factor(Dmrt)、Feminization-1(fem-1)和5-hydroxytryptamine receptor gene(5-HT), 以及参与原始生殖细胞的迁移和维持的Nanos和调节成熟促进因子亚基的周期蛋白cyclin B等6个候选生殖相关基因并进行引物设计(表1), 以β-actin作为内参基因对孤雌卤虫卵巢不同繁殖模式的mRNA表达量进行qRT-PCR检测, 采用2–∆∆Ct法[16]计算目的基因在生殖腺不同繁殖模式的相对表达水平。

表1 qRT-PCR引物Tab.1 qRT-PCR primers

1.6 候选生殖相关基因的蛋白结构域和系统发育分析

通过在线软件对候选生殖相关基因进行保守结构域预测(https : //www.ncbi.nlm.nih.gov/Structure/cdd/docs/cdd_search.html); 采用NCBI(http://www.ncbi.nlm.nih.gov)数据库对候选生殖相关基因的氨基酸序列进行BLAST比对, 分析与其他物种的相似性; 利用MEGA7.0软件将不同物种的氨基酸序列进行多重比对, 构建系统进化树。

2 结果

2.1 Unigene的组装和功能注释

从原始数据中将低质量序列剔除后, 共得到了94901316条序列, 从头组装后得到97146条unigene(表2), 长度大于2000 bp的有12338条, 平均长度为1023 bp, N50 unigene有1558个, N90 unigene有450个。

表2 Unigene序列在公共数据库中的注释结果Tab.2 Comment results of the unigene sequence in the public database

将获得的97146条Unigene 数据信息分别在蛋白质家族域数据库(Pfam)、蛋白质序列数据库(Swiss Prot)、非冗余蛋白数据(NR)、真核生物蛋白相邻类的聚簇(KOG)、基因本体论(GO)数据库和基因组百科全书(KEGG)七大数据库进行比对注释。结果显示, 注释率最高的是GO数据库, 达到37.39%; 注释率最低的为NT数据库中, 仅为13.01%(表2)。转录组原始数据已上传至基因组数据中心,登录号为: PRJNA788349。

2.2 差异表达基因分析

转录组数据分析获得卵生和卵胎生卤虫卵巢的1452个差异表达基因, 上调基因与下调基因分别为601个和851个(图1)。

图1 卵胎生与卵生卤虫卵巢中上调和下调的差异表达基因Fig.1 Up-regulated and down-regulated differentially expressed genes in ovoviviparous vs oviparous Artemia ovary

GO数据库涵盖3个方面(图2): 生物过程(Biological process)、细胞组成(Cellular component)和分子功能(Molecular function)。统计结果显示: 注释到生物过程数据库中的差异表达基因有1243个,其中22.84%富集在有机氮化合物代谢, 碳水化合物衍生物代谢过程占11.19%; 注释富集至分子功能数据库的有530个差异表达基因, 17.13%为结构分子活性, 几丁质结合占6.06%; 富集至细胞组分的有306个, 细胞组分占11.54%, 核蛋白复合体占7.11%,核糖体蛋白占6.41%。

图2 GO富集分析图Fig.2 Enrichment scatter diagram of GO pathway

在本研究中, 对卵生、卵胎生的差异表达基因进行KEGG富集分析, 共富集到220条代谢调节通路, 挑选了富集最显著的20条信号通路条目在图3中进行展示。表3为卵生和卵胎生差异表达上调和下调基因富集的前10个信号通路, 其中上调富集基因数目最多的通路为蛋白质消化吸收, 下调富集基因数目最多的为核糖体通路。

表3 卤虫表达上调和下调基因富集的10个信号通路Tab.3 The top -ten enriched pathways of up -regulated and down-regulated unigene in Artemia

图3 KEGG pathway富集散点图Fig.3 Enrichment scatter diagram of KEGG pathway

从KEGG富集通路中筛选后利用NCBI比对,挑选可能对卤虫不同生殖模式起到调控作用的基因, 并结合文献查找获得了Serine/threonine-like kinase(SLK)、Double sex and mab-3 related transcription factor(Dmrt3)、Feminization-1(fem-1)和5-hydroxytryptamine receptor gene(5-HT), 以及参与原始生殖细胞的迁移和维持的Nanos和调节成熟促进因子亚基的周期蛋白cyclin B共6个候选生殖相关基因。

2.3 候选生殖相关基因编码蛋白的保守结构域预测

对SLK蛋白的保守结构域预测发现该蛋白共包含1个结构域, 其中第334—581位氨基酸为 MPP_PPP_family结构域(图4A)。对Dmrt3蛋白的保守结构域预测发现该蛋白共包含1个结构域, 其中第68—208位氨基酸为 DM superfamily结构域(图4B)。对Fem-1蛋白的保守结构域预测发现该蛋白共包含5个结构域, 其中第92—183位氨基酸、第158—240位氨基酸、第44—116位氨基酸都为Ank_2结构域, 第527—571位氨基酸为 Ank_2超家族结构域, 第486—582位氨基酸为ANKYR结构域(图4C)。对CyclinB蛋白的保守结构域预测发现该蛋白共包含2个结构域, 其中第61—190位氨基酸为CYCLIN_CCNB1-like_rpt1结构域, 第195—315位氨基酸为CYCLIN_SF 超家族结构域(图4D)。对5-HT蛋白的保守结构域预测发现该蛋白共包含1个结构域,其中第32—371位氨基酸为7tmA_5-HT7结构域(图4E)。对Nanos蛋白的保守结构域预测发现该蛋白共包含1个结构域, 其中第232—279位氨基酸为zf-nanos结构域(图4F)。

2.4 候选生殖相关基因的相似性分析及系统进化树构建

分别对候选生殖相关基因构建系统进化树分析可知,SLK基因氨基酸序列与昆虫纲中的马铃薯的SLK序列相似性最高(图5A)。Dmrt3基因氨基酸序列与两栖纲、鳃足纲(节肢动物门)相关物种的Dmrt3基因氨基酸序列相对较近, 其中Dmrt3基因与两性卤虫和盐湖卤虫的Dmrt3序列相似性最高(图5B)。fem-1基因氨基酸序列与鳃足纲(节肢动物门)和桡足纲相关物种的fem-1基因氨基酸序列相对较近,其中fem-1基因与两性卤虫的fem-1序列相似性最高(图5C)。cyclinB基因氨基酸序列与辐鳍鱼纲(脂鲤科)相关物种的cyclinB基因氨基酸序列相对较近(图5D)。5-HT基因氨基酸序列与介形虫的5-HT序列相似性最高(图5E)。Nanos基因氨基酸序列与昆虫纲的相关物种的Nanos基因氨基酸序列相对较近, 其中Nanos基因与多胚跳小蜂和竹节虫的Nanos序列相似性最高(图5F)。

图5 SLK、Dmrt3、Fem-1、CyclinB、5-HT和Nanos蛋白系统进化树Fig.5 Phylogenetic tree of SLK, DMRT3, Fem-1, CyclinB, 5-HT and Nanos proteins

2.5 qRT-PCR验证和卤虫候选生殖相关基因的表达分析

qRT-PCR的结果表明, 所选的与生殖相关候选基因在卵生和卵胎生两种繁殖模式中的表达量具有显著性差异。SLK、Dmrt3、fem-1、cyclin B、5-HT和Nanos在卤虫卵巢中的表达量均为卵生高于卵胎生(P<0.05; 图6)。

图6 SLK、Dmrt3、fem-1、cyclinB、5-HT和Nanos基因在不同繁殖模式下的表达情况Fig.6 Expression of SLK, Dmrt3, fem-1, cyclinB, 5-HT and Nanos genes in different reproductive modes of Artemia oocysts

3 讨论

3.1 转录组测序分析

本研究以卵生和卵胎生生殖方式的孤雌卤虫为材料进行转录组测序, 获得了97146条非冗余Unigene序列, 其中注释到已知基因的Unigene为37.39%, 注释率相对较低。在其他关于甲壳动物的转录组测序比对注释的报道中也出现了许多未被注释到的情况, 如日本沼虾(Macrobrachium nipponensis)转录组中只有23.89%为被注释的序列[17], 银色小长臂虾(Palaemonetes argentinusfan)中注释到的序列只有24%[18], 而凡纳滨对虾(Litopenaeus vannamei)的注释率也只有37.28%[19]。出现转录组注释率低的原因可能是在目前公共基因数据库中相关物种基因注释信息有限。而随着卤虫基因组测序结果的完成与公布, 这些潜在的与生殖相关功能基因有望得到进一步挖掘。

已有研究证实, 甲壳动物的卵巢发育受到内部因子如类固醇激素、神经多肽、神经递质以及外源环境因子的综合调控[14,20—24]。本研究从GO富集数据库中筛选到一些参与卵泡细胞发育、卵母细胞分化、生殖细胞发育、卵壳形成和性腺发育等过程的基因, 这些基因可能参与卤虫不同繁殖模式的生殖过程。Unigene在KOG真核直系同源数据库中注释到参与信号转导机制、细胞周期控制、细胞分裂、染色体分裂及核苷酸运输和代谢等相关基因, 可能与卤虫生殖有关。另外从KEGG代谢通路注释分析中, 筛选到一些代谢通路如甾体激素生物合成[25]、PI3K-Akt信号通路[26]、TGF-beta信号通路[27—29]及Wnt信号通路[30], 这些都与甲壳动物卵巢的发育和成熟相关, 因此我们推断这些代谢通路可能参与卤虫不同生殖模式的发生过程。

3.2 生殖相关候选基因表达分析

卤虫产卵和产幼的机制复杂, 可能与p26基因[31]、卵壳[32]、卵黄蛋白和外部环境因素有关[3,33]。本研究的不同生殖模式中的319个差异表达基因显著富集至抗原处理和呈递、内质网蛋白加工过程、PI3K-Akt信号通路、核糖体、吞噬体等KEGG通路中, 这些代谢通路对于卤虫生长发育相关蛋白的合成和转运起到调控作用。为进一步探究这些过程是否与卤虫不同繁殖模式的生殖过程有关, 从这些通路中筛选到了可能与卤虫不同生殖模式相关的候选基因SLK、Dmrt3、cyclin B、Nanos、fem-1和5-HT, 对卤虫在不同繁殖模式下生殖相关基因进行了蛋白结构域预测和表达情况分析。研究发现上述6个基因预测得到的蛋白结构域与之前报道的相应基因保守结构域一致, 且均与生殖相关。这预示着本研究所选6个生殖候选基因在卤虫中可能同样参与生殖发育过程。

在关于卵胎生卤虫研究中发现,SLK主要参与调控细胞周期、细胞凋亡等过程[34], 并且胚胎发育时期SLK蛋白表达量比卵母细胞发育时期高[35]。对Artemia franciscana的雌虫及雄虫做转录组分析发现,Dmrt相关基因调控性腺分化过程中蜕皮素的生物合成, 并参与性别决定过程, 表明卤虫中Dmrt基因参与生殖腺的早期发育[36]。对中华绒螯蟹(Eriocheir sinensis)fem-1基因研究发现, 该基因的表达量在卵巢组织和睾丸组织中最高, 并且在早期胚胎发育时高度表达, 说明该基因对于性腺发育后期有潜在作用, 同时与中华绒螯蟹的生殖密切相关[37]。在克氏原螯虾(Procarabarus clarkii)中发现cyclin B基因在细胞分化、脂质转运体活性和减数分裂周期过程中都发挥着重要作用, 说明cyclin B在甲壳类卵母细胞成熟调控中起着不可或缺的翻译激活作用[38]。5-HT既是神经递质又是神经激素, 存在于卵巢和雄激素腺中, 并对性腺成熟产生影响。5-HT基因在远洋梭子蟹(Portunus pelagicus)的中枢神经系统和卵巢中存在, 并可能通过作用于中枢神经系统的某些生殖激素对远洋梭子蟹的繁殖起着关键作用[39]。Nanos基因在成年雌性斑马鱼(Brachydanio rerio)生殖系的初期卵母细胞中表达, 对于斑马鱼卵母细胞的持续产生起着重要作用[40]; 在多倍体银鲫(Carassius gibelio)中发现,Nanos2在生殖干细胞中呈现出偏向表达的特征, 调控卵巢早期胞囊的发育过程[41]。在本研究中, 通过实时荧光定量PCR得到SLK、Dmrt、fem-1、cyclin B、5-HT和Nanos的表达情况, 在卤虫不同繁殖模式下, 以上基因在卤虫卵生繁殖模式中的表达量均显著高于卵胎生表达量, 因此推测这些基因可能更多参与卤虫休眠卵形成过程。

4 结论

利用高通量测序技术分析了卤虫卵生和卵胎生差异转录组文库, 结合KEGG富集通路及GO功能分类对两种繁殖模式间的差异表达基因进行了分析。挑选6个生殖相关候选基因, 分别对其蛋白结构域进行预测并分析其在不同生殖模式下的表达情况, 发现SLK、Dmrt3、cyclin B、fem-1、5-HT和Nanos在卵生卤虫卵囊的表达量显著高于卵胎生卤虫, 说明这些基因可能更多地参与卤虫休眠卵形成过程, 影响卤虫卵生繁殖模式形成。在后续研究中将进一步利用RNAi技术精确地靶向敲降目的基因,检测出下游基因表达量的变化, 并结合蛋白印迹技术对目标基因功能进行深入研究。