CN111243665A - 一种核糖体印记测序数据分析方法及系统 - Google Patents
一种核糖体印记测序数据分析方法及系统 Download PDFInfo
- Publication number
- CN111243665A CN111243665A CN202010014759.0A CN202010014759A CN111243665A CN 111243665 A CN111243665 A CN 111243665A CN 202010014759 A CN202010014759 A CN 202010014759A CN 111243665 A CN111243665 A CN 111243665A
- Authority
- CN
- China
- Prior art keywords
- analysis
- ribosome
- sequencing data
- reads
- data
- 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.)
- Pending
Links
- 238000004458 analytical method Methods 0.000 title claims abstract description 81
- 238000012163 sequencing technique Methods 0.000 title claims abstract description 74
- 108090000623 proteins and genes Proteins 0.000 claims abstract description 88
- 230000014509 gene expression Effects 0.000 claims abstract description 52
- 238000013519 translation Methods 0.000 claims abstract description 47
- 238000000034 method Methods 0.000 claims abstract description 39
- 238000012216 screening Methods 0.000 claims abstract description 22
- 230000000007 visual effect Effects 0.000 claims abstract description 3
- 108020004705 Codon Proteins 0.000 claims description 25
- 238000012545 processing Methods 0.000 claims description 23
- 108091032973 (ribonucleotides)n+m Proteins 0.000 claims description 21
- 210000003705 ribosome Anatomy 0.000 claims description 16
- 108020004566 Transfer RNA Proteins 0.000 claims description 12
- 238000001914 filtration Methods 0.000 claims description 12
- 108091070501 miRNA Proteins 0.000 claims description 12
- 238000013518 transcription Methods 0.000 claims description 12
- 230000035897 transcription Effects 0.000 claims description 12
- 239000002679 microRNA Substances 0.000 claims description 10
- 102000042773 Small Nucleolar RNA Human genes 0.000 claims description 8
- 108020003224 Small Nucleolar RNA Proteins 0.000 claims description 8
- 238000010586 diagram Methods 0.000 claims description 7
- 238000010201 enrichment analysis Methods 0.000 claims description 7
- 238000010835 comparative analysis Methods 0.000 claims description 6
- 230000037361 pathway Effects 0.000 claims description 6
- 238000003559 RNA-seq method Methods 0.000 claims description 5
- 238000012098 association analyses Methods 0.000 claims description 5
- 238000010219 correlation analysis Methods 0.000 claims description 5
- 238000007619 statistical method Methods 0.000 claims description 5
- 238000007621 cluster analysis Methods 0.000 claims description 4
- 238000001303 quality assessment method Methods 0.000 claims description 3
- 238000013441 quality evaluation Methods 0.000 claims description 3
- 238000003672 processing method Methods 0.000 claims description 2
- 238000005516 engineering process Methods 0.000 abstract description 5
- 238000003908 quality control method Methods 0.000 abstract description 5
- 238000007405 data analysis Methods 0.000 abstract description 4
- 238000012165 high-throughput sequencing Methods 0.000 abstract description 3
- 238000010276 construction Methods 0.000 abstract description 2
- 230000009897 systematic effect Effects 0.000 abstract description 2
- 230000014616 translation Effects 0.000 description 42
- 239000012634 fragment Substances 0.000 description 10
- 102000039471 Small Nuclear RNA Human genes 0.000 description 7
- 108020004688 Small Nuclear RNA Proteins 0.000 description 7
- 108091081024 Start codon Proteins 0.000 description 6
- 238000011144 upstream manufacturing Methods 0.000 description 6
- 108090000765 processed proteins & peptides Proteins 0.000 description 5
- 108020004418 ribosomal RNA Proteins 0.000 description 5
- 108091026890 Coding region Proteins 0.000 description 4
- 150000001413 amino acids Chemical class 0.000 description 3
- 230000008859 change Effects 0.000 description 3
- 230000007547 defect Effects 0.000 description 3
- 230000008569 process Effects 0.000 description 3
- 102000004169 proteins and genes Human genes 0.000 description 3
- 108020005345 3' Untranslated Regions Proteins 0.000 description 2
- 108091026898 Leader sequence (mRNA) Proteins 0.000 description 2
- 108091030146 MiRBase Proteins 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 2
- 108020004999 messenger RNA Proteins 0.000 description 2
- 238000001556 precipitation Methods 0.000 description 2
- -1 snorRNA Proteins 0.000 description 2
- 239000000126 substance Substances 0.000 description 2
- 238000011282 treatment Methods 0.000 description 2
- 102000006382 Ribonucleases Human genes 0.000 description 1
- 108010083644 Ribonucleases Proteins 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000004422 calculation algorithm Methods 0.000 description 1
- 210000004027 cell Anatomy 0.000 description 1
- 239000003153 chemical reaction reagent Substances 0.000 description 1
- 150000001875 compounds Chemical class 0.000 description 1
- 238000004590 computer program Methods 0.000 description 1
- 238000011109 contamination Methods 0.000 description 1
- 230000000593 degrading effect Effects 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 230000006870 function Effects 0.000 description 1
- 238000007417 hierarchical cluster analysis Methods 0.000 description 1
- 230000002401 inhibitory effect Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000001124 posttranscriptional effect Effects 0.000 description 1
- 239000002243 precursor Substances 0.000 description 1
- 238000001243 protein synthesis Methods 0.000 description 1
- 230000001105 regulatory effect Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 210000004708 ribosome subunit Anatomy 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
- 230000014621 translational initiation Effects 0.000 description 1
Images
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
- G16B25/00—ICT specially adapted for hybridisation; ICT specially adapted for gene or protein expression
- G16B25/10—Gene or protein expression profiling; Expression-ratio estimation or normalisation
-
- 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
- G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
-
- 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
- G16B45/00—ICT specially adapted for bioinformatics-related data visualisation, e.g. displaying of maps or networks
Landscapes
- Life Sciences & Earth Sciences (AREA)
- Physics & Mathematics (AREA)
- Health & Medical Sciences (AREA)
- Engineering & Computer Science (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Medical Informatics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Biotechnology (AREA)
- Evolutionary Biology (AREA)
- General Health & Medical Sciences (AREA)
- Theoretical Computer Science (AREA)
- Biophysics (AREA)
- Genetics & Genomics (AREA)
- Data Mining & Analysis (AREA)
- Molecular Biology (AREA)
- Chemical & Material Sciences (AREA)
- Analytical Chemistry (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Artificial Intelligence (AREA)
- Bioethics (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Databases & Information Systems (AREA)
- Epidemiology (AREA)
- Evolutionary Computation (AREA)
- Public Health (AREA)
- Software Systems (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明涉及高通量测序和生物信息技术分析领域，特别是涉及翻译组和核糖体印记测序数据分析方法及系统。本发明的核糖体印记测序数据分析方法，包括：将测序数据进行筛选；将数据筛选结果进行比对，保留目标reads；将所述比对的结果进行RFs分布分析；将获得的目标reads进行编码基因的统计；将获得的目标reads的编码基因进行表达差异分析。本发明分析方法更为系统、详细和全面，从下机数据的质量控制到可视化分析以及表达差异分析的一套全面分析方法，既可以评估核糖体印记建库数据特征是否满足要求，判断其是否准确，又可以进行后续更丰富的数据生物信息分析内容，挖掘实验数据生物学意义。
Description
技术领域
本发明涉及高通量测序和生物信息技术分析领域，特别是涉及翻译组和核糖体印记测序数据分析方法及系统。
背景技术
核糖体印记测序(Ribosome profiling sequencing，简称Ribo-seq)是一种基于高通量测序的检测全基因组水平RNA翻译的技术。Ribo-seq也是目前翻译组学(研究RNA到蛋白的翻译过程的组学)研究的主流方法，其具体方法是用低浓度RNase处理核糖体-新生肽链复合物，降解掉无核糖体覆盖的RNA片段，再去除核糖体，最后利用二代测序技术检测被核糖体保护的约～30bp的正在翻译RNA小片段。这些被核糖体保护的RNA片段，准确指示了核糖体正在进行翻译的“足迹”/“印记”，因此这些被核糖体保护的RNA片段，又被称为核糖体足迹(ribosome footprints，简称RFs)。
发明专利CN 201910407961.7一种从低质量核糖体印记数据预测基因编码框的方法和系统仅公开了如何从multitaper算法和复杂度上对低质量核糖体数据的分析。现有技术对此项测序数据生物信息分析的方法不够全面、对数据质控不够全面、不能评估核糖体建库特征是否满足需求、且未有与其他组学关联分析的方法内容。
发明内容
鉴于此，有必要针对上述问题提供一种核糖体印记测序数据分析方法及系统，以解决现有技术中分析不够全面具体的问题。
本发明是通过以下技术方案实现的：
一种核糖体印记测序数据的处理方法，包括如下操作：
将测序数据进行筛选；将数据筛选结果进行比对，保留目标reads；将所述比对的结果进行RFs分布分析；将获得的目标reads进行编码基因的统计；将获得的目标reads的编码基因进行表达差异分析。
进一步的，所述测序数据的筛选包括测序数据统计、过滤和质量评估。
优选的，测序数据统计是指对每个样本进行总测序数据量统计；由于Ribo-seq测序的目标RNA片段约为30bp，单端reads就已经可以将RNA片段测通，因此仅保留双末端测序的一端reads用于后序分析。同时，为了保证数据质量，需要在信息分析前对下机数据进行质控与过滤。本发明对下机后经过初步过滤得到的reads碱基总数进行进一步更严格的过滤，得到clean reads，用于后续的信息分析，将过滤后的reads结果以饼状图或柱状图的形式可视化展示。
优选的，所述过滤的包括：
1)去除含有接头的reads(截去接头及接头后面的部分)；
2)去除全部都是A碱基的reads；
3)去除含N比例大于10％的reads；
4)去除低质量的reads(质量值Q≤20的碱基数占整条reads的50％以上)。
优选的，所述质量评估操作包括：将各类过滤的reads以及高质量reads占总reads的比例通过图示化(饼状图)展示，以及对过滤后的数据通过绘制碱基组成和质量值分布图，来直观地看数据质量情况。
进一步的，对测序数据进行筛选后，还可包括以下操作：将数据筛选结果可视化展示。
进一步的，所述“将数据筛选结果进行比对，保留目标reads”的操作包括：将数据筛选结果进行比对，通过reads长度分布统计，保留长度为20bp～40bp的reads。
优选的，所述“将数据筛选结果进行比对，保留目标reads”包括如下操作：
比对去除核糖体RNA：本发明使用reads比对工具bowtie2将Clean Reads比对到该物种的核糖体，去除比对上核糖体RNA的reads，保留下来的数据用于后续进一步的分析。
比对去除转运RNA：本发明通过比对数据库来去除样本数据中的tRNA。具体操作为：将比对去除核糖体RNA后的reads通过blast比对到GenBank及Rfam数据库，尽可能地发现并去除样本中的tRNA。
比对去除snoRNA、snRNA、miRNA：通过比对数据库来去除样本数据中的snoRNA、snRNA、miRNA。具体为：将比对去除核糖体RNA和tRNA后的reads通过blast比对到GenBank及Rfam数据库，尽可能地发现并去除样本中的snoRNA和snRNA。然后比对miRBase中该物种的miRNA序列，去除比对上的miRNA序列。
比对参考基因组，得到核糖体印记。
本发明使用短reads比对软件Bowtie2将以上比对去除了核糖体RNA、tRNA、snoRNA、snRNA、miRNA，长度过滤后得到的reads比对到参考基因序列上；能够比对到参考基因组的符合预期长度的reads即为核糖体印记(ribosome footprints，RFs)。
进一步的，所述“将所述比对结果进行RFs分布分析”包括如下操作：
RFs在编辑基因上的分布：本发明方法根据RFs在编码基因上的比对位置，本发明将RFs分为四类：CDS、5’UTR、3’UTR、Intron。一般来说RFs多分布在CDS区，在UTR区则数量较少；将RFs在编码基因上的位置分布情况进行统计并做饼图；
RFs在起始和终止密码子周边的分布：本发明方法根据RFs 5’端在基因组的比对位置，统计位于编码基因CDS起始密码子和终止密码子周边的RFs的数量。
一般而言，起始密码子上游或终止密码子下游不会被翻译或翻译丰度极低，因此对应的RFs信号整体弱于编码区。理论上RFs 5’端比对位置开始于起始密码子上游12～13np，停止于终止密码子上游18nt的位置。得到实验各样本所有RFs在起始和终止密码子周边的分布丰度图(图2)。
RFs比对密码子位置的分布：由于核糖体在转录本滑动翻译蛋白的过程中，会每隔3个碱基(1个密码子)产生一个停顿，完成一个氨基酸的肽段延伸。
本发明将比对到CDS区的RFs，按照RFs 5’端比对位置对应的密码子位置归为三类(密码第1～3个碱基)；然后计算每个基因中三类RFs比例，并绘制柱状图(图3)。由于核糖体在密码子第一碱基位置停留时间最长，因此RFs比对位置对应密码子第一个碱基的比例通常最高。
进一步的，所述“将获得的目标reads进行编码基因的统计”包括但不限于基因表达结果统计、基因覆盖度统计、基因表达量统计。
优选的，所述基因覆盖度指每个基因被reads覆盖的百分比，其值等于基因中被比对的reads覆盖的碱基数跟基因编码区所有碱基数的比值。
优选的，所述基因表达量统计方法具体为：使用Rsem软件计算编码基因ORF区内Ribo-seq水平的reads数，并换算为FPKM值，从而得到基因在翻译水平的表达量。
进一步的，所述“将获得的目标reads的编码基因进行表达差异分析“包括但不限于分组间差异翻译统计分析、表达模式聚类分析、差异基因GO/Pathway富集分析。
优选的，使用edgeR软件对组间基因进行差异翻译分析，利用FDR与log2FC来筛选差异翻译基因(differentially translated gene,DTG)，筛选条件为FDR<0.05且|log2FC|>1。
优选的，表达模式聚类分析为基于基因表达量，对样本和基因间的关系进行层级聚类，并使用热图来呈现聚类结果。
进一步的，若进行Ribo-seq的样品，同时进行了RNA-seq或lncRNA-seq，进行了以上所述的分析处理外，还进一步进行两组学关联分析，其分析包括但不限于：基因翻译表达量与转录表达量相关性分析、翻译差异和转录差异比较分析、翻译效率分析(translational efficiency，TE)、TE差异与转录差异比较分析。
本发明还公开了一种核糖体印记测序数据分析的系统，所述的系统具体包括：
测序数据筛选单元，用于筛选原始测序数据；
对比分析单元，用于对筛选所得的测序数据进行比对分析，以找出目标reads；
RFs分布分析单元，用于将获得的比对分析数据结果进行RFs分布分析；
编码基因统计单元；所述编码基因统计单元实现的统计任务包括：基因表达结果统计、基因覆盖度统计、基因表达量统计；
表达差异分析单元；所述表达差异分析单元实现的份吸任务包括：分组间差异翻译统计分析、表达模式聚类分析、差异基因GO/Pathway富集分析。
进一步的，所述系统还可以包括：两组学关联分析单元。
本发明还公开了一种核糖体印记测序数据处理平台，所述平台包括：
处理器、存储器以及核糖体印记测序数据处理控制程序；
其中在所述的处理器执行所述的核糖体印记测序数据处理平台控制程序，所述的核糖体印记测序数据处理平台控制程序被存储在所述的存储器中，所述的核糖体印记测序数据处理平台控制程序，实现上述的核糖体印记测序数据处理方法步骤。
本发明还公开了一种计算机可读取存储介质，所述计算机可读取存储介质存储有核糖体印记测序数据处理平台控制程序，所述的核糖体印记测序数据处理平台控制程序，实现上述的核糖体印记测序数据处理方法步骤。
本发明有益效果：
本发明的核糖体印记测序数据分析方法及系统通过生物信息技术，解决现有技术中核糖体印记测序数据分析方法不全面的缺点，可以准确、系统、全面的对核糖体印记测序数据进行评估和生物信息分析。同时，可以对核糖体印记数据特征进行可视化展示。
本发明分析方法更为系统、详细和全面，从下机数据的质量控制到可视化分析以及表达差异分析的一套全面分析方法，既可以评估核糖体印记建库数据特征是否满足要求，判断其是否准确，又可以进行后续更丰富的数据生物信息分析内容，挖掘实验数据生物学意义。
附图说明
图1是本发明一种核糖体印记测序数据分析方法的流程图。
图2是本发明实施例1中的核糖体印记测序数据分析方法流程图。
图3本发明实施例1中样本所有RFs在起始和终止密码子周边的分布丰度。
图4是本发明实施例1中样本RFs比对密码子位置的分布图。
图5是本发明实施例2中基因翻译表达量与转录本表达量相关性分析散点图。
图6是本发明实施例2中转录组与翻译组差异的比较散点图。
图7是本发明实施例2中组TE丰度分布密度图。
图8是本发明实施例2中组TE分布盒型图。
具体实施方式
为了更好地说明本发明所解决的问题、所采用的技术方案和所达到的效果，现结合具体实施例和相关资料进一步阐述。需要说明的是，本发明内容包含但不限于以下实施例及其组合实施方式。
本发明实施例中未注明具体技术或条件的，按照本领域内的文献所描述的技术或条件或者按照产品说明书进行。所用试剂或仪器未注明生产厂商者，均为可以通过市购等途径获得的常规产品。
经分析发现，本发明的核糖体印记测序数据分析方法及系统(图1)，可克服目前核糖体印记测序数据分析方法不全面的缺点，可以准确、系统、全面的对核糖体印记测序数据进行评估和生物信息分析。
实施例1
本发明提出一种核糖体印记测序数据分析方法(图2)，具体内容如下：
步骤S1：对测序数据进行筛选，包括统计、过滤和质量评估，并可视化展示筛选结果；
优选的，测序数据统计是指对每个样本进行总测序数据量统计；由于Ribo-seq测序的目标RNA片段约为30bp，单端reads就已经可以将RNA片段测通，因此仅保留双末端测序的一端reads用于后序分析。同时，为了保证数据质量，需要在信息分析前对下机数据进行质控与过滤。本发明对下机后经过初步过滤得到的reads碱基总数进行进一步更严格的过滤，得到clean reads，用于后续的信息分析，将过滤后的reads结果以饼状图或柱状图的形式可视化展示。
在一些实施例中，所述过滤的步骤如下：
1)去除含有接头的reads(截去接头及接头后面的部分)；
2)去除全部都是A碱基的reads；
3)去除含N比例大于10％的reads；
4)去除低质量的reads(质量值Q≤20的碱基数占整条reads的50％以上)；
在一些实施例中，所述质量评估操作包括：将各类过滤的reads以及高质量reads占总reads的比例通过图示化(饼状图)展示，以及对过滤后的数据通过绘制碱基组成和质量值分布图，来直观地看数据质量情况。
步骤S2：将S1得到的数据结果进行比对分析，通过reads长度分布统计，保留长度为20bp～40bp的reads用于RFs分布分析(步骤S3分析所用数据)，包括如下步骤：
步骤S2.1：比对去除核糖体RNA：
受样本质量和物种的影响，实验方法去除核糖体RNA的效率可能不太稳定，而核糖体RNA的污染会影响后续的分析。因此本发明使用reads比对工具bowtie2将Clean Reads比对到该物种的核糖体，去除比对上核糖体RNA的reads，保留下来的数据用于后续进一步的分析；
步骤S2.2：比对去除转运RNA：
tRNA(转运RNA)在翻译过程中负责识别mRNA上的密码子并转运相应的氨基酸，在翻译过程中起着重要的作用，是翻译起始复合物的一个重要组成部分。tRNA的种类和浓度高度影响着蛋白质合成的种类和翻译速度。由于沉淀核糖体-新生肽链复合物的过程中也有很大可能把tRNA一起沉降下来，因此本发明通过比对数据库来去除样本数据中的tRNA。将比对去除核糖体RNA后的reads(步骤S2.1所得结果)通过blast比对到GenBank及Rfam数据库，尽可能地发现并去除样本中的tRNA；
步骤S2.3：比对去除snoRNA,snRNA,miRNA：
细胞中的snoRNA和snRNA，由于参与了基因转录后加工过程的RNA剪接、rRNA前体加工及核糖体亚基的组装，所以也有可能在沉淀核糖体-新生肽链复合物的过程中一起沉降下来。另外，miRNA可以结合mRNA，通过降解或抑制mRNA的翻译来调控基因的表达。因此本发明通过比对数据库来去除样本数据中的snoRNA,snRNA,miRNA。将比对去除核糖体RNA和tRNA后的reads(步骤S2.2所得结果)通过blast比对到GenBank及Rfam数据库，尽可能地发现并去除样本中的snoRNA和snRNA。然后比对miRBase中该物种的miRNA序列，去除比对上的miRNA序列；
步骤S2.4：比对参考基因组，得到核糖体印记：本发明使用短reads比对软件Bowtie2将以上比对去除了核糖体RNA、tRNA、snoRNA、snRNA、miRNA，长度过滤后得到的reads(步骤S2.3所得结果)比对到参考基因序列上。能够比对到参考基因组的符合预期长度的reads即为核糖体印记(ribosome footprints，RFs)。
步骤S3：将S2步骤获得的比对分析数据结果进行RFs分布分析，包括如下内容：
S3.1：RFs在编辑基因上的分布：本发明方法根据RFs在编码基因上的比对位置，本发明将RFs分为四类：CDS,5’UTR,3’UTR,Intron。一般来说RFs多分布在CDS区，在UTR区则数量较少。将RFs在编码基因上的位置分布情况进行统计并做饼图；
S3.2：RFs在起始和终止密码子周边的分布：本发明方法根据RFs 5’端在基因组的比对位置，统计位于编码基因CDS起始密码子和终止密码子周边的RFs的数量。一般而言，起始密码子上游或终止密码子下游不会被翻译或翻译丰度极低，因此对应的RFs信号整体弱于编码区。理论上RFs 5’端比对位置开始于起始密码子上游12～13np，停止于终止密码子上游18nt的位置。得到实验各样本所有RFs在起始和终止密码子周边的分布丰度图(图3)。
S3.3：RFs比对密码子位置的分布：由于核糖体在转录本滑动翻译蛋白的过程中，会每隔3个碱基(1个密码子)产生一个停顿，完成一个氨基酸的肽段延伸。本发明将比对到CDS区的RFs，按照RFs 5’端比对位置对应的密码子位置归为三类(密码第1～3个碱基)。然后计算每个基因中三类RFs比例，并绘制柱状图(图4)。由于核糖体在密码子第一碱基位置停留时间最长，因此RFs比对位置对应密码子第一个碱基的比例通常最高。
步骤S4：进一步进行编码基因的统计，包括但不限于基因表达结果统计、基因覆盖度统计、基因表达量统计。
在一些实施例中，所述基因覆盖度指每个基因被reads覆盖的百分比，其值等于基因中被比对的reads覆盖的碱基数跟基因编码区所有碱基数的比值。
在一些实施例中，基因表达量统计方法具体为：使用Rsem软件计算编码基因ORF区内Ribo-seq水平的reads数，并换算为FPKM值，从而得到基因在翻译水平的表达量。
在一些实施例中，基因表达量的计算使用FPKM(Fragments Per Kilobase oftranscript per Million mapped reads)法，其计算公式为：
设FPKM(A)为基因A的表达量，则C为比对到基因A的测序片段数，N为比对到参考基因的总测序片段数，L为基因A的碱基数。FPKM法能消除基因长度和测序量差异对计算基因表达的影响，计算得到的基因表达量可直接用于比较不同样品间的基因表达差异。
步骤S5：表达差异分析，包括但不限于分组间差异翻译统计、表达模式聚类分析、差异基因GO/Pathway富集分析。
在一些实施例中，使用edgeR软件对组间基因进行差异翻译分析，利用FDR与log2FC来筛选差异翻译基因(differentially translated gene,DTG)，筛选条件为FDR<0.05且|log2FC|>1。
在一些实施例中，表达模式聚类分析为基于基因表达量，对样本和基因间的关系进行层级聚类，并使用热图来呈现聚类结果。
实施例2
6例实验样本，分为3组(标识A/B/C)，进行Ribo-seq和RNA-seq，进行S1至S5的分析外，还可以进一步进行步骤S6两组学关联分析，包括但不限于：基因翻译表达量与转录本表达量相关性分析、翻译差异和转录差异比较分析、翻译效率分析(translationalefficiency，TE)、TE差异与转录差异比较分析。
S6.1：基因翻译表达量与转录本表达量相关性分析，为计算组内的基因翻译表达量与转录本丰度的皮尔森相关系数，并绘制散点图，可分析翻译组和转录组两组学之间的相关性高低(图5)。
S6.2：翻译差异和转录差异比较分析，具体为根据基因在两个组学中的变化规律，将基因分为5类，对应分别为：1)Transcription：仅在转录组中差异显著的基因2)Translation：仅在翻译组中差异显著的基因3)Homodirection:在两个组学都有显著差异，且上下调方向相同4)Opposite:在两个组学都有显著差异，且上下调方向相反5)Unchanged：在两个组学中都没有显著差异我们统计各类基因的数量，并绘制散点图(图6)。散点图图中X轴为转录水平差异倍数的log2值，Y轴为翻译水平的差异倍数的log2值。从图可以得出处理组间转录水平的变化与翻译水平的变化是否一致。进一步可进行各分类基因的GO/Pathway功能富集分析。
S6.3：翻译效率分析代表样本中某个基因的总RNA分子(通常指mRNA)与核糖体结合并进行翻译的比例。利用Ribo-seq和RNA-seq的数据，计算TE，其计算公式为：TE＝(FPKMin Ribo-seq/(FPKM in RNA-seq)，并将结果通过组TE丰度分布密度图(图7)和组TE分布盒型图(图8)可视化展示。
通过计算翻译效率TE，可以进一步进行更为丰富的生物信息分析，本案例中主要进行翻译效率与转录丰度相关性分析、差异翻译效率分析、差异翻译效率基因的聚类分析、以及对TE差异显著的基因集分别进行GO和KEGG富集分析。
在一些实施例中，翻译效率与转录丰度相关性分析，为计算组内的基因翻译效率与转录丰度的皮尔森相关系数，并绘制散点图，可分析基因翻译效率与转录丰度之间的相关性高低
在一些实施例中，差异翻译效率分析，为使用Ribodiff软件对组间基因翻译效率进行差异分析，利用FDR与log2FC来筛选差异翻译效率基因DTEG，筛选条件为FDR<0.05且|log2FC|>1。
在一些实施例中，差异翻译效率基因的聚类分析，为基于每个比较组合都会得到的一个差异TE基因集，将所有比较组合的差异基因集的并集在每个实验组中的的TE值，用于层次聚类分析图。
实施例3
一种核糖体印记测序数据分析的系统，包括但不限于计算机可读介质或云平台，所述计算机可读介质或云平台存储有能进行核糖体印记测序数据分析的计算机程序，能够实现以上方法的步骤。
综上所述，本发明一种核糖体印记测序数据分析方法及系统，通过生物信息方法结合，进行核糖体印记测序数据的分析，解决现有技术中可克服目前核糖体印记测序数据分析方法不全面的缺点，可以准确、系统、全面的对核糖体印记测序数据进行评估和生物信息分析。
以上所述实施例仅表达了本发明的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对本发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进，这些都属于本发明的保护范围。因此，本发明专利的保护范围应以所附权利要求为准。
Claims (11)
1.一种核糖体印记测序数据的处理方法，其特征在于，包括如下操作：
对测序数据进行筛选，并可视化展示结果；
将数据筛选结果进行比对，保留目标reads；
将所述比对的结果进行RFs分布分析；
将获得的目标reads进行编码基因的统计；
将获得的目标reads编码基因进行表达差异分析。
2.根据权利要求1所述的核糖体印记测序数据的处理方法，其特征在于，所述“对测序数据进行筛选”包括测序数据的过滤和质量评估；
所述过滤包括：去除含有接头的reads、去除全部都是A碱基的reads、去除含N比例大于10％的reads和去除低质量的reads；
所述质量评估包括：将各类过滤的reads以及高质量reads占总reads的比例通过图示化展示，以及对过滤后的数据通过绘制碱基组成和质量值分布图，来直观地看数据质量情况。
3.根据权利要求1所述的核糖体印记测序数据的处理方法，其特征在于，所述“将数据筛选结果进行比对，保留目标reads”的操作中，所述比对包括：比对去除核糖体RNA、比对去除转运RNA、比对去除snoRNA、snRNA、miRNA，最后保留比对至参考基因组的reads。
4.根据权利要求1所述的核糖体印记测序数据的处理方法，其特征在于，所述“将所述比对的结果进行RFs分布分析”的操作中，所述RFs分布分析包括：核糖体印记在编辑基因上的分布分析、核糖体印记在起始和终止密码子周边的分布分析、核糖体印记比对密码子位置的分布分析。
5.根据权利要求1所述的核糖体印记测序数据的处理方法，其特征在于，所述的“将获得的目标reads进行编码基因的统计”包括：基因表达结果统计、基因覆盖度统计和基因表达量统计。
6.根据权利要求1所述的核糖体印记测序数据的处理方法，其特征在于，所述的“将获得的目标reads编码基因进行表达差异分析”包括：分组间差异翻译统计、表达模式聚类分析和差异基因GO/Pathway富集分析。
7.根据权利要求1所述的核糖体印记测序数据的处理方法，其特征在于，所述核糖体印记测序数据的处理方法还包括：两组学关联分析；
所述两组学关联分析包括：结合样本RNA-seq或lncRNA-seq，进行了如权利要求1所述的分析后，再进行基因翻译表达量与转录表达量相关性分析、翻译差异和转录差异比较分析、翻译效率分析和TE差异与转录差异比较分析。
8.一种核糖体印记测序数据处理系统，其特征在于，所述的系统具体包括：
测序数据筛选单元，用于筛选原始测序数据；
对比分析单元，用于对筛选所得的测序数据进行比对分析，以找出目标reads；
RFs分布分析单元，用于将获得的比对分析数据结果进行RFs分布分析；
编码基因统计单元；所述编码基因统计单元实现的统计任务包括：基因表达结果统计、基因覆盖度统计、基因表达量统计；
表达差异分析单元；所述表达差异分析单元实现的份吸任务包括：分组间差异翻译统计分析、表达模式聚类分析、差异基因GO/Pathway富集分析。
9.根据权利要求8所述的核糖体印记测序数据处理系统，其特征在于，所述系统还包括：两组学关联分析单元。
10.一种核糖体印记测序数据处理平台，其特征在于，包括：
处理器、存储器以及核糖体印记测序数据处理控制程序；
其中在所述的处理器执行所述的核糖体印记测序数据处理平台控制程序，所述的核糖体印记测序数据处理平台控制程序被存储在所述的存储器中，所述的核糖体印记测序数据处理平台控制程序，实现如权利要求1至9中任一项所述的核糖体印记测序数据处理方法步骤。
11.一种计算机可读取存储介质，其特征在于，所述计算机可读取存储介质存储有核糖体印记测序数据处理平台控制程序，所述的核糖体印记测序数据处理平台控制程序，实现如权利要求1至7中任一项所述的核糖体印记测序数据处理方法步骤。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010014759.0A CN111243665A (zh) | 2020-01-07 | 2020-01-07 | 一种核糖体印记测序数据分析方法及系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010014759.0A CN111243665A (zh) | 2020-01-07 | 2020-01-07 | 一种核糖体印记测序数据分析方法及系统 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN111243665A true CN111243665A (zh) | 2020-06-05 |
Family
ID=70866079
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010014759.0A Pending CN111243665A (zh) | 2020-01-07 | 2020-01-07 | 一种核糖体印记测序数据分析方法及系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN111243665A (zh) |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111635935A (zh) * | 2020-06-10 | 2020-09-08 | 华中农业大学 | 一种低成本、高通量间接量化蛋白表达水平的方法及其应用 |
| CN112908407A (zh) * | 2021-02-02 | 2021-06-04 | 北京大学 | 一种用tRNA组学来质控蛋白生物合成体系的方法 |
| CN118038991A (zh) * | 2024-04-12 | 2024-05-14 | 宁波甬恒瑶瑶智能科技有限公司 | 基因序列处理方法、系统、电子设备及存储介质 |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20160237487A1 (en) * | 2015-02-10 | 2016-08-18 | The Texas A&M University System | Modeling and Predicting Differential Alternative Splicing Events and Applications Thereof |
| CN108624651A (zh) * | 2018-05-14 | 2018-10-09 | 深圳承启生物科技有限公司 | 一种构建Ribo-seq测序文库的方法 |
| US20180334679A1 (en) * | 2017-05-21 | 2018-11-22 | Ramot At Tel-Aviv University Ltd. | Methods for modifying the growth rate of a cell |
| CN110136776A (zh) * | 2019-05-15 | 2019-08-16 | 深圳大学 | 一种从低质量核糖体印迹数据预测基因编码框的方法和系统 |
| CN110317771A (zh) * | 2019-05-07 | 2019-10-11 | 深圳大学 | 一种高质量水稻核糖体印记文库的构建方法 |
| CN110556163A (zh) * | 2019-09-04 | 2019-12-10 | 广州基迪奥生物科技有限公司 | 一种基于翻译组的长链非编码rna翻译小肽的分析方法 |
-
2020
- 2020-01-07 CN CN202010014759.0A patent/CN111243665A/zh active Pending
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20160237487A1 (en) * | 2015-02-10 | 2016-08-18 | The Texas A&M University System | Modeling and Predicting Differential Alternative Splicing Events and Applications Thereof |
| US20180334679A1 (en) * | 2017-05-21 | 2018-11-22 | Ramot At Tel-Aviv University Ltd. | Methods for modifying the growth rate of a cell |
| CN108624651A (zh) * | 2018-05-14 | 2018-10-09 | 深圳承启生物科技有限公司 | 一种构建Ribo-seq测序文库的方法 |
| CN110317771A (zh) * | 2019-05-07 | 2019-10-11 | 深圳大学 | 一种高质量水稻核糖体印记文库的构建方法 |
| CN110136776A (zh) * | 2019-05-15 | 2019-08-16 | 深圳大学 | 一种从低质量核糖体印迹数据预测基因编码框的方法和系统 |
| CN110556163A (zh) * | 2019-09-04 | 2019-12-10 | 广州基迪奥生物科技有限公司 | 一种基于翻译组的长链非编码rna翻译小肽的分析方法 |
Non-Patent Citations (2)
| Title |
|---|
| HSIN-YEN LARRY WU，等: "The Tomato Translational Landscape Revealed by Transcriptome Assembly and Ribosome Profiling", 《PLANT PHYSIOLOGY》, vol. 181, no. 1, pages 367 - 380, XP093047348, DOI: 10.1104/pp.19.00541 * |
| SONIA CHOTHANI，等: "deltaTE: Detection of TranslationallyRegulated Genes by Integrative Analysis of Ribo-seq and RNA-seq Data", 《CURRENT PROTOCOLS IN MOLECULAR BIOLOGY》, vol. 129, pages 1 - 22 * |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111635935A (zh) * | 2020-06-10 | 2020-09-08 | 华中农业大学 | 一种低成本、高通量间接量化蛋白表达水平的方法及其应用 |
| CN112908407A (zh) * | 2021-02-02 | 2021-06-04 | 北京大学 | 一种用tRNA组学来质控蛋白生物合成体系的方法 |
| CN118038991A (zh) * | 2024-04-12 | 2024-05-14 | 宁波甬恒瑶瑶智能科技有限公司 | 基因序列处理方法、系统、电子设备及存储介质 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN111243665A (zh) | 一种核糖体印记测序数据分析方法及系统 | |
| Müller et al. | Capturing the dynamics of genome replication on individual ultra-long nanopore sequence reads | |
| CN111261229B (zh) | 一种MeRIP-seq高通量测序数据的生物分析流程 | |
| CN111354418B (zh) | 基于参考基因组注释文件的高通量测序技术动物tRFs数据分析方法 | |
| CN114657238B (zh) | 一种枸杞40k液相芯片及应用 | |
| CN110021347B (zh) | 一种基于miRBase数据库的动物有参的miRNA数据分析方法 | |
| CN112133368A (zh) | 一种基于三代测序技术的宏基因组测序数据的自动化分析方法 | |
| CN110544509B (zh) | 一种单细胞ATAC-seq数据分析方法 | |
| CN112270953A (zh) | 基于bd单细胞转录组测序数据的分析方法、装置及设备 | |
| CN111192637B (zh) | 一种lncRNA鉴定和表达定量的分析方法 | |
| CN112349346A (zh) | 检测基因组区域中的结构变异的方法 | |
| CN115101128B (zh) | 一种杂交捕获探针脱靶危险性评估的方法 | |
| CN113066532A (zh) | 基于高通量测序技术的宿主中病毒来源sRNA数据分析方法 | |
| CN108192893B (zh) | 基于转录组测序开发艾纳香ssr引物的方法 | |
| CN108165620A (zh) | 标签及其制备方法和应用 | |
| CN110970091B (zh) | 标签质控的方法及装置 | |
| CN114708910A (zh) | 一种利用单细胞测序数据计算池测序中细胞亚群富集分数的方法 | |
| CN109461473B (zh) | 胎儿游离dna浓度获取方法和装置 | |
| CN108728515A (zh) | 一种使用duplex方法检测ctDNA低频突变的文库构建和测序数据的分析方法 | |
| CN111192636B (zh) | 一种适用于oligodT富集的mRNA二代测序结果分析方法 | |
| KR102347463B1 (ko) | 핵산 서열 분석에서 위양성 변이를 검출하는 방법 및 장치 | |
| WO2023184330A1 (zh) | 基因组甲基化测序数据的处理方法、装置、设备和介质 | |
| Eché et al. | A Bos taurus sequencing methods benchmark for assembly, haplotyping, and variant calling | |
| CN114530200A (zh) | 基于计算snp熵值的混合样本鉴定方法 | |
| WO2022006443A1 (en) | Systems and methods for detecting cell-associated barcodes from single-cell partitions |
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 |