文|爱叨叨的小胖

编辑 | 爱叨叨的小胖

背景

性别二型性的分子机制仅仅只有一部分有具体的阐明。基因表达分析已经通过转录组差异分析来回答相关的问题。虽然这些研究中有许多都是根据来自牛群、啮齿动物和果蝇的数据，但其中进行描述牛的性别二型性的研究还是只有少数。

此外，大多数牛转录组学旨在识别与组织生长和发育相关的性别二型性机制的研究，只在植入前胚胎上进行；对其他组织的探索性分析可能有助于进一步了解这些机制。

在本研究中，我们将鉴定有助于牛的性别二型性的性别二型基因，采用两种统计方法分析，从多个不同组织收集的样本RNA-seq数据在多个组织中鉴定性二型基因的优缺点，以及确定牛性二型性的不同机制。目前我们将展开本次研究结果的论述。

结果：RNA-seq分析管道的描述

为了提取RNA-seq基因表达数据，我们在40个样本中的平均存活读取率为98.75%，这些序列使用Bowtade2，与牛参考基因组对齐平均定位率为81.91%。

我们使用UCSC基因组浏览器，通用转移格式（GTF）文件计算基因表达，并在HTSeq GTF包上绘制映射结果。共使用13,148个基因进行DEG分析。

为了鉴定多个组织中的两性二型基因，我们使用M1对每个组织进行了两组测试，结果显示24、14、86和57个基因，分别在肝脏、脂肪、肌肉和垂体组织中，存在显著差异表达（FDR调整后的P值<0.05）。

在所有四种组织中，均有3个常见的基因：DDX3Y、USP9Y和ZFY，显著检测到的deg具有较高的组织特异性；只有少数基因，同时在多个组织中发现显著。在肝组织显著检测到的24个DEGs中，有21个基因仅在肝组织中被发现，包括CUX2（FDR调整P值：4.41E-04）、CYP7A1（4.41E-04）、AK4（1.05E-03）、COL27A1（1.18E-03）、TNC（1.32E-03）。

这种检测到的deg的强组织特异性模式也在其他组织中观察到，在脂肪组织里14个总DEGs中，有8个被鉴定为脂肪组织特异性的性别二型基因包括IGFBP1（2.02E-04）、TECTB（6.77E-03）和ACR（1.24E-02）。

在肌肉组织中，检测到最多的性别差异基因，86个基因被显著检测，只有7个deg在其他组织中被发现。MYH1（7.51E-16）、MMP12（1.02E-07）、MCHR1（4.15E-05）和SH3KBP1（2.89E-04）在肌肉组织中极为显著，仅限于肌肉组织。

最后在垂体组织检测到的57个显著的DEGs中，有51个具有组织特异性，包括GRP（1.79E-14）、LOC781146（2.70E-09）和LYSB（9.12E-08），为了进行更详细的组织调查特异性，我们进行了层次聚类，以揭示使用所有基因的组织之间的关系。

图1显示了雄性和雌性聚集的样本结果显示性别间的表达差异，小于不同组织间的差异。在此图中观察到脂肪组织样本的组内类似性低于其他样本。

不过这种现象不仅可以在必定的FDR截止处观察到；还通过FDR调整的P值截止值0.05至0.2来计算组织特异性。

图1检测到的性别相关基因的组织特异性

使用基于统计模型的方法同时思考整个组织的影响，对性别二型基因的鉴定

为了同时思考组织，对统计模型的影响并识别性别二型基因，在edgeR中实现的ANODEV允许在数据结构更复杂的情况下识别deg，当调整组织对模型的影响时M2用于检测性别二型基因。

结果，采用似然比检验（LRT）确定了655个显著基因（FDR调整P值<0.05），包括USP9Y（5.47E-72）、DDX3Y（6.56E-55）、ZFY（8.79E-55）、XIST基因（1.37E-17）、KDM6A（9.01E-05）和常染色体基因如LOC780876（3.71E-12）、LYSB（5.71E-06）等。

我们使用维恩图比较了显著基因的数量，改变了fdr调整后的P值截断值，如0.05、0.1和0.2。M2比M1鉴定出更多的性二型基因，使用M1得到的显著检测到的deg中约有一半与集成模型（M2）重叠；重叠基因的比例为40~92%。

不过有大量的非重叠基因只有在使用M1或M2时，才被显著检测到，我们研究了导致M1和M2之间有明显区别的潜在因素。

在FDR调整后的P值截断值为0.05的情况下，使用M1和M2显著检测到752个基因，在这些基因中，有三个具有代表性的基因在M1和M2之间显示出明显不同的模式，它们被可视化为盒图和线线图。

从这些图中可以观察到，从M2中显著检测到的deg与从M1中显著鉴定出的基因相比，在组织中显示出相对类似的斜率，AK4在M1中显著检测到，而在M2中没有显著检测到，同样ANXA9和STRA6基因在肌肉和垂体腺的M1中分别具有显著性。

M2报告的FDR调整后的肌肉腺，和垂体腺的P值分别为0.782和0.891，相比之下，AGPHD1、KDM6A和SRSF2基因仅使用M2具有显著性（FDR调整后的P值分别为0.001、7.99E-05和0.0007），M1并没有在任何组织中显著检测到任何这些基因。

利用M2中显著的两性二态性基因，进行染色体富集分析

众所周知，性二态基因不仅在染色体上广泛存在，而且在不同物种的性染色体中，表现出高富集，以M2显著检测到的655个DEGs为基础进行染色体富集分析，以检测染色体富集的程度。调查程度利用M2鉴定的deg，对牛染色体的富集采用Fisher准确检验（表1）。

表1 染色体富集试验结果（Fisher准确检验）P值<为0.05时的显著基因

正如之前的研究结果所预期的那样，两性二型基因在牛的X和Y染色体上高度富集。这一结果表明，有性二态基因不仅包括性染色体还包括常染色体，而曼哈顿图的创建显示，M2比M1检测到更多的deg。

基于M2中检测到的655个基因，利用染色体信息计算其性别偏倚性（图2），除在Y染色体中发现的基因外，所检测到的大量基因，均为雌性偏倚基因。

可以理解的是，在性染色体上，检测到的两性二型基因表现出特异性的性别偏倚表达，无一例外。

图2区分显著的性别二态性基因（FDR<0.05）与性别偏差、组织和染色体位置等多个类别

在常染色体和未知染色体中，虽然观察到少数雄性偏倚基因，但大多数检测到的基由于雌性偏倚基因，相比之下X和Y连锁基因，分别只在雌性和雄性的偏倚类别中被发现。

在这个图中，可以观察到这种模式与使用M1的肌肉组织结果超级类似，为了证实这一结果，我们将显著观察到的655个检测基因的基于logFC的曼哈顿图（图3）。

图3 M2的性二态基因模式总结（FDR<0.05）

利用qRT-PCR，检测两个不同物种的两性二型基因的技术验证

采用qRT-PCR进行技术验证，1.通过比较牛的RNA-seq和qRT-PCR表达基因，技术验证检测到的性别二态性基因；2，比较牛和大鼠，以调查已验证的性别二态性基因，在这两种不同哺乳动物之间是否共享。

利用M1和M2两种统计方法，我们鉴定了752个显著的（FDR<.05）性别二态基因，在这些检测到的DEGs里从牛和大鼠物种的已知引物序列中，随机选择40个基因，并进行qRT-PCR，不仅比较牛的RNA-seq和qRT-PCR，还比较牛和大鼠的qRT-PCR（图4）。

图4 鉴定了752个显著的（FDR<.05）性别二态基因

为了检测三个实验中的整体基因表达，我们使用分位数归一化值可视化了一个比较热图（图5）。

第一我们观察到qRT-PCR和RNA-seq之间的大多数基因表达模式是类似的，为了定量评估这一观察结果，我们使用三个实验中，雄雌之间的对数倍变化比值，计算了皮尔逊相关系数，如图5-(d)，所示结果表明，牛RNA-seq和qRT-PCR之间的相关系数较高（0.43~0.87）。

图5 从牛RNA-seq分析中随机选择40个性二态基因绘制相对热图

另一方面，牛和大鼠之间的相关系数较小（-0.18~0.13），为了在另一个方面比较这个验证结果，我们对我们的qRT-PCR结果进行了统计检验（表2）。

由于40个性别二型基因，是使用广义线性模型（GLM）进行RNA-seq分析，如M1（两组比较）和M2（使用ANODEV），因此，采用t检验和方差分析（ANOVA），来应用类似的方法qRT-PCR分析的统计学检验。qRT-PCR检测到牛和大鼠中性二型基因。

表2 对40个性别二型基因进行40个性别二型基因统计检验

这一统计结果表明虽然两个物种的基因表达模式并不一样，但在两个物种中普遍可以发现性别二型基因，为了检验这一现象我们使用聚类R包，对40个基因的基因表达数据进行了层次聚类分析。

采用欧氏距离来确定样本之间的距离，并计算聚类各分支的代表距离。此外，计算轮廓分数，以检测三种类型的数据集的最佳聚类数量，结果表明，牛的大鼠RNA-seq和qRT-PCR的最佳簇数k估计为4，而牛的qRT-PCR估计为7。

基于这些估计的最优聚类数量，分层切割树，为每个数据集进行可视化，这三棵树，不仅揭示了样本可以通过组织清楚地区分，而且样本也可以根据每个组织的性别清楚地划分。

这些经过技术验证的40个基因，在牛和大鼠物种中都是性别二型基因。最后，在随机选择的40个deg中有33个在M2中显著检测到。

为了测量m2的准确性，进行统计学我们比较了qRT-PCR和RNA-seq的结果（表3），结果有21个基因在qRT-PCR中进行了技术验证。目前我们将展开本次实验研究的讨论。

表3 4用0个基因进行qRT-PCR和RNA-seq统计学

讨论：通过rna序列分析检测到牛的性染色体基因

由于我们在本研究中，排除了表达水平超级低的基因我们预计所有y连锁基因，都表现出显著的雄性倾向性表达。

不过根据我们的分析只有少数y连锁的非deg，基于bosTau7参考基因组的8个注释基因中，有5个基因（EIF1AY、PRAME、SRY、TBL1X、ZNF280B）在M1和M2中检测到不显著。3个基因（EIF1AY、PRAME、SRY）的TMM归一化值极低。

为了研究每个基因与其类似基因之间的交叉配对的可能性，我们使用一个基本的局部比对搜索工具，使用ZNF280Bmrna序列，对参考基因组（BLAST）进行了局部比对。

结果我们不仅观察到ZNF280B基因的查询覆盖率为98%，而且还观察到了17号染色体上一个未注释基因的表达，查询覆盖率为92%。

此外，先前的研究发现，ZNF280B是一个常染色体衍生的y染色体基因，因此y-连锁非deg的出现，可能是由于它们的低表达或重复之间的差异表达，也可能是由内参基因中两个序列类似的基因引起的。

在FDR调整后的P值<为0.05时，我们在X染色体上没有观察到雄性偏倚基因，不过在FDR调整后的P值<为0.1时，三个x连锁基因，在雄性中显著上调三个基因；ZIC3、CXHXorf34和TMEM35在雄雌间的表达差异相对较小，但与其他752个显著的性二型基因相比，个体差异也较高。

在325个x染色体注释基因中，有26个在M1或M2中，检测到雌性偏倚的deg，这些基因可能为牛物种的x染色体失活逃逸基因，提供了一个很好的候选基因。

在肌肉组织中，共鉴定出5个这样的基因：CA5B、RBM3、SH3KBP1、XPNPEP2和XIST。此外，NXF3、XIST和ZFX在垂体组织中检测为DEGs。

在这些基因中只有一个x-连锁基因XIST，一般在一些组织中被发现使用M1和M2中显著检测到XIST。

不过XIST在肝组织中，没有明显检测到；这一结果令人惊讶，由于我们预计这个雌性特异性x失活相关基因，将包含在整个组织中。

对原始数据的检查显示，雌性和雄性的肝脏中，均有表达异常值，结果没有统计学意义，最终将XIST从DEG列表中排除。

在代谢组织和脑垂体中，表现出性别二型性

我们研究了为什么雄性偏倚基因，只在垂体组织中富集的缘由，DAVID显著报道了大量与肌肉发育相关的基因集，包括横纹肌组织发育（1.30E-03）、肌肉组织发育（1.60E-03）、心肌组织形态发育（1.10E-02）、肌肉组织形态发育（1.10E-02）、心肌组织发育（5.10E-02）、骨骼肌组织发育（7.40E-02）和骨骼肌器官发育（7.40E-02）等。

之前的研究结果表明，雄激素在肌肉和骨骼的发育中起着超级重大的作用，特别是在雄性中，似乎在垂体组织中显著检测到的性别二型基因显示出雄性偏倚的富集，由于许多在该组织中发现的显著基因在肌肉发育中发挥了重大作用。

结果表明，肌肉组织对基因表达的性别差异的影响高于其他组织，这些结果与早期使用M1进行的研究的结果一致，这些结果和检测到的基因数量表明肌肉组织中的基因表达，比其他组织更受性别二型性的影响。

脑下垂体，是代谢和生殖节律性控制的关键结构；这种功能是生长激素（GH）等激素和促性腺激素以脉动方式分泌的结果，虽然性别二型性、代谢和节律过程之间的关系复杂且不清楚，但对重大的节律过程相关基因，和通路的鉴定，可能表明垂体的生理作用及其对外周组织的影响。

例如，众所周知肝脏基因表达的大量性别差异是由中枢昼夜节律系统控制的，以及由循环GH的脉动性造成的，不过不能排除核心昼夜节律振荡器的非昼夜节律作用，昼夜节律相关基因，可能是这些信号级联下游的DE调控。

综上所述本研究结果表明，垂体表达基因可能至少部分参与了从激素分泌、外周组织基因表达开始，最后在代谢过程中建立性别二型性的连锁事件。

腹部脂肪组织基因表达的性别二型性程度

使用M1进行的RNA-seq分析结果显示，从腹部脂肪组织中，检测到的DEGs不仅与肝脏、肌肉和垂体相比显示出很小的组织特异性，而且检测到的性二型基因数量也最低，这表明在脂肪组织中，发现的大多数显著的基因对性别二型性的影响小于其他组织。

脂肪酸代谢和区域脂肪分布的性别偏倚差异，在以往的研究中得到了广泛的认识，不过关于腹部脂肪基因表达的性别差异的研究数量有限，不足以阐明脂肪组织中，性别二型性的程度。

与此同时一些研究报告称，内脏脂肪组织功能的某些方面似乎是独立于性别行为的，例如，之前的一项研究报道，在人类脂肪组织中没有观察到雄激素结合和雌激素受体表达的性别差异。

此外在脂肪组织中，检测到比其他组织中更少的另一个可能缘由，可能是韩国本土牛的育种策略，由于高大理石纹牛肉是韩国牛肉行业的消费者和制造商的首选，牛的育种系统主要聚焦在增加脂肪大理石纹，这很可能导致雄雌基因表达差异逐渐减少。

需要在转录组学、蛋白质组学和生理水平上进行进一步的研究，以得出结论脂肪库是否与其他代谢组织相比具有较小程度的性别二型性。

我们观察到雄性中大量的脑垂体表达上调，而雌性中大多数的肝脏、脂肪和肌肉deg表达上调。之前有几项研究比较了组织之间基因表达水平上的性别二型性，包括垂体。

虽然这些结果，并不能解释为什么垂体基因似乎表现出性别特异性的表达模式，但为未来的研究确定垂体中雄性偏倚的基因表达模式，在人类和啮齿动物等其他物种中可能有利于未来的研究。

结语

为了研究牛的性别二型性，我们使用两种不同的统计方法进行了RNA-seq分析。

结果，在不同组织中成功鉴定了许多性二型基因和通路；表达表现出较强的组织特异性模式。

为了验证鉴定出的牛两性二型基因，我们分别对牛和大鼠物种进行了qRT-PCR实验。结果显示，虽然这两种哺乳动物的性别二型基因是共享的但基因表达模式在不同的组织中有所不同。

我们的研究结果表明，许多生物过程可能参与了牛的代谢组织的性别二型性。特别是，脑垂体中性别二型性基因的表达模式，不仅表明了大脑本身的性别二型性的影响，而且也表明了大脑本身的影响中枢神经系统的外周组织在决定性别差异。最后，我们得出结论，在思考多个组织的RNA-seq研究中，两种统计方法都有其优缺点。