I have been using a custom 'strandedBamImport' function suggested by Florian Gviz - show plus and minus strand from a bam to plot plus and minus strand coverage from bam files with gviz and it has been really useful for my work. However I want to plot strand-specific library-size normalized RPM coverage with gviz instead of raw coverage, and I am not sure how to easily do this for multiple libraries. I want to provide the library sizes for each libraries rather than using the no. of mapped reads (we performed some filtering on the mappings), but I couldn't figure out how to incorporate this within the custom stranded bam import function. Could anyone show me how to adapt the custom function code or suggest other function which can be used to calculate normalized coverage based on user provided library sizes? Thanks!
Here is a sample code I am using now to plot raw coverage for plus and minus strand from a bam file for an arabidopsis gene:
biomTrack = BiomartGeneRegionTrack(genome="TAIR10",biomart=mart, chromosome="chr2", start=16531688, end=16533931,showId=T, geneSymbols=T, rotate.title=TRUE, filters=list(biotype="protein_coding"),collapseTranscripts=FALSE) bamfile1<-"lib1.bam" dataTrack1 = DataTrack(bamfile1, genome = "TAIR10",chromosome="chr2", name="Coverage",importFunction=strandedBamImport, stream=TRUE, legend=TRUE) plotTracks(list(biomTrack,dataTrackin1), from = 16531688, to = 16533931, type="hist", col.histogram=NA, cex.title=1, cex.axis=1, title.width=1.2, col=c("#1d91c0","#54278f"),groups=c("+", "-"))