CN110021352B - 一种基于miRBase数据库的植物有参的miRNA数据分析方法 - Google Patents
一种基于miRBase数据库的植物有参的miRNA数据分析方法 Download PDFInfo
- Publication number
- CN110021352B CN110021352B CN201811639517.XA CN201811639517A CN110021352B CN 110021352 B CN110021352 B CN 110021352B CN 201811639517 A CN201811639517 A CN 201811639517A CN 110021352 B CN110021352 B CN 110021352B
- Authority
- CN
- China
- Prior art keywords
- mirna
- analysis
- sequence
- data
- species
- 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.)
- Active
Links
Images
Classifications
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Landscapes
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明公开了一种基于miRBase数据库的植物有参的miRNA数据分析方法,其特征在于,包括如下步骤:文件准备步骤;下机数据过滤步骤;sRNA分类注释步骤;miRNA差异分析步骤;miRNA功能和通路分析步骤;miRNA序列特征分析步骤;结果整理步骤。本发明的有益效果在于:针对植物miRNA的特点采用合适的分析方法。结果全面,包含涉及到的miRNA分析内容以及其他测到的小RNA信息注释。自动整理所有分析结果,完成各个部分分析之后,自动对结果进行统计,可视化,以及归类整理,使结果排布一目了然,直接用于报告生成。所有操作步骤可见,方便错误查询。
Description
技术领域
本发明涉及转录组测序领域,具体涉及一种在miRBase数据库中有参考数据的植物miRNA测序的数据分析方法。
背景技术
miRNA是一类由内源基因编码非编码单链RNA分子,在动植物中参与转录后基因表达调控。多数miRNA以单拷贝、多拷贝或基因簇的形式存在于基因组中。miRNA在很多物种中被广泛发现,且在进化进程中高度保守,因此研究miRNA的确切功能、目的靶基因、以及其作用机制,是转录组学数据分析中的重要一环,对于了解生物体内基因的表达调控机制有重要意义。
miRNA的作用机制在动物和植物之间存在明显差异,且有的物种有丰富的miRNA参考数据,但有的物种缺乏参考数据,甚至有些物种没有参考基因组信息,这些情况下的miRNA测序的数据分析方法十分不同。目前还没有针对植物小RNA分析的工具。植物miRNA和动物miRNA在生物体内的作用机制不同,保守性程度也不同。目前还没有现成的流程分析植物小RNA测序数据;尤其是没有自动化的分析平台实现植物小RNA测序结果的流程化分析工具,包括后续的sRNA注释,miRNA序列的特征分析,表达量分析和差异分析,靶基因位点分析,等各个步骤的自动化整合。
发明内容
为了克服现有技术所存在的上述缺陷,本发明的目的在于提供一种基于miRBase数据库的植物有参的miRNA数据分析方法。
为了实现本发明的目的之一,所采用的技术方案是:
一种基于miRBase数据库的植物有参的miRNA数据分析方法,包括如下步骤:
步骤一、文件准备步骤:准备并读取config文件,软件读取相关信息后,会生成进行以下列出的所有分析步骤对应的shell脚本,按顺序运行即可,在运行同时每一步都会有运行日志,方便结果检查;
步骤二、下机数据过滤步骤:
下机后的原始数据,去除接头,然后过滤低质量序列,即:以5个碱基长度为窗口对原始序列进行搜索,当窗口中碱基的平均测序质量低于20时,将从窗口最前端开始的部分截断并舍弃。将过滤后的数据进行去重,获得无重复的序列,并标记所有序列数量并统计,过滤序列用于后续分析;
步骤三:sRNA分类注释步骤:
将去重后的序列与Rfam数据库进行blast比对,筛选出碱基错配数小于2的结果,注释出其中的非编码RNA序列,
将其余的小RNA序列与miRBase数据库中该物种的miRNA成熟体序列进行比对,筛选出碱基错配数小于2的结果,注释为已知的miRNA序列,同时计算测到的miRNA表达量,进行表达模式分析;
步骤四、miRNA差异分析步骤:
根据上一步注释到的miRNA信息以及表达量结果,使用DESeq进行差异表达分析,并按照差异倍数(FoldChange>2)和显著性(Pvalue<0.05)筛选差异表达的miRNA并绘制图像;
步骤五、miRNA功能和通路分析步骤:
以目标物种的mRNA序列为目标序列,使用psRNATarget或者psRobot软件对差异表达的miRNA序列,进行靶基因位点搜索;
对上一步预测到的miRNA靶基因进行GO功能和KEGG通路的富集分析,获得差异miRNA可能参与的功能和代谢通路;
步骤六、miRNA序列特征分析步骤:
miRNA碱基偏好性分析;
将该物种成熟miRNA序列与近缘物种进行blast比对,筛选出物种间保守的miRNA,并标记其相似度;
对检测到的已知miRNA进行家族归类,并查找相应miRNA家族在其他物种中的存在情况;
步骤七、结果整理步骤:
将所有用于生成miRNA结题报告的统计分析结果进行整理。
在本发明的一个优选实施例中,所述文件准备步骤中,所述的文件包括:下机数据位置以及对应的样本名和分组名、用于差异分析的分组、分析结果保存路径、任务名称、物种简称、测序接头序列、该物种名miRNA的成熟体序列、基因组序列及其index文件的位置、用于功能注释的基因注释文件,mRNA序列、GTF文件中的任意一种或多种。
在本发明的一个优选实施例中,所述下机数据过滤步骤当中,所述的统计为同时对原始数据和过滤数据量进行统计,并以柱状图展示不同长度的序列的数量分布特征。
在本发明的一个优选实施例中,所述sRNA分类注释步骤当中,还包括新的miRNA预测:使用mapper.pl将剩余的序列与基因组进行比对,并使用mireap.pl对比对上的序列进行新的miRNA预测,并使用RNAfold获得结构信息,最后对所有的小RNA序列的注释结果进行统计。
在本发明的一个优选实施例中,所述miRNA差异分析步骤当中,所述绘制图像包括:采用R语言的ggplot2软件包绘制差异表达miRNA的火山图(直观了解差异miRNA的分布情况)和、或MA图(评估文库标准化的好坏)、采用Pheatmap包对差异表达miRNA的表达量绘制热图。
在本发明的一个优选实施例中,所述miRNA碱基偏好性分析为不同长度的miRNA的首位碱基的偏好性,以及所有miRNA每个位置上的碱基偏好性。
本发明的主要创新点在于:
针对植物miRNA的特点采用合适的分析方法。
结果全面,包含涉及到的miRNA分析内容以及其他测到的小RNA信息注释。
自动整理所有分析结果,完成各个部分分析之后,自动对结果进行统计,可视化,以及归类整理,使结果排布一目了然,直接用于报告生成。
所有操作步骤可见,方便错误查询,在进行每一步分析时,都会记录所用到的命令行和参数,以及运行中产生的日志结果,一旦程序运行出错,可以快速检查错误。
附图说明
图1为本发明的流程示意图。
图2为运行日志示意图。
图3为本发明的MA示意图。
图4为本发明的火山图示意图。
图5为本发明的psRobot结果示意图。
图6为本发明的blast结果示意图。
图7为结果目录示意图。
图8为原始数据以及相关统计图表示意图。
图9为去冗余序列以及相关统计图表示意图。
图10为各种小RNA注释结果示意图。
图11为miRNA特征分析结果示意图。
图12为miRNA表达量相关分析结果示意图。
图13为miRNA靶基因预测相关分析结果示意图。
具体实施方式
以下通过实施例对本发明作进一步的说明,但这些实施例不得用于解释对本发明的限制。
在步骤S1)中接受用户的小RNA测序数据,以及相关的数据库信息,然后对所有的数据进行相关的分析,得到每个样本中所有小RNA的注释信息,并对miRNA进行序列特征分析和表达量分析,以及样本间差异表达分析,功能和通路富集分析。
首先是对下机数据进行过滤和数量统计。本发明实施例中,对下机数据进行去除接头和低质量序列的过滤处理,得到高质量的测序结果。作为示例地,采用perl语言脚本去除接头序列(filter_data.py脚本),并通过5bp的滑动窗口,对原始序列进行搜索,当窗口中碱基的平均测序质量低于20时,将从窗口最前端开始的部分截断并舍弃过滤低质量序列。然后过滤掉长度小于18或者大于36bp的序列。然后对高质量数据的重复序列进行归纳,得到所有的无冗余序列。并对原始数据和高质量进行数量统计。
接下来先通过比对注释出ncRNA序列。作为示例的,使用Blast将这些序列与Rfam数据库比对,注释其他如rRNA,tRNA,snRNA,snoRNA等非编码RNA信息。然后使用perl脚本对结果筛选出碱基错配数小于2的结果,注释出其中的非编码RNA序列。
然后注释出miRNA序列。作为示例的,将其余的小RNA序列与miRBase数据库中该物种的miRNA成熟体序列进行Blast比对,筛选出碱基错配数小于2的结果,注释为已知的miRNA序列,同时计算测到的miRNA表达量,进行表达模式分析。
然后从剩余的序列预测新的miRNA信息。作为示例的,使用mapper.pl将剩余的序列与基因组进行比对,并使用mireap.pl对比对上的序列进行新的miRNA预测,并使用RNAfold获得结构信息。最后对所有的小RNA序列的注释结果进行统计。
对于之前检测到的保守miRNA序列根据其表达量,进行差异表达分析。作为示例的,使用DESeq进行差异表达分析,并按照差异倍数(FoldChange>2)和显著性(Pvalue<0.05)筛选差异表达的miRNA。同时采用R语言的ggplot2软件包绘制差异表达miRNA的火山图(直观了解差异miRNA的分布情况)和MA图(评估文库标准化的好坏)。采用Pheatmap包对差异表达miRNA的表达量绘制热图。
根据序列相似性,对筛选到的显著差异表达的miRNA进行靶基因预测。作为示例的,以本物种的mRNA的3’UTR序列为目标序列,使用psRNATarget或者psRob软件对差异表达的miRNA序列,进行靶基因位点搜索。然后使用R语言通过超几何检验计算靶基因富集到哪些GO功能和KEGG代谢通路上,从而了解这些差异miRNA所发挥的功能。
psRobot结果参见图5,为miRNA:
vvi-miR156a与mRNA:VIT_03s0097g00240.t01的关系预测结果以及位点详情。
blast结果示例参见图6,第一列为本物种miRNA,第二列为与之最匹配的其他物种的miRNA。
本发明还对预测到的保守的miRNA序列进行序列特征分析,包括碱基偏好性分析,保守性分析和家族分析。
作为示例的,采用perl脚本,先对不同长度的miRNA序列,分别统计第一位碱基的种类分布数量;以及所有miRNA每个位置上的碱基种类分布数量,并使用R语言画图展示结果。然后将该物种的miRNA序列与近缘物种进行比对,找出物种间存在的保守性miRNA,并标记之间的相似度。根据每个miRNA的家族信息,找出在近缘物种中是否包含对应家族的miRNA信息。
最终整理所有的分析结果,所所有分析内容按类别排放在不同的目录下。作为示例的,将原始数据单独存放;将数据过滤的统计结果,序列长度分布图形单独存放;将所有小RNA的注释结果及注释结果统计都单独存放;将miRNA序列特征分析结果单独存放;将miRNA表达量以及差异表达相关的分析内容单独存放;将差异表达的miRNA对应的靶基因预测结果,以及功能和通路富集分析结果单独存放,结果目录排布参见7-13。
Claims (6)
1.一种基于miRBase数据库的植物有参的miRNA数据分析方法,其特征在于,包括如下步骤:
步骤一、文件准备步骤:准备并读取config文件,软件读取相关信息后,会生成进行以下列出的所有分析步骤对应的shell脚本,按顺序运行即可,在运行同时每一步都会有运行日志,方便结果检查;
步骤二、下机数据过滤步骤:
下机后的原始数据,去除接头,然后过滤低质量序列,即:以5个碱基长度为窗口对原始序列进行搜索,当窗口中碱基的平均测序质量低于20时,将从窗口最前端开始的部分截断并舍弃;
将过滤后的数据进行去重,获得无重复的序列,并标记所有序列数量并统计,过滤序列用于后续分析;
步骤三:sRNA分类注释步骤:
将去重后的序列与Rfam数据库进行blast比对,筛选出碱基错配数小于2的结果,注释出其中的非编码RNA序列,
将其余的小RNA序列与miRBase数据库中所述其余的小RNA序列的物种的miRNA成熟体序列进行比对,筛选出碱基错配数小于2的结果,注释为已知的miRNA序列,同时计算测到的miRNA表达量,进行表达模式分析;
步骤四、miRNA差异分析步骤:
根据上一步注释到的miRNA信息以及表达量结果,使用DESeq进行差异表达分析,并按照差异倍数:FoldChange>2和显著性:Pvalue<0.05筛选差异表达的miRNA并绘制图像;
步骤五、miRNA功能和通路分析步骤:
以目标物种的mRNA序列为目标序列,使用psRNATarget或者psRobot软件对差异表达的miRNA序列,进行靶基因位点搜索;
对上一步预测到的miRNA靶基因进行GO功能和KEGG通路的富集分析,获得差异miRNA参与的功能和代谢通路;
步骤六、miRNA序列特征分析步骤:
miRNA碱基偏好性分析;
将该物种成熟miRNA序列与近缘物种进行blast比对,筛选出物种间保守的miRNA,并标记其相似度;
对检测到的已知miRNA进行家族归类,并查找相应miRNA家族在其他物种中的存在情况;
步骤七、结果整理步骤:
将所有用于生成miRNA结题报告的统计分析结果进行整理。
2.如权利要求1所述的一种基于miRBase数据库的植物有参的miRNA数据分析方法,其特征在于,所述文件准备步骤中,所述的文件包括:下机数据位置以及对应的样本名和分组名、用于差异分析的分组、分析结果保存路径、任务名称、物种简称、测序接头序列、该物种miRNA的成熟体序列、基因组序列及其index文件的位置、用于功能注释的基因注释文件,mRNA序列、GTF文件中的任意一种或多种。
3.如权利要求1所述的一种基于miRBase数据库的植物有参的miRNA数据分析方法,其特征在于,所述下机数据过滤步骤当中,所述的统计为同时对原始数据和过滤数据量进行统计,并以柱状图展示不同长度的序列的数量分布特征。
4.如权利要求1所述的一种基于miRBase数据库的植物有参的miRNA数据分析方法,其特征在于,所述sRNA分类注释步骤当中,还包括新的miRNA预测:使用mapper.pl将剩余的序列与基因组进行比对,并使用mireap.pl对比对上的序列进行新的miRNA预测,并使用RNAfold获得结构信息,最后对所有的小RNA序列的注释结果进行统计。
5.如权利要求1所述的一种基于miRBase数据库的植物有参的miRNA数据分析方法,其特征在于,所述miRNA差异分析步骤当中,所述绘制图像包括:采用R语言的ggplot2软件包绘制差异表达miRNA的火山图和/或MA图、采用Pheatmap包对差异表达miRNA的表达量绘制热图。
6.如权利要求1所述的一种基于miRBase数据库的植物有参的miRNA数据分析方法,其特征在于,所述miRNA碱基偏好性分析为不同长度的miRNA的首位碱基的偏好性,以及所有miRNA每个位置上的碱基偏好性。
Applications Claiming Priority (2)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811159378 | 2018-09-30 | ||
CN2018111593780 | 2018-09-30 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110021352A CN110021352A (zh) | 2019-07-16 |
CN110021352B true CN110021352B (zh) | 2023-06-16 |
Family
ID=65849002
Family Applications (3)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201811394656.0A Active CN109545281B (zh) | 2018-09-30 | 2018-11-20 | 一种基于二代高通量测序的trio家系遗传突变模式的分析方法 |
CN201811644819.6A Active CN110021347B (zh) | 2018-09-30 | 2018-12-29 | 一种基于miRBase数据库的动物有参的miRNA数据分析方法 |
CN201811639517.XA Active CN110021352B (zh) | 2018-09-30 | 2018-12-29 | 一种基于miRBase数据库的植物有参的miRNA数据分析方法 |
Family Applications Before (2)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201811394656.0A Active CN109545281B (zh) | 2018-09-30 | 2018-11-20 | 一种基于二代高通量测序的trio家系遗传突变模式的分析方法 |
CN201811644819.6A Active CN110021347B (zh) | 2018-09-30 | 2018-12-29 | 一种基于miRBase数据库的动物有参的miRNA数据分析方法 |
Country Status (1)
Country | Link |
---|---|
CN (3) | CN109545281B (zh) |
Families Citing this family (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110648722B (zh) * | 2019-09-19 | 2022-05-31 | 首都医科大学附属北京儿童医院 | 新生儿遗传病患病风险评估的装置 |
CN111354418B (zh) * | 2020-01-19 | 2023-02-10 | 上海欧易生物医学科技有限公司 | 基于参考基因组注释文件的高通量测序技术动物tRFs数据分析方法 |
CN111653316A (zh) * | 2020-05-27 | 2020-09-11 | 上海寻因生物科技有限公司 | 基于二代测序的可视化分析方法、系统及存储介质 |
CN112908412A (zh) * | 2021-02-10 | 2021-06-04 | 北京贝瑞和康生物技术有限公司 | 用于复合杂合变异致病证据适用性的方法、设备和介质 |
CN112967756B (zh) * | 2021-03-30 | 2022-07-26 | 上海欧易生物医学科技有限公司 | 基于snakemake语言快速批量可自动邮件反馈结果的高通量测序质控分析方法 |
Family Cites Families (11)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102222175A (zh) * | 2011-05-06 | 2011-10-19 | 西南大学 | 一种与miRNA相关的功能性分子标记的开发方法 |
CN105349617A (zh) * | 2014-08-19 | 2016-02-24 | 复旦大学 | 一种对高通量rna测序数据的质量控制方法及装置 |
CN105528532B (zh) * | 2014-09-30 | 2019-08-16 | 深圳华大基因科技有限公司 | 一种rna编辑位点的特征分析方法 |
CN105112569B (zh) * | 2015-09-14 | 2017-11-21 | 中国医学科学院病原生物学研究所 | 基于宏基因组学的病毒感染检测及鉴定方法 |
CN105925685A (zh) * | 2016-05-13 | 2016-09-07 | 万康源(天津)基因科技有限公司 | 一种基于家系的外显子组潜在致病变异检测方法 |
CN106119353A (zh) * | 2016-06-25 | 2016-11-16 | 广州泰因生物科技有限公司 | 一种显性家系遗传病致病位点的快速筛查方法 |
GB201615837D0 (en) * | 2016-09-16 | 2016-11-02 | Genomics Plc | Methods and apparatus for identifying one or more genetic variants associated with disease in an individual or group of related individuals |
CN107103205A (zh) * | 2017-05-27 | 2017-08-29 | 湖北普罗金科技有限公司 | 一种基于蛋白质质谱数据注释真核生物基因组的生物信息学方法 |
CN107577919A (zh) * | 2017-08-21 | 2018-01-12 | 上海派森诺生物科技股份有限公司 | 一种基于高通量测序技术的宏基因组数据分析方法 |
CN107828857A (zh) * | 2017-11-23 | 2018-03-23 | 南宁科城汇信息科技有限公司 | 一种转录组测序及RNAseq数据分析方法 |
CN108220403B (zh) * | 2017-12-26 | 2021-07-06 | 北京科迅生物技术有限公司 | 特定突变位点的检测方法、检测装置、存储介质及处理器 |
-
2018
- 2018-11-20 CN CN201811394656.0A patent/CN109545281B/zh active Active
- 2018-12-29 CN CN201811644819.6A patent/CN110021347B/zh active Active
- 2018-12-29 CN CN201811639517.XA patent/CN110021352B/zh active Active
Also Published As
Publication number | Publication date |
---|---|
CN110021347B (zh) | 2023-08-08 |
CN110021347A (zh) | 2019-07-16 |
CN109545281B (zh) | 2023-06-20 |
CN109545281A (zh) | 2019-03-29 |
CN110021352A (zh) | 2019-07-16 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN110021352B (zh) | 一种基于miRBase数据库的植物有参的miRNA数据分析方法 | |
CN109859797B (zh) | 一种基于miRBase数据库的无参的miRNA数据分析方法 | |
CN106446600B (zh) | 一种基于CRISPR/Cas9的sgRNA的设计方法 | |
Fei et al. | Tomato Functional Genomics Database: a comprehensive resource and analysis package for tomato functional genomics | |
CN111354418B (zh) | 基于参考基因组注释文件的高通量测序技术动物tRFs数据分析方法 | |
CN106951731B (zh) | 一种大片段插入或缺失的预测方法及系统 | |
CN106021985B (zh) | 一种基因组数据压缩方法 | |
CN113066532B (zh) | 基于高通量测序技术的宿主中病毒来源sRNA数据分析方法 | |
CN115631789A (zh) | 一种基于泛基因组的群体联合变异检测方法 | |
Stokes et al. | Transcriptomics for clinical and experimental biology research: hang on a seq | |
CN114107444A (zh) | 一种发掘调控植物发育关键调控因子的方法及其应用 | |
Grones et al. | Best practices for the execution, analysis, and data storage of plant single-cell/nucleus transcriptomics | |
KR101506916B1 (ko) | miRNA 탐색 자동화 시스템을 이용하여 시료로부터 miRNA를 자동으로 동정하는 방법 | |
WO2012155296A1 (zh) | 获得基因组大小和误差的方法 | |
CN112233722A (zh) | 品种鉴定的方法、其预测模型的构建方法和装置 | |
CN109754844B (zh) | 一种在全基因组水平上预测植物内源siRNA的方法 | |
Kuang et al. | Machine Learning Approaches for Plant miRNA Prediction: Challenges, Advancements, and Future Directions | |
Liu et al. | Gene coexpression and evolutionary conservation analysis of the human preimplantation embryos | |
CN110684830A (zh) | 一种石蜡切片组织rna分析方法 | |
US20190172553A1 (en) | Using k-mers for rapid quality control of sequencing data without alignment | |
CN112102880A (zh) | 品种鉴定的方法、其预测模型的构建方法和装置 | |
CN110867212A (zh) | 一种猪的品种溯源方法及装置 | |
Salas Zamora | Comparative analysis of piRNA precursor annotation across different primates | |
CN116590435B (zh) | 一种与猪背膘厚相关的因果候选基因及其鉴定方法和应用 | |
US20230093253A1 (en) | Automatically identifying failure sources in nucleotide sequencing from base-call-error patterns |
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 |