Search
Question: .wig files for strand-specific paired-end RNA-Seq
0
gravatar for Guest User
5.4 years ago by
Guest User12k
Guest User12k wrote:
Hi, Is there a simple way to make strand-specific .wig file (i.e., a separate track for + and - strand) from paired-end data (where the second read maps to the other strand)? I've tried using this: library(Rsamtools) library(rtracklayer) myReads <- readGappedAlignments("RNAseqMapping.bam") coveragePlus <- coverage(myReads[strand(myReads) == '+']) export(coveragePlus, "RNAplus.wig") coverageMinus <- coverage(myReads[strand(myReads) == '-']) export(coverageMinus, "RNAminus.wig") But it appears that the second read in the pair contributes to the other strand, generating similar tracks for the + and the - strands. Is there a way to deal with this better? Thanks! Igor. -- output of sessionInfo(): R version 2.13.1 (2011-07-08) Platform: i386-pc-mingw32/i386 (32-bit) locale: [1] LC_COLLATE=English_United States.1252 [2] LC_CTYPE=English_United States.1252 [3] LC_MONETARY=English_United States.1252 [4] LC_NUMERIC=C [5] LC_TIME=English_United States.1252 attached base packages: [1] stats graphics grDevices utils datasets methods base loaded via a namespace (and not attached): [1] tools_2.13.1 -- Sent via the guest posting facility at bioconductor.org.
ADD COMMENTlink modified 5.4 years ago by Michael Lawrence9.8k • written 5.4 years ago by Guest User12k
0
gravatar for Michael Lawrence
5.4 years ago by
United States
Michael Lawrence9.8k wrote:
On Wed, Jun 20, 2012 at 4:53 AM, Igor Ulitsky [guest] < guest@bioconductor.org> wrote: > > Hi, > > Is there a simple way to make strand-specific .wig file (i.e., a separate > track for + and - strand) from paired-end data (where the second read maps > to the other strand)? I've tried using this: > > library(Rsamtools) > library(rtracklayer) > myReads <- readGappedAlignments("RNAseqMapping.bam") > coveragePlus <- coverage(myReads[strand(myReads) == '+']) > export(coveragePlus, "RNAplus.wig") > coverageMinus <- coverage(myReads[strand(myReads) == '-']) > export(coverageMinus, "RNAminus.wig") > > But it appears that the second read in the pair contributes to the other > strand, generating similar tracks for the + and the - strands. > Is there a way to deal with this better? > > Are you trying to generate coverages for the actual strand of transcription? If so, you would probably get that information from the XS tag and set it as your strand prior to export, but unless you used a special protocol the XS information would be incomplete. Btw, I would recommend BigWig export of those coverage tracks. Michael > Thanks! > > Igor. > > -- output of sessionInfo(): > > R version 2.13.1 (2011-07-08) > Platform: i386-pc-mingw32/i386 (32-bit) > > locale: > [1] LC_COLLATE=English_United States.1252 > [2] LC_CTYPE=English_United States.1252 > [3] LC_MONETARY=English_United States.1252 > [4] LC_NUMERIC=C > [5] LC_TIME=English_United States.1252 > > attached base packages: > [1] stats graphics grDevices utils datasets methods base > > loaded via a namespace (and not attached): > [1] tools_2.13.1 > > -- > Sent via the guest posting facility at bioconductor.org. > [[alternative HTML version deleted]]
ADD COMMENTlink written 5.4 years ago by Michael Lawrence9.8k
Yes, I'm trying to get the strand of stranscription. If the XS flags are correct (I assume they will be if I run tophat with the appropriate library-type?) then the export will also be on the right strand? Thanks! Igor. On Wed, Jun 20, 2012 at 3:59 PM, Michael Lawrence <lawrence.michael at="" gene.com=""> wrote: > > > On Wed, Jun 20, 2012 at 4:53 AM, Igor Ulitsky [guest] > <guest at="" bioconductor.org=""> wrote: >> >> >> Hi, >> >> Is there a simple way to make strand-specific .wig file (i.e., a separate >> track for + and - strand) from paired-end data (where the second read maps >> to the other strand)? I've tried using this: >> >> library(Rsamtools) >> library(rtracklayer) >> myReads <- readGappedAlignments("RNAseqMapping.bam") >> coveragePlus <- coverage(myReads[strand(myReads) == ?'+']) >> export(coveragePlus, "RNAplus.wig") >> coverageMinus <- coverage(myReads[strand(myReads) == ?'-']) >> export(coverageMinus, "RNAminus.wig") >> >> But it appears that the second read in the pair contributes to the other >> strand, generating similar tracks for the + and the - strands. >> Is there a way to deal with this better? >> > > Are you trying to generate coverages for the actual strand of transcription? > If so, you would probably get that information from the XS tag and set it as > your strand prior to export, but unless you used a special protocol the XS > information would be incomplete. Btw, I would recommend BigWig export of > those coverage tracks. > > Michael > >> >> Thanks! >> >> Igor. >> >> ?-- output of sessionInfo(): >> >> R version 2.13.1 (2011-07-08) >> Platform: i386-pc-mingw32/i386 (32-bit) >> >> locale: >> [1] LC_COLLATE=English_United States.1252 >> [2] LC_CTYPE=English_United States.1252 >> [3] LC_MONETARY=English_United States.1252 >> [4] LC_NUMERIC=C >> [5] LC_TIME=English_United States.1252 >> >> attached base packages: >> [1] stats ? ? graphics ?grDevices utils ? ? datasets ?methods ? base >> >> loaded via a namespace (and not attached): >> [1] tools_2.13.1 >> >> -- >> Sent via the guest posting facility at bioconductor.org. > >
ADD REPLYlink written 5.4 years ago by Igor Ulitsky10
0
gravatar for Michael Lawrence
5.4 years ago by
United States
Michael Lawrence9.8k wrote:
As long as you change the strand assignment to the XS values, it should be what you want. Michael On Wed, Jun 20, 2012 at 6:06 AM, Igor Ulitsky <ulitskyi@gmail.com> wrote: > Yes, I'm trying to get the strand of stranscription. If the XS flags > are correct (I assume they will be if I run tophat with the > appropriate library-type?) then the export will also be on the right > strand? > > Thanks! > > Igor. > > On Wed, Jun 20, 2012 at 3:59 PM, Michael Lawrence > <lawrence.michael@gene.com> wrote: > > > > > > On Wed, Jun 20, 2012 at 4:53 AM, Igor Ulitsky [guest] > > <guest@bioconductor.org> wrote: > >> > >> > >> Hi, > >> > >> Is there a simple way to make strand-specific .wig file (i.e., a > separate > >> track for + and - strand) from paired-end data (where the second read > maps > >> to the other strand)? I've tried using this: > >> > >> library(Rsamtools) > >> library(rtracklayer) > >> myReads <- readGappedAlignments("RNAseqMapping.bam") > >> coveragePlus <- coverage(myReads[strand(myReads) == '+']) > >> export(coveragePlus, "RNAplus.wig") > >> coverageMinus <- coverage(myReads[strand(myReads) == '-']) > >> export(coverageMinus, "RNAminus.wig") > >> > >> But it appears that the second read in the pair contributes to the other > >> strand, generating similar tracks for the + and the - strands. > >> Is there a way to deal with this better? > >> > > > > Are you trying to generate coverages for the actual strand of > transcription? > > If so, you would probably get that information from the XS tag and set > it as > > your strand prior to export, but unless you used a special protocol the > XS > > information would be incomplete. Btw, I would recommend BigWig export of > > those coverage tracks. > > > > Michael > > > >> > >> Thanks! > >> > >> Igor. > >> > >> -- output of sessionInfo(): > >> > >> R version 2.13.1 (2011-07-08) > >> Platform: i386-pc-mingw32/i386 (32-bit) > >> > >> locale: > >> [1] LC_COLLATE=English_United States.1252 > >> [2] LC_CTYPE=English_United States.1252 > >> [3] LC_MONETARY=English_United States.1252 > >> [4] LC_NUMERIC=C > >> [5] LC_TIME=English_United States.1252 > >> > >> attached base packages: > >> [1] stats graphics grDevices utils datasets methods base > >> > >> loaded via a namespace (and not attached): > >> [1] tools_2.13.1 > >> > >> -- > >> Sent via the guest posting facility at bioconductor.org. > > > > > [[alternative HTML version deleted]]
ADD COMMENTlink written 5.4 years ago by Michael Lawrence9.8k
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 2.2.0
Traffic: 120 users visited in the last hour