用gffread从基因组中提取基因、转录本、蛋白、启动子、非编码序列

https://jishuin.proginn.com/p/763bfbd5e36e
1.获取转录本序列
gffread GRCh38.gtf -g GRCh38.fa -w GRCh38.transcripts.fa

获取CDS序列
gffread GRCh38.gtf -g GRCh38.fa -x GRCh38.cds.fa
获取蛋白序列
gffread GRCh38.gtf -g GRCh38.fa -y GRCh38.protein.fa
4.提取基因启动子序列
首先确定启动子区域，这里定义转录起始位点上游1000 bp和下游500 bp为启动子区域。
sed 's/"/\t/g' GRCh38.gtf | awk 'BEGIN{OFS=FS="\t"}{if( $3=="gene") {if($ 7=="+") {start= $4-1000; end=$ 4+500;} else {if( $7=="-") start=$ 5-500; end= $5+1000; } if(start<0) start=0; print$ 1,start,end, $14,$ 10,$7;}}' >GRCh38.promoter.bed

然后提取序列。这里用到了bedtools工具，官方有提供编译好的二进制文件，下载下来即可使用。
-name: 输出基因名字（bed文件的第四列）
-s: 考虑到正反链（对于启动子区域，是否考虑链的信息关系不太大）
bedtools getfasta -name -s -fi GRCh38.fa -bed GRCh38.promoter.bed >GRCh38.promoter.fa
5.提取基因序列
提取基因序列的操作也类似于提取启动子序列。这里要注意GFF文件的序列位置是从1开始，而bed文件的位置是从0开始，前闭后开，所以要对序列的起始位置进行-1的操作。
type="gene"
sed 's/"/\t/g' GRCh38.gtf | awk -v type=" ${type}" 'BEGIN{OFS=FS="\t"}{if($ 3==type) {print $1,$ 4-1, $5,$ 14,".", $7}}' >GRCh38.gene.bed 提取基因序列 bedtools getfasta -name -s -fi GRCh38.fa -bed GRCh38.gene.bed >GRCh38.gene.fa 6.提取非编码RNA的序列在GTF文件中有转录本类型的注释，包含下面这些注释类型 ntisense_RNA；lincRNA；miRNA；misc_RNA；processed_pseudogene；processed_transcript；protein_coding；rRNA；scaRNA；sense_intronic；sense_overlapping；snoRNA；snRNA；TEC；transcribed_processed_pseudogene；transcribed_unitary_pseudogene；transcribed_unprocessed_pseudogene；unitary_pseudogene；unprocessed_pseudogene 我们只筛选lincRNA grep 'transcript_biotype "lincRNA"' GRCh38.gtf >GRCh38.lincRNA.gtf gffread GRCh38.lincRNA.gtf -g GRCh38.fa -w GRCh38.lincRNA.fa 7.提取一个个外显子序列获取外显子的坐标 type="exon" sed 's/"/\t/g' GRCh38.gtf | awk -v type="$ {type}" 'BEGIN{OFS=FS="\t"}{if( $3==type) {print$ 1, $4-1,$ 5, $14,$ 20, $7}}' >GRCh38.exon.bed 提取序列 -name: 输出基因名字（bed文件的第四列） -s: 考虑到正反链（对于启动子区域，是否考虑链的信息关系不太大） bedtools getfasta -name -s -fi GRCh38.fa -bed GRCh38.exon.bed >GRCh38.exon.fa 8.提取一个个内含子序列确定内含子区域 sed 's/"/\t/g' GRCh38.gtf | awk 'BEGIN{OFS=FS="\t";oldtr="";}{if($ 3=="exon") {tr= $14; if(oldtr!=tr) {start=$ 5; oldtr=tr;} else {print $1,start,$ 4-1,tr, $20,$ 7; start=$5;} } }' >GRCh38.intron.bed

最后编辑于：2022.04.24 16:33:28