【3.5.1.3】GTF与GFF的处理--gffread

gffread 不仅可以实现GTF与GFF的互相转换，而且还可以对GFF文件进行过滤处理。可以直接读取GTF文件。

一、安装

cd /data4/software
git clone https://github.com/gpertea/gffread
cd gffread
make


#下载gtf/gff文档及hg19文件
wget -c ftp://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_19/gencode.v19.annotation.gtf.gz
wget -c ftp://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_19/gencode.v19.annotation.gff3.gz
wget -c http://hgdownload.cse.ucsc.edu/goldenPath/hg19/bigZips/chromFa.tar.gz

gunzip gencode.v19.annotation.gtf.gz
gunzip gencode.v19.annotation.gff3.gz
tar  xzf chromFa.tar.gz
#使用gunzip及tar分别解压，并合并下载内容获得hg19.fa
cat *.fa >hg19.fa

二、示例

示例：

Usage:
gffread <input_gff> [-g <genomic_seqs_fasta> | <dir>][-s <seq_info.fsize>]
[-o <outfile.gff>] [-t <tname>] [-r [[<strand>]<chr>:]<start>..<end> [-R]]
[-CTVNJMKQAFGUBHZWTOLE] [-w <exons.fa>] [-x <cds.fa>] [-y <tr_cds.fa>]
[-i <maxintron>]
Filters and/or converts GFF3/GTF2 records.
<input_gff> is a GFF file, use '-' if the GFF records will be given at stdin

常用参数介绍：

-g   序列文件，即GFF/GTF文件第一列ID对应的序列文件。
-i   丢弃掉内含子大于的转录本（mRNA/transcript）
-r   起始和终止位置，填写示例100.10000即为输出与100到10000有重叠的所有转录组，也可以限制序列ID及链，填写示例：+Chr1:100..10000。
-R   丢弃掉此范围的转录本，与-r相反。
-U   丢弃掉 single-exon的转录本
-C   丢低调无CDS的转录本。
-V   丢弃掉含有移码突变的转录本。
-H   如果使用了-V，则重新检查并调整内含子相位，避免由于翻译起始位点选择的位置不对导致移码突变的产生。
-B   如果使用了-V， 对于单外显子基因，则重新检查相反的链，是否存在移码突变。
-N   丢弃掉多外显子基因剪接位点不是常见的 GT-AG, GC-AG or AT-AC序列。
-J   丢弃掉没有起始密码子或者终止密码子的转录本，仅保留有完整编码框的转录本。
--no-pseudo:  过滤掉含有 'pseudo' 的注释信息
-M/--merge :  合并完全相同的或者存在包含关系的转录本。
-d:使用 -M  将合并信息输出到文件中
--cluster-only: 类似于 --merge 但是不合并转录本
-K   对于-M 选项：also collapse shorter, fully contained transcripts
      with fewer introns than the container
-Q   对于-M 选项：移除包含关系的转录本的限制条件：多外显子转录本将会合并，如果他们内含子位置完全一样，单外显子转录本只需要有80%一样即可合并。
--force-exons:  使GFF features的最小层级为exon
-E   对于重复的 ID或者 GFF/GTF 其他潜在的格式问题给出警告信息。
-Z   将内含子小于4 bp的邻近的两个外显子合并为一个。
-w   输出每个转录本的外显子序列
-x   输出CDS序列
-W   对于 -w 和 -x 选项，输出外显子位置坐标到FASTA序列的ID中
-y   输出蛋白质序列
-L   将Ensembl GTF 转换为 GFF3 conversion (implies -F; should be used with -m)
-o   输出"filtered" 后的GFF文件 。
-T -o  参数将输出 GTF格式。

常用示例：

GFF转换GTF
gffread input.gff3 -T -o out.gtf

GTF转换GFF3
gffread input.gtf -o out.gff3

根据GFF或者GTF提取CDS，蛋白质和外显子序列
gffread gene.gff3 -g genome.fa -x cds.fa -y pep.fa -w cdna.fa

只提取翻译后蛋白序列
gffread gencode.v19.annotation.gff3 -g hg19.fa -y tr_pep.fa

根据reference提取CDS序列
gffread gencode.v19.annotation.gff3 -g hg19.fa -x cds.fa

只提取外显子序列
gffread gencode.v19.annotation.gff3 -g hg19.fa -w exons.fa

三、注意

1）如果genome.fa中有N、或非ATCG字符，会在转成CDS的时候变为H。：

提取后的cds file：

2）gff 注意正负链问题：

genome：XXXXXXXXXXXXTaGXXXXXXXXXXXX
cds: XXXXXXXXXXXXXCtAXXXXXXXXXXXX

当genome中负链上的碱基有改变：某个A改变成a，这时候提取CDS序列后，会发现是T变t。

那么后面插入一个T就相当于前面插入一个a；

3）gffread是提取序列，不是生成bed文件。

三、我的案例

/data4/software/gffread/gffread /data4/database/human/gencode.v19.annotation.gff3 -g /data4/database/human/hg19.fa  -r chr3:38182631..38182651  -y tr_pep.fa -x tr_nucl.fa -W

参考资料

• https://github.com/gpertea/gffread
• https://www.jianshu.com/p/cadfbdcf938f