10月31 【成功案例】湖羊骨骼肌发育过程中全基因组范围LncRNAs的显著变化

Genome-Wide Analysis Reveals Extensive Changes in LncRNAs during Skeletal Muscle Development in Hu Sheep

湖羊骨骼肌发育过程中全基因组范围LncRNAs的显著变化

 

杂志： genes 

影响因子：3.600

PMID: 28763026

 

 

 

研究背景：

羊肉具有蛋白含量高、脂肪和胆固醇含量低等优点，促进羊的肌肉生长可提高羊肉产量。过去关于羊骨骼肌生长的研究主要集中在蛋白质编码基因，然而基因组绝大部分是非编码序列。长非编码RNA（lncRNAs）是一种重要的非编码RNA，越来越多研究发现某些功能型的lncRNAs是新的肌肉调节元件，发挥重要作用。关于羊肌肉相关lncRNA的转录组学研究较少，尤其关于湖羊的研究更少， lncRNA在羊骨骼肌生长过程中的表达类型和潜在作用大部分都是未知的，因此了解羊不同生长时期肌肉转录组的动态变化具有重要意义。

本文拟研究湖羊三个关键生长阶段（110天胎儿期，5天幼年期，2岁成年期）肌肉lncRNA的表达模式和潜在作用，筛选差异表达lncRNA的和DEGs（差异表达基因）的靶基因。同时通过lncRNA-gene网络预测湖羊肌肉生长的lncRNA潜在调控因子，获得湖羊肌纤维生长相关的转录水平的候选调节因子。

材料和方法

材料：相同的自然光环境下饲养的胎儿期、幼年期、成年期湖羊，每个时期选取3个随机样本，取既定部位（12和13胸椎之间）的LD肌肉样本，采样后即刻液氮冷冻提取RNA。

测序平台：Illumina platform

分析平台：百迈客云平台（BMKCloud），分析内容如下：

在百迈客云平台上，对胎儿期、幼年期、成年期湖羊肌肉细胞进行RNA-seq测序后转录组拼接。进行注释后，先排除＜200个核苷酸或单外显子的转录子,再用CPC,CNCI,PFAM,CPAT等工具从未知转录组中筛选候选lncRNA。四种方法FPKM＞0.1的转录子交集定义为lncRNA转录子。

应用DEGseq packages (1.10.1) 来分析组间的表达差异，FDR值＜5%，log2(倍数变化)的绝对值被定义为差异化表达。

预测差异表达的lncRNA的顺式靶基因和反式靶基因。

为分析lncRNA的主要功能，通过NCBI的Nr、GO、KEGG、COG数据库比对进行靶基因和DEGs注释，KS≤0.05，信号通路矫正p值≤0.05的GO term被定义为显著富集的。

为进一步研究lncRNA和他们的靶基因的相互作用，建立LncRNA-Gene共表达网络。

结果：

1.测序拼接和转录组分析

胎儿期、幼年期、成年期3组RNA-seq序列质控后，经过去接头、多聚尾、低质量序列后，每组平均获得了65578070，65591958,71241551的拼接序列，9个文库平均GC含量为54.39%，每个样本Q30值≥91.71%， 92%以上与O. aries参考基因组特异性匹配，分析显示胎儿期有45.1%的序列与外显子区域匹配，后期阶段匹配度更高（幼年期53.8%、成年期53.7%），在幼年、成年期组的内含子和基因间序列比例低于胎儿期。（表1）

表1.矫正后序列与不同阶段湖羊参考基因组比对

2.lncRNA筛选和特征性描述
为研究湖羊肌肉lncRNA的基本特征，筛选lncRNA并与mRNA比对。CPC、CNCI、PFKM、CPAT交集部分有6924个lncRNA转录子，包括已知的保守lncRNA，肌肉分化相关的lncMD（图1A）

图1 A 韦恩图  图1B. lncRNA和mRNA在每个染色体上的分布比例。红色和黑色线条分布代表lncRNA和mRNA，蓝色线条代表相应染色体在基因组中的大小比例

lncRNA转录子分为4606个lincRNA（66.5%），1131个内含子lncRNA（16.3%）,1187个反义lncRNA（17.1%）。与mRNA类似，lncRNA转录子在26对染色体和X染色体中分布广泛，但线粒体中不存在。而且lncRNA在几个染色体中的比例与染色体大小比例一致，尤其是18号染色体，表明在此研究中相应lncRNA可能体现重要功能（图1B）

比较lncRNA和mRNA的外显子特征，结果显示大部分lncRNA每个转录子包含2-5个外显子（均值3.3），比mRNA少（均值7.9，图1C）

图1C 湖羊lncRNA每个转录子外显子数      图1D 湖羊lncRNA外显子大小分布

另外大部分lncRNA包含2个外显子，而包含2个外显子的mRNA仅占比3.9%，明显低于lncRNA。lncRNA中外显子平均长度相对长于mRNA，大部分小于200bp（图1D）。与蛋白编码基因一致，在肌肉生长过程中同一阶段lncRNA表达趋势相似，且平均表达水平比与蛋白编码基因低（图2A,B）

每组重复间的相关性分析显示相关性高（图3A-C），而且在6924个表达的lncRNA转录子中，于某一生长阶段特异性表达的占比33.03%，这个比例在蛋白编码基因中较低（4%），提示lncRNA的特殊意义和动态特性。胎儿期特异性lncRNA为1042个，远高于幼年期（626）和成年期（619），表明lncRNA在早期发展阶段的重要意义。

图2A胎儿期、幼年期、成年期湖羊肌肉lncRNA的FPKM分布

图2B胎儿期、幼年期、成年期湖羊肌肉蛋白编码基因的FPKM分布

图3 每组重复间的相关系数分析（A胎儿期组，B幼年期组，C成年期组）

3.差异表达分析及靶点基因预测

3个时期两两比较分析显示在胎儿期vs幼年期、幼年期vs成年期、胎儿期vs成年期（对照组vs实验组），分别有27、14、92个lncRNA，239 270 1437个基因是特异性表达的。(|log2FC| > 1, FDR < 0.05).值得一提的是，幼年期vs成年期组差异表达的lncRNA量最低，在胎儿期vs幼年期组和幼年期vs成年期组，DEGs的数量几乎一样多，表明产前和产后尽管在时间只相差1一个月，但是转录水平的差异较大（图4A-D）。

图4 对比组差异表达lncRNA和基因的数量（A 不同对比组差异表达lncRNA数量的韦恩图 B 不同对比组DEGs数量的韦恩图C不同对比组差异表达lncRNA总数 D不同对比组DEGs总数）

差异表达的lncRNA中，在胎儿期vs幼年期有36个上调lncRNA，42个下调lncRNA，幼年期vs成年期有13个上调lncRNA，28个下调lncRNA，胎儿期vs成年期有68个上调lncRNA，78个下调lncRNA。DEGs中，胎儿期vs幼年期组有1028个上调基因和487个下调基因、幼年期vs成年期组有659个上调基因和900个下调基因、胎儿期vs成年期有1862个上调基因和1749个下调基因。差异表达的lncRNA和DEGs的分层集群显示幼年期和成年期的表达模式相似而与胎儿期有所差异（图4E,F）。

图4 E 差异表达的lncRNA的分层集群  F DEGs的分层集群

为进一步评估RNA测序的结果，选择MYOG(从胎儿期富集的与晚期肌肉细胞分化相关的基因)，MYH7(肌肉结构基因)，5种差异表达的lncRNA，4种DEGs进行qRT-PCR分析。表达量与RNA-Seq结果一致。结果表明表达与测序结果有较好相关性，表明测序结果可行度高。

lncRNA可以通过与靶基因顺式和反式作用发挥功能，通过两两对比分析，预测相邻上下游100kb的和或差异表达lncRNA的互补蛋白编码基因。得到共201个靶基因。

图5 qRT-PCR和RNA-Seq分别验证不同时期湖羊肌肉差异表达lncRNA和基因的表达水平

4.差异表达的lncRNA和mRNA靶基因的生信分析

为进一步研究差异表达的lncRNA,通过与NR / GO/ COG和KEGG数据库比对，对lncRNA和DEGs的靶基因做注释。

表2 差异表达lncRNA和DEGs的功能注释

基于GO数据库，靶基因和DEGs被分别归到生物学过程，细胞构成，分子功能三个功能类别，而且在所有比对中，在肌肉生长过程中发挥重要作用的器官形成，骨骼肌系统发展与应激被归类为显著富集的GO terms。

图6不同对比组的差异表达lncRNA和基因的靶基因的top GO分析（A不同对比组的差异表达lncRNA top20个GO terms B不同对比组的DEGs top20个GO terms）

COG功能聚类分析把差异表达的lncRNA和基因分别归为17和24个类别。在胎儿vs幼年期组，氨基酸转运合成，碳水化合物转运合成，离子转录合成的靶基因和DEGs的比例比其他两个对比组高，提示这类基因在产前阶段肌肉发育过程中的重要性。此外，对靶基因和DEGs分析显示，在幼年期vs成年期组信号转导机制和细胞骨架比其他对比组比例更高，表明这类terms和相应基因在产后肌肉生长中的重要意义（图7）。

图7 差异表达lncRNA和DEGs的靶基因的COG分类（A 差异表达lncRNA的靶基因的COG分析 BDEGs的靶基因的COG分析）

 

根据KEGG分析，胎儿组vs幼年组，幼年组vs成年组，胎儿组vs成年组的21、6、51个lncRNA的靶基因分别与40、16、48个通路对应。虽然每个对比组的差异表达lncRNA的靶基因数量少，通路（MAPK信号通路，间隙结合，钙信号通路，胰岛素信号通路，激动蛋白细胞骨架调节）是可被诱导的。这些通路中，MAPK信号通路在胎儿组vs幼儿组，胎儿组vs成年组中最常见，提示此通路和相关lncRNA可能参与了湖羊肌肉的生长调控。

在胎儿期vs幼年期，胎儿期vs成年期对比组，氧化磷酸化，碳代谢和心肌收缩是DEGs数量最高的前三个通路，在幼年期vs成年期组却不是。值得注意的是，途径如糖酵解和生成，脂肪酸降解，氨基酸合成的和过氧物酶体增殖激活受体 (PPAR)信号通路特异性存在于胎儿期vs幼年期组，提示这些途径中DEGs在早期肌纤维生长过程中的重要作用。此外磷酸肌醇激酶-3-激酶蛋白激酶B(PI3k-Akt)信号通路和类固醇合成在幼年期vs成年期组中特异性存在，强化了它们在产后肌肉生长中的作用。总之，差异表达的lncRNA和DEGs在肌肉生长调节方面显示了巨大潜力。

表3 不同对比组的DEGs的KEGG通路

5.肌纤维增长相关的潜在功能性lncRNA筛选

为研究lncRNA如何与靶基因来调节肌肉生长，依据FPKM进行差异表达lncRNA和对应的差异表达靶基因的共表达分析,共获得了15（胎儿组vs幼年组），7（幼年期vs成年期），37（胎儿期vs成年期）个共调节的lncRNA-基因对。三个网络都显示了与肌肉生长相关的候选lncRNA（图8）。差异表达的lncRNA和相应的差异表达的顺式和反式作用的靶基因来构建lncRNA-基因互作网络。差异表达的靶基因直接与肌肉生长过程相关，对靶基因对应的差异表达lncRNA也进行了分类。

图8 不同对比组lncRNA-基因关系网络

（A胎儿期/幼年期组，B幼年期/成年期组，C胎儿期/成年期组；基因以紫色表示，lncRNA以绿色表示，上调用三角形表示，下调用圆形表示，顺式作用实心线表示，反式用虚线表示）

随后通过GO分析得到了8个生物学过程和相关20个mRNA和41个lncRNA，如RTL1和JPH2(表S5)。

qRT-PCR分析显示在肌肉生长过程中lncRNA表达和对应靶基因变化一致，包括TCONS_00606329和它的靶基因elastin，TCONS_00758916、 TCONS_00685981 和他们的靶基因phosphofructokinase muscle and ankyrin repeat SOCS box containing 8; TCONS_00297401和靶基因泛素特异性蛋白酶2编码基因，
TCONS_00377352, TCONS_00381991、TCONS_00381994和他们的靶基因RTL1 (图 9) 。这些结果进一步揭示了lncRNA和他们靶基因的互作关系。

图9胎儿期、幼年期、成年期差异表达lncRNA和共表达的靶基因的表达水平的qRT-PCR验证

 

研究结论

本研究基于Illumina平台，筛选出一系列与湖羊三个关键生长时期相关的lncRNA和基因，筛选得到的6924个差异表达的lncRNA转录子中，特异性表达于某一生长阶段的占比高于蛋白编码基因（33.03%：4%），提示lncRNA的特殊意义和动态特性。且在胎儿期阶段特异性lncRNA数量较多，提示lncRNA在早期发展阶段的重要意义。

通过两两比较分析筛选特异性表达的lncRNA和基因，在胎儿期vs幼年期、幼年期vs成年期、胎儿期vs成年期（对照组vs实验组）分别得到特异性表达的27、14、92个lncRNA，239 270 1437个基因。差异表达的lncRNA和DEGs的分层集群显示幼年期和成年期的表达模式相似，与胎儿期有差异。

选择MYOG，MYH7，5种差异表达的lncRNA，4种DEGs进行qRT-PCR分析。表达量与RNA-Seq结果一致，测序结果可信度高。

通过与NR，GO,COG,KEGG数据库进行比对，对差异表达lncRNA和基因的靶基因进行了多种肌肉相关生物学过程的注释。且构建了lncRNA-gene 共调控网络，提供了候选lncRNA的有价值的信息， qRT-PCR分析显示在肌肉生长过程中lncRNA表达和对应靶基因变化一致。

 

创新点

首次对三个关键生长阶段的湖羊肌肉生长相关lncRNA进行系统性描述。从转录子结构和表达模式上研究lncRNA，然后进一步预测lncRNA的靶基因来进行功能研究。为进一步研究湖羊lncRNA功能提供了依据，为促进湖羊肌肉生长提供了有价值的信息。