[BioC] Gviz - Plot genes and data from - strand (3'-5') in 5'-3' direction

[Robert Ivanek]

I realized it is better to change the coordinates than the actual values:

reverseCounts <- function(reads, genes) {
    genes <- genes[strand(genes)=="-"]
    strand(genes) <- "*"
    ##
    ov <- findOverlaps(ranges(reads), ranges(genes))
    ww <- width(reads)[queryHits(ov)]
    start(reads)[queryHits(ov)] <- start(genes)[subjectHits(ov)] +
       (end(genes)[subjectHits(ov)] - start(reads)[queryHits(ov)])
    width(reads)[queryHits(ov)] <- ww
    ##
    return(reads)
}

Best
Robert


On 12/09/13 17:01, Robert Ivanek wrote:
> I see, sorry I got it wrong.
> It cannot be done automatically in Gviz but maybe following can help you:
> 
> library(GenomicRanges)
> library(Gviz)
> 
> genes.gr <- GRanges("chr1", IRanges(c(10,20,30,40), width=8),
>                     strand=c("+","-","+","-"))
> aT <- AnnotationTrack(genes.gr)
> 
> reads.gr <- GRanges("chr1", IRanges(seq(1,50,2), width=1),
>                     strand="+", count=seq(1,50,2))
> dT1 <- DataTrack(reads.gr)
> 
> ## this function expects DataTrack in which you would like to flip the
> points and AnnotationTrack with genes which would define regions you
> would like to flip, those regions (genes) have to be on minus strand
> 
> reverseCounts <- function(reads, genes) {
>     genes <- genes[strand(genes)=="-"]
>     strand(genes) <- "*"
>     ##
>     ov <- findOverlaps(ranges(reads), ranges(genes))
>     indx <- split(queryHits(ov), subjectHits(ov))
>     indx <- lapply(indx, rev)
>     indx <- unlist(indx)
>     values(reads)[,queryHits(ov)] <- values(reads)[,indx]
>     ##
>     return(reads)
> }
> ##
> dT2 <- reverseCounts(dT1, aT)
> ##
> plotTracks(list(aT, dT1, dT2))
> 
> Is that what you want?
> Best
> Robert
> 
> 
> On 12/09/13 16:21, Jager, Dominik wrote:
>> ...yes, thats correct.
>>
>> So what I really want is to plot my data (see attached PNG) flipped by 180 degrees. So now the genes and the data were plotted in 3' to 5' direction, but I want in 5' to 3'/
>>
>> Here some code:
>>
>>> setwd("~/WORKING/")
>>> library(Gviz)
>>> options(ucscChromosomeNames = FALSE)
>>> tkgff <- import.gff3("~/WORKING/NC_006624.gff", format = "gff3", genome= "NC_06624.1", asRangedData = FALSE)
>>> atrack <- AnnotationTrack(tkgff, strand = "*", chromosome = chr, genome = gen, name = TK)
>>> gtrack <- GenomeAxisTrack()
>>> WIG <- import.wig("~/S0_20min_reverse.wig", asRangedData = FALSE)
>>> dtrack <- DataTrack(WIG, name = "S0+")
>> >from <- 1049512
>>> to <- 1052176
>>> plotTracks(list(gtrack,atrack, dtrack), from = from, to = to, cex = 1, type = "h")
>>
>> I also attached the plot as PNG.
>>
>>
>> Dominik
>>
>> Dr. Dominik Jäger
>> Biological Sciences Department of Microbiology
>> 405 Biological Science Building | 484 West 12th Avenue Columbus, OH 43210
>> (614)-682-6890 Office
>> jager.9 at osu.edu osu.edu
>>
>> ________________________________________
>> Von: bioconductor-bounces at r-project.org [bioconductor-bounces at r-project.org]" im Auftrag von "Steve Lianoglou [lianoglou.steve at gene.com]
>> Gesendet: Donnerstag, 12. September 2013 09:57
>> An: Robert Ivanek
>> Cc: bioconductor at r-project.org list
>> Betreff: Re: [BioC] Gviz - Plot genes and data from - strand (3'-5') in 5'-3' direction
>>
>> Hi,
>>
>> On Thu, Sep 12, 2013 at 6:14 AM, Robert Ivanek <robert.ivanek at unibas.ch> wrote:
>>> Hi Dominik,
>>>
>>> You can achieve that by changing the strand in corresponding track
>>> objects.
>> [snip]
>>
>> While this might change the strand of a gene, I do not think this is
>> what the OP is really after -- as I have previously been keen to try
>> to implement the approach that (I think) the OP wants.
>>
>> I believe the OP would like to plot "the mirror image" of the plot
>> that Gviz would "normally" show given the start/stop coordinates that
>> the plotting functionality infers (or is given).
>>
>> The approach you suggest wouldn't properly "flip" my gene models and
>> the data associated with them -- it would only show them on the
>> positive strand. The 5'utr of the gene that was originally on the
>> negative strand would look as if it were the 3'utr of the gene on the
>> (now annotated) positive strand, and vice versa, and if I were (say)
>> plotting ChIP-seq data over this locus, the TF would appear to be
>> binding the 3'utr of a gene on the positive strand (or downstream of
>> it), when that's not really the case.
>>
>> -steve
>>
>> --
>> Steve Lianoglou
>> Computational Biologist
>> Bioinformatics and Computational Biology
>> Genentech
>>
>> _______________________________________________
>> 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
>>
>>
>> _______________________________________________
>> 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
>>
> 
> _______________________________________________
> 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
>