Coverage Question: GenomicRanges
1
0
Entering edit mode
@lakshmanan-iyer-1829
Last seen 9.2 years ago
United States
I am calculating the coverage on a subset of a gapped aliment that obtained using "subsetByOverlaps" between a "GappedAlignments" and a "GRanges" object #> class (Dendaligns) #[1] "GappedAlignments" #attr(,"package") #[1] "GenomicRanges" #> class (p3UTRsChr) #[1] "GRanges" #attr(,"package") #[1] "GenomicRanges" DendfiltData2 <- subsetByOverlaps (Dendaligns, p3UTRsChr); #> class (DendfiltData2) #[1] "GappedAlignments" #attr(,"package") #[1] "GenomicRanges" # # Find coverage Dendcov <- coverage(DendfiltData2) # # Find islands of coverage with min value of 1 Dendislands <- slice (Dendcov, lower=1); # # convert the island ranges in to a data frame DendislandsData <- as.data.frame (ranges(Dendislands[[chr]]) ); # # extract the scores of the coverage in a list Dendisland_scores <- viewApply (Dendislands[[chr]], function (x) as.integer(x)) When I look at the scores, the first island contains a series of 1 (45 of them in my case). However, when I visualize the same regions in IGV with the raw bam file (see attached picture), I see that the counts are > 1 for bases after first 10 bases. I am wondering where the counting becomes different from summing up the reads at any position! Thanks for any thoughts -best -Lax -------------- next part -------------- A non-text attachment was scrubbed... Name: igv_panel.png Type: image/png Size: 6317 bytes Desc: not available URL: <https: stat.ethz.ch="" pipermail="" bioconductor="" attachments="" 20120309="" 93fcaa7e="" attachment.png="">
Coverage convert Coverage convert • 1.4k views
ADD COMMENT
0
Entering edit mode
@michael-lawrence-3846
Last seen 3.0 years ago
United States
One possibility is that you need to use ignore.strand = TRUE in your subsetByOverlaps call, i.e., you are filtering out all of those negative strand reads. Michael On Fri, Mar 9, 2012 at 9:29 AM, Lakshmanan Iyer <laxvid@gmail.com> wrote: > I am calculating the coverage on a subset of a gapped aliment that obtained > using "subsetByOverlaps" between a "GappedAlignments" and a "GRanges" > object > #> class (Dendaligns) > #[1] "GappedAlignments" > #attr(,"package") > #[1] "GenomicRanges" > #> class (p3UTRsChr) > #[1] "GRanges" > #attr(,"package") > #[1] "GenomicRanges" > DendfiltData2 <- subsetByOverlaps (Dendaligns, p3UTRsChr); > #> class (DendfiltData2) > #[1] "GappedAlignments" > #attr(,"package") > #[1] "GenomicRanges" > # > # Find coverage > Dendcov <- coverage(DendfiltData2) > # > # Find islands of coverage with min value of 1 > Dendislands <- slice (Dendcov, lower=1); > # > # convert the island ranges in to a data frame > DendislandsData <- as.data.frame (ranges(Dendislands[[chr]]) ); > # > # extract the scores of the coverage in a list > Dendisland_scores <- viewApply (Dendislands[[chr]], function (x) > as.integer(x)) > > When I look at the scores, the first island contains a series of 1 (45 of > them in my case). > However, when I visualize the same regions in IGV with the raw bam file > (see attached picture), I see that the counts are > 1 for bases after first > 10 bases. > > I am wondering where the counting becomes different from summing up the > reads at any position! > > Thanks for any thoughts > -best > -Lax > > _______________________________________________ > Bioconductor mailing list > Bioconductor@r-project.org > https://stat.ethz.ch/mailman/listinfo/bioconductor > Search the archives: > http://news.gmane.org/gmane.science.biology.informatics.conductor > [[alternative HTML version deleted]]
ADD COMMENT
0
Entering edit mode
You are right. From the picture attached, I could see that only the first read is on + strand, the other three are - strand and are being ignored. Let me rerun and see Best Lax Sent from my iPhone On Mar 9, 2012, at 1:05 PM, Michael Lawrence <lawrence.michael@gene.com> wrote: > One possibility is that you need to use ignore.strand = TRUE in your subsetByOverlaps call, i.e., you are filtering out all of those negative strand reads. > > Michael > > On Fri, Mar 9, 2012 at 9:29 AM, Lakshmanan Iyer <laxvid@gmail.com> wrote: > I am calculating the coverage on a subset of a gapped aliment that obtained > using "subsetByOverlaps" between a "GappedAlignments" and a "GRanges" object > #> class (Dendaligns) > #[1] "GappedAlignments" > #attr(,"package") > #[1] "GenomicRanges" > #> class (p3UTRsChr) > #[1] "GRanges" > #attr(,"package") > #[1] "GenomicRanges" > DendfiltData2 <- subsetByOverlaps (Dendaligns, p3UTRsChr); > #> class (DendfiltData2) > #[1] "GappedAlignments" > #attr(,"package") > #[1] "GenomicRanges" > # > # Find coverage > Dendcov <- coverage(DendfiltData2) > # > # Find islands of coverage with min value of 1 > Dendislands <- slice (Dendcov, lower=1); > # > # convert the island ranges in to a data frame > DendislandsData <- as.data.frame (ranges(Dendislands[[chr]]) ); > # > # extract the scores of the coverage in a list > Dendisland_scores <- viewApply (Dendislands[[chr]], function (x) > as.integer(x)) > > When I look at the scores, the first island contains a series of 1 (45 of > them in my case). > However, when I visualize the same regions in IGV with the raw bam file > (see attached picture), I see that the counts are > 1 for bases after first > 10 bases. > > I am wondering where the counting becomes different from summing up the > reads at any position! > > Thanks for any thoughts > -best > -Lax > > _______________________________________________ > Bioconductor mailing list > Bioconductor@r-project.org > https://stat.ethz.ch/mailman/listinfo/bioconductor > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor > [[alternative HTML version deleted]]
ADD REPLY

Login before adding your answer.

Traffic: 823 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6