泥蚶(Tegillarca granosa)属软体动物门(Mollusca)、双壳纲(Bivalvia)、列齿目(Taxodonta)、蚶科(Arcidae), 富含营养, 且含有特有的维生素B₁₂和血红蛋白, 味道鲜美, 具有补血养胃和除淤散积的功效, 是我国传统养殖贝类(李太武等, 2003)。血细胞是机体抵抗病原微生物的重要载体, 贝类血细胞能吞噬病原微生物或形成包囊包裹异物, 还在伤口修复、炎症反应过程中发挥重要作用(孙敬锋等, 2006)。

转录组测序(RNA sequencing)是利用大规模测序技术直接对cDNA序列进行测序的技术, 能获得细胞在某一状态下所能转录出来的所有mRNA。RNA-seq具有通量高、成本低、灵敏度高等优点, 是一种有效的基因组研究和功能基因鉴定方法, 利用转录组测序研究转录谱特征近年来呈爆发式增长(陈雪峰等, 2019; Wang et al, 2019)。目前有关泥蚶研究主要集中在养殖技术、病害防控以及单一基因克隆(李太武等, 2003; 董迎辉等, 2013; Liu et al, 2017; 许家辉等, 2018)等方面。本研究采用RNA-Seq技术测序了鳗弧菌刺激前后泥蚶血细胞转录组, 通过转录组序列组装获得了7个免疫相关基因蛋白序列, 并对7个免疫相关基因蛋白序列理化性质、二级结构、三级结构、保守结构域、进化发生和mRNA表达水平进行了分析, 探索其在细菌感染过程中所发挥的作用, 以期为泥蚶抗病机理研究提供理论基础。

1 材料与方法 1.1 样本准备

实验用泥蚶来自浙江省舟山市沈家门国际水产城, 选择贝壳无损伤、对外界刺激反应迅速的泥蚶, 带回实验室暂养5d, 暂养盐度为20±1、温度(18±2)℃、溶解氧高于5mg/L。暂养结束后, 随机挑选5个大小相一致的泥蚶置于含有5×10⁶cfu/mL鳗弧菌悬液中, 5个泥蚶置于无菌海水中。24h后, 采用无菌注射器从两组泥蚶外套腔中取血, 获得的血淋巴立即与Trizol (Invitrogen, CA, USA)试剂采按1 : 5混合, 将混合好的血淋巴冻存于液氮中备用。

1.2 RNA提取、cDNA文库构建与测序

按照TRIzol^® Reagent试剂盒操作步骤提取总RNA, 提取的总RNA采用Bioanalyzer 2100和RNA 1000 Nano LabChip Kit (Agilent, CA, USA)定性定量试剂盒进行质检, 获得的RNA满足总量≥5μg、浓度≥200ng/μL、OD₂₆₀/OD₂₈₀为1.8—2.2、RIN > 8.0, 进行文库构建。先使用连接有Oligo(dT)的磁珠富集mRNA, 富集的mRNA被片段化试剂(Fragmentation Buffer)随机打断成短片段, 以片段化的mRNA为模板, 用六碱基随机引物(Random hexamers)合成一链cDNA, 随后加入缓冲液、dNTPs、RNaseH和DNA Polymerase I进行二链cDNA合成。AMPure XP beads纯化双链产物, 利用T4 DNA聚合酶和Klenow DNA聚合酶将DNA的粘性末端修复为平末端, 3′末端加碱基A并加接头, AMPureXP beads进行片段选择, 之后用USER酶降解含有U的cDNA第二链, 然后进行PCR扩增获得最终测序文库。文库质检合格后采用Illumina Hiseq4000进行测序, 测序读长为双端2*150bp (PE150)。

1.3 转录组de novo组装

测序产生的原始数据(Raw Date)采用使用Cutadapt去除测序reads中接头后, 使用Fqtrim过滤掉不合格的序列后得到有效数据。具体为先对测序reads进行窗口法质量扫描, 扫描窗口默认为6bp, 当窗口内平均质量值低于20时, 将read从窗口起始到3′终止的部分截掉, 去除poly-A/T、去除截断后长度小于100bp以及N的含量在5%以上的序列, 余下的序列为有效数据。过滤后的数据(Clean reads), 使用Trinity软件进行从头组装, reads首先被拼接成长片段contigs, 去除冗余contigs后将contigs进行分组, 最后连接成两端不能再延长的转录片段(Unigene)。

1.4 基因注释和结构分析

采用软件DIAMOND, 设定期望值E Value < 1e-5, 在NCBI-NR (http://www.ncbi.nlm.nih.gov/)、GO (http://www.geneontology.org)、KEGG (http://www.genome.jp/kegg/)、Pfam (http://pfam.xfam.org/)和Swiss-Prot (https://www.uniprot.org/)数据库中进行比对, 获得Unigenes注释信息。基于注释信息, 利用NCBI数据库SmartBlast工具对7个候选基因进一步进行比对, 以保证注释信息正确, 7个候选基因编码区由NCBI网站获得(http://www.ncbi.nlm.nih.gov/)。

利用ExPaSy中的ProtParam (http://web.expasy.org/protparam/)对蛋白序列进行理化性质分析; 利用SMART (http://smart.embl-heidelberg.de/)分析序列保守结构域(包括信号肽)特征; 利用Prabi程序(https://npsa-prabi.ibcp.fr/cgi-bin/secpred_hnn.pl)对蛋白质二级结构进行预测; 利用在线工具Swiss-model (http://swissmodel.expasy.org)进行三级结构分析, 选取同源性大于30% (Seq Identity)、GMQE和QMEAN值大的模板生成蛋白质的三级结构, 然后使用PyMOL Viwer软件进行编辑, 最终绘制出蛋白质的三级结构模型图。

1.5 差异表达基因分析

采用Salmon软件(Patro et al, 2017)将测序获得的Clean reads比对到Trinity组装的序列, 从而获得每个unigene上的read数目, 然后计算TPM (Transcripts Per Kilobase Million), 从而获得Unigene的表达量。采用R软件包中的edgeR筛选差异表达基因, 筛选阈值为P值< 0.05且|log₂Fold Change| > 1。

1.6 系统发生分析

以人类Homo sapiens、斑马鱼Danio rerio及其他水生生物作为参考物种, 利用MEGA-X对7种免疫相关基因进行Clustal W序列比对, 用最大似然法(ML法)构建系统发育树, 构建系统发育树前, 先计算最优模型, 采用最优模型进行建树, 同时进行步长循环1000检验, 所用模型及序列见表 1。

2 结果 2.1 测序数据统计、组装和注释

通过高通量测序泥蚶血细胞中共获得9.01G raw base和60001612 raw reads, 过滤后获得8.75G valid bases和59199480 valid reads, valid reads的GC (guanine-cytosine)含量为42.41%, 组装后, 共获得21823条Unigenes, 原始数据已上传至国家基因组科学数据中心(https://bigd.big.ac.cn/, Project:PRJCA002006, 提交编号: subCRA002389)。通过与Pfam、NR、GO、Swissprot和KEGG数据库比对, 共有10691条Unigenes与已知基因同源, GO数据库中注释到8301 (38.04%)条, KEGG中注释到5376 (24.63%)条, Swissprot中注释到7686 (35.22%)条, Pfam中注释到9868 (36.05%)条, NR中注释到10636条(48.74%)条(图 1)。物种同源性分析表明, 与NR数据库中太平洋牡蛎Crassostrea gigas同源序列比例最高, 其次是栉孔扇贝Lottia gigantea(图 1)。

2.2 溶酶体酶

泥蚶血细胞转录组序列中鉴定出4种溶酶体酶, 包括2种糖苷酶、硫酸酯酶B和天冬酰胺氨基葡萄糖酶。溶酶体α-葡萄糖苷酶(GAA)长度为3168bp, 开放阅读框长度为2347bp, 编码779个氨基酸, 预测分子量约为87.9kDa、理论等电点6.68;丝氨酸和苏氨酸含量最高, 占7.7%;不稳定系数为34.03, 属于稳定蛋白。其N端(第1—96位氨基酸)是一个β桶状结构, 一般位于真核生物麦芽糖酶-葡糖淀粉酶催化域的N端, 属于糖基水解酶31的N-末端保守结构域; 第97—163氨基酸属于糖苷水解酶家族31的保守结构域; 第183—651氨基酸为O-糖基水解酶31结构域(图 2a)。二级结构分析其含有α螺旋、β折叠及无规则卷曲。SWISS同源性分析显示泥蚶溶酶体α-葡萄糖苷酶与模板5kzw.1.A、5nn3.1.A同源性均超过50%, 自动建模只有一种三级结构被建立(图 6a)。多重序列比对显示泥蚶α-葡萄糖苷酶与人、小鼠和斑马鱼同源性为50%左右, 与太平洋牡蛎Crassostrea gigas、美洲巨牡蛎Crassostrea virginica和虾夷扇贝Mizuhopecten yessoensis同源性超过60%;系统进化树表明泥蚶α-葡萄糖苷酶与虾夷扇贝优先聚为一支, 然后与太平洋牡蛎和美洲巨牡蛎相聚(图 2b)。

溶酶体α-甘露糖苷酶(LAMAN)长度为3473bp, 开放阅读框长度为3180bp, 编码1060个氨基酸, 预测分子量约为120.9kDa、理论等电点8.77;含量最高的氨基酸为亮氨酸, 占9.4%;不稳定系数为32.94, 属于稳定蛋白; 二级结构由α螺旋、β折叠及无规则卷曲组成。其N端(第1—24位氨基酸)具有信号肽, 第43—360氨基酸属于糖苷水解酶家族38的保守结构域, 第365—442为α-甘露糖苷酶中间保守结构域, 第491—1054氨基酸为糖苷水解酶家族38的C端保守结构域(图 3a); 三级结构显示其与Swiss-Model数据库中α-甘露糖苷酶模板1o7d.1.A和1o7d.1.B同源性超过60%, 与模板1o7d.1.D同源性超过55%, 自动建模获得地三级结构见图 6b。多重序列比对显示其与人、斑马鱼、美洲巨牡蛎和虾夷扇贝同源性超过50%。系统进化树分析显示泥蚶α-甘露糖苷酶与美洲巨牡蛎亲缘关系最近(图 3b)。

芳基硫酸酯酶B(ASB)长度为1743bp, 开放阅读框为1536bp, 编码512个氨基酸, 预测分子量约为59.1kDa、理论等电点6.32;含量最高的氨基酸为甘氨酸, 占8.2%;不稳定系数为47.08, 属于不稳定蛋白。二级结构分析其含有α螺旋、β折叠及无规则卷曲。其N端具有信号肽和硫酸酯酶结构域(图 4a)。其与Swiss-model数据库中模板B1fsu.1.A同源性高达47.01%, 自动建模只有一种三级结构被建立(图 6c)。多重序列比对显示其与太平洋牡蛎、美洲巨牡蛎和虾夷扇贝同源性高达55%以上, 与人类芳基硫酸酯酶B同源性超过40%;系统进化分析显示泥蚶芳基硫酸酯酶B优先与虾夷扇贝聚为一支, 然后与太平洋牡蛎和美洲巨牡蛎相聚(图 4b)。

天冬酰胺氨基葡萄糖酶(AGA)长度为1256bp, 开放阅读框为1044bp, 编码348个氨基酸, 预测分子量约为120.9kDa、理论等电点8.77;含量最高的氨基酸为亮氨酸, 占9.4%;不稳定系数为32.94, 属于稳定蛋白; 具有α螺旋、β折叠及无规则卷曲组成的二级结构。第16—342位氨基酸残基是天冬酰胺酶2结构域(图 5a)。三级结构分析显示其与模板1p4k.1.A、4r4y.1.A、同源性超过40%, 与模板1apz.1.A同源性为超过60%, 以模板1apz.1.A建立的三级结构见图 6d。多重序列比对显示其与太平洋牡蛎、美洲巨牡蛎和虾夷扇贝相似度高达55%以上, 与人类、斑马鱼、大黄鱼等相似度超过50%。系统进化分析显示泥蚶天冬酰胺氨基葡萄糖酶优先与虾夷扇贝聚为一支(图 5b)。

2.3 整合素

泥蚶血细胞转录组序列中, 鉴定出整合素β3亚基(ITB3), β3亚基长度为3455bp, 编码793个氨基酸, 预测分子量约为87.9kDa、理论等电点4.94;苏氨酸含量最高, 占7.9%;不稳定系数为44.51, 属不稳定蛋白。SMART保守结构域分析其具有人整合素β3相似结构域, 胞外序列N端有信号肽和整合素亚基保守结构域(integrin beta subunits domain, INB)(第27—452位氨基酸), 第631—714位氨基酸是整合素β尾部保守区域, 第714—736位氨基酸为跨膜区域, 胞内有整合素β-cyt结构域, 泥蚶整合素第461—622位氨基酸位置上还有3个EGF样(epidermal growth factor)重复(图 7b)。具有α螺旋、β折叠及无规则卷曲组成的二级结构。通过Swiss-model进行蛋白质三级结构同源建模, 结果显示泥蚶整合素β3与模板3ije.1.B、6avr.1.B、6bxj.1.A、4mmy.1.B、4mmx.1.B等同源性均超过35%, 模板6avr.1.B建模效果最优, 其GMQE(0.63)和QMEAN(-3.31)值均最大, 泥蚶整合素β3三级结构见图 6e。NCBI结果显示泥蚶整合素β3序列与美洲牡蛎整合素β3同源性超过45%, 与人和斑马鱼同源性超过35%。系统进化树表明:泥蚶整合素β3亚基与太平洋牡蛎和美洲巨牡蛎亲缘关系最近(图 7c)。

2.4 热休克蛋白

泥蚶血细胞转录组中鉴定出1种新型热休克蛋白, 为热休克蛋白9 (strss-70 protein, mitochondrial, HSPA9)。热休克蛋白9长度为2189bp, 开放阅读框长度为2067bp, 编码688个氨基酸, 编码蛋白分子量大小约为75.6kDa、理论等电点6.32;谷氨酸含量最高, 占8.4%;不稳定系数为41.96, 属于稳定不蛋白; 具有α螺旋、β折叠及无规则卷曲组成的二级结构; 三级结构分析显示其与模板2kho.1.A、1yuw.1.A和4fl9.1.A同源性超过55%, 以同源性高的模板2kho.1.A生成三级结构见图 8a。NCBI多重序列比对结果显示泥蚶热休克蛋白9氨基酸序列与哺乳类(人Homo sapiens, 小鼠Mus musculus)、斑马鱼、太平洋牡蛎和虾夷扇贝等序列相似度都在70%以上。系统进化树表明:泥蚶HSPA9与虾夷扇贝优先聚为一支, 然后与太平洋牡蛎和美洲巨牡蛎相聚, 表明它们之间亲缘关系比较近(图 8b)。

2.5 Toll相互作用蛋白

Toll相互作用蛋白(TOLLIP)长度为995bp, 开放阅读框为867bp, 编码289个氨基酸, 分子量约为32.4kDa、理论等电点4.77;含量最高的氨基酸为谷氨酰胺, 占10.03%;不稳定系数为57.98, 属于稳定不蛋白; 二级结构由α螺旋、β折叠及无规则卷曲组成, α螺旋占比最低。第68—165位置上具有C2(蛋白激酶C保守结构域2)结构域, C端第245—287具有toll相互作用蛋白CUE保守结构域(图 9b)。Swiss-model数据库中暂时没有toll相互作用蛋白模板。多重序列比对显示其与太平洋牡蛎和虾夷扇贝相似度高达60%以上, 与人类和斑马鱼等相似度超过50%。系统进化分析显示泥蚶toll相互作用蛋白与太平洋牡蛎、美洲巨牡蛎和虾夷扇贝亲缘关系比较近(图 9c)。

2.6 感染前后7种免疫相关基因表达差异

基因丰度反映了基因表达水平, 鳗弧菌刺激后, 7种免疫相关基因丰度TPM见图 10, 整合素β3以及Toll相互作用蛋白刺激前后丰度变化不大, 四种溶酶体酶和热休克蛋白9刺激后表达丰度显著增加, 以P < 0.05且|log₂Fold Change| > 1进行差异基因表达分析, 四种溶酶体酶和热休克蛋白9表达上调。

3 讨论 3.1 免疫相关基因的进化特点

本文采用转录组测序获得了泥蚶7个免疫相关基因序列, 并分析了其一级、二级和三级结构, 采用Swiss-model进行同源建模, 一般要求序列同源性大于30% (Arnold et al, 2006; Marco et al, 2014), 本研究中6个蛋白序列与Swiss-model数据库中的模板同源性均超过35%, 表明三级结构可靠。利用最大似然法(ML法)并采用步长循环1000进行进化树检验, 构建了7种免疫相关基因的系统进化树, 系统进化树显示泥蚶与虾夷扇贝Mizuhopecten yessoensis、太平洋牡蛎Crassostrea gigas和美洲巨牡蛎Crassostrea virginica亲缘关系比较近, 四者常聚为一个紧密的簇; 福寿螺Pomacea canaliculata与泥蚶亲缘关系次于三种双壳贝类; 斑马鱼Danio rerio、大黄鱼Larimichthys crocea和人Homo sapiens与泥蚶亲缘关系较远。与多数双壳贝类不同, 泥蚶血液呈红色, 以血红蛋白作为呼吸色素, 系统进化关系显示免疫相关基因与无脊

椎动物同源性更高些, 这与血细胞形态学研究结果相一致, 朱泽闻等(2011)通过研究发现泥蚶血细胞由红色颗粒细胞、嗜碱性颗粒细胞和透明细胞为主, 血细胞组成上不同于脊椎动物。

3.2 溶酶体酶

泥蚶属低等无脊椎动物, 血细胞在免疫防御中起重要作用, 血细胞直接参与异物的免疫黏附、包囊、吞噬等过程, 细胞内溶酶体酶的水解作用是杀伤和清除吞噬病原微生物主要途径, 部分溶酶体酶还能释放到血淋巴中, 通过改变病原菌表面的分子结构, 有利于吞噬细胞识别病原微生物(Cheng, 1983)。Carballal等用半定量光密度法(API ZYM试剂盒)测定了贻贝血细胞中与免疫防御功能相关的酶类, 贻贝血细胞具有α-葡萄糖苷酶活性(Carballal et al, 1997); Xue等用API ZYM试剂盒在健康欧洲扁牡蛎血细胞中发现14种酶, 没有发现α-葡萄糖苷酶活性, 原生动物Bonamia ostreae感染后血细胞具有α-葡萄糖苷酶活性(Xue et al, 2000)。通过转录组测序, 健康泥蚶血细胞中有α-葡萄糖苷酶表达, 细菌感染后泥蚶血细胞中α-葡萄糖苷酶表达量显著增加。生理状态下, 血细胞中α-葡萄糖苷酶是否具有活性, 不同贝类研究结果不一样, 表明其存在形式具有种属特异性, 欧洲扁牡蛎血细胞中α-葡萄糖苷酶可能以无活性前体形式存在, 而贻贝是以活性形式存在, 泥蚶α-葡萄糖苷酶分子特征显示其没有信号肽, 长牡蛎和斑马鱼溶酶体α-葡萄糖苷酶二级结构分析显示其没有信号肽, 而人类α-葡萄糖苷酶前体N端含有信号肽, 暗示泥蚶血细胞中溶酶体α-葡萄糖苷酶以活性形式存在。

溶酶体α-甘露糖苷酶也属于糖苷酶, 参与甘露糖苷键水解, 贻贝、欧洲扁牡蛎、太平洋牡蛎血细胞中均发现α-甘露糖苷酶, 鳗弧菌感染后泥蚶血细胞中溶酶体α-甘露糖苷酶表达量增加, 与太平洋牡蛎研究结果相似(Xue et al, 2000)。泥蚶溶酶体α-甘露糖苷酶其序列N端具有信号肽, 人类和斑马鱼溶酶体α-甘露糖苷酶前体二级结构显示其N端具有信号肽(Paciotti et al, 2017), 表明泥蚶血细胞中溶酶体α-甘露糖苷酶以前体的形式存在。

芳基硫酸酯酶的作用是水解体内的硫酸盐, 芳基硫酸酯酶B主要水解体内的糖胺聚糖(gags), 特别是水解硫酸皮肤素和硫酸软骨素, 二级结构分析显示其与人芳基硫酸酯酶前体一样, N端具有信号肽, 以前体存在。天冬酰胺氨基葡萄糖酶能水解与天冬酰胺残基相连的N-糖苷键(GlcNAc-Asn), 二级结构显示其没有信号肽, 表明其以活性形式存在溶酶体中。

鳗弧菌感染后, 泥蚶血细胞中四种酶表达量均增加, 表明其在水解细菌过程中发挥作用, 溶酶体酶水解作用是血细胞杀伤和清除细菌的重要途径, 与魁蚶血细胞超微结构研究相一致(Zhou et al, 2019)。

3.3 整合素

整合素是一种跨膜蛋白, 由α和β两个亚基通过非共价结合形成二聚体, 介导细胞与细胞以及细胞与细胞外基质相互识别和粘附(马进, 2004)。在人类研究中发现, 整合素β3有两种存在形式:与αIIb形成二聚体以及与αV形成二聚体, αIIbβ3主要表达于巨核细胞和血小板上, αVβ3存在于多种细胞上(姜泓等, 2007)。泥蚶没有巨核细胞和血小板, β3应以αVβ3形式存在于血细胞上。整合素β3的天然配体都含有一段精氨酸-甘氨酸-天冬氨酸(RGD)序列, 包括纤连蛋白、玻连蛋白和层粘连蛋白等配体, 整合素β3与配体识别后会引起细胞内信号转导过程, 使下游激酶发生磷酸化修饰, 如丝裂原活化蛋白激酶(Mitogen-activated protein kinases, MAPK)、磷酸肌醇激酶(Phosphoinositide 3-kinase, PI3K)和细胞外信号调节激酶(Extracellular signal regulated kinase, ERK)磷酸化, 从而发挥生理学功能, 如细胞增殖、细胞迁移等(Sun et al, 2016;Zoppi et al, 2018 )。转录组测序获得的泥蚶整合素β3序列具有典型的整合素β结构域, 包括胞外区域、跨膜区域和胞内区域(Xiong et al, 2001), 泥蚶整合素β3具有3个EGF样结构域, 与人类整合素3保守结构域稍有差异。鳗弧菌刺激后, 血细胞中整合素β3表达没有明显变化, 暗示泥蚶整合素β3不直接参与鳗弧菌识别与吞噬, 而Lv等(2019)研究发现长牡蛎血细胞中整合素β5参与结合和吞噬大肠埃希杆菌。

3.4 热休克蛋白和Toll相互作用蛋白

HSPA9是热休克蛋白70家族成员之一, HSPA9存在于线粒体, 是线粒体主要蛋白质之一, 与细胞多种生理活动密切相关, 如负责精确调控线粒体蛋白的输入和输出、线粒体产能、分子伴侣以及应激反应等(张国彬等, 2009)。鳗弧菌感染后, 泥蚶血细胞中HSPA9表达显著增高, 表明其参与机体应激反应, 与斑点叉尾鮰研究结果相一致(Song et al, 2016)。

Toll相互作用蛋白是一种接头蛋白, 在天然免疫中发挥重要作用, 其是Toll样受体信号通路的组成部分。未受刺激的泥蚶血细胞中, Toll相互作用蛋白表达量较高, 表明其也是一种负调控因子, Tollip能下调MyD88依赖性的TLRS信号通路转导从而抑制细胞活化和炎症反应, 使NF-κB的转录活性降低, 以避免免疫反应对机体本身的损伤(陈维英等, 2011)。受鳗弧菌刺激24h后, 泥蚶血细胞中Tollip表达量变化不显著, 鳗弧菌刺激后虾夷扇贝血细胞中Tollip表达量24h时显著上升(李若佼, 2015), 而日本鳗鲡感染嗜水气单胞菌24h后肝脏、脾脏和肾中Tollip表达量没有明显变化(Feng et al, 2019), 表明Tollip响应细菌感染的表达模式有种属差异。

4 结论

本研究通过对泥蚶血细胞进行RNA-Seq转录组测序, 对测序序列进行基因信息注释和生物信息学分析, 获得了7个免疫相关基因的蛋白序列, 包括溶酶体酶、整合素、热休克蛋白和Toll相互作用蛋白, 分析了7个基因蛋白序列理化性质、保守结构域、二级结构和三级结构, 鳗弧菌感染后4种溶酶体酶和热休克蛋白9表达量显著升高, 获得地这些基因序列有利于更深入研究泥蚶的抗病机理。