How to get exon number from position from script?

Question

[sacha]

Hello,

Is there a way from genomic position to get the exon number ? Would be greate if It can be provided from an API REST.

For example, get the gene name and exon or intron id from:

chr1    6500665    6500888

Using hg19.

Answer 1

[Juke34]

Hi,

I think the easier way is to download the whole annotation file in gtf or gff, and typed in your shell the following command lines (It doesn't work on Windows :) ):

To get all the features at this location you type:

awk '{if($1=="chr1" && $3>6500665 && $4<6500888) print $0}' your_annotation.gtf

If you want filtered by feature type, as example gene and exon type:

awk '{if($1=="chr1" && $4>10000 && $5<20000 && ($3=="gene" || $3=="exon")) print $0}' your_annotation.gtf

Finally to know the exon number between these locations type:

awk '{if($1=="chr1" && $4>10000 && $5<20000 &&  $3=="exon") arr++}END{print arr}' your_annotation.gtf

Answer 2

[Ram]

You can get the exon starts and ends from UCSC Genome Browser's MySQl tables, and whip up a quick script to figure out which exons the co-ordinates fall in. I'd suggest bedtools for the second step, but converting the exonStarts and exonEnds fields to BED format would be an added burden in your case.