全转录组测序分析精子发生中RNA结合蛋白质的动态表达

李 凯，邙新雨，邹定峰，李梦真，缪时英，王琳芳，宋 伟

(中国医学科学院基础医学研究所 北京协和医学院基础学院 生物化学与分子生物学系医学分子生物学国家重点实验室，北京 100005)

精子发生是从精原细胞发育成为成熟精子的一个复杂有序的连续细胞分化过程，主要分为3个时期：精原细胞有丝分裂期、精母细胞减数分裂期和精子形成期[1]。该过程众多阶段特异性基因的表达与调控是维持精子发生正常进行的分子基础。目前，大量单细胞组学研究高精度解析了精子发生中生精细胞有序性发生的转录图谱[2-6]，然而，在转录后水平如何操控这些转录生成物的命运和功能进而影响生精细胞的增殖或分化仍未详细阐明。

RNA结合蛋白质(RNA-binding proteins，RBPs)是一类通过其功能结构域与RNA互作并操控RNA命运和功能的蛋白质[7]。RBPs在转录后水平以多种方式参与调控RNA命运，例如mRNA选择性剪接、运输、编辑和翻译等，这些方式均会引起相应的基因表达变化[8]。目前，在精子发生中已发现部分RBPs在转录后水平发挥关键的基因表达调控作用[9-10]，然而有关RBPs在精子发生全程的动态表达图谱仍缺乏完整认识。本研究整合小鼠6种类型生精细胞的全转录组测序数据，系统分析RBPs在精子发生中的动态表达全貌、阶段特异性及协同表达模式，并对其潜在功能进行预测，为阐释精子发生的分子机制及诊治男性不育相关疾病提供新的科学依据。

1 材料与方法

1.1 材料

小鼠睾丸组织6种类型生精细胞的全转录组测序数据(Bulk RNA-seq)来源于本实验室前期的研究成果[11]。该数据可从美国国立生物信息中心的基因表达综合数据库GEO下载(https://www.ncbi. nlm.nih.gov/geo/)，其登录号为GSE145130。6种类型生精细胞分别为原始A型精原细胞(primitive type A spermatogonia, priSG-A)、B型精原细胞(type B spermatogonia, SG-B)、前细线期精母细胞(preleptotene spermatocytes, plpSC)、粗线期精母细胞(pachy-tene spermatocytes, pacSC)、圆形精子(round spermatids, rST)和长形精子(elongating spermatids, elST)。

1.2 方法

1.2.1 RNA-seq数据处理：6种类型生精细胞RNA-seq的文库构建、测序数据质控、基因组比对及基因表达分析方法见参考文献[11]。

1.2.2 基因表达热图：利用R语言pheatmap工具包展示RBPs在6种类型生精细胞中的基因表达水平(Fragments Per Kilobase of transcript sequence per Millions base pairs mapped, FPKM≥1)。

1.2.3 RBPs差异表达分析及功能预测：利用R语言DESeq2工具包分析6种类型生精细胞差异表达的RBPs，筛选标准为：|log2FoldChange|≥1.5且P.adjust<0.05。利用时间序列分析软件(Short Time-series Expression Miner, STEM)分析差异表达RBPs动态表达模式，基因簇最大数目设置为50。利用ClusterProfiler工具包对6种类型生精细胞中差异表达的RBPs分别进行GO(Gene Ontology)功能富集分析(P.adjust<0.05)。

1.2.4 RBPs共表达调控网络分析：WGCNA(Weighted Gene Co-Expression Network Analysis)称为加权基因共表达网络分析，通过计算基因间表达关系鉴定表达模式相似的基因模块(Module，ME)，位于同一模块的基因共表达程度较高并且具有相似的调控作用。利用该方法分析RBPs在精子发生中的共表达网络，算法软阈值设置为22，其他为默认参数。

2 结果

2.1 RBPs在精子发生中的阶段特异性及动态表达模式

目前小鼠物种RBPs的数目预计为1 913个[7]。基于6种类型生精细胞的全转录组测序数据分析RBPs在精子发生过程的全局转录图谱，结果显示，在这6种类型生精细胞中共检测到1 835个RBPs(FPKM≥1)(图1A)。根据差异基因筛选标准进一步在精子发生中鉴定了519个阶段相对特异表达的RBPs，其在priSG-A、SG-B、plpSC、pacSC、rST和elST中的数目分别为71、102、97、134、80和35个，其中减数分裂时期(plpSC与pacSC)的RBPs比例最高(44.5%)，有丝分裂时期(priSG-A与SG-B)RBPs比例次之(33.3%)，精子形成时期(rST与elST)的RBPs比例最低(22.2%)。STEM软件分析结果显示阶段特异表达的RBPs在精子发生中主要具有7种动态表达模式(Cluster 1-7)(图1B)。

2.2 RBPs在精子发生中的潜在调控作用

为了进一步预测RBPs在精子发生中的潜在调控作用，分别对阶段特异表达的RBPs进行GO功能富集分析(P.adjust<0.01)。结果发现priSG-A特异表达RBPs主要富集在RNA剪接(RNA splicing)和mRNA加工(mRNA processing)等(图2A)；SG-B特异表达RBPs主要富集在核糖核蛋白复合物生成(ribonucleoprotein complex biogenesis)和核糖体生成(ribosome biogenesis)等(图2B)；plpSC特异表达RBPs主要富集在RNA剪接(RNA splicing)和mRNA加工(mRNA processing)等(图2C)；pacSC特异表达RBPs主要富集在核糖核蛋白复合物生成(ribonucleoprotein complex biogenesis)和核糖体生成(ribosome biogenesis)等(图2D)；rST特异表达RBPs主要富集在RNA代谢(mRNA metabolic process)和翻译(regulation of translation)等(图2E)；elST特异表达RBPs主要富集在RNA翻译(regulation of translation)和RNA代谢(mRNA metabolic process)(图2F)。

2.3 RBPs在精子发生中的共表达网络分析

RBPs在细胞发育过程中通常会形成互作复合体与RNA互作并操控RNA命运和功能，因此鉴定RBPs共表达模式将有利于进一步发现其在精子发生中的重要调控作用。利用WGCNA分析阶段特异表达RBPs在精子发生中共表达调控网络，结果共获得5个共表达基因模块(ME 1-5)，并计算各模块间的关联性(相关系数>0.9)(图3A)。根据RBPs之间的表达量进行聚类绘制得到RBPs共表达调控网络, 这5个共表达基因模块内的RBPs数量分别为285、145、219、77和76个(图3B)。将关联性得分最高的RBPs作为某一RBPs的潜在共表达对象，取交集去掉重复后共获得246个共表达RBPs。STEM软件分析结果显示这些共表达RBPs在精子发生中主要具有7种动态表达模式(Cluster 1-7)，每一种表达模式内的RBPs数目分别为19、28、45、36、21、30和26个，其中减数分裂时期(Cluster 3-5)共表达RBPs比例最高(49.8%)，精子形成时期(Cluster 6-7)共表达RBPs比例次之(27.3%)，有丝分裂时期(Cluster 1-2)共表达RBPs比例最低(22.9%)(图3C)。

A.heat map showing global transcriptional profile of RBPs in spermatogenesis; B.STEM showing the stage specificity and dynamic expression pattern of RBPs in spermatogenesis

A-F.enriched GO terms of stage-specific RBPs in priSG-A, SG-B, plpSC, pacSC, rST and elST图2 阶段特异表达的RBPs在精子发生中的GO功能富集分析Fig 2 GO enrichment analysis of stage-specific RBPs in spermatogenesis

A.heatmap showing the correspondence between co-expression modules; B.heatmap showing WGCNA analysis of RBPs; C.heatmap showing the dynamic expression pattern of RBPs co-expression in spermatogenesis

2.4 共表达RBPs在精子发生中的潜在调控作用

利用GO功能富集分析预测共表达RBPs在精子发生中的潜在调控作用(P.adjust<0.05)。priSG-A共表达RBPs数量较少，并未富集到明显的GO功能条目。SG-B共表达RBPs主要富集在核糖核蛋白复合物生成(ribonucleoprotein complex biogenesis)和核糖体RNA代谢(rRNA metabolic process)等(图4A)；plpSC共表达RBPs主要富集在RNA剪接(RNA splicing)和mRNA加工(mRNA processing)等(图4B)；pacSC共表达RBPs主要富集在蛋白质定位(protein localization)和核糖核蛋白复合物生成(ribonucleoprotein complex biogenesis)等(图4C)；rST和elST共表达RBPs表达模式较为相似，作为一个基因集进行GO富集分析，主要富集在核糖体生成(ribosome biogenesis)和翻译调控(regulation of translation)等(图4D)。

3 讨论

根据分子生物学中心法则，以往被视为遗传信息传递中间站的信使RNA (mRNA)其实具有远超人们所理解的转录后调控方式，例如mRNA选择性剪接、编辑、运输和翻译等，而这些调控方式主要是由RBPs介导完成，最终引起相应的基因表达变化。本研究系统描绘了RBPs在精子发生中的动态表达全貌、阶段特异性及协同表达模式，并对其潜在功能进行了预测。

基因表达的转录后调控对于维持精子发生的正常进行至关重要[12]。精子发生早期阶段基因转录异常活跃，大量mRNA被转录生成后，可与RBPs相互作用形成mRNP复合物进行存储，在精子发生后期阶段由于染色质高度压缩，基因转录活性逐渐降低。然而，为了维持生精细胞的正常发育，早期转录并存储的mRNA在此时开始进行翻译，该现象称为“转录-翻译”解偶联[13]。目前，在精子发生中已发现部分RBPs在转录后水平发挥关键的基因调控作用[14-16]。中科院生化与细胞研究所刘默芳教授团队发现MIWI/piRNA通过与翻译起始因子eIF3f及AU-rich元件结合蛋白HuR等相互作用，激活生精细胞mRNA的翻译并调控生精细胞发育[9]。本研究发现RBPs在精子发生过程呈现阶段特异性表达，并且在早期阶段(例如priSG-A、SG-B及plpSC)主要参与调控mRNA加工、选择性剪接或稳定等过程，而在后期阶段(例如rST及elST)主要参与调控核糖体组装或mRNA翻译等过程，提示这些阶段特异表达的RBPs在生精细胞发育过程发挥了重要调控作用，并为精子发生中的“转录-翻译”解偶联提供了潜在分子基础。此外，RBPs在细胞内经常以共表达形式存在并能驱动靶标mRNA的协同表达，进而调控细胞分化或组织发育。本研究发现RBPs在精子发生过程具有明显的共表达特征，并主要参与mRNA代谢、mRNP组装或翻译等过程，提示精子发生中共表达的RBPs是操控RNA命运和功能的潜在重要参与者。然而，这些共表达的RBPs需要在细胞和动物水平进行分子生物学验证，深入解析RBPs与RNA的相互作用机制、生化特征及生理功能可为理解、诊断和治疗男性不育相关疾病提供新的线索。

A-D.enriched GO terms of RBPs co-expression in SG-B, plpSC, pacSC, rST and elST图4 共表达RBPs在精子发生中的GO功能富集分析Fig 4 GO enrichment analysis of RBPs co-expression in spermatogenesis