利用长链非编码RNAs与mRNAs联合分析半番鸭及番鸭胚胎期性腺发育差异
1.
2.
3.
4.
Combined analysis of embryonic gonadal development differences of mulard duck and muscovy duck using long non-coding RNAs and mRNAs
1.
2.
3.
4.
通讯作者:
收稿日期: 2022-07-07 接受日期: 2022-10-31
基金资助: |
|
Received: 2022-07-07 Accepted: 2022-10-31
作者简介 About authors
李丽(https://orcid.org/0000-0002-9953-7785),E-mail:
关键词:
Keywords:
本文引用格式
李丽, 章琳俐, 辛清武, 缪中纬, 朱志明, 邱君志, 郝晓娜, 黄勤楼, 郑嫩珠.
LI Li, ZHANG Linli, XIN Qingwu, MIAO Zhongwei, ZHU Zhiming, QIU Junzhi, HAO Xiaona, HUANG Qinlou, ZHENG Nenzhu.
半番鸭是番鸭(Cairina moschata)与家鸭(Anas platyrhynchos var. domestica)属间远缘杂交的后代,兼有双亲的优点并具有很强的杂种优势,如生长速度快、肉质细嫩和瘦肉率高等[1]。但在实际生产中,半番鸭仍存在一个显著问题,即没有实际种用价值,生产中往往采用父母本杂交繁殖,自然交配率低,公鸭利用率差,主要依赖人工授精技术,其过程烦琐且培育成本较高,极大地限制了半番鸭的生产。从性成熟期和胚胎期系统对比分析半番鸭与繁殖力正常的父母本性腺的差异基因和转录本,更易于发现半番鸭不育性状的遗传基础。长链非编码RNAs(long non-coding RNAs, lncRNAs)是一类长度大于200 nt并且一般不编码蛋白质的RNA分子,长期被当作无用的转录物而被忽视[2]。远缘杂交不育的遗传基础较复杂。目前,学者们先后从解剖形态学、染色体学说、生理生化(包括激素含量)和差异基因等方面对半番鸭的不育机制进行了分析。许玉德等[3]认为半番鸭性腺中激素分泌不足导致其不育。檀俊秩等[4]发现家鸭和番鸭的染色体核型不一致,从染色体组学角度解析了半番鸭不育的原因。随着功能基因组时代的到来,学者们运用高通量测序技术研究畜禽性腺发育的分子机制,大量基因和转录本得以挖掘。PENG等[5]发现lncRNAs通过调控FOLX2、Wnt6等靶基因参与卵母细胞减数分裂等过程,进而影响鸡卵泡发育。AYERS等[6]发现钙蛋白酶5基因(calpain 5 gene, CAPN5)、G蛋白偶联受体56基因(G protein coupled receptor 56 gene, GPR56)和成纤维细胞生长因子受体3基因(fibroblast growth factor receptor 3 gene, FGFR3)等与鸡性别分化相关。REN等[7]筛选到6个与鸭卵巢在不同阶段转化显著相关的lncRNAs。随着高通量测序技术的快速发展,尽管众多研究都强调了差异转录本、微RNAs(microRNAs, miRNAs)及基因在畜禽性腺发育方面的重要作用,但关于半番鸭不育性状的生物学意义还鲜为人知,尤其是lncRNAs在半番鸭不育方面的作用还鲜有报道。近年来,我们在半番鸭不育性状方面进行了一系列的研究,完成了半番鸭和家鸭染色体差异分析,以及半番鸭和亲本在性成熟期精巢与卵巢组织形态学差异及基因表达差异分析。结果发现,与父母本北京鸭或番鸭相比,半番鸭在性成熟期生殖系统发育不全,多表现为精巢或输卵管缺陷,或缺乏性激素分泌关键基因和在减数分裂中发挥重要作用的着丝粒蛋白基因,从而直接导致半番鸭繁殖力低下,且其成年后的表型取决于胚胎期组织器官的分化发育[8-9]。为进一步探究半番鸭不育性状分化发育差异的根源,本研究对半番鸭和番鸭雄性性腺组织进行全转录组高通量测序,探究半番鸭及其亲本在胚胎期性腺发育差异的关键lncRNAs和mRNAs,阐释半番鸭不育的分子机制,深入了解半番鸭及番鸭胚胎期性腺发育差异,为禽类性腺发育机制和生殖发生研究提供基础资料,为解决半番鸭繁育难题提供科学的参考依据。
1 材料与方法
1.1 试验材料
选取半番鸭(福建千里洋种禽有限公司)和番鸭(漳州昌龙农牧有限公司)的受精蛋各50枚进行预试验,在相同条件下孵化,番鸭孵化期为35 d,半番鸭孵化期为30 d左右。在预试验中,分别取3~22 d(7、8、9、10、11、12、13 d为重要时间节点)的鸭胚,于无菌超净台中提取其胚胎性腺,对各鸭胚的外观形态和性腺进行观察和拍照。其中,7 d时可见半番鸭性腺,而由于番鸭孵化期较长,7~9 d时其性腺不可见,10 d时可在显微镜下隐约分辨其性腺,12.5 d时可完整采集2个鸭品种的性腺。在前期观察胚胎发育预试验基础上,选择入孵后第12.5 d的胚蛋,每个品种各取10枚。在无菌超净台中对胚蛋进行消毒后,将鸭胚置于显微镜下,提取胚体并置于装有磷酸盐缓冲液(phosphate buffered saline, PBS)(武汉赛维尔生物科技有限公司)的培养皿中,借助SMZ745T体视显微镜(日本Nikon公司)、显微剪和显微镊,对各组性腺组织的外观形态进行观察,结果发现,雄性性腺为两侧发育,雌性性腺因右侧卵巢退化而左右大小不同;剥离性腺组织并迅速置于液氮中,之后转移至-80 ℃冰箱中,待用;剪下腿部并置于-20 ℃冰箱中,待用。利用聚合酶链反应(polymerase chain reaction, PCR)再次鉴定性别后,每个品种各取3个雄性性腺组织进行lncRNAs和mRNAs转录组测序。试验于福建省农业科学院畜牧兽医研究所动物房内进行,并参照国家实验动物处理行为准则,由福建省农业科学院畜牧兽医研究所动物实验伦理委员会批准(批准号:FAAS-IAHV-AEC-2020-0310)。
1.2 雌雄鉴定
参考美国国家生物技术信息中心基因和基因组信息数据库(NCBI GenBank)中鸟类CHD1基因序列,提取供试鸭胚腿部DNA用于性别鉴定。引物序列如下:CHD1-F,5
1.3 总RNA提取及cDNA文库构建和测序
以TRIzol试剂(美国Thermo Fisher Scientific公司)提取供试鸭胚性腺组织总RNA,用DNA酶Ⅰ去除残余的基因组DNA,并对RNA样品进行严格质控:1)采用琼脂糖凝胶电泳分析样品RNA完整性及是否存在DNA污染;2)利用Nanodrop微量分光光度计(美国Thermo Fisher Scientific公司)进行初步定量,检测RNA浓度及纯度;3)利用2100生物分析仪(美国Agilent公司)精确检测RNA的完整度。
采用去除核糖体RNA的方法构建链特异性文库。先用Qubit荧光定量仪(美国Thermo Fisher Scientific公司)对文库进行初步定量,并稀释至1 ng/μL;然后用2100生物分析仪对文库的插入片段长度进行检测,发现片段长度为250~300 bp,符合预期;接着,使用实时荧光定量聚合酶链反应(real-time fluorescence quantitative polymerase chain reaction, qRT-PCR)方法对文库的有效浓度进行准确定量,结果表明,文库的有效浓度>2 nmol/L,保证了样本文库的质量;最后,根据得到的文库有效浓度及数据产出需求,采用Illumina平台进行PE150测序。建库和测序由北京诺禾致源科技股份有限公司完成。
1.4 数据过滤及转录水平分析
为了保证数据分析的质量及可靠性,利用FastQC软件过滤掉原始数据中少量带有测序接头或测序质量较低的片段。利用TopHat2软件将测序得到的高质量读长(clean reads)比对到参考基因组[
1.5 差异基因功能富集分析
采用GOseq软件进行基因本体(gene ontology, GO)富集分析,分为分子功能(molecular function)、生物过程(biological process)和细胞组分(cellular component)3部分;采用KOBAS 2.0软件的京都基因和基因组数据库(Kyoto Encyclopedia of Genes and Genomes, KEGG)进行通路富集分析,从而确定候选靶基因参与的最主要的生化代谢途径和信号转导途径。
1.6 lncRNAs靶基因预测与功能分析
通过lncRNAs与蛋白编码基因的位置关系(co-location)和表达相关性(co-expression)预测lncRNAs的靶基因:根据lncRNAs与mRNAs的位置关系进行cis靶基因预测,筛选范围为100 kb以内;根据lncRNAs与mRNAs的表达相关性进行靶基因预测,筛选条件是相关系数大于0.95(当样品数<5时不进行该项分析)。然后对差异lncRNAs的靶基因进行GO/KEGG功能富集分析,预测lncRNAs的主要功能。
1.7 测序结果可靠性验证
为了验证测序结果的可靠性,对差异表达的mRNAs和lncRNAs进行验证。用Primer Premier 6.0和Beacon Designer 7.8软件进行qRT-PCR引物设计,然后由生工生物工程(上海)股份有限公司负责合成,引物序列如表1所示。在基因表达差异统计分析中,每个样品重复3次,各基因的相对表达水平以2(CT,内参基因-CT,目的基因)值进行统计分析。
表1 qRT-PCR引物序列信息
Table 1
名称 Name | GenBank登录号 GenBank accession No. | 引物序列(5 Primer sequence (5 | 扩增长度 Size/bp | 退火温度 Annealing temperature/℃ |
---|---|---|---|---|
GAPDH | XM_005016745.2 | F: GGAGCTGCCCAGAACATTATC R: GCAGGTCAGGTCCACGACA | 141 | 60 |
AKR1D1 | 101802488 | F: CGTGTCCTCCCCTCCTTTA R: CCTCGCTGGATGCTGAAAC | 105 | 60 |
17β-HSD3 | 101799136 | F: GGGACTGATCTGTTTCTGTGCT R: CCATCTCCTGCTCCTGTGA | 138 | 60 |
LOC101794158 | 101794158 | F: CCTTTGGCATCGGGATTGGA R: GCCGGATGACTCTGTTGGA | 138 | 60 |
PTGFRN | 101794936 | F: CACAATGGGAGCTGGCATAAAGT R: GGTGGGGTCATCTGAAAACTTG | 132 | 60 |
TCONS_00189165 | XLOC_112996 | F: GCTATGAAGAGGAAAAGGGATGTC R: GCGCTTAAACGAGAGGAACCA | 79 | 60 |
XR_003494977.1 | 106018179 | F: GGAACCATAATGTGCCAGCTGAA R: CCTTGCAGTTGATAAGACGATGGTA | 96 | 60 |
TCONS_00087264 | 101791443 | F: GAGCCACAGACACGAAGAGAA R: GGGGAATCAGTAACCCATAGC | 96 | 60 |
XR_003498542.1 | 113844190 | F: GAAGAGGACACCGATCAAGAG R: GTAACAACCTTCAGCCACTCAAAC | 87 | 60 |
1.8 重组载体17β -HSD3-pcDNA3.1 3×Flag的构建
1.8.1 17β-羟基类固醇脱氢酶3基因(17β-hydroxy-steroid dehydrogenase 3 gene, 17β-HSD3)的PCR扩增
以cDNA为模板,设计引物17β-HSD3-F和17β-HSD3-R并进行PCR扩增。引物序列如下:17β-HSD3-F,5
1.8.2 载体和目的片段重组(In-Fusion法)
PCR扩增产物采用QIAquick PCR & Gel Cleanup试剂盒进行切胶回收(具体操作参照试剂盒说明书),然后采用GeneArt® Seamless Cloning & Assembly Enzyme Mix进行重组质粒的无缝克隆。17β-HSD3蛋白和pcDNA3.1融合反应体系:灭菌去离子水2 μL,17β-HSD3蛋白纯化产物2 μL,pcDNA3.1(BamHⅠ/EcoRⅠ 50 ng,美国Thermo Fisher Scientific公司)1 μL,GeneArt® 2×Enzyme Mix 5 μL。室温放置20 min后,在冰上放置2~3 min,然后转化大肠埃希菌DH5α感受态细胞。
1.9 靶向关系验证
1.9.1 17β-HSD3-pcDNA3.1 3×Flag免疫共沉淀
鸭胚成纤维细胞采用α-MEM培养液[含10%胎牛血清(fetal bovine serum, FBS)和青霉素-链霉素双抗]进行复苏,在37 ℃、5% CO2条件下培养至汇合度为90%左右;加入0.25%胰酶-乙二胺四乙酸进行消化,以5×105 mL-1细胞密度接种于6 cm细胞培养皿中进行传代。按LipofectamineTM 3000说明书进行转染,待细胞裂解后,进行优质预偶联蛋白磁珠A/G和抗体复合物制备,RNA结合蛋白免疫沉淀(RNA binding protein immunoprecipitation, RIP)试验按EZ-Magna RIP试剂盒(美国Millipore公司)说明书进行。
1.9.2 反转录试验
TCONS_00246198 cDNA模板合成采用SuperScriptTM Ⅲ First Strand Synthesis SuperMix系统(美国Thermo Fisher Scientific公司)。第1链cDNA合成反应体系:分离纯化的RNA 4 μL,2×RT Reaction Mix 5 μL,RT Enzyme Mix 1 μL。反应条件:25 ℃孵育10 min,50 ℃孵育30 min,85 ℃终止反应5 min。将反转录产物保存在-20 ℃冰箱中,备用。
1.9.3 qRT-PCR扩增
采用PowerUpTM SYBRTM Green Master Mix系统(美国Thermo Fisher Scientific公司)进行qRT-PCR扩增。引物序列:TCONS_00246198-F,5
2 结果与分析
2.1 鸭胚雌雄鉴定
利用CHD1和PCR进行雌雄鉴定,结果如图1所示,呈现明显双条带的为雌性(ZW型),呈现单一条带的为雄性(ZZ型)。在本试验中,选取单一且较亮条带的半番鸭和番鸭样本各3个进行测序,其中半番鸭标记为BF1、BF2、BF3,番鸭标记为F1、F2、F3。
图1
图1
半番鸭和番鸭雌雄鉴定
bf/BF:半番鸭;b:北京鸭;f/F:番鸭;M:DNA标志物。
Fig. 1
Female and male identification of the mulard duck and the muscovy duck
bf/BF: Mulard duck; b: Pekin duck; f/F: Muscovy duck; M: DNA marker.
2.2 测序数据的质量及比对分析
本研究通过Illumina PE150对半番鸭和番鸭的性腺组织进行测序,获得13.25~15.94 Gb的数据量(原始碱基数),Q20(碱基错误识别概率小于1%)均≥97.97%,Q30(碱基错误识别概率小于0.1%)均≥94.25%(表2)。测序数据经过滤后,平均69.57%的高质量读长比对到了鸭参考基因组上,说明测序数据的质量良好,可用于后续分析。
表2 测序数据的质量及比对分析
Table 2
样品 Sample | 原始读长 Raw read | 高质量读长 Clean read | 原始碱基数 Raw base number/Gb | 高质量碱基数 Clean base number/Gb | Q20/% | Q30/% | GC含量 GC content/% | 比对读长(比对率) Mapped read (mapped rate/%) |
---|---|---|---|---|---|---|---|---|
BF1 | 88 347 436 | 87 273 712 | 13.25 | 13.09 | 98.01 | 94.25 | 49.53 | 66 096 582(75.73) |
BF2 | 106 246 148 | 104 905 270 | 15.94 | 15.74 | 97.97 | 94.25 | 49.20 | 77 624 844(74.00) |
BF3 | 88 259 882 | 87 107 522 | 13.24 | 13.07 | 98.16 | 94.74 | 48.71 | 66 689 926(76.56) |
F1 | 106 867 756 | 105 328 090 | 15.59 | 15.80 | 98.15 | 94.74 | 48.78 | 67 180 376(63.78) |
F2 | 92 794 782 | 91 434 968 | 15.27 | 13.72 | 98.29 | 94.99 | 50.08 | 57 816 779(63.23) |
F3 | 96 607 494 | 95 297 670 | 15.70 | 14.29 | 98.29 | 95.09 | 49.06 | 61 078 250(64.09) |
2.3 mRNAs及lncRNAs差异表达分析
以padj<0.05且|log2差异倍数[log2 (fold change), log2 (FC)]|>1作为显著差异基因的筛选条件,在此条件下,从半番鸭和番鸭中共筛选出1 109个差异基因。相对于番鸭,半番鸭中上调表达基因共857个,下调表达基因共252个。其中,与繁殖相关的基因有醛酮还原酶家族1成员D1基因(aldo-keto reductase family 1 member D1 gene, AKR1D1)、17β-HSD3、前列腺素F2受体负调节因子基因(prostaglandin F2 receptor negative regulator gene, PTGFRN)和胆固醇侧链裂解酶基因(cholesterol side-chain cleavage enzyme gene, P450scc)等,提示这些基因可能影响半番鸭性腺组织的分化发育。
为了分析lncRNAs是否参与胚胎期鸭的性腺发育,本研究筛选了半番鸭和番鸭雄性性腺组织中差异表达的lncRNAs。结果(图2)显示,半番鸭和番鸭性腺组织中共有733个显著差异表达的lncRNAs,相对于番鸭,半番鸭中显著上调的lncRNAs共660个,显著下调的lncRNAs共73个。
图2
图2
半番鸭中相对于番鸭的差异表达的lncRNAs
Fig. 2
lncRNAs differentially expressed in the mulard duck compared with the muscovy duck
表3 差异倍数前10位的差异基因
Table 3
GenBank 登录号 GenBank accession No. | 基因名称 Gene name | FPKM值 FPKM value | log2 (FC) | p值 p-value | padj | |
---|---|---|---|---|---|---|
半番鸭 Mulard duck | 番鸭 Muscovy duck | |||||
101792001 | THSD4 | 30.576 36 | 0 | 21.438 5 | 2.38×10-23 | 8.19×10-19 |
101805128 | MACROD2 | 23.394 99 | 0 | 20.716 1 | 5.52×10-24 | 3.80×10-19 |
101792143 | TSPAN4 | 10.613 14 | 0 | 20.461 3 | 1.06×10-21 | 2.44×10-17 |
101791830 | RPL7A | 1 032.393 00 | 0 | 20.430 9 | 1.73×10-21 | 2.98×10-17 |
101794364 | PTBP3 | 41.773 16 | 0 | 19.867 6 | 4.27×10-20 | 4.91×10-16 |
101794113 | TCF4 | 14.043 13 | 0 | 19.581 0 | 8.10×10-20 | 6.97×10-16 |
101790148 | PRKAR1B | 21.604 88 | 0 | 19.434 5 | 7.84×10-20 | 6.97×10-16 |
101 798636 | BMPR1A | 30.705 36 | 0 | 19.327 0 | 1.58×10-19 | 1.09×10-15 |
101799292 | MXRA7 | 83.668 07 | 0 | 19.209 6 | 4.22×10-19 | 2.64×10-15 |
101803865 | LOC101803865 | 190.221 80 | 0 | 19.199 0 | 1.21×10-18 | 5.13×10-15 |
101796808 | LRP2 | 0 | 37.974 58 | -20.637 3 | 1.44×10-7 | 1.92×10-5 |
101790849 | BBS9 | 0 | 7.721 36 | -19.545 7 | 2.49×10-7 | 3.06×10-5 |
101803047 | AUTS2 | 0 | 2.290 90 | -19.088 0 | 1.13×10-6 | 1.06×10-4 |
101793582 | TEAD1 | 0 | 9.605 26 | -18.878 1 | 1.53×10-6 | 1.36×10-4 |
101795217 | ACTN1 | 0 | 12.786 68 | -18.458 2 | 1.95×10-4 | 5.97×10-3 |
101798600 | SLC4A4 | 0 | 4.550 25 | -18.312 7 | 7.69×10-6 | 5.08×10-4 |
101803985 | TBC1D5 | 0 | 2.929 80 | -18.191 8 | 2.24×10-4 | 6.58×10-3 |
101790849 | BBS9 | 0 | 2.916 37 | -18.165 6 | 2.27×10-4 | 6.64×10-3 |
101791925 | GAK | 0 | 12.974 59 | -18.117 0 | 5.60×10-16 | 9.41×10-13 |
101794936 | PTGFRN | 0 | 21.933 99 | -17.955 5 | 6.66×10-15 | 8.34×10-12 |
表4 差异倍数前10位的差异表达的lncRNAs
Table 4
转录本号 Transcript ID | 转录本名称 Transcript name | FPKM值 FPKM value | log2 (FC) | p值 p-value | padj | |
---|---|---|---|---|---|---|
半番鸭 Mulard duck | 番鸭 Muscovy duck | |||||
TCONS_00036148 | LOC101796976-OT2 | 57.680 6 | 0 | 20.229 8 | 1.18×10-4 | 4.07×10-3 |
TCONS_00258829 | LINC4586 | 131.759 7 | 0 | 19.484 8 | 1.56×10-19 | 1.09×10-15 |
TCONS_00260238 | LINC4642 | 94.465 7 | 0 | 18.437 0 | 3.71×10-17 | 9.82×10-14 |
TCONS_00258290 | LOC101805096-OT1 | 5.334 5 | 0 | 18.284 4 | 3.62×10-6 | 2.73×10-4 |
TCONS_00203531 | LOC113845374-OT5 | 7.646 8 | 0 | 17.936 7 | 9.07×10-15 | 1.10×10-11 |
TCONS_00258785 | LINC4564 | 19.433 6 | 0 | 16.873 3 | 4.84×10-12 | 2.03×10-9 |
TCONS_00262674 | LOC113842089-OT1 | 14.538 8 | 0 | 16.850 8 | 3.31×10-13 | 2.17×10-10 |
TCONS_00036199 | LOC101796281-OT2 | 9.120 5 | 0 | 16.749 4 | 4.36×10-4 | 1.07×10-2 |
XR_003494977.1 | LOC106018179 | 22.879 7 | 0 | 16.701 0 | 6.73×10-13 | 3.89×10-10 |
TCONS_00258797 | LINC4579 | 17.896 5 | 0 | 16.600 1 | 1.55×10-12 | 7.73×10-10 |
TCONS_00032957 | PHF21B-OT1 | 0 | 2.490 1 | -17.482 1 | 2.30×10-14 | 2.37×10-11 |
TCONS_00087264 | ESRRG-OT2 | 0 | 0.921 4 | -15.881 6 | 1.61×10-10 | 4.31×10-8 |
TCONS_00250319 | LOC101789623-OT1 | 0 | 1.653 1 | -14.860 8 | 1.80×10-3 | 2.95×10-2 |
TCONS_00054873 | TRHR-OT1 | 0 | 1.634 7 | -14.786 3 | 8.36×10-5 | 3.16×10-3 |
TCONS_00254443 | LOC113841128-OT6 | 0 | 1.830 1 | -14.273 0 | 2.08×10-4 | 6.24×10-3 |
TCONS_00254482 | LINC4324 | 0 | 10.620 7 | -14.202 6 | 3.15×10-3 | 4.43×10-2 |
TCONS_00110367 | APELA-OT1 | 0 | 3.733 0 | -13.559 6 | 5.12×10-3 | 6.24×10-2 |
TCONS_00200982 | LINC3452 | 0 | 0.893 3 | -13.375 1 | 5.94×10-6 | 4.11×10-4 |
TCONS_00153334 | NR5A2-OT1 | 0 | 0.515 0 | -13.217 4 | 6.14×10-3 | 7.04×10-2 |
TCONS_00190897 | LINC3338 | 0 | 0.230 9 | -13.158 0 | 7.03×10-4 | 1.50×10-2 |
2.4 lncRNAs靶基因预测
lncRNAs生物学功能丰富,广泛参与生物体各类重要生理过程,可在转录水平和转录后水平上对靶基因的表达进行调控,并找出不同处理产生的差异功能通路。lncRNAs主要通过位置关系或表达相关性调控靶基因的表达。针对位置关系调控机制共发现1 587个靶基因,针对表达相关性调控机制共发现4 507个靶基因。
选取靶基因并与差异表达的mRNAs进行分析。通过靶基因预测和差异基因维恩图(图3)可知,共有136个下调的lncRNAs和893个上调的lncRNAs参与调控差异基因表达,并可能存在潜在的调控关系。
图3
图3
差异lncRNAs靶向差异mRNAs的维恩图
Fig. 3
Venn diagram of differential lncRNAs targeting differential mRNAs
在靶基因预测结果中,已知与繁殖相关的lncRNAs TCONS_00246198靶向17β-HSD3,TCONS_00229529靶向四次穿膜蛋白2基因(tetraspanin-2 gene, TSPAN2),XR_217588.4靶向TSPAN5,TCONS_00087816靶向5α-类固醇还原酶2型基因(5α-steroid reductase type 2 gene, SRD5A2),TCONS_00223199靶向类固醇激素合成急性调节蛋白基因(steroidogenic acute regulatory protein gene, STAR),XR_003497885.1、XR_001193797.2、XR_003497886.1、TCONS_00139307靶向骨形态发生蛋白受体1A基因(bone morpho-genetic protein receptor 1A gene, BMPR1A),XR_003494073.1、TCONS_00048332、TCONS_00048329、TCONS_00240705靶向辅肌动蛋白α2基因(actinin alpha 2 gene, ACTN2),XR_003494639.1、XR_003494505.1靶向AKR1D1。推测以上lncRNAs或差异基因参与了鸭胚胎期性腺发育。
2.5 差异基因和差异lncRNAs靶基因功能富集分析
GO功能注释结果(图4)显示,有43个生物过程显著富集,主要参与细胞蛋白定位、DNA代谢过程、细胞大分子定位等。在差异lncRNAs靶基因表达相关性的GO分析中,蛋白质结合、S-腺苷甲硫氨酸依赖性甲基转移酶活性、自噬、O-甲基转移酶活性、蛋白质糖基化、蛋白质修饰过程等显著富集。在差异lncRNAs靶基因位置关系的GO分析中,蛋白质ADP-核糖基化、蛋白质糖基化、糖基化、糖蛋白生物合成过程、糖蛋白代谢过程、转录因子TFIID复合物等显著富集。
图4
图4
差异基因和差异lncRNAs靶基因的GO富集过程
A.差异基因;B.差异lncRNAs(表达相关性);C.差异lncRNAs(位置关系)。
Fig. 4
GO enrichment process of differential genes and differential lncRNAs targeting genes
A. Differential genes; B. Differential lncRNAs (co-expression); C. Differential lncRNAs (co-location).
为明确半番鸭的不育性状是否存在某些特殊的过程,对差异基因和差异lncRNAs靶基因进行KEGG通路注释,结果共鉴定出13条显著通路(图5),主要包括内质网中蛋白质加工、类固醇激素生物合成、谷胱甘肽代谢、自噬调节、细胞外基质(extracellular matrix, ECM)-受体相互作用和三磷酸腺苷结合盒转运体等。在差异lncRNAs靶基因表达相关性的KEGG分析中,显著富集的有真核生物中核糖体生物发生和ECM-受体相互作用。在差异lncRNAs靶基因位置关系的KEGG分析中,鞘糖脂生物合成——乳糖和新乳糖系列显著富集。ECM-受体相互作用在差异基因和差异lncRNAs靶基因中均显著富集。
图5
图5
差异基因和差异lncRNAs靶基因的KEGG富集过程
Fig.5
KEGG enrichment process of differential genes and differential lncRNAs targeting genes
由差异基因和差异lncRNAs靶基因功能富集注释结果可知,类固醇激素生物合成通路可能与性腺发育有关,该通路上包含AKR1D1、17β-HSD3和睾丸细胞色素P450家族1亚家族B1等显著差异基因。另外,ECM-受体相互作用也可能在鸭性腺发育中起到了重要作用,该通路上包含差异基因突触囊泡糖蛋白2B基因。
2.6 测序数据验证
为了验证测序结果,分别挑选4个差异基因和4个lncRNAs进行qRT-PCR验证,其中,AKR1D1、17β-HSD3和PTGFRN可能与性腺分化发育相关,而LOC101794158、TCONS_00189165、XR_003494977.1、TCONS_00087264和XR_003498542.1为随机挑选的不同差异倍数的差异基因或lncRNAs,以便全面验证测序结果。如图6所示,在qRT-PCR检测中,差异基因或lncRNAs的表达趋势均与lncRNAs测序结果中基因的表达趋势一致,即上调的AKR1D1、TCONS_00189165和XR_003494977.1在qRT-PCR中也为上调,下调的17β-HSD3、TCONS_00087264、XR_003498542.1等在qRT-PCR中均为下调,表明测序结果较为准确,可用于后续分析。
图6
图6
qRT-PCR验证半番鸭和番鸭胚胎期性腺组织中差异表达基因和lncRNAs
Fig. 6
Differential genes and lncRNAs in embryonic gonadal tissues of the mulard duck and the muscovy duck verified by qRT-PCR
2.7 17β -HSD3-pcDNA3.1 3×Flag免疫共沉淀检测
为研究lncRNAs TCONS_00246198与17β-HSD3的靶向关系,通过转染17β-HSD3-pcDNA3.1 3×Flag重组质粒,采用Flag抗体进行RIP检测。结果(图7)表明,与IgG相比,TCONS_00246198富集水平达71.51倍,说明TCONS_00246198与17β-HSD3蛋白存在直接或间接的相互作用,即两者可能存在靶向关系。
图7
图7
17β -HSD3蛋白对lncRNAs TCONS_00246198的富集分析
Fig. 7
Enrichment analysis of lncRNAs TCONS_00246198 by 17β-HSD3 protein
3 讨论
近年来研究发现,lncRNAs可与mRNAs或miRNAs共同参与诸多性腺分化发育生物学过程,如lncPGCR靶向Btrc基因在原始生殖细胞(primordial germ cells, PGCs)发育中具有生长促进作用[10],lnc_13814通过apla-miR-145-4促进鸭颗粒细胞的凋亡[11],lncPGCAT-1通过MAPK信号转导调节原始生殖细胞的形成[12],lncRNAs可调控鸡胚胎干细胞向精原干细胞的分化过程[13]。LI等[14]通过测序探究了鸡胚胎干细胞向精原干细胞分化过程中与lncRNAs相关的调控机制。近些年,笔者对半番鸭与亲本性成熟期精巢和卵巢组织形态学差异和基因差异进行了系统研究,为了进一步探究半番鸭不育性状分化发育差异的根源,本研究以不育半番鸭和繁殖性状正常番鸭的胚蛋为研究对象,对胚胎期性腺组织进行高通量测序,筛选出可能与胚胎期性腺发育相关的差异基因和差异lncRNAs,通过功能注释推测其参与繁殖过程的代谢通路。
在雄性动物中,睾酮是主要的雄激素[15],由睾丸间质细胞分泌,可影响精巢组织生长发育并维持精子发生[16],其合成相关酶主要来自羟基类固醇脱氢酶家族(如17β-HSD3)和细胞色素P450家族(如P450scc)[17]。本研究中,TCONS_00246198的靶基因为17β-HSD3,推测两者为半番鸭性腺发育的lncRNAs和关键基因。17β-HSD3蛋白可将由胆固醇转化而成的雄烯二酮最终氧化为睾酮[18],对睾丸激素合成和雄性繁育性状至关重要[19]。在人类医学研究中,46,XY性别分化障碍患者缺乏17β-HSD3且体内雄激素水平较低[20]。RHOUMA等[21]和ENGELI等[22]研究发现,17β-HSD3错义突变导致睾酮含量过低,并使46,XY性别分化障碍。以上研究提示,17β-HSD3可能与性腺发育相关。因此,可将TCONS_00246198和靶基因17β-HSD3作为下一步研究半番鸭性腺发育的关键因子。此外,本研究发现,P450scc可能参与类固醇激素合成通路并影响半番鸭性腺发育。NURLIANI等[23]研究发现,胆固醇侧链裂解酶(cholesterol side-chain cleavage enzyme, P450scc/Cyp11a)在睾丸间质细胞中合成雄激素;MENG等[24]的研究进一步证实了P450scc影响睾酮合成并可能与性腺发育相关;LIANG等[25]通过免疫荧光分析发现,Cyp11a位于间质细胞中,催化鱼类中所有类固醇激素的合成,并参与比目鱼的性腺分化和发育。
本研究中,lncRNAs TCONS_00229529靶向TSPAN2,XR_217588.4靶向TSPAN5,即可能参与鸭性腺发育,而KASHEF等[26]研究发现,四次穿膜蛋白家族参与爪蟾胚胎协调发育过程,但未发现其参与性腺发育过程。因此,TSPAN2和TSPAN5在性腺发育方面的功能还需要进一步探究。
RENAULT等[27]认为,BMPR1A和BMPR1B在卵巢功能中起关键作用,其新序列发生变异可导致原发性卵巢功能不全。EDSON等[28]研究发现BMPR1A和BMPR1B可调节卵巢生育力。本研究中,BMPR1A是lncRNAs XR_003497885.1、XR_001193797.2、XR_003497886.1、TCONS_00139307的靶基因,并且在半番鸭睾丸组织中高表达,在番鸭睾丸组织中未表达,与以上研究结果不尽相同。因此,lncRNAs XR_003497885.1、XR_001193797.2、XR_003497886.1、TCONS_00139307靶向BMPR1A是促进还是抑制半番鸭性腺发育,还有待进一步研究。
PTGFRN是促性腺激素抑制剂。AHMED等[29]研究发现,PTGFRN为miRNAs miR-7的直接靶标,通过调节前列腺素和骨形态生成蛋白信号通路来负调节垂体激素转录和分泌,进而调节动物性成熟和生殖功能。本研究中,PTGFRN在半番鸭睾丸组织中未表达,在番鸭睾丸组织中高表达,还需要通过定量试验来进一步验证PTGFRN在半番鸭和番鸭睾丸组织中的表达情况。
综上所述,本研究筛选到可能影响胚胎期性腺发育的差异基因和差异lncRNAs,如TCONS_00246198靶向17β-HSD3,TCONS_00229529靶向TSPAN2,TCONS_00087816靶向SRD5A2。其中,17β-HSD3是类固醇激素合成途径中的关键差异基因,已有相关报道将17β-HSD3视为性腺发育的靶基因[21-22],表明该关键基因可能影响半番鸭胚胎期性腺发育过程。另外,关键lncRNAs TCONS_00246198对生殖细胞发生的作用也还需要进一步验证。未来将继续研究上述关键基因和关键lncRNAs对细胞增殖或凋亡的作用,从分子水平科学解释属间杂交后代半番鸭不育的分子遗传机制,为半番鸭不育研究提供新的途径。
4 结论
本研究以福建省特色优势杂交品种半番鸭为试验对象,依据父本番鸭正常性腺分化发育规律,获得了一批可能影响鸭胚胎期性腺发育的关键mRNAs和lncRNAs,同时,确定了关键lncRNAs TCONS_00246198和17β-HSD3的靶向关系。本研究结果有助于了解半番鸭不育分子机制、鸭胚胎期性腺发育差异及禽类性腺发育机制。
参考文献
福建省白羽半番鸭及其母本选育
[J].,
White mulard duck and female parent breeding in Fujian Province
[J]. ,DOI:10.16372/j.issn.1004-6364.2010.14.002 [本文引用: 1]
Long noncoding RNAs and exosomal lncRNAs: classification, and mechanisms in breast cancer metastasis and drug resistance
[J]. ,
半番鸭(番鸭 ×家鸭 )雄性不育的辨析
[J].,
Studies on male sterility of hybrid duck
[J]. ,
半番鸭繁殖性状的研究
[J].,
Studies on reproductive character of mule duck
[J]. ,
Genome-wide differential expression of long noncoding RNAs and mRNAs in ovarian follicles of two different chicken breeds
[J]. ,
Identification of candidate gonadal sex differentiation genes in the chicken embryo using RNA-seq
[J]. ,
Divergently expressed gene identification and interaction prediction of long noncoding RNA and mRNA involved in duck reproduction
[J]. ,
半番鸭与番鸭精巢组织差异表达转录组测序分析
[J].,
Transcriptome analysis of differential gene expression associated with testis tissue in mule duck and muscovy duck
[J]. ,DOI:10.3864/j.issn.0578-1752.2017.18.016 [本文引用: 1]
Comparative transcriptome and histomorphology analysis of testis tissues from mulard and Pekin ducks
[J]. ,
Long noncoding RNA LncPGCR mediated by TCF7L2 regulates primordial germ cell formation in chickens
[J]. ,
LncRNA lnc_13814 promotes the cells apoptosis in granulosa cells of duck by acting as apla-miR-145-4 sponge
[J]. ,
P53 and H3K4me2 activate N6-methylated LncPGCAT-1 to regulate primordial germ cell formation via MAPK signaling
[J]. ,
Analysis of lncRNA expression profile during the formation of male germ cells in chickens
[J]. ,
Regulation of crucial lncRNAs in differentiation of chicken embryonic stem cells to sperma-togonia stem cells
[J]. ,
Regulation of seminiferous tubule-associated stem Leydig cells in adult rat testes
[J]. ,
De novo testosterone production in luteinizing hormone receptor knockout mice after transplantation of Leydig stem cells
[J]. ,
Leydig cells: formation, function, and regulation
[J]. ,
Pathways and genes involved in steroid hormone metabolism in male pigs: a review and update
[J]. ,
Pig Hsdl7b3: alternative splice variants expression, insertion/deletion (indel) in promoter region and their associations with male repro-ductive traits
[J]. ,
46, XY disorder of sex development (DSD) due to 17β-hydroxy-steroid dehydrogenase type 3 deficiency
[J]. ,
Novel cases of Tunisian patients with mutations in the gene encoding 17β-hydroxysteroid dehydrogenase type 3 and a founder effect
[J]. ,
Biochemical analyses and molecular modeling explain the functional loss of 17β-hydroxysteroid dehydrogenase 3 mutant G133R in three Tunisian patients with 46, XY disorders of sex development
[J]. ,
An immunohistochemical study on testicular steroidogenesis in the Sunda porcupine (Hystrix javanica)
[J]. ,
Roles of two cyp11 genes in sex hormone biosynthesis in Japanese flounder (Paralichthys olivaceus)
[J]. ,
Characteristics of Cyp11a during gonad differentiation of the olive flounder Paralichthys olivaceus
[J]. ,
Expression of the tetraspanin family members Tspan3, Tspan4, Tspan5 and Tspan7 during Xenopus laevis embryonic develop-ment
[J]. ,
BMPR1A and BMPR1B missense mutations cause primary ovarian insufficiency
[J]. ,
Granulosa cell-expressed BMPR1A and BMPR1B have unique functions in regulating fertility but act redundantly to suppress ovarian tumor development
[J]. ,
Loss of microRNA-7a2 induces hypogonadotropic hypogonadism and infertility
[J]. ,
Identification of genome-wide SNP-SNP interactions associated with important traits in chicken
[J]. ,
CDC14B acts through FZR1 (CDH1) to prevent meiotic maturation of mouse oocytes
[J]. ,
46, XY disorder of sexual development resulting from a novel monoallelic mutation (p.Ser31Phe) in the steroid 5α‐reductase type‐2 (SRD5A2) gene
[J]. ,
Novel mutations of the SRD5A2 and AR genes in Thai patients with 46, XY disorders of sex development
[J]. ,
Spectrum of pathogenic variants in SRD5A2 in Indian children with 46, XY disorders of sex development and clinically suspected steroid 5α-reductase 2 deficiency
[J]. ,
A large cohort of disorders of sex development and their genetic characteristics: 6 novel mutations in known genes
[J]. ,
Pig StAR: mRNA expression and alternative splicing in testis and Leydig cells, and association analyses with testicular morphology traits
[J]. ,
Chemerin impairs in vitro testosterone production, sperm motility, and fertility in chicken: possible involvement of its receptor CMKLR1
[J]. ,
Molecular and morphological sex differentiation in sablefish (Anoplopoma fimbria), a marine teleost with XX/XY sex determination
[J]. ,
/
〈 | 〉 |