朱淑斌, 徐 盼, 周春宝, 许琴瑟, 吴嘉韵
(1.江苏农牧科技职业学院,江苏 泰州 225300; 2.姜曲海猪保种场,江苏 泰州 225300)
母猪高产和健康的保持一直是养猪业面临的最大挑战之一。母猪繁殖障碍主要表现为不发情或发情异常、不孕或受孕率低、早产、流产以及产胎率低等临床症状。产仔数、产活仔数和窝质量是评估母猪繁殖性能的重要指标,母猪的繁殖能力低会直接影响生猪生产效益[1-2]。母猪的繁殖器官结构复杂,任何一个组成部分出现异常都会导致繁殖障碍,其中子宫和卵巢的发育情况备受重视,它们在母猪繁殖过程中具有重要的作用[3-4]。因此,围绕母猪子宫和卵巢发育过程中的转录变化进行深入研究,挖掘影响母猪生殖器官发育的重要功能基因,可以进一步通过分子选育对功能基因进行筛选和鉴定,有望减少母猪繁殖障碍现象的发生,进而提高母猪的繁殖性能。
RNA-Seq是目前分子生物学研究领域使用最为广泛的组学技术之一,具有通量高、效率高和成本低的优点。RNA-Seq可以全面提供样本基因转录本的表达信息,让研究者快速了解样本转录差异的整体规律。近十年来,RNA-Seq技术的快速发展,使其在畜牧研究中得到广泛应用。段修军等[5]通过对黑羽番鸭胸肌、腿肌中挥发性风味物质进行测定,并利用RNA-Seq技术筛选得到可能通过形成肌内挥发性风味物质进而影响肌肉风味的候选基因。安清明等[6]利用RNA-Seq技术筛选贵州白山羊肌肉生长发育的候选基因。苏比奴尔·艾力等[7]基于RNA-Seq技术对塔里木马鹿毛色相关基因进行筛选及分析,筛得7个候选基因可能在塔里木马鹿毛色形成过程中发挥重要作用。梁鹏等[8]利用RNA-Seq技术分析滩羊、杜泊羊和小尾寒羊肉品质变化和肌肉内基因的表达情况,挖掘它们肉质性状差异的相关调控基因,为其他绵羊品种选育与改良提供理论依据。此外,在猪养殖中,目前也有许多研究人员应用RNA-Seq技术挖掘猪种优良性状的功能基因,相关研究结果对增加猪产肉量、增强猪抗病力和提高猪群繁殖力具有重要的参考价值[9]。Fraser等[10]利用RNA-Seq技术比较不同冷冻公猪精子的转录谱,为进一步阐明精子与其低温存活相关基因的研究提供了基础,鉴定到的关键候选基因具有作为预测猪精液冷冻性标记的潜力。徐盼等[11]基于RNA-Seq技术鉴定90日龄与180日龄苏姜猪背最长肌之间以及180日龄苏姜猪背最长肌与腿肌之间的差异表达基因,为后续苏姜猪肌肉生长发育的分子改良提供了一定的研究基础。
姜曲海猪的主要产区位于泰州市姜堰区、南通市曲塘镇、南通市海安市一带,此类猪种全身被毛黑色,其中有部分猪在鼻吻处偶有白斑,外号“花鼻子”,其具有早熟易肥、繁殖性能强的优势[12-13]。母猪1月龄生殖系统发育还未完全,未达到性成熟,8月龄姜曲海母猪已达到性成熟,具备完整的繁育能力。因此,本研究拟利用RNA-Seq技术对1月龄和8月龄的姜曲海猪子宫和卵巢组织进行研究,筛选和分析差异表达基因,挖掘发育性状相关的有效基因,以期为研究姜曲海猪繁育的分子机理水平提供基础数据和理论支撑。
本次试验动物材料来源于江苏姜曲海猪猪种有限公司。随机选取6头同窝体质量相近[体质量(6.87±0.21) kg]的健康姜曲海母猪进行标准化饲养,在1月龄和8月龄时期各屠宰3头,于屠宰后30 min内采集子宫(子宫角部位)和卵巢组织样品置于冻存管,分别命名为:1月龄-子宫-编号(1M-ute-1、1M-ute-2、1M-ute-3)、8月龄-子宫-编号(8M-ute-1、8M-ute-2、8M-ute-3)、1月龄-卵巢-编号(1M-ova-1、1M-ova-2、1M-ova-3)、8月龄-卵巢-编号(8M-ova-1、8M-ova-2、8M-ova-3),放入液氮罐带回实验室,再转移到-80 ℃冰箱内长期保存,用于总RNA的提取。
TRIzol(货号:15596018)、Oligo磁珠(货号:B-3002)购自广东润科生物工程股份有限公司;琼脂糖(货号:C042)购自上海如吉生物科技有限公司;反转录试剂盒(货号:B300537)、实时荧光定量PCR试剂盒(货号:B639271)购自生工生物工程(上海)股份有限公司;分光光度计(型号:K5500)购自北京凯奥科技发展有限公司;细胞分析仪(型号:2100Bio-analyzer)购自Agilent公司;定量PCR仪(型号:CG-02)购自HealForce公司;测序仪(型号:HiseqTM2500)购自Illumina公司。
利用TRIzol法提取各子宫和卵巢组织样本的总RNA。总RNA质量经检测后分离纯化mRNA,合成cDNA并加测序接头后回收目的片段,进行PCR扩增。构建好的文库利用PE150测序策略进行测序,以此得到原始数据raw read[7]。
通过软件FastQC(v0.11.2)对测序的原始数据进行质量评估,通过软件Trimmomatic(v0.36)进行质量剪切,去掉测序接头及低质量序列后,得到相对准确的高质量的有效数据(clean read)。利用HiSAT2(v2.1.0)软件将样本有效数据与猪的参考基因组[Scrofa10.2 (fp://ftp.ensembl.org pub/release 87/fasta/sus_ scrofa/dna/)]进行比对,统计Mapping信息。将唯一比对到基因组的基因用于差异基因表达分析,使用每千个碱基转录每百万映射读取的碎片(TPM)法对每个样本中基因的表达量进行标准化。以|log2FoldChange|>2(FoldChange表示表达量差异倍数)和P<0.05为标准,统计筛选差异表达基因。使用DESeq2(v1.26.0)软件进行基因差异表达分析。利用GO功能富集分析和KEGG信号通路富集分析将差异表达基因注释到相应的GO和KEGG信号通路上,进一步将差异表达基因与猪QTL数据库进行比对,以鉴定姜曲海猪子宫和卵巢组织中与发育性状有关的候选基因。
随机挑选RNA-Seq筛选的10个差异表达基因,使用实时荧光定量PCR试剂盒进行实时荧光定量PCR(Quantitative real-time PCR, qRT-PCR)检测。20.0 μl qRT-PCR反应体系为:10.0 μl实时荧光定量试剂盒预混液,0.4 μl上游引物(10 μmol/L),0.4 μl下游引物(10 μmol/L),2.0 μl cDNA,7.2 μl双蒸水(RNase free ddH2O)。扩增程序为:95 ℃ 3 min;95 ℃ 5 s,60 ℃ 30 s,45个循环。利用Primer 6.0软件,跨外显子设计所筛选基因的qRT-PCR引物,以GAPDH基因作为内参基因。所有引物均由上海生工生物工程技术服务有限公司合成,引物详细信息见表1。每组qRT-PCR检测3个重复,利用SPSS 20.0软件进行t检验[14],本研究对qRT-PCR和RNA-esq结果进行了Pearson相关性分析。
表1 qRT-PCR引物信息
通过对RNA-Seq测序数据进行质控分析发现,姜曲海猪的子宫和卵巢组织样本平均read数如表2所示,clean read的数量和比对结果符合要求。同时,子宫和卵巢的3D主成分分析(PCA)结果、基因表达量密度以及组间的相关性分析结果表明,组间样本差异较小(图1A~图1C)。上述结果均表明测序质量较高,为下一步数据分析的可靠性奠定了基础。此外,由图1D可知,1月龄和8月龄的子宫和卵巢组织样本共有的基因表达数为19 435个,后续生物信息学分析将基于该共有表达基因开展。
表2 子宫和卵巢组织样本比对参考基因组结果
A:子宫和卵巢的3D主成分分析(PCA)结果;B:基因表达量密度曲线图,TPM:每千个碱基转录每百万映射读取的碎片;C:组间的相关性分析,数值表示组间相关性;D:组间共表达基因韦恩图。重叠区数字代表不同组共有的基因表达数,非重叠区数字代表特异的基因表达数。1M-ute-1、1M-ute-2、1M-ute-3、8M-ute-1、8M-ute-2、8M-ute-3、1M-ova-1、1M-ova-2、1M-ova-3、8M-ova-1、8M-ova-2、8M-ova-3见表2注。
1M-ute-1、1M-ute-2、1M-ute-3分别表示1月龄子宫组织样本的3组重复;8M-ute-1、8M-ute-2、8M-ute-3分别表示8月龄子宫组织样本的3组重复;1M-ova-1、1M-ova-2、1M-ova-3分别表示1月龄卵巢组织样本的3组重复;8M-ova-1、8M-ova-2、8M-ova-3分别表示8月龄卵巢组织样本的3组重复。
本研究设置阈值为|log2FoldChange|>2、P<0.05的标准来筛选差异表达基因。结果(图2A)显示,1月龄姜曲海猪子宫和卵巢组织的基因表达水平相近;8月龄姜曲海猪子宫组织的基因表达水平低于卵巢组织;姜曲海猪子宫和卵巢组织的基因表达水平都表现为8月龄组高于1月龄组。图2B显示,1月龄和8月龄姜曲海猪卵巢组织上调差异表达基因个数明显多于下调差异表达基因个数。由图3可知,1月龄和8月龄姜曲海猪子宫组织中有1 688个差异表达基因,其中表达上调的基因有844个,表达下调的基因有844个(图3A);1月龄和8月龄姜曲海猪卵巢组织中存在3 833个差异表达基因,有2 831个表达基因上调,1 002个表达基因下调(图3B)。此外,1月龄姜曲海猪子宫和卵巢组织中的差异表达基因有978个基因表达上调,765个基因表达下调(图3C);8月龄姜曲海猪子宫和卵巢组织中的差异表达基因有1 132个基因表达上调,2 919个基因表达下调(图3D)。
A:比较组表达量箱型图,纵轴为每千个碱基的转录每百万映射的读数(TPM)的对数;B:表达差异分析统计柱状图,纵轴为上调和下调的差异基因数量。1M-ute-1、1M-ute-2、1M-ute-3、8M-ute-1、8M-ute-2、8M-ute-3、1M-ova-1、1M-ova-2、1M-ova-3、8M-ova-1、8M-ova-2、8M-ova-3见表2注。
A:8月龄和1月龄子宫组织样本间比较的差异表达基因数量;B:8月龄和1月龄卵巢组织样本间比较的差异表达基因数量;C:1月龄子宫组织样本和卵巢组织样本间比较的差异表达基因数量;D:8月龄子宫组织样本和卵巢组织样本间比较的差异表达基因数量。图中每个点代表一个基因,不同颜色代表不同基因类型,红色、绿色和黑色分别表示上调基因、下调基因和无差异表达基因。1M-ute-1、1M-ute-2、1M-ute-3、8M-ute-1、8M-ute-2、8M-ute-3、1M-ova-1、1M-ova-2、1M-ova-3、8M-ova-1、8M-ova-2、8M-ova-3见表2注。
GO功能富集分析结果显示,1月龄和8月龄卵巢组织的差异表达基因显著富集于多细胞生物的发育、系统发育和细胞发育等过程(图4A);1月龄和8月龄子宫组织的差异表达基因显著富集于动物器官发育、组织发育和细胞分化等过程(图4B);1月龄子宫和卵巢组织的差异表达基因极显著富集于动物器官发育、组织发育、多细胞生物的发育等生物学过程(图4C);8月龄子宫和卵巢组织的差异表达基因显著富集于动物器官发育、组织发育、繁殖和生殖等过程(图4D)。
通过KEGG信号通路数据库对不同组的差异表达基因进行通路注释发现,1月龄和8月龄卵巢组织的差异表达基因显著富集于碱基切除修复和DNA复制(图5A);1月龄和8月龄姜曲海猪子宫组差异表达基因无显著通路富集,主要涉及的通路为Wnt信号通路和细胞色素P450对异种生物的代谢作用(图5B);1月龄姜曲海猪子宫和卵巢组织的差异表达基因显著富集于紧密连接信号通路(图5C);8月龄姜曲海猪子宫和卵巢组织的差异表达基因无显著通路富集,主要涉及的通路为调节干细胞多能性的信号通路和轴突导向等(图5D)。
图中点的大小表示差异表达基因的数量,点越大表示包含差异表达基因越多。1M-ute-1、1M-ute-2、1M-ute-3、8M-ute-1、8M-ute-2、8M-ute-3、1M-ova-1、1M-ova-2、1M-ova-3、8M-ova-1、8M-ova-2、8M-ova-3见表2注。
图中点的大小表示差异表达基因的数量,点越大表示包含差异表达基因越多。1M-ute-1、1M-ute-2、1M-ute-3、8M-ute-1、8M-ute-2、8M-ute-3、1M-ova-1、1M-ova-2、1M-ova-3、8M-ova-1、8M-ova-2、8M-ova-3见表2注。
结合GO功能和KEGG信号通路显著性富集分析结果,将差异表达基因与猪数量性状基因座(QTL)数据库(https://www.animalgenome.org/cgi-bin/QTLdb/SS/index)进行比对,一共筛选出11个子宫中与发育性状有关的候选基因(ADCY2、DIO2、MAGEL2、LHFPL1、CDH13、ITIH3、MTUS2、LIF、RBP4、EPHB2、TESC)以及31个卵巢中与发育性状有关的候选基因(VRTN、FRAS1、LHFPL1、SDCCAG8、ITIH4、CPNE5、ARL4C、MTUS2、FGF14、SMTNL2、LIF、KIF5C、LEF1、PIWIL1、ZFYVE9、OSBPL10、RBP4、EPHB2、SPEF2、ESR2、ELOVL2、DCC、AREG、RIMKLA、WNT10B、HECW1、FUT2、TESC、WDHD1、DLGAP5、WNT3)。其中,2个组织中共同筛选得到的发育性状差异表达基因有6个,分别为LHFPL1(LHFPL tetraspan subfamily member 1)、MTUS2(Microtubule associated scaffold protein 2)、LIF(Leukemia inhibitory factor)、RBP4(Retinol binding protein 4)、EPHB2(Erythropoietin-producing hepatocellular B2)和TESC(Tescalcin)。此外,我们将这6个共有的发育性状差异表达基因在不同组间的表达绘制热图,发现除LIF和RBP4外的其他基因随着母猪月龄的增大其表达水平在2个组织中均呈上升趋势(图6)。
A:卵巢和子宫共有的发育性差异表达基因Venn图;B:共有的发育性差异表达基因在不同组间的表达热图。a:1M-ova-1;b:1M-ova-2;c:1M-ova-3;d:1M-ute-1;e:1M-ute-2;f:1M-ute-3;g:8M-ova-1;h:8M-ova-2;i:8M-ova-3;j:8M-ute-1;k:8M-ute-2;l:8M-ute-3。1M-ute-1、1M-ute-2、1M-ute-3、8M-ute-1、8M-ute-2、8M-ute-3、1M-ova-1、1M-ova-2、1M-ova-3、8M-ova-1、8M-ova-2、8M-ova-3见表2注。
为进一步验证差异表达基因的表达模式,随机选取10个基因(EGR1、ACKR3、DUSP1、CCN2、SAT1、ASPN、ADGRA2、CCDC3、SFRP2、BGN)进行qRT-PCR验证,Pearson相关性分析结果表明,qRT-PCR和RNA-Seq相关系数为0.811(P<0.01),这些基因的表达模式在RNA-Seq和qRT-PCR之间保持一致(图7),表明转录组数据分析的可靠性和准确性高。
a:EGR1;b:ACKR3;c:DUSP1;d:CCN2;e:SAT1;f:ASPN;g:ADGRA2;h:CCDC3;i:SFRP2;j:BGN。
RNA-Seq能够检测物种特定组织或器官在某一状态下的转录本序列信息,从而通过生物学分析手段筛选和鉴定该物种优良性状的关键候选基因[15-19]。因此,目前RNA-seq已经广泛应用于畜禽经济性状的研究当中,并且许多研究结果也进一步表明RNA-Seq技术在畜禽重要经济性状选育中的重要作用,基于RNA-Seq可以高效筛选特定性状的关键候选基因[5-11]。因此,本研究利用RNA-Seq技术鉴定1月龄和8月龄姜曲海猪子宫和卵巢组织差异基因的表达谱,结果发现姜曲海猪卵巢组织不同月龄的差异表达基因多于子宫,表明在发育过程中卵巢组织发生了更多的转录变化,并且2个组织间的差异表达基因比较结果也表明不同组织的发育存在一定的差异性。
本研究对不同月龄不同组织中的差异表达基因进行富集分析,发现这些差异表达基因显著富集于多细胞生物的发育生物学过程,相较于1月龄时期,8月龄的子宫和卵巢组织的差异表达基因在繁殖和生殖过程中显著富集。1月龄和8月龄卵巢组的差异表达基因显著富集于碱基切除修复和DNA复制信号通路。值得注意的是,在子宫和卵巢的发育过程中,差异表达基因显著富集于Wnt信号通路,Wnt信号通路在女性生殖系统的发育中起着重要的作用, 是其分化形成的关键因子, 在胚胎干细胞定向分化为子宫内膜样细胞的过程中,Wnt信号通路同样起着关键的调控作用[20]。由此推测,在母猪子宫和卵巢的发育过程中Wnt信号通路发挥着调控细胞分化的关键作用。根据上述结果,本研究推测在姜曲海母猪子宫和卵巢的发育过程中以上信号通路发挥了潜在的调控作用,Wnt、碱基切除修复和DNA复制等信号通路及通路中的差异表达基因可能对姜曲海猪生殖器官的发育具有一定的调控作用,进而影响姜曲海母猪的繁殖性能。
为了进一步筛选与姜曲海猪子宫和卵巢发育相关的关键候选基因,本研究将差异表达基因与猪QTL数据库进行比对后鉴定到6个潜在的重要候选基因。其中,LHFPL1基因位于X染色体上,是脂肪瘤HMGIC融合伙伴(LHFP)基因家族的一个成员,LHFPL1在人的大部分组织中广泛表达,特别是在肺、胸腺、骨骼肌、结肠和卵巢中表达较高[21],其在卵巢中的高表达可能与其参与卵巢的发育相关。MTUS2是一种蛋白质编码基因,研究发现MTUS2与肺非鳞状非小细胞癌相关[22]。LIF基因位于第14号染色体上,参与造血、神经发育和内分泌调控等多种生命活动,在母体和胎儿界面的免疫耐受性方面也可能发挥作用[23-25]。RBP4基因位于第14号染色体上,属于脂蛋白家族,可将视黄醇从肝脏储存到外围组织[26]。另有研究结果表明,RBP4可能通过抑制IRS/PI3K/Akt的磷酸化参与高尿酸血症诱导的胰岛素抵抗[27]。EPHB2属于受体酪氨酸激酶跨膜糖蛋白的Eph受体家族,其在结直肠癌、胃癌、肝细胞癌等多种癌症中均被证实异常表达,导致肿瘤发生和发展[28]。TESC亦称为钙调磷酸酶B同源蛋白3(CHP3),其在细胞生长和分化中起重要作用[29]。子宫和卵巢的发育过程中细胞会不断发生分化,而本研究筛选到的候选基因中有很多参与细胞生长和分化相关的生命活动,因此,本研究初步推测这些候选基因在姜曲海猪的子宫和卵巢发育过程中发挥重要作用,今后有必要对这些候选基因进行进一步的功能验证和探究。
本研究基于姜曲海这一地方猪品种繁殖性能强的优势,利用RNA-Seq技术对影响姜曲海母猪子宫和卵巢组织发育相关的功能基因进行了深入的挖掘,揭示了姜曲海母猪发育过程中子宫和卵巢组织的基因表达图谱,并成功鉴定到LHFPL1、MTUS2、LIF、RBP4、EPHB2和TESC等6个潜在的发育候选基因,有望为姜曲海猪繁殖性状的分子改良提供良好的试验依据和理论基础,同时也为姜曲海猪子宫和卵巢发育分子层面的研究提供了新的研究方向。