一种基因注释文件格式及针对该基因注释文件格式的解析工具

技术领域

本发明涉及一种基因注释文件格式，具体涉及一种基因注释文件格式及针对该基因注释文件格式的解析工具，属于生物信息领域。

背景技术

一个物种的基因组测序完成后，需要对这些数据进行解读，首先要先找到这些序列中转录起始/终止位点、基因、外显子、内含子等组成元件在染色体中的位置信息(即注释)后才能再进行深入的分析。GFF/GTF是目前贮存这些注释信息的两种文件格式，其中，GFF(general feature format)主要是用来存储基因组基因注释结果，并且现大部分利用的是第三版，即GFF3；GTF(gene transfer format)主要是用来存储基因组转录本注释结果。当前广泛使用的GTF格式为第二版，即GTF2。其中，GFF3允许使用#作为注释符号，除去注释行外，主体部分共有9列，并且GFF3中每一列的含义：序列编号(seqid)、来源(source)、类型(type)、起始位置(start)、终止位置(end)、得分(score)、链的方向(strand)、相位(phase)和属性(attributes)，其中

第一列：基因组序列编号

第二列：注释信息的来源

第三列：注释信息的类型，比如基因(gene)、信使RNA(mRNA)、外显子(exon)、编码序列(CDS)、非编码序列(UTR)等

第四列：第三列的注释类型在基因组序列上的起始位置

第五列：第三列的注释类型在基因组序列上的终止位置

第六列：得分，是注释信息可靠性的说明，“.”表示为空

第七列：该基因或转录本位于基因组序列的正链(+)或负链(-)上

第八列：相位，第三列为“CDS"时的起始编码的相对位置，有效值为0、1、2

第九列：包含众多注释信息，以多个键值对组成的注释信息描述，不同属性之间以分号相隔，如：

ctg123.mRNA 1300 9000.+.ID＝mrna0001；Name＝sonichedgehog

ctg123.exon 1300 1500.+.ID＝exon00001；Parent＝mrna0001

其中：

ID，注释信息的编号，在一个GFF文件中必须唯一。

Parent，指明ID所从属的父序列ID。用于表述外显子(exons)从属于转录本(transcript)的关系，转录本(transripts)从属于基因(gene)的关系。

GTF是GFF2的扩展，且GTF的前八列GTF与GFF相同，不同之处在于：

(1)GTF还包含可选字段:如5’端非编码区(5’UTR),3’端非编码区(3’UTR),基因间区(inter),基因间保守的非编码序列区(inter_CNS),以及内含子中的保守非编码序列区(intron_CNS)。

(2)GTF文件的第9列，以分号分隔，且每一组信息内部以空格分隔。且基因名和转录本名为gene_id和transcript_id。

从基因注释格式来说，随着基因组研究和三代测序PacBio(Wang,et al.,2016)和Nanopore(Garalde,et al.,2018)等技术的不断发展，基因的鉴定及后续研究更加深入，而基于以前基因研究设计的GFF3和GTF基因注释格式缺乏对某些基因或结构的表述。存在以下缺陷：(1)无法准确描述一些新发现的或复杂的基因类型及其转录本，比如原始转录本(primary transcript),融合基因(fusion gene),环状RNA(circleRNA)等；(2)对调控元件的表述单一，而在不同生物过程中调控元件是多样的，比如在原始转录本(primarytranscript)形成加工转录本(processed transcript)的过程中，编辑位点(editedsite),修饰位点(modified site),polyA信号序列(polyA signal sequence),弱化子(attenuator)可能行使调控功能(Roundtree,et al.,2017)；(3)缺乏相对坐标体系，即以序列自身为参考的坐标系统，这种特征对于真核生物多外显子序列的分析研究尤其重要；(4)缺乏基于序列相似性或基于功能的基因关系描述体系：相似性较高的序列可能在生物体内行使相同或相近的生物学功能，如由于基因复制导致的旁系同源基因彼此之间有较高的序列相似性。因此，基因注释格式从序列相似程度描述组织和继承关系(通过序列相似性程度，将序列信息聚集在一起，形成一个集合)有利于生物学功能的研究。

从解析基因注释工具来说，目前的解析工具功能较为单一，往往只有单独针对某种注释格式文件(GFF3或GTF)的写入和读取结构信息的功能，或者某两种格式的转化功能，或者是结合基因组序列文件和注释格式文件得到基因序列的功能。缺乏绝对位置(基于基因组序列的坐标)和相对位置(基于序列自身为参考的坐标)的转化功能，而该功能有助于对记录和直观分析非线性组合的序列，如多外显子形成的CDS的坐标的快速提取。

综上所述，更新目前的基因注释文件格式以及开发一种更为全面的解析工具来满足日益增加的分析研究需求是至关重要的。

发明内容

为了克服现有技术所存在的上述缺陷，本发明提供了一种新型的基因注释文件格式(GFF4)及针对该基因注释格式文件的解析工具，兼容了现有注释格式和新型注释格式的各项信息提取和转换。本发明的解析工具包含的各注释格式间的转化(GFF4、GTF、GFF3)、检索(如合并、提取和查询)、序列提取和位置转换多项功能，满足了生物信息下游分析的多样化需求。本发明的基因注释文件格式能更适应于目前生物领域遗传基因相关的更加多样、精细的研究。

为了实现本发明的目的，本发明的技术方案如下：

本发明提供了一种基因注释文件格式，所述文件格式共9列，分别为：序列编号、起始、终止、链的方向、层级、特征类型、编号、父序列编号和属性，其中：

按照序列合成的顺序，所述层级从上到下的序列性质为：基因座位->原始转录本->加工转录本->产物。

在本发明的一优选实施方式中，所述属性，按照不同的所述层级，有其特定的生物类型值以及调控元件生物类型值。

在本发明的一优选实施方式中，在所述基因注释文件格式中每一个层级下有共享的特征类型。

在本发明的一更优选实施方式中，当所述特征类型为调控元件时,在第九列的属性中包含调控元件生物类型的值。

在本发明的一优选实施方式中，通过所述属性中的特定属性来实现生物学功能上的逻辑组织，其序列相似程度从大到小为：产物->产物簇->基因座位名->基因。

在本发明的一优选实施方式中，所述基因注释文件格式增加了相似程度和交叉关系的分析结果。

本发明还提供了一种针对基因注释文件格式的解析工具，所述解析工具包括：

格式转化模块，用于GTF、GFF3和GFF4格式之间的相互转换；

查询模块，用于查询基因注释文件中目标信息；

注释合并模块，用于将多种注释格式、信息丰富度不同的注释文件有机地合并成一个完整的注释文件；

单元提取模块，用于提取目标序列的子/父序列方向的注释信息；

位置转化模块，用于相对位置与绝对位置间的相互转化；

序列提取模块，用于根据注释信息和对应的基因组序列提取相应的基因和调控元件序列。

在本发明的一优选实施方式中，所述查询模块用于查询第一至九列的信息。在本发明的一优选实施方式中，所述单元提取模块用于提取给定ID的子/父序列方向的注释信息，例如当给定一个加工转录本ID时，可提取形成该转录本的原始转录本和该转录本的产物信息。

本发明的基因注释文件格式(GFF4)增加了相似程度和交叉关系的分析结果，相对于现有的基因注释文件格式更加精细，同时也减少了下游分析的繁杂程度。

本发明的基因注释文件格式的解析工具不仅顺应了生命科学研究的发展，还综合了多种功能，具有更好的兼容性，使得下游的信息提取、分析更加便捷，同时有利于注释信息的及时补充更新。

附图说明

图1为本发明的各层级和生物类型间的关系(富含生物过程)示意图；

图2为本发明的解析工具的格式转化模块示意图；

图3为本发明的解析工具的查询模块示意图；

图4为本发明的解析工具的注释合并模块示意图；

图5为本发明的解析工具的单元提取模块示意图；

图6为本发明的解析工具的位置转化模块示意图；

图7为本发明的解析工具的序列提取模块示意图；

图8显示了本发明GFF4格式的示例；

图9为本发明的解析工具的示意图。

具体实施方式

1、本发明的基因注释格式(GFF4)：

(1)共有9列，分别为：第一列：序列编号(seqId)、第二列：起始(start)、第三列：终止(end)、第四列：链的方向(strand)、第五列：层级(level)、第六列：特征类型(featureType)、第七列：编号(id)、第八列：父序列编号(parentId)和第九列：属性(attributes)。

a)按照序列合成的顺序，层级(level)从上到下的序列性质为：基因座位(locus)->原始转录本(primary)->加工转录本(processed)->产物(product)。

b)属性(attributes)按照不同的层级(level)，会有其特定的生物类型值(Biotype value)以及调控元件生物类型值(regulator_biotype value)，具体每个层级下的特征，如下：

注：

在GFF4格式中每一个层级下会有一些共享的特征类型。

调控元件是起调控作用一段序列，当特征类型为调控元件时,在第九列的属性(attributes)中必须包含调控元件生物类型值。

层级由高向低依次为：基因座位(locus)->原始转录本(primary)->加工转录本(processed)->产物(product)。

表1各层级下调控元件的生物类型

1.1.属性(attributes)

1.1.1.原始转录本生物类型(primary_biotype),加工转录本生物类型(processed_biotype),产物生物类型(product_biotype)

原始转录本生物类型(primary_biotype),加工转录本生物类型(processed_biotype),产物生物类型(product_biotype)分别是层级(level)为原始转录本(primary),加工转录本(processed),产物(product)下的一个属性标签。

该属性标签值的类型：字符串

该属性标签是否必需:否

该属性标签是否允许多个值:是

表2原始转录本生物类型(primary_biotype),加工转录本生物类型(processed_biotype),产物生物类型(product_biotype)可选的属性值

1.1.2.产物簇(product_cluster)/基因座位名(gene)/基因(gene_name)

该属性标签值的类型：字符串

该属性标签是否必需:否

该属性标签是否允许多个值:否

1.1.3.序列来源(seqSource)

该属性标签值的类型：字符串，可选：核染色体(major)，质粒(plasmid)，B染色体(B)，线粒体(mitochondrion)，动基体(kinetoplast)，纺锤体(mitosome)，质体(plastid)，叶绿体(Chloroplast)，色质体(Chromoplast)，质体样细胞器(apicoplast)

该属性标签默认值：未知(unknown)

该属性标签是否必需:否

该属性标签是否允许多个值:否

1.1.4.环状(circular)

该序列是否成环或未知。

该属性标签值：否(FALSE)，是(TRUE)，不确定(ND)

该属性标签默认值：不确定(ND)

该属性标签是否必需:否

该属性标签是否允许多个值：否

1.1.5.相位(phase)

该属性标签值：0,1,2

该属性标签是否必需：是(对于编码序列(CDS))

该属性标签是否允许多个值：否

1.1.6.序列(sequence)

该属性标签值的类型：字符串

该属性标签是否必需：否

该属性标签默认值：“”

该属性标签是否允许多个值：否

1.1.7.移码(frameshift)

该属性标签值的类型：字符串

该属性标签默认值：“”

该属性标签是否必需：否

该属性标签是否允许多个值：否

标签值的格式为“X:Y”，X是位置，Y是偏移量，可能的偏移量包括：

1.-2:核糖体后退两个碱基(5端的方向)当移动到位置X时；

2.-1:核糖体后退一个碱基(5端的方向)当移动到位置X时；

3.+1:核糖体前进一个碱基(3端的方向)当移动到位置X时；；

4.+2:核糖体前进两个碱基(3端的方向)当移动到位置X时；；

1.1.8.特征编号(featureUid)

该属性标签值的类型：字符串

该属性标签默认值：序列编号_起始_终止_链的方向(seqId_start_end_strand)

该属性标签是否必需：否

该属性标签是否允许多个值：否

1.1.9.名字(name)

该属性标签值的类型：字符串

该属性标签默认值：“”

该属性标签是否必需：否

该属性标签是否允许多个值：否

1.1.10.别称(alias)

该属性标签值的类型：字符串

该属性标签默认值：“”

该属性标签是否必需：否

该属性标签是否允许多个值：是

1.1.11.注释来源(annotation Source)

该属性标签值的类型：字符串

该属性标签默认值：“”

该属性标签是否必需：否

该属性标签是否允许多个值：是

注释来源用于描述得到这个特征的算法或者处理流程，特别的可能是一些软件或者数据库的名称，比如"Genescan"，"Genbank"。

1.1.12.打分(score)

该属性标签值的类型：浮点型或整数型

该属性标签默认值：NA

该属性标签是否必需：否

该属性标签是否允许多个值：否

得分表明这行特征或坐标的可靠程度。

1.1.13.注解(note)

该属性标签值的类型：字符串

该属性标签默认值：“”

该属性标签是否必需：否

该属性标签是否允许多个值：是

1.1.14.序列类型(seqType)

该属性标签值：DNA,RNA,蛋白(Protein),不确定(ND)

该属性标签默认值：不确定(ND)

该属性标签是否必需：否

该属性标签是否允许多个值：否

1.1.15.基因座位类型(locus_type)

该属性标签值：无基因(gene_free),包含基因(gene_containing),不确定(ND)

该属性标签默认值：不确定(ND)

该属性标签是否必需：否

该属性标签是否允许多个值：否

描述最终的产物是否包含基因

1.1.16.基因数目(gene_number)

该属性标签值的类型：整数型

该属性标签默认值：NA

该属性标签是否必需：否

该属性标签是否允许多个值：否

在该基因座位(locus)下基因的数目

1.1.17.假基因化(pseudo)

该属性标签值：否(FALSE)，是(TRUE)，不确定(ND)

该属性标签默认值：不确定(ND)

该属性标签是否必需：否

该属性标签是否允许多个值：否

表明在突变之后这条序列是否经历假基因化，成为丧失功能的序列

1.1.18.表达量(expression)

该属性标签值的类型：浮点型，整数型

该属性标签是否必需：否

该属性标签是否允许多个值：否

(2)按照生物学过程的组织与继承：

编号(id)-->父序列编号(parentId)，箭头左侧承自箭头右侧，代表的是图8中第7列和第8列间的关系

(3)按照生物学功能的组织和继承：

通过属性(attributes)中特定属性来实现生物学功能上的逻辑组织，其序列相似程度从大到小依次为：产物(product)->产物簇(product_cluster)->基因座位名(gene)->基因(gene_name)。

即：不同的产物可能聚集成一个产物簇，不同的产物簇可能同属一个基因，不同基因座位(locus)的基因可能同属一个基因(包含了旁系同源基因)。而它们聚集的条件是相似度。

其中，第七列相同的编号(id)的各行代表组成同一产物的各个片段。

对于同一个基因座位(locus)中的产物(product),不同加工转录本(processed)序列来源的产物(product)可能在序列上完全一样，这些序列相同的产物(product)，在逻辑上被归为一个产物簇(product_cluster)(如图8中L01.t4.p1和L01.t1.p4这两个产物(product)，序列一模一样，因此具有相同的产物簇编号(product_cluster id)：L01.g1.pc3)；共享序列的片段往往具有相关的生物学功能(比如同一个蛋白编码基因的不同可变剪接形式)，因此还将同一基因座位(locus)上同种产物类型(product_biotype)中具有共享片段的不同产物簇(product_cluster)归为一个基因座位名(gene)(如图8中，L01.g4.pc1和L01.g4.pc2这2个产物簇(product_cluster)形成基因座位名(gene):L01.g4，但是由于L01.g8.pc1序列来源完全独立，因此为独立的另外一个基因座位名(gene)：L01.g8)。

不同基因座位(locus)的基因座位名(gene)，其序列可能高度相似(旁系同源基因)。因此，为序列高度相似(仅具有个别位点的变异)或已被实验验证的不同基因座位(locus)来源的基因座位名(gene)指定相同的基因(gene_name)。

当特征类型(featureType)为调控元件(regulator)时，能描述不同层级的下调控元件，继承方式同上。

根据编号(id)，父序列编号(parentId)以及属性(attributes)中加工转录本生物类型(processed_biotype)；产物生物类型(product_biotype)的关系，可以描述同一个原始转录本下生成不同的加工转录本最后形成不同功能的产物的情况(如图8中的原始转录本:L01.pt4)。

2、针对本发明的基因注释格式的解析工具功能丰富，包括以下模块：格式转化模块(form)、查询模块(grep)、注释合并模块(merge)、单元提取模块(unit)、位置转化模块(abs2rel和rel2abs)和序列提取模块(seq)。其中，

格式转化模块用于GTF、GFF3和GFF4格式之间的相互转换，兼容性强；

查询模块用于提取满足一定条件的注释信息，特别是对第九列属性的信息(用一般linux操作系统上的awk命令很难进行快速提取)，使检索更加快捷；

注释合并模块用于将多种注释格式、信息丰富度不同的注释文件有机地合并成一个完整的注释文件，有利于整合不同版本的注释信息，例如对于某一个基因组，已有一个基因注释版本，现新注释了某些基因，可以使用这个模块选择已有基因注释文件为参考，根据位置关系，将新注释的基因合并，提高了基因组信息的完整性；

单元提取模块用于提取某个产物、转录本等的子/父序列方向的注释信息，有利于对单个基因的全面研究，即便于了解一段序列所经历的详细的生物过程；

位置转化模块包含绝对位置转化相对位置模块(abs2rel)(绝对->相对)和相对位置转化绝对位置模块(rel2abs)(相对->绝对)，用于相对位置与绝对位置间的相互转化，有利于分析研究非线性组合的序列，如多外显子形成的CDS的坐标的快速提取。对基因突变方面的研究也同样有益，例如，从文章中知晓某个基因转录本的编码序列(CDS)上的第50个碱基发生了突变，可以通过位置转化模块，得到该突变碱基在基因组上的位置，便于提取、分析子序列的信息，如突变分析。其中，绝对位置指基于基因组上的位置；相对位置指基于序列自身为参考的位置。此外，通过在起始列和终止列两列中，数值前面加负号来表示该数值为相对位置，这样更好描述了转录本类型为环状RNA(circleRNA)的序列(如图8中的产物：L01.t6.p1)；

序列提取模块用于根据注释信息和对应的基因组序列(核酸、蛋白)提取相应的基因、调控元件序列，提高了序列提取的效率。

本发明的基因注释文件解析工具采用perl计算机程序语言编写。

3、针对本发明的基因注释文件格式的解析方法包括：

a)格式转化模块：

格式检验：从参数中获取输入文件，并传递给解析步骤。

格式解析：根据不同格式解析文件，并将解析后的数据以JSON格式存放。

格式输出：根据参数得到的输出格式，将数据按照所需格式输出。

上述输入输出均支持三种注释格式(GTF、GFF3、GFF4)。

b)查询模块：

给定一个模式匹配字符串，该模块搜索文件中能匹配给定模式列表的行，旨在快速定位满足一定条件的基因注释条目，目前输入只支持GFF4注释格式。

步骤：格式解析->grep模块–>格式输出

格式解析：从参数中获取输入文件(GFF4格式注释文件)，并传递给解析步骤，解析的结果数据以JSON的格式存放(JSON1)；

grep模块：利用模式匹配，从JSON1中获取满足条件的注释数据，并以JSON的格式存放(JSON2)；

其中，模式匹配公式为：[标签操作符值]：

标签：第1-8列的列名(seqId,start,end,strand,level,featureType,id,parentId)，以及属性(attribute)列中的所有标签(见1.1中所有标签)

操作符：(数字类比较：>,<,＝＝,<＝,<＝,！＝)、(字符类比较：＝＝)、(模式匹配：＝～,！～)

值：数字类、字符串(双引号括起来)、模式匹配(perl运行的所有正则表达式，//)

示例：

1.'seqId！～/ChrM/'

2.'level＝＝"processed"&&product_biotype＝＝"ORF"'

格式输出：以GFF4格式输出

c)注释合并模块：

对多个注释文件进行解析，存储在JSON格式的数据结构中，以一个注释文件作为参考，其他文件依次与参考注释文件合并，位置有交集的基因会合并成为一个基因。最后按统一的格式输出。

上述输入输出均支持三种注释格式(GTF、GFF3、GFF4)。

d)单元提取模块：

对输入文件(一个任意注释文件)进行解析，以JSON格式存放(JSON1)；之后，提取目标序列编号(id)的父序列(parent)或子序列(child)(也可以同时提取子父序列)方向的注释信息，以JSON格式存放(JSON2)；最后以GFF4格式输出。该模块便于得到目标序列编号(id)的序列来源或/和序列产物。

注：

父序列方向：根据父序列编号递归向上提取；

子序列方向：根据子序列编号递归向下提取；

命令格式：[目标序列编号1][:][parent/child][目标序列编号2][:][parent/child]

示例：

父序列方向提取：L01.t1:parent；

子序列方向提取：L01.t1:child；

子父双向提取：L01.t1:parent:L01.t1:child。

e)位置转化模块：

示例1：

一条加工转录本t1，它的一条产物是p1，从文章中知道p1(产物类型为编码序列(CDS))第100个碱基A突变为T，通过该条编码序列在基因组上的位置，利用相对位置转化绝对位置模块，研究者可以得到该突变碱基在基因组上的位置，即绝对位置。

该模块只支持GFF4格式。

输入为一个目标序列注释记录(相对位置形式，GFF4格式)，和一个包含目标序列的注释文件(绝对位置形式，任意注释格式)。

示例2：

已知一条产物p1(产物类型为编码序列(CDS))的绝对位置，为获得突变位点在编码序列上的相对位置，以便分析突变对蛋白编码的影响，则需通过绝对位置转相对位置模块，得到相对产物p1的突变位置，即可分析相应的蛋白序列发生了怎样的变化，是否影响到了蛋白的功能。

该模块同时支持GFF4、GFF3、GTF格式输入，但输出只支持GFF4格式。

只需要输入目标注释记录(以注释文件形式输入)。

f)序列提取模块：

根据注释文件和基因组序列文件，生成注释文件中每个编号(id)的核苷酸序列文件(*_nucleotide.fasta)；如果包含编码序列的转录本，会有编码序列文件(CDS)(*_CDS.fasta)以及蛋白序列文件(*_pep.fasta)；如果有调控元件序列，比如终止子，还会生成调控元件的序列文件(*_regulator.fasta)。

具体步骤：对输入的注释文件进行解析，并以JSON格式存放(JSON1)；之后，对输入的基因组序列文件依照每个条目的位置、产物类型，对序列进行提取、翻译、输出，生成各条目相应的序列文件。

该模块同时支持GFF4、GFF3、GTF格式。

以上显示和描述了本发明的基本原理和主要特征以及本发明的优点。本行业的技术人员应该了解，本发明不受上述实施例的限制，上述实施例和说明书中描述的只是说明本发明的原理，在不脱离本发明精神和范围的前提下，本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明范围内。本发明要求保护范围由所附的权利要求书及其等效物界定。