The support.bioconductor.org editor has been updated to markdown! Please see more info at: Tutorial: Updated Support Site Editor

Question: Gviz color alignments by read strand
0
gravatar for Chris Stubben
4 months ago by
Salt Lake City, Utah
Chris Stubben40 wrote:

Is there a way to use fill.reads to color alignments by read strand?  I would like to plot these on a single track

library(Gviz)
alTrack <- AlignmentsTrack(system.file(package="Gviz", "extdata", 
  "gapped.bam"), isPaired=FALSE)
plotTracks(alTrack, from=2997000, to=2998200, chromosome="chr12", 
  fill.reads=c("red", "blue"))

I tried using the strandedBamImport function from a previous question (Gviz - show plus and minus strand from a bam) but get an Error in .Call2("cigar_ranges" ...  

 

gviz • 197 views
ADD COMMENTlink modified 6 weeks ago by Robert Ivanek590 • written 4 months ago by Chris Stubben40
Answer: Gviz color alignments by read strand
2
gravatar for Robert Ivanek
6 weeks ago by
Robert Ivanek590
Switzerland
Robert Ivanek590 wrote:

The way how Gviz handles colouring for AlignmentsTrack is a bit special. It recycles colours specified in the fill.reads along the length of unique read ids (uid column). On top of it, if you plot reads as arrows, it also sorts the reads before plotting according to the strand.

Maybe this piece of code does what you want. It is similar to what you suggested but still retain all plotting functionalities for the AlignmentsTrack class. I tested it on few regions (even with gapped reads) and it seems to work.

library(Gviz)

bf <- system.file(package="Gviz", "extdata", "gapped.bam")

chrom <- "chr12"
from <- 2997500
to <- 2998200

alTrack <- Gviz:::.import.bam.alignments(bf, GRanges(chrom, IRanges(from, to)))
alTrack <- AlignmentsTrack(sort(alTrack))

displayPars(alTrack)$fill.reads <- c(plus="#ECAFAE", minus="#AFAEED")[as.numeric(sort(ranges(alTrack)$readStrand))]

plotTracks(alTrack, from=from, to=to, chromosome=chrom)

Let me know if that works for you.

ADD COMMENTlink written 6 weeks ago by Robert Ivanek590

That works great, thanks.  I just had to add options(ucscChromosomeNames=FALSE)

ADD REPLYlink modified 6 weeks ago • written 6 weeks ago by Chris Stubben40
Answer: Gviz color alignments by read strand
0
gravatar for Chris Stubben
4 months ago by
Salt Lake City, Utah
Chris Stubben40 wrote:

I can color by read strand if I subset a data.frame before creating the alignments track.  

bam1 <- readGAlignments(BamFile(system.file(package="Gviz", "extdata", "gapped.bam")))
bamdf <- data.frame(bam1)
names(bamdf)[1] <- "chromosome"
bamdf <- subset(bamdf, start > 2997000 & end < 2998200)
clrs <- ifelse(bamdf$strand == "+", "#ECAFAE", "#AFAEED")
alTrack <- AlignmentsTrack(bamdf, name="Alignments", fill = clrs)
plotTracks(alTrack)
ADD COMMENTlink written 4 months ago by Chris Stubben40

This only works if ALL reads are plotted in the display.

ADD REPLYlink written 4 months ago by Chris Stubben40

Hi Chris,

what exactly do you want to achieve? Do you want to plot individual reads (coloured by the read strand) or coverage (histogram) for both strands separately?

Robert

ADD REPLYlink written 6 weeks ago by Robert Ivanek590

I have some Digenome-seq data and want to plot potential cut sites that should have reads stacked on opposite strands (so individual reads colored by read strand).  My code works if there's only a few hundred reads that all fit in the plot.

ADD REPLYlink written 6 weeks ago by Chris Stubben40
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 16.09
Traffic: 333 users visited in the last hour