CN110619926A - 一种识别全部rna剪切位点的分析方法及分析系统 - Google Patents
一种识别全部rna剪切位点的分析方法及分析系统 Download PDFInfo
- Publication number
- CN110619926A CN110619926A CN201910726790.4A CN201910726790A CN110619926A CN 110619926 A CN110619926 A CN 110619926A CN 201910726790 A CN201910726790 A CN 201910726790A CN 110619926 A CN110619926 A CN 110619926A
- Authority
- CN
- China
- Prior art keywords
- splicing
- junctions
- asja
- linear
- read
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Classifications
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
- G16B20/30—Detection of binding sites or motifs
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B35/00—ICT specially adapted for in silico combinatorial libraries of nucleic acids, proteins or peptides
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Chemical & Material Sciences (AREA)
- Theoretical Computer Science (AREA)
- Molecular Biology (AREA)
- Biophysics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Biotechnology (AREA)
- Evolutionary Biology (AREA)
- General Health & Medical Sciences (AREA)
- Medical Informatics (AREA)
- Biochemistry (AREA)
- Library & Information Science (AREA)
- Analytical Chemistry (AREA)
- Genetics & Genomics (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明提出了一种识别全部RNA剪切位点的分析方法,英文名称为Assembling Splice Junctions Analysis,简称ASJA,包括以下步骤:步骤A:ASJA算法识别全部的剪接事件的连接点,包括以下:步骤A1:对RNA‑seq数据进行比对以及产生拼接后的转录本;步骤A2:提取所有的剪接连接点,包括线性剪切位点、反式剪切位点和融合剪切位点;步骤A3:注释和整合所述剪接连接点;步骤B:评估ASJA的效能。本发明方法提供唯一识别剪接事件的位点以及对每种连接点的标准化定量。本发明的ASJA方法比现有方法运行速度快,且准确性高。本发明通过已经发表过的RNA‑seq数据进行进行评估,通过ASJA可以有效地分析检测和定量剪接连接点,同时也发现了新的剪接连接点。本发明还提出了一种识别全部RNA剪切位点的分析系统。
Description
技术领域
本发明涉及RNA剪切技术领域,尤其涉及一种RNA剪切位点的分析方法及分析系统。
背景技术
RNA剪接是指细胞核内基因信息在转录过程中或是在转录过后的一种修饰,是内含子从一个新合成的前体RNA(pre-RNA)中被移除,外显子被合并,产生一个成熟的含有剪接位点的RNA的过程。人类多外显子基因的精确可变剪接形式造就了基因组的多样性。异常的剪接与很多疾病相关,例如:癌症、神经衰退和肌肉萎缩。通过高通量测序技术(RNA-seq),大量的非线性的剪接事件连接点包括反式剪接连接点、融合剪接连接点被识别。研究者发现融合转录本与肉瘤和恶性血液疾病有关,例如,BCR-ABL1在成人急性髓系白血病中扮演着重要的角色,同时被用于慢性髓性白血病的生物学标记物。反式剪切形成的环状RNA(circRNA)具有单链、共价闭合环结构的特点。从高通量测序计算的结果中已经得知circRNA在转录组中广泛存在。尽管大多数circRNA的功能基本上还未被阐明,但已知的功能包括结合miRNA或蛋白质,调节其亲本基因和生成蛋白质。此外,目前的研究已经揭示一些circRNA在神经系统、癌症发展和先天免疫反应中的关键作用。
尽管目前RNA-seq分析工具在预测剪切位点方面有较好的性能,但是在剪接连接点的检测以及其定量上还是有局限性的。Nellore等人通过比对GEO数据库中大约2万个RNA-seq样本最终得到大量的线性的外显子-外显子连接点。然而,由于错误的读段位置、样本特异性变异和基因组错误组装的区域,他们直接从映射读段中提取的剪接连接点中可能具有较高的假阳性率。当前工具只倾向于检测一、两种类型剪接连接点。为了在RNA测序数据的基础上对整个转录组进行全面的评价,有必要开发一个综合的工具来识别和描述不同类型的剪接连接点。
发明内容
本发明提供了一个能够从RNA-seq测序数据中检测和注释所有剪接事件连接点的拼接剪接连接点分析(Assembling Splice Junctions Analysis,ASJA)的方法。
本发明通过对GSE77661数据集合的分析,发现了大量的位于基因间隔区域的、未被注释过的线性剪接连接点(10,353),其功能有待深度探索。
本发明提出了一种识别全部RNA剪切位点的分析方法,包括以下步骤:
步骤A:ASJA算法识别全部的剪接事件的连接点,包括以下:
步骤A1:对RNA-seq数据进行比对以及产生拼接后的转录本;
步骤A2:提取所有的剪接连接点,包括线性剪切位点、反式剪切位点和融合剪切位点;
步骤A3:注释和整合所述剪接连接点;
步骤B:评估ASJA的效能:
步骤B1:从GEO中下载GSE77661数据集以及SRR934794SRR934744SRR934930SRR2976716和SRR2976715的SRA文件并用FastQC将其转化fastq 文件;
步骤B2:在Linux环境中运行对每个样本执行ASJA命令,检测和注释数据集中全部的剪接连接点;同时检查ASJA输出的日志文件以及观察linux的load average、Cpu(s)、Mem参数;
步骤B3:统计分析数据集检测到的剪接连接点并与其他方法进行比较。统计ASJA的运行时间和占用Linux的内存。
相比于剪接位点提取的另一个软件Mapsplice2,ASJA是从拼接的转录本中提取的线性剪接位点,而且是使用剪接位点两侧的外显子最小覆盖度来定义剪接位点的表达量而不是简单的位点的横跨读段数目。Mapsplice2使用带有"ZF:Z" 标签的sam文件来提取融合基因,而ASJA是基于STAR的嵌合比对的结果。
本发明所述步骤A1中对RNA-seq数据进行比对包括:使用FastQC软件对 RNA-seq数据进行质量控制和数据过滤;过滤后的读段使用STAR比对软件进行两轮基因组比对。
本发明所述步骤A1中产生拼接后的转录本包括:将STAR比对的结果输入到StringTie中得到用于线性剪切连接点检测的拼接后的转录本;已知的转录本注释文件被用于转录本拼接模型的参考文件,以此来指导装配过程,其结果得到已知的转录本和未在GENCODE中注释的转录本。
本发明所述步骤A2提取剪接连接点包括:线性剪接连接点、反式剪接连接点、融合连接点。
本发明提取所述线性剪接连接点包括:在比对和拼接结束后,使用perl脚本提取来自拼接后转录本的线性剪接连接点,同时用下式对连接点进行定量:
COV(AB)=min(∑cov(A),∑cov(B));
其中,cov(A)cov(B)表示每个外显子的覆盖度;
所述连接点的表达量需要根据已知注释的连接点的覆盖度进行标准化:
CPT(AB)=COV(AB)(*10e7)/TC;
其中,TC表示已知全部连接点的覆盖度。
本发明提取所述反式剪接连接点包括:使用perl脚本从嵌合比对的结果中提取反式剪接连接点;所述反式剪接连接点满足:
1.嵌合读段比对到同一条染色体和同一条链上;
2.3’末端和5’末端位点之间的距离小于3,000,000bp;
3.剪切方式符合GT/AG法则;
4.基因组中反式剪接连接的起始和结束位置不合理;比对到线粒体和未知染色体的读段要被剔除。
本发明提取所述融合连接点包括:从嵌合比对的结果中提取所述融合连接点,并应用以下步骤降低假阳性率:
1.比对到线粒体和其他未标记的染色体或未映射到重叠群的读段被过滤掉;
2.反式剪接的读段被剔除;
3.读段数目大与1的融合连接点被保留;ASJA算法同时计算了类似 STARChip中的SpanningReads变量验证融合连接点的有效性。
本发明所述步骤A3包括:准备剪接连接点的注释文件、计算剪接率、整合分析三种类型的剪接事件和过滤剪接事件;其中,
所述准备剪接连接点的注释文件包括:ASJA算法提供的初级线性连接点注释文件、从UCSC table browser中下载的BED格式的外显子注释文件;
所述计算剪接率包括:根据每个线性连接点注释到的基因来计算该位置发生剪接的概率,对未注释的连接点,ASJA算法将计算该连接点在对应的转录本中的剪接概率;Weight ratioi=CPTi/CPTm,其中CPTm表示连接点注释到基因的最大的剪接值;所述ASJA算法利用如下计算3’末端和5’末端的反式剪切率:
5′_ratio=(5′back_splicedread)/linerread;
3′_ratio=(3′back_splicedread)/linerread;
其中,5’back_splicedread代表5’末端的读段数目,linearread表示对应该5’端位置对应的线性连接的读段数目;
ASJA使用类似反式剪切率的方法来计算发生融合连接点的概率:
donorratio=donorread/sum(donorread,linerread);
acceptorratio=acceptorread/sum(acceptorread,linerread);
其中,linearread表示在线性剪接中对应donor、acceptor剪接点的读段数目;
所述整合分析三种类型的剪接事件包括:基于基因注释和剪接位点的信息, ASJA算法将三种不同的剪接事件整合到一个文件中,以便发现剪接事件之间的相互关系;该部分输出内容包括:线性剪接的识别号,基因名称,反式剪接识别号以及融合剪接识别号;
所述过滤剪接事件包括:ASJA算法通过对读段数目和剪接率的过滤得到可信度较高的剪接事件。
本发明所述步骤B中,使用来自验证集合的数据来评估ASJA的效能,包括 3个脑胶质瘤样本(已验证包含9个融合基因)和核糖体RNA-/RNaseR处理的人类PA1样品。使用所述GSE77661数据集的12个正常组织,7个癌组织和7 个匹配的癌旁组织的RNA-seq数据实际应用ASJA;其中,所述正常组织包括:脑,结肠,心脏,肝脏,肺和胃;所述癌组织包括:膀胱尿路上皮癌,乳腺癌,结直肠癌,肝细胞癌,胃癌,肾透明细胞癌和前列腺腺癌。
基于以上方法,本发明还提出了一种识别全部RNA剪切位点的分析系统,所述系统包括:
连接点识别模块,其通过ASJA算法识别全部的剪接事件的连接点;
能效评估模块,其用于评估ASJA的效能。
本发明提出的RNA剪切位点的分析方法包括使用STAR比对软件中的嵌合比对以及StringTie拼接软件中的拼接转录本部分。ASJA提供唯一识别剪接事件的位点以及对每种连接点的标准化定量。
本发明通过已经发表过的RNA-seq数据进行评估,以及与现有软件进行比较。结果表明,ASJA都显示出良好的精确度,尤其是对于融合剪切位点的分析。与ASJA相比,大多数融合基因检测软件的召回率较高,但假阳性也很高,需要进一步的复杂筛选才能获得真正的融合基因。此外,ASJA的运行速度比其他的软件更快。通过ASJA可以有效地分析检测和定量剪接连接点,同时也发现了新的剪接连接点。
附图说明
图1为ASJA工作流程概要图。ASJA的结构自上而下分为三层分别是:使用STAR对RNA-seq的嵌合比对;根据剪接事件的不同特征识别剪接事件连接点;对剪接连接点的注释和整合。
图2为ASJA与其他软件在结果准确性和运行效率上的比较,其中(A)韦恩图显示了使用三种circRNA预测工具在12个正常组织中预测的circRNA的数量;(B)韦恩图展示了ASJA,ACFS和circRNA_finder在poly(A)-/RNaseR 和poly(A)-样品中预测circRNA的结果;图(C)软件分析验证数据集合所花费的时间。
图3为ASJA识别到的线性连接点的特征图,其中图(A)饼图展示了原始的线性连接点数目和高可信度的线性连接点;图(B)线性连接点注释到已知基因上的分布,非注释到基因上的认为是新型连接点;图(C)饼图展示了新型连接点在同性异构体和基因间隔区域的分布;图(D)新型连接点在不同组织中的个数分布;图(E)使用线性连接点数据对癌症和癌旁样本进行无监督聚类。热图展示的结果是基于log2(fold-change)>1and p<0.01(Wilcoxon test)。
图4为ASJA检测到的反式剪接连接点(circRNA)的特征图,其中(A) circRNA注释和未注释的分布;图(B)柱状图显示了基因中发生反式剪接事件的数目;图(C)circRNA在基因组中的位置分布;图(D)多维度筛选法用于鉴定高度丰富的circRNA。
具体实施方式
结合以下具体实施例和附图,对发明作进一步的详细说明。实施本发明的过程、条件、实验方法等,除以下专门提及的内容之外,均为本领域的普遍知识和公知常识,本发明没有特别限制内容。
本发明提供了一种识别全部RNA剪切位点的分析方法,具体包括以下步骤:
步骤A:ASJA算法,ASJA算法可通过下面的3个步骤识别全部的剪接事件的连接点:
步骤A1:对RNA-seq数据进行比对以及产生拼接后的转录本;
步骤A2:提取所有的剪接连接点,包括线性剪切位点、反式剪切位点和融合剪切位点;
步骤A3:注释和整合剪接连接点。
目前,大量检测剪接连接点的方法是根据已有的转录本信息出发,未考虑到转录本拼接的实际情况,ASJA算法检测到的线性剪接连接点是来源于从拼接后的转录本,且能检测到未知的剪切位点。在线性连接点的定量方面,ASJA算法不同于以往直接给出连接点的读段数目(read count),而是根据其注释和整体读段数目数目进行了标准化,这有利于连接点在不同样本之间的相互比较。对于非经典的剪接连接点的检测和注释,ASJA算法格外增加了其剪接率的指标,该指标是描述非经典在经典剪接方式(线性)中的占比情况,为非经典剪接方式的真实性提供了参考。
(A1)比对RNA-seq数据和产生拼接转录本
设置STAR比对软件的参数
使用FastQC软件对RNA-seq数据进行质量控制和数据过滤。所有过滤后的读段都使用STAR比对软件进行两轮基因组比对。第一轮,使用STAR在默认参数下形成基因组index作为比对的检索index。在第一轮拼接结束后,将其输出的含有剪接位点的文件重新利用生成用于第二比对的新检索index。在第二轮比对过程中,需要添加嵌合比对的参数-chimSegmentMin。
a、产生拼接后的转录本
将STAR比对的结果输入到StringTie中得到用于线性剪切连接点检测的拼接后的转录本。已知的转录本注释文件(GENCODE)被用于转录本拼接模型的参考文件,以此来指导装配过程。其结果不仅可以得到已知的转录本,未在 GENCODE中注释的转录本也会被记录。
StringTie的设置如下:stringtie<input mapped_bam>-f 0.1-o<out file>-p 4-G<GTF>.
(A2)提取剪接事件的连接点
b、提取线性剪接连接点
在比对和拼接结束后,使用perl脚本提取来自拼接后转录本的线性剪接连接点,同时用公式对连接点进行定量。
COV(AB)=min(∑cov(A),∑cov(B)),
其中,cov(A)cov(B)表示每个外显子的覆盖度。
连接点的表达量还需要根据已知注释的连接点的覆盖度进行标准化:
CPT(AB)=COV(AB)(*10e7)/TC,
其中TC表示已知全部连接点的覆盖度。
使用perl脚本计算剪接连接点的表达量相比于R语言可节约时间和内存。但同时需要对perl中哈希表和数据集合灵活运用。本发明中将连接点的每个外显子覆盖度存储到数组中且将该数组作为二维哈希表的值(value)。
c、提取反式剪接连接点
使用perl脚本从嵌合比对的结果中提取反式剪接连接点。反式剪接连接点需要满足以下几点:
1.嵌合读段需要比对到同一条染色体和同一条链上;
2.3’末端和5’末端位点之间的距离要小于3,000,000bp;
3.剪切方式符合GT/AG法则;
4.基因组中反式剪接连接的起始和结束位置不合理,即起止点的位置在终止点位置的后面。比对到线粒体和未知染色体的读段也要被剔除。
d、提取融合连接点
融合连接点也是从嵌合比对的结果中提取的。应用以下步骤来降低假阳性率:
1.比对到线粒体和其他未标记的染色体或未映射到重叠群(contigs)的读段被过滤掉;
2.反式剪接的读段被剔除;
3.读段数目大与1的融合连接点被保留。ASJA同时计算了类似STARChip 中的SpanningReads变量来验证融合连接点的有效性。
(A3)注释和整合剪接连接点
a、准备剪接连接点的注释文件
ASJA提供初级的线性连接点注释文件(直接从GENCODE中提取线性连接点信息),同时还提供剔除了通读基因(read-through)和同源基因(paralog gene) 的后注释文件。为了注释反式剪切和融合连接点,ASJA提供了从UCSC table browser中下载的BED格式的外显子注释文件。
b、计算剪接率
根据每个线性连接点注释到的基因来计算该位置发生剪接的概率,对未注释的连接点,ASJA将计算该连接点在对应的转录本中的剪接概率。Weight ratioi=CPTi/CPTm,其中CPTm表示连接点注释到基因的最大的剪接值。ASJA 利用如下计算3’末端和5’末端的反式剪切率:
5′_ratio=(5′back_splicedread)/linerread
3′_ratio=(3′back_splicedread)/linerread
其中5’back_splicedread代表5’末端的读段数目,linearread表示对应该5’端位置对应的线性连接的读段数目。
ASJA使用类似反式剪切率的方法来计算发生融合连接点的概率。
donorratio=donorread/sum(donorread,linerread)
acceptorratio=acceptorread/sum(acceptorread,linerread)
其中,linearread表示在线性剪接中对应donor、acceptor剪接点的读段数目。
c、整合分析三种类型的剪接事件
基于基因注释和剪接位点的信息,ASJA将三种不同的剪接事件整合到一个文件中,以便发现剪接事件之间的相互关系。该部分输出内容包括:线性剪接的识别号,基因名称,反式剪接识别号以及融合剪接识别号。
d、过滤剪接事件
ASJA通过对读段数目和剪接率的过滤得到可信度较高的剪接事件。ASJA 默认设置如下:线性剪接事件的读段数目大于2和剪接率大于0.01;反式剪接事件的读段数目大于1和反式剪接率大于0.01;融合剪接事件的读段数目大于1 和融合位点的一端与外显子剪切位点一致。
步骤B:评估ASJA的效能包括:
步骤B1,从GEO中下载GSE77661数据集以及SRR934794
SRR934744SRR934930SRR2976716和SRR2976715的SRA文件并用FastQC将其转化fastq文件。
步骤B2:在Linux环境中运行对每个样本执行ASJA命令,检测和注释数据集中全部的剪接连接点。同时检查ASJA输出的日志文件以及观察linux的load average、Cpu(s)、Mem等参数。
步骤B3:统计分析数据集检测到的剪接连接点并与其他方法进行比较。统计ASJA的运行时间和占用Linux的内存等。
本发明使用来自GEO(GSE77661)的12个正常组织,7个癌组织和7个匹配的癌旁组织(NT)的RNA-seq数据实际应用ASJA。
其中,正常组织包括脑,结肠,心脏,肝脏,肺和胃,7种癌组织包含膀胱尿路上皮癌(BLCA),乳腺癌(BRCA),结直肠癌(CRC),肝细胞癌(HCC),胃癌(GC),肾透明细胞癌(KCA)和前列腺腺癌(PRAD)。
本发明用来自NCBI SRA的数据构建验证集合来评估ASJA的效能,包括脑胶质瘤样本SRR934794SRR934744和SRR934930,以及poly(A)-/RNaseR和 poly(A)-处理样品RR2976716SRR2976715。
本发明还提出了一种识别全部RNA剪切位点的分析系统,系统包括:
连接点识别模块,其通过ASJA算法识别全部的剪接事件的连接点;
能效评估模块,其用于评估ASJA的效能。
实施例
1.ASJA流程概述
ASJA软件包是基于参考基因组的拼接(线性)和嵌合基因组(反式剪接、融合)的比对的结果检测三种不同形式的剪接,并对其进行注释和整合(图1)。 ASJA的工作流程如下:1.ASJA利用STAR和StringTie的优势得到嵌合基因组比对的结果以及拼接后的转录本;2.ASJA基于拼接后的转录本检测线性剪接事件,反式剪接和融合连接点的检测是基于嵌合基因组比对的结果;3.ASJA不仅对不同类型的剪接事件进行定量和注释,还可以根据剪接事件之间的关系对其进行文件整合。
2.对ASJA效能的评估
2.1评估线性剪接位点
本发明引用了已发表文章中定义金标准的方法。对样本进行有参考基因组注释的一轮STAR比对,把满足线性剪接位点的表达量(读段深度)大于1同时其对应转录本表达(FPKM)大于1的剪接位点视为金标准(正常肠癌样本01,总计:20,618)。ASJA检测到已知线性剪接连接点的灵敏度为97.3%。对于新的线性连接,将两轮无参考基因组注释的STAR结果中的已知的线性剪接位点与金标准进行比较,灵敏度为89.8%。不仅如此,以相同的金标准为参考,MapSplice2灵敏度为91.5%(1,750/20,618),其灵敏度要低于ASJA。
2.2评估反式剪接位点
对于反式剪接位点(环状RNA)的检测,本发明将ASJA与另外两种算法进行了比较,circRNA_finder和ACFS。本发明使用了来自12个正常组织(GSE77661) 的RNAseq数据集。ASJA在三种工具中观察到相同环状RNA的比例较高 (75.5%),而且ASJA检测到了最多的环状RNA(图2A)。为了评估假阳性率,本发明使用了RNase R(用于验证circRNA)处理的RNA-seq和相应的 poly(A)-RNA-seq数据集(GSE75733)探讨环状RNA的假阳性。本发明观察到的ASJA假阳性率(31.2%)也与其他算法包括ACFS(43%)和circRNA_finder(31.5%) 的结果相似(图2B)。
2.3评估融合剪接位点
为了评估ASJA融合剪接位点的准确性,本发明将ASJA与其他两个检测融合基因的软件进行了比较,包括MapSplice2和deFuse。本发明使用GBM的三个样本作为阳性集合,分别为SRR934794、SRR934744和SRR934930。表1的结果表明,尽管测试软件对融合基因都有较高的召回率,但ASJA的精确度明显高于其他软件。
2.4运行时间
本发明将ASJA与其他软件的运行时间进行比较。ASJA在运行速度方面的表现要优于其他的软件,而且其速度比其他方法快2-10倍。统计的运行时间都是从FASTQ开始计时,直到其产生结果为止。
3.使用ASJA检测到正常和癌症组织的剪接连接点
对近期发布的RNA-seq数据(GSE77661)使用ASJA软件包检测其中的剪接事件,并对其进行量化。本发明总共检测到了唯一的233,675个线性连接点, 81,484个反式剪接连接点和33个融合连接点(table2)。平均每个样本中检测到 165,997个线性连接点和5,668个反式剪接连接点和1个融合基因。线性连接点的数目要远大于反式剪接和融合连接点。除了脑组织占有4%,大部分组织中反式剪接读段的覆盖度仅占剪接事件的1%。该结果与circRNA在有丝分裂组织中数目较少的概念一致。此外,在正常组织中几乎检测不到融合连接点,但是在癌症中存在一些融合连接,如BRCA。
3.1ASJA检测到的线性连接点的特征
本发明使用ASJA默认的程序参数得到线性连接点。在所有组织中共发现 322,675个不同的线性连接点,其中284,287个连接点至少含有两个读段数目且剪接率大于0.08(图3A)。将注释到已知100,7444个转录本的剪接位点根据基因组位置其分为三类,其中78.5%的线性连接点是位于蛋白编码基因,还有少部分是位于非编码基因和假基因区域(图3B)。而且,本发明还发现每个基因平均发生10.4个已知的线性剪接事件,240,453中的10,870个已知的线性剪接位点注释到编码基因的5’非编码区域。同时检测到43,834(15.6%)个全新的剪接位点(图3B)。虽然大多数全新的连接(33,480)与已知基因重叠,但其余10,353 个连接(23.6%)起源于基因间区域(图3C)。本发明的研究结果表明,许多新型连接点似乎在各种组织中特异性表达。脑组织中有7,107个新的线性连接点,数目是远高于其他组织的(图3D)。本发明还比较了癌组织和匹配的非癌组织 (NCT)的连接点表达谱,并鉴定了癌症中109个下调和765个上调的线性连接点(图3E)。
3.2ASJA检测到的反式剪接连接点的特征
在ASJA检测的26个样本中总共鉴定的31,346个circRNA,但是20,475个尚未在circBase中注释(图4A)。每个基因仅有2(中位数,反式剪接事件的数目排序是从1到72)个反式剪接事件的发生(图4B)。本发明使用GENCODE 注释文件研究了这些circRNA候选物的基因组起源。超过90%的circRNA由蛋白质编码外显子组成,而较小的部分与长链非编码RNA和已知转录本的反义区域(图3C)。ASJA通过估计5'末端或3'末端的反式剪接比率来量化每个circRNA 相对于其替代线性同种型的丰度。虽然在某些情况下不存在线性剪接产物,但这些位点的背拼接比变化很大。当使用严格的反式剪接比率和读段数目来剔除低质量的circRNA(平均反式剪接比率>0.15;log2(read_count))>-1)时,本发明得到404个高丰度circRNA(图4D)。
表1:不同软件检测胶质瘤样本中融合基因的准确性。
注:*样本编号中括号里面的数字表示该样本中真阳性融合基因的数目.n。
表2:每个样本中的连接数。
注:BLCA:膀胱尿路上皮癌,BRCA:乳腺癌,CRC:结直肠癌,HCC:肝癌,GC:胃癌,KCA:肾透明细胞癌,PRAD:前列腺腺癌。
本发明的保护内容不局限于以上实施例。在不背离发明构思的精神和范围下,本领域技术人员能够想到的变化和优点都被包括在本发明中,并且以所附的权利要求书为保护范围。
Claims (9)
1.一种识别全部RNA剪切位点的分析方法,其特征在于,包括以下步骤:
步骤A:ASJA算法识别全部的剪接事件的连接点,包括以下:
步骤A1:对RNA-seq数据进行比对以及产生拼接后的转录本;
步骤A2:提取所有的剪接连接点,包括线性剪切位点、反式剪切位点和融合剪切位点;
步骤A3:注释和整合所述剪接连接点;
步骤B:评估ASJA的效能。
2.如权利要求1所述的RNA剪切位点的分析方法,其特征在于,所述步骤A1中对RNA-seq数据进行比对包括:使用FastQC软件对RNA-seq数据进行质量控制和数据过滤;过滤后的读段使用STAR比对软件进行两轮基因组比对。
3.如权利要求2所述的RNA剪切位点的分析方法,其特征在于,所述步骤A1中产生拼接后的转录本包括:将STAR比对的结果输入到StringTie中得到用于线性剪切连接点检测的拼接后的转录本;已知的转录本注释文件被用于转录本拼接模型的参考文件,以此来指导装配过程,其结果得到已知的转录本和未在GENCODE中注释的转录本。
4.如权利要求1所述的RNA剪切位点的分析方法,其特征在于,所述步骤A2提取剪接连接点包括:线性剪接连接点、反式剪接连接点、融合连接点。
5.如权利要求4所述的RNA剪切位点的分析方法,其特征在于,提取所述线性剪接连接点包括:在比对和拼接结束后,使用perl脚本提取来自拼接后转录本的线性剪接连接点,同时用下式对连接点进行定量:
COV(AB)=min(∑cov(A),∑cov(B));
其中,cov(A)cov(B)表示每个外显子的覆盖度;
所述连接点的表达量需要根据已知注释的连接点的覆盖度进行标准化:
CPT(AB)=COV(AB)(*10e7)/TC;
其中,TC表示已知全部连接点的覆盖度。
6.如权利要求4所述的RNA剪切位点的分析方法,其特征在于,提取所述反式剪接连接点包括:使用perl脚本从嵌合比对的结果中提取反式剪接连接点;所述反式剪接连接点满足:
1.嵌合读段比对到同一条染色体和同一条链上;
2.3’末端和5’末端位点之间的距离小于3,000,000bp;
3.剪切方式符合GT/AG法则;
4.基因组中反式剪接连接的起始和结束位置不合理;比对到线粒体和未知染色体的读段要被剔除。
7.如权利要求4所述的RNA剪切位点的分析方法,其特征在于,提取所述融合连接点包括:从嵌合比对的结果中提取所述融合连接点,并应用以下步骤降低假阳性率:
1.比对到线粒体和其他未标记的染色体或未映射到重叠群的读段被过滤掉;
2.反式剪接的读段被剔除;
3.读段数目大与1的融合连接点被保留;ASJA算法同时计算了类似STARChip中的SpanningReads变量验证融合连接点的有效性。
8.如权利要求4所述的RNA剪切位点的分析方法,其特征在于,所述步骤A3包括:准备剪接连接点的注释文件、计算剪接率、整合分析三种类型的剪接事件和过滤剪接事件;其中,
所述准备剪接连接点的注释文件包括:ASJA算法提供的初级线性连接点注释文件、从UCSC table browser中下载的BED格式的外显子注释文件;
所述计算剪接率包括:根据每个线性连接点注释到的基因来计算该位置发生剪接的概率,对未注释的连接点,ASJA算法将计算该连接点在对应的转录本中的剪接概率;Weightratioi=CPTi/CPTm,其中CPTm表示连接点注释到基因的最大的剪接值;所述ASJA算法利用如下计算3’末端和5’末端的反式剪切率:
5′_ratio=(5′back_splicedread)/linerread;
3′_ratio=(3′back_splicedread)/linerread;
其中,5’back_splicedread代表5’末端的读段数目,linearread表示对应该5’端位置对应的线性连接的读段数目;
ASJA使用类似反式剪切率的方法来计算发生融合连接点的概率:
donorratio=donorread/sum(donorread,linerread);
acceptorratio=acceptorread/sum(acceptorread,linerread);
其中,linearread表示在线性剪接中对应donor、acceptor剪接点的读段数目;
所述整合分析三种类型的剪接事件包括:基于基因注释和剪接位点的信息,ASJA算法将三种不同的剪接事件整合到一个文件中,以便发现剪接事件之间的相互关系;该部分输出内容包括:线性剪接的识别号,基因名称,反式剪接识别号以及融合剪接识别号;
所述过滤剪接事件包括:ASJA算法通过对读段数目和剪接率的过滤得到可信度较高的剪接事件。
9.一种识别全部RNA剪切位点的分析系统,其特征在于,采用如权利要求1-8之任一项所述的识别全部RNA剪切位点的分析方法,所述系统包括:
连接点识别模块,其通过ASJA算法识别全部的剪接事件的连接点;
能效评估模块,其用于评估ASJA的效能。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910726790.4A CN110619926B (zh) | 2019-08-07 | 2019-08-07 | 一种识别全部rna剪切位点的分析方法及分析系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910726790.4A CN110619926B (zh) | 2019-08-07 | 2019-08-07 | 一种识别全部rna剪切位点的分析方法及分析系统 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110619926A true CN110619926A (zh) | 2019-12-27 |
CN110619926B CN110619926B (zh) | 2023-03-31 |
Family
ID=68921536
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910726790.4A Active CN110619926B (zh) | 2019-08-07 | 2019-08-07 | 一种识别全部rna剪切位点的分析方法及分析系统 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110619926B (zh) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111508563A (zh) * | 2020-05-22 | 2020-08-07 | 四川大学华西医院 | 一种长非编码rna的癌症相关可变剪接数据库系统 |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20150254397A1 (en) * | 2014-01-11 | 2015-09-10 | Cytognomix Inc | Method of Validating mRNA Splciing Mutations in Complete Transcriptomes |
CN109390037A (zh) * | 2018-10-08 | 2019-02-26 | 齐齐哈尔大学 | 基于SVM-AdaBoost的成熟miRNA全位点识别方法 |
CN110010201A (zh) * | 2019-04-16 | 2019-07-12 | 山东农业大学 | 一种rna选择性剪接位点识别方法及系统 |
-
2019
- 2019-08-07 CN CN201910726790.4A patent/CN110619926B/zh active Active
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20150254397A1 (en) * | 2014-01-11 | 2015-09-10 | Cytognomix Inc | Method of Validating mRNA Splciing Mutations in Complete Transcriptomes |
CN109390037A (zh) * | 2018-10-08 | 2019-02-26 | 齐齐哈尔大学 | 基于SVM-AdaBoost的成熟miRNA全位点识别方法 |
CN110010201A (zh) * | 2019-04-16 | 2019-07-12 | 山东农业大学 | 一种rna选择性剪接位点识别方法及系统 |
Non-Patent Citations (1)
Title |
---|
王端青: "《基于转录组测序数据计算识别RNA编辑位点和可变剪接事件》", 《万方学位论文数据库》 * |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111508563A (zh) * | 2020-05-22 | 2020-08-07 | 四川大学华西医院 | 一种长非编码rna的癌症相关可变剪接数据库系统 |
Also Published As
Publication number | Publication date |
---|---|
CN110619926B (zh) | 2023-03-31 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Lun et al. | csaw: a Bioconductor package for differential binding analysis of ChIP-seq data using sliding windows | |
CN112951418B (zh) | 基于液体活检的连锁区域甲基化评估方法和装置、终端设备及存储介质 | |
CN106599616B (zh) | 基于duplex-seq的超低频突变位点检测分析方法 | |
CN106021984A (zh) | 一种全外显子组测序数据分析系统 | |
CN106951731B (zh) | 一种大片段插入或缺失的预测方法及系统 | |
CN110021347B (zh) | 一种基于miRBase数据库的动物有参的miRNA数据分析方法 | |
CN113257350A (zh) | 基于液体活检的ctDNA突变程度分析方法和装置、ctDNA性能分析装置 | |
CN105930690A (zh) | 一种全外显子组测序数据分析方法 | |
CN111584006B (zh) | 基于机器学习策略的环形rna识别方法 | |
CN109686439A (zh) | 遗传病基因检测的数据分析方法、系统及存储介质 | |
CN113838533B (zh) | 一种癌症检测模型及其构建方法和试剂盒 | |
CN106295246A (zh) | 找到与肿瘤相关的lncRNA并预测其功能 | |
CN110910950A (zh) | 一种联合分析单细胞scRNA-seq和scATAC-seq的流程方法 | |
CN112599188B (zh) | 一种融合驱动基因单端锚定的dna融合断点注释方法 | |
CN115132274B (zh) | 循环无细胞dna转录因子结合位点的甲基化水平分析方法及装置 | |
CN110021346A (zh) | 基于RNAseq数据的基因融合与突变检测方法及系统 | |
CN110619926B (zh) | 一种识别全部rna剪切位点的分析方法及分析系统 | |
CN109859797A (zh) | 一种基于miRBase数据库的无参的miRNA数据分析方法 | |
Miller et al. | Quality-controlled R-loop meta-analysis reveals the characteristics of R-loop consensus regions | |
CN113096737B (zh) | 一种用于对病原体类型进行自动分析的方法及系统 | |
CN105483210A (zh) | 一种rna编辑位点的检测方法 | |
CN117275585A (zh) | 基于lp-wgs和dna甲基化的肺癌早筛模型构建方法及电子设备 | |
CN110164504B (zh) | 二代测序数据的处理方法、装置及电子设备 | |
CN111192632A (zh) | 整合dna和rna的深度测序数据提取基因融合免疫治疗新抗原的方法和装置 | |
CN116230109A (zh) | 一种基于深度学习的手性分离预测方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |