I have a batch of 'noteworthy positions' (50,000+) and I'm trying to determine the best way to batch annotate them with my annotated genome (GTF format or anything which it can be converted to).

Does anyone know a package that already does this? I have considered writing a basic script but perhaps such a thing exists already.

My intention is to provide the coordinates and return columns such as:

```
Nearest    Is within    Intron/exon    Distance    positional annotation
gene       gene/UTR/                   from TSS    (if exists)
           Promotor?
```

Thanks
To get the nearest gene, you could use BEDOPS <a href="http://bedops.readthedocs.org/en/latest/content/reference/set-operations/closest-features.html">*closest-features*</a>:

    $ closest-features --closest positions.bed genes.bed > answer.bed

To answer if within gene/UTR/promoter/intron/exon, you could use BEDOPS <a href="http://bedops.readthedocs.org/en/latest/content/reference/statistics/bedmap.html#indicator">*bedmap --indicator*</a> (or similar options, depending on what level of detail you need):

```
$ bedmap --echo --indicator positions.bed genes.bed > answer_genes.bed
$ bedmap --echo --indicator positions.bed UTR.bed > answer_UTR.bed
$ bedmap --echo --indicator positions.bed promoter.bed > answer_promoter.bed
```

etc.

To get the distance from the nearest TSS, you could use BEDOPS *closest-features* again:

    $ closest-features --closest --dist positions.bed TSS.bed > answer.bed

To get annotations of any overlapping elements, you could use BEDOPS *bedmap --echo-map*:

<pre>
$ bedmap --echo --echo-map positions.bed annotations.bed > answer.bed</pre>

If your annotations are in a non-BED format, you could use BEDOPS <a href="http://bedops.readthedocs.org/en/latest/content/reference/file-management/conversion/convert2bed.html">*convert2bed*</a> to turn them into BED:

```
$ gtf2bed < annotations.gtf > annotations.bed
$ gff2bed < annotations.gff > annotations.bed
$ vcf2bed < annotations.vcf > annotations.bed
```

etc.

Because BEDOPS tools use Unix I/O streams, you can convert inside a larger BEDOPS operation pipeline:

    $ gff2bed < annotations.gff | closest-features --closest --dist positions.bed - > answer.bed

Or you can use process substitution, if you use `bash`:

    $ closest-features --closest --dist positions.bed <(gff2bed < annotations.gff) > answer.bed

To turn BED-formatted annotations into BED-formatted lists of UTRs, introns, exons and promoters, etc. is easy and a cursory search of Biostars or seqanswers can help locate small scripts in *awk* or Perl etc. that can do this for various annotation subcategories.

For instance, to turn a BED-formatted list of stranded gene annotations into 1k proximal promoter regions:

```
$ awk '{ \
        if ($6=="+") { \
            print $1"\t"($2 - 1000)"\t"$2"\t"$4"\t"$5"\t"$6; \
        } \
        else { \
            print $1"\t"$3"\t"($3 + 1000)"\t"$4"\t"$5"\t"$6; \
        } \
    }' genes.bed \
    > promoters.bed
```

Then it is a matter of using `promoters.bed` in downstream BEDOPS operations.