I find that I have various uses for a function that generates "gene model" GRanges on the basis of a symbol.  Even more detailed gene model computations are performed in Gviz.

Do other users find such a function useful?  Are there other aspects of model structure that should be available?  Should the GenomicFeatures package include a function of this nature?

​> genemodel

function (sym, genome = "hg19", annoResource = Homo.sapiens, 

    getter = exonsBy, byattr = "gene") 

{

    stopifnot(is.atomic(sym) && (length(sym) == 1))

    if (!exists(dsa <- deparse(substitute(annoResource)))) 

        require(dsa, character.only = TRUE)

    num = select(annoResource, keys = sym, keytype = "SYMBOL", 

        columns = c("ENTREZID", "SYMBOL"))$ENTREZID

    getter(annoResource, by = byattr)[[num]]

}

<environment: namespace:erma>

> genemodel("CD28")

GRanges object with 6 ranges and 2 metadata columns:

      seqnames                 ranges strand |   exon_id   exon_name

         <Rle>              <IRanges>  <Rle> | <integer> <character>

  [1]     chr2 [204571198, 204571471]      + |     35150        <NA>

  [2]     chr2 [204591356, 204591421]      + |     35151        <NA>

  [3]     chr2 [204591356, 204591712]      + |     35152        <NA>

  [4]     chr2 [204591674, 204591712]      + |     35153        <NA>

  [5]     chr2 [204594371, 204594495]      + |     35154        <NA>

  [6]     chr2 [204599507, 204603636]      + |     35155        <NA>

  -------

  seqinfo: 93 sequences (1 circular) from hg19 genome

I like this suggestion Vince.  There is another high level function that I recently added to OrganismDbi and which I think could be a good complement to something quite close to this.  I plan to explore this some more when I work on enhancing the other function.

Here is a version of genemodel that uses current OrganismDb infrastructure and autopropagates seqinfo.  Code is available in devel branch of erma package.

genemodel <-
  function (key, keytype = "SYMBOL", annoResource = Homo.sapiens) 
{
    oblig = c("EXONCHROM", "EXONSTART", "EXONEND", "EXONSTRAND", 
        "EXONID")
    addrs = select(annoResource, keys = key, keytype = keytype, 
        columns = oblig)
    ans = GRanges(addrs$EXONCHROM,
        IRanges(addrs$EXONSTART, addrs$EXONEND), 
        strand = addrs$EXONSTRAND, exon_id = addrs$EXONID)
    mcols(ans)[[keytype]] = key
    useq = unique(as.character(seqnames(ans)))
    si = seqinfo(annoResource)
    seqinfo(ans) = si[useq, ]
    ans
}