GFF文件和基因组文件提取mRNA,cds,protein序列

首先是用gffread提取cds序列,蛋白序列,转录本序列

gffread genome.gff3 -g genome.fa -x  cds.fa
gffread genome.gff3 -g genome.fa -y  protein.fa
gffread genome.gff3 -g genome.fa -w  transcripts.fa

接下来我们利用组合工具来提取mRNA,和gene序列

python -m jcvi.formats.gff bed --type=mRNA --key=ID a.gff3 > mRNA.bed
这个需要借助Python 包jcvi包
抽取 GFF文件中mRNA, key是ID
基因也是如此
python -m jcvi.formats.gff bed --type=gene --key=ID  a.chr.gff3  > gene.bed

image.png

注意用bedtools提取序列，要加-s 参数（能区分正反链）

bedtools getfasta -fi ./a.fasta  -bed gene.bed   -nameOnly  -s -fo a.gene.fa

a.gene.fa就是我们需要的文件，同理也可以得到a.mRNA.fa。
再利用提取最长转录本脚本，获得基因中最长可变剪切的序列。