我们知道,无参考基因组的物种进行全长转录组测序最宝贵的就是得到物种最完善的全长转录本(isoform)序列信息,基于这些isoform序列信息方可进行转录本的结构和功能分析。另外,结合二代数据还可进行基因及转录本定量和差异分析(如图1)。
图1无参全长转录组基本分析流程
为了打造优质、高效、实用的Iso-seq服务,菲沙基因特地为大家带来无参全长转录组分析结果快速解读指南。结合实际案例和关键问题进行全方位剖析答疑,让您一拿到分析结果就能迅速抓到重点进行生物学故事阐述。
正题之前,我们先来了解下无参全长转录组三代和二代联合分析的结果文件组成,主要包括两部分内容,分别是“Iso-Seq_Result”和“RNA-Seq_Result”。其中的RNA-Seq展示的是全长isoform的表达定量、表达差异及功能富集信息,Iso-Seq展示的是全长isoform的功能及结构信息。今天小编将以Iso-Seq部分为例,结合实际疑问为您揭晓结果文件的具体含义。
首先,结果的第一部分,我们展示的是数据产出和全长isoform的统计。数据产出指标根据PacBio官方Iso-Seq pipeline指定为Polymerase read(酶读长)、Subread(原始序列的每个子序列,也叫插入序列)、CCS(ROI,subreads校正后的cDNA序列)和FLNC(全长非嵌合序列)。分别统计相应指标的数据量(Total bases)、测序深度(# of Reads)、平均长度(Mean length)和Read N50。Read N50是将得到的相应指标read按照长度从长到短排序,依次累加相应指标read的长度直至达到总长50%时的相应指标read的长度,一定程度上能够反映相应指标的优质度。同时进行各个指标的长度分布分析(图2)。
图2 reads长度分布图
(B、C、D,x轴分别表示Subreads、CCS、FLNC reads的长度,左边Y轴为柱形图坐标,表示长度在一定范围内(X轴)的reads数量。右边Y轴为曲线图坐标,表示长度大于一定数值(X轴)的reads数量)
图3 Iso-Seq流程图
图4 ICE聚类和自我纠错后的isoform序列
(cb: 序列编号;c:根据序列相似性将FLNC分为多个聚类簇(c120618代表第120619个聚类簇);f3p0:f代表全长非嵌合序列条数,p代表非全长非嵌合序列条数;2124:代表该序列碱基数)
全长转录本序列功能预测
全长转录本基础功能注释包括NR、GO、KO、KOG、Swiss-Prot等数据库的注释,注释情况统计如图5。
图5功能注释结果统计
(NR、GO、KO、KOG及Swiss-Prot:分别代表在这五个数据库中注释成功的isoforms数目及其占isoforms总数的比例;Unannotated:在以上五个数据库中都没有任何功能注释的isoforms数目及其占isoforms总数的比例)
图6 GO分类图
(纵坐标为GO三个大类的下一层级的GO term,横坐标为注释到该term下(包括该term的子term)的转录本个数。3种不同分类表示Go term的三种基本分类(从上往下依次为生物学过程,细胞成分,分子功能)
图7 GO分类图
(横坐标为KOG的26个group,纵坐标为注释到该group下的转录本个数)
图8 KEGG分类图
【纵坐标为KEGG代谢通路的名称,横坐标为注释到该通路下的转录本个数。将转录本根据参与的KEGG代谢通路分为5个分支:细胞过程(Cellular Processes),环境信息处理(Environmental Information Processing),遗传信息处理(Genetic Information Processing),代谢(Metabolism),有机系统(Organismal Systems)】
基于之前的注释结果,对在蛋白数据库没有注释信息的转录本,用CPAT进一步评估其编码潜能;再通过一定的编码潜能筛选过滤,留下来的作为最终的非编码RNA预测的结果,预测结果如图9。
图9 lncRNA预测结果展示
(ID:转录本ID;mRNA_size:序列长度;ORF_size:开放阅读框长度;Fickett_score,Hexamer_score:用于识别编码区域、评估编码潜能的值;coding_prob:编码潜能。coding_prob越小代表lncRNA序列更有依据)
过滤掉被预测为lncRNA的部分,利用transDecoder预测序列的ORF基因结构注释及功能预测,如图10。其中的transdecoder.cds文件为预测到的CDS序列,transdecoder.pep文件为预测到的氨基酸序列,secretory_protein文件为预测到的分泌蛋白。
图10 orf预测结果
可变剪切分析
利用特定软件及pipeline对去冗余后的全长转录本序列预测可变剪切,见图11。
图11可变剪切示意图
(比较同一聚类下的任意两条转录本发生可变剪切的位置,红色区域的gap即为可变剪切发生的位置)
由于重复单位的核苷酸不同以及重复次数不完全相同,造成了SSR长度的高度变异性,其中最常见的是双核苷酸重复类型,如(CA)n。我们的分析中采用MISA软件(1.0版,默认参数);对应的各个unit size的最少重复次数分别为:1-10,2-6, 3-5,4-5,5-5,6-5(如:1-10,以单核苷酸为重复单位时,其重复数至少为10才可被检测到;2-6,以双核甘酸为重复单位时,其最少重复数为6)对转录本进行SSR检测,如图12。
图12 SSR预测结果
(ID:转录本ID;SSR nr.:转录本上SSR编号;SSR:SSR信息;SSR type:SSR类型;size:SSR位点长度;start:SSR起始位点;end:SSR终止位点)
随着PacBio sequelⅡ的引入,Iso-Seq成本大幅下降,周期大大缩减,进行多个样本的三代全长转录组比较研究越发便(bian)宜。例如不同时期或组织,我们可以通过比较组间的转录本序列及功能差异从而深入挖掘生长发育、细胞分化与疾病发生机制。
小编针对广大客户对结果文件的疑惑做了相应的答案总结,相信以上的详细解读能够cover到您关注的point,如仍有疑问可直接拨打官网电话或email菲沙技术support@frasergen.com。
PacBio ISO-seq技术的长度长可以解决二代测序技术由于读长过短引发的一系列问题(例如拼接组装,GC偏好等),能够识别二代转录组无法识别的同源异构体 (isoform)、同源基因、超家族基因或等位基因表达的转录本。菲沙基因是国内首家拥有PacBio Sequel及SequelⅡ三代测序仪的高新技术企业,已与广大科研学者合作完成各种三代测序分析科研项目累积数百项,拥有成熟的实验技术、强大的生信分析团队和丰富的项目经验,期待与您的合作,助您的科研快速起飞!