[BioC] get promoter regions (bed format) of genes hg19

[Harris A. Jaffee]

Thanks, Val, for the details.  Now I have what is either a shameless plug or
a request for explanation.  matchGenes() in _bumphunter_ is similar, although
definitely not the same.  I am not sure how to translate between the results,
except to note that matchGenes uses nearest() with either select="arbitrary",
by default, or select="all" -- rather than findOverlaps().  It is much faster
than what you did, in effect because your transcriptsBy call is "compiled in".
At any rate, the user just does
	
	library(bumphunter)
	
	# this uses nearest(select="arbitrary")
	matchGenes(htr_genes_gr, promoterDist=2000)
or
	# this uses nearest(select="all")
	matchGenes(htr_genes_gr, promoterDist=2000, all=TRUE)

The results are annotated based on the matching transcript's TSS, but neither
answer has length 6 in this example.  Maybe someone can elucidate the other
differences (and perhaps some similarities).

On Jan 22, 2014, at 11:56 AM, Valerie Obenchain wrote:

> Hi Ninni,
> 
> A more complete description of the promoters() function is on the man page in the GenomicRanges package. You can find this by specifying the full method name:
> 
>  ?'promoters,GenomicRanges-method'
> 
> The GRanges or GRangesList you supply as 'x' is assumed to contain gene ranges. Here is a snip from the man page stating that the assumption is that start(x) is the TSS:
> 
>>        • ‘promoters’ returns an object of the same type and length as
>>          ‘x’ containing promoter ranges. Promoter ranges extend around
>>          the transcription start site (TSS) which is defined as
>>          ‘start(x)’. The ‘upsteam’ and ‘downstream’ arguments define
>>          the number of nucleotides in the 5' and 3' direction,
>>          respectively. The full range is defined as,
>> 
>>          (start(x) - upstream) to (start(x) + downstream - 1).
> 
> So promoters() simply returns the region upstream and/or downstream of the range you supply as 'x'. If you have a generic range and want to determine if it falls within a gene you can use findOverlaps() on a TranscriptDb and your ranges.
> 
>  library(TxDb.Hsapiens.UCSC.hg19.knownGene)
>  txdb <- TxDb.Hsapiens.UCSC.hg19.knownGene
>  genes <- transcriptsBy(txdb, "gene")
>  hits <- findOverlaps(htr_genes_gr, genes)
> 
> 'queryHits' is an index of 'htr_genes_gr' ranges that maps to the 'subjectHits' index of 'genes' ranges. Here we see that both range 2 and 4 in 'htr_genes_gr' each map to 2 different genes.
> 
>>> hits
>> Hits of length 6
>> queryLength: 4
>> subjectLength: 23459
>>  queryHits subjectHits
>>   <integer>   <integer>
>> 1         1       19185
>> 2         2       17448
>> 3         2       20105
>> 4         3       17448
>> 5         4        5026
>> 6         4        8270
> 
> Isolate the genes that the ranges in 'htr_genes_gr' hit:
> 
>  mygenes <- genes[subjectHits(hits)]
> 
> Call promoters() on the gene ranges of interest.
> 
>  promoters(mygenes, upstream=2000)
> 
> 
> Valerie
> 
> On 01/22/2014 01:12 AM, Ninni Nahm [guest] wrote:
>> 
>> Dear all,
>> 
>> I'm trying to obtain promoter regions of a list of genes.
>> This is my granges example:
>> 
>> htr_genes_gr = GRanges( seqnames = rep("chr1", 4),
>> ranges = IRanges(start = c(131025, 761586, 762988, 860260), end = c(134836, 794826, 794826, 879955)))
>> 
>> I wanted to use promoters(htr_genes_gr, upstream=2000)
>> to get the upstream regions. However, I get something, but Im afraid that it does not match with the USCS genome browser. So, I take the output and view it in IGV with hg19 loaded. But the promoter seems shifted, i.e. it does not start where the gene begins. When I look at the IRanges web page, there is an indication that the genome version might be hg18. I am not sure if this is the problem. Anyways, does anyone have a suggestion how I can get the bed file of the promoter regions. I need the promoter regions of all genes in the end.
>> Thanks a lot!
>> 
>> 
>> 
>>  -- output of sessionInfo():
>> 
>>  sessionInfo()
>> R version 3.0.2 (2013-09-25)
>> Platform: x86_64-pc-linux-gnu (64-bit)
>> 
>> locale:
>>  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
>>  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
>>  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8
>>  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C
>>  [9] LC_ADDRESS=C               LC_TELEPHONE=C
>> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
>> 
>> attached base packages:
>> [1] stats     graphics  grDevices utils     datasets  methods   base
>> 
>> --
>> Sent via the guest posting facility at bioconductor.org.
>> 
>> _______________________________________________
>> Bioconductor mailing list
>> Bioconductor at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/bioconductor
>> Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor
>> 
> 
> 
> -- 
> Valerie Obenchain
> 
> Program in Computational Biology
> Division of Public Health Sciences
> Fred Hutchinson Cancer Research Center
> 1100 Fairview Ave. N, M1-B155
> P.O. Box 19024
> Seattle, WA 98109-1024
> 
> E-mail: vobencha at fhcrc.org
> Phone:  (206) 667-3158
> Fax:    (206) 667-1319
> 
> _______________________________________________
> Bioconductor mailing list
> Bioconductor at r-project.org
> https://stat.ethz.ch/mailman/listinfo/bioconductor
> Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor