从gff3中提取gene坐标等转为bed

13325312-c089d3fbed4d1f2e.png
gff3文件

看下第三列:

13325312-0e29975df06f720e.png
主要分成这八种

awk '$3 == "gene"' gencode.vM19.chr_patch_hapl_scaff.annotation.gff3 | awk 'BEGIN{FS="\t|=|;";OFS="\t"}{print $1,$4,$5,$16}' > genes.bed

13325312-6091b632b7815b14.png
提取基因

去除非chr染色体的基因:

sed -n '/^chr/p' genes.bed > genes.fine.bed

查看了下有: 在chr染色体上有54446个基因

#######################################

其他方法:

ref:RSeQC使用笔记 – 生信笔记

1.##gtf转化为bed:

cat reference/genome/hg19/gencode.v26lift37.annotation.gtf |awk 'OFS="\t" {if($3=="transcript") {print $1,$4-1,$5-1,$12,$6,$7}}' |tr -d '";' >hg19.bed

2.用ucsc的gtfToGenePred配合shell来将gtf转化为bed12:

gtfToGenePred -genePredExt -geneNameAsName2 ~/reference/genome/hg19/gencode.v26lift37.annotation.gtf gene.tmp

awk '{print $2"\t"$4"\t"$5"\t"$1"\t0\t"$3"\t"$6"\t"$7"\t0\t"$8"\t"$9"\t"$10}' gene.tmp >  hg19.bed12

展开阅读全文

没有更多推荐了,返回首页