how to get Homo_sapiens.GRCh38.84.bed file

Question

[yuxinghai]

I'm using geneBody_coverage2.py in RseQC package to calculate Read coverage over gene body, which need a bed12 file. I did alignment based on Homo_sapiens.GRCh38.84.gtf , so i need a Homo_sapiens.GRCh38.84.bed. But gtfToBed can't produce a BED12 file, how can I get it? anyone can help me?

Answer 1

[Devon Ryan]

awk '{if($3 != "gene") print $0}' Homo_sapiens.GRCh38.84.gtf | grep -v "^#" | gtfToGenePred /dev/stdin /dev/stdout | genePredToBed stdin Homo_sapiens.GRCh38.84.bed

You can download the two tools from UCSC.

If you're curious, this is from here and is essentially the code I use to generate the reference data for our institute.

[RNAseqer]

Hi there, I was hoping you mightclarify a problem I had attempting your suggestion. I tried implementing this with the ensemble .gtf file Homo_sapiens.GRCh38.95.abinitio.gtf with the usage:

myapple$ awk '{if($3 != "gene") print $0}' Homo_sapiens.GRCh38.95.abinitio.gtf | grep -v "^#" | gtfToGenePred /dev/stdin /dev/stdout | genePredToBed stdin Homo_sapiens.GRCh38.95.bed

but got the error message:

/dev/stdin doesn't appear to be a GTF file (GFF not supported by this program)

Do you know why this may be?

[Devon Ryan]

I've never looked at the ab initio GTF files, perhaps they're a bit wonky. Is there a reason you're using that?

[RNAseqer]

Hmm. There must be something different about abinito files, since I just tried again with another gtf file from the same source (there were several available and I just happened to grab abinito one, im testing and RNAseq pipeline and learning as I go, being new to the field). With Homo_sapiens.GRCh38.95.gtf the process goes to completion and generates a bed file.

$ awk '{if($3 != "gene") print $0}' Homo_sapiens.GRCh38.95.gtf | grep -v "^#" | gtfToGenePred /dev/stdin /dev/stdout | genePredToBed stdin Homo_sapiens.GRCh38.95_not_abinito_file.bed