[Gmod-ajax] tracks for bam coverage

classic Classic list List threaded Threaded
3 messages Options
Reply | Threaded
Open this post in threaded view
|

[Gmod-ajax] tracks for bam coverage

Ke Jiang
Hi,

I'm trying to visualize bam coverage on JBrowse. I noticed this potential problem (see the attached picture): when doing this for mRNA-Seq data, the reads that spanning more than one exons could be correctly visualized using 'Alignments2', the partial reads are connected by a thin black line. But if using 'SNPCoverage' or 'FeatureCoverage', the coverage between the partial reads (introns) are the same as the parital reads, instead of zero, which may lead the confusion of certain coverage in the intron. So on JBrowse with mRNA-Seq data, we don't see the 'pulse' of coverage on exons with zero coverage in introns, instead, we see flat coverage in introns, reflecting the coverage at the exon/intron boundaries, which is misleading. I tried to convert bam to BigWig format for coverage visualization, but it seems this problem carries on with the bedGraph and BigWig format. I wonder is there a way to work around this?

Thanks!

Ke

Inline image 1

--
Ke Jiang, Ph.D. 
Post-doctoral Fellow
Delbruck Laboratory
Cold Spring Harbor Laboratory
Cold Spring Harbor, NY 11724


------------------------------------------------------------------------------
Meet PCI DSS 3.0 Compliance Requirements with EventLog Analyzer
Achieve PCI DSS 3.0 Compliant Status with Out-of-the-box PCI DSS Reports
Are you Audit-Ready for PCI DSS 3.0 Compliance? Download White paper
Comply to PCI DSS 3.0 Requirement 10 and 11.5 with EventLog Analyzer
http://pubads.g.doubleclick.net/gampad/clk?id=154622311&iu=/4140/ostg.clktrk
_______________________________________________
Gmod-ajax mailing list
[hidden email]
https://lists.sourceforge.net/lists/listinfo/gmod-ajax
Reply | Threaded
Open this post in threaded view
|

Re: tracks for bam coverage

vkrishna
Hi Ke,

Not sure why you’re seeing problems when visualizing RNA-seq coverage. On our Medicago JBrowse, we have several tracks which display read coverage, where the underlying data is in the form of BigWig files generated from BAM alignment files.

Here is an example region:
http://www.jcvi.org/medicago/jbrowse/?data=data%2Fjson%2Fmedicago&loc=chr1%3A26628677..26632655&tracks=gene_models%2CNodule%20RNA-seq%20Coverage&highlight=

Here is the protocol I use for the conversion from BAM to BigWig:
files=Mt*fastq
for file in $files; do
    outdir=`basename $file .fastq`
    if [ ! -d "$outdir" ]; then
        mkdir $outdir
    fi
    if [ ! -e "$outdir/accepted_hits.bam" ]; then
        qsub -N "$outdir" -P 0372 -cwd -pe threaded 16 -l medium -o $outdir/log -e $outdir/err \
            tophat2 -p 16 --output-dir ./$outdir \
                --min-intron-length 50 --max-intron-length 10000 \
                --library-type fr-unstranded \
                ./JCVI.Medtr.v4 $file
    else
        trackopts="name=\"$outdir RNA-seq coverage\""
        bedtools genomecov -ibam $outdir/accepted_hits.bam -g JCVI.Medtr.v4.fasta.sizes -bg -split > $outdir/$outdir.bedGraph
        bedGraphToBigWig $outdir/$outdir.bedGraph JCVI.Medtr.v4.fasta.sizes $outdir/$outdir.bw
    fi
done


I recall that it is important to specify the `-split` parameter for the bedtools genome coverage calculation step.

Hope this helps!
Vivek


On Sep 24, 2014, at 2:16 PM, Ke Jiang <[hidden email]> wrote:

Hi, 

I'm trying to visualize bam coverage on JBrowse. I noticed this potential problem (see the attached picture): when doing this for mRNA-Seq data, the reads that spanning more than one exons could be correctly visualized using 'Alignments2', the partial reads are connected by a thin black line. But if using 'SNPCoverage' or 'FeatureCoverage', the coverage between the partial reads (introns) are the same as the parital reads, instead of zero, which may lead the confusion of certain coverage in the intron. So on JBrowse with mRNA-Seq data, we don't see the 'pulse' of coverage on exons with zero coverage in introns, instead, we see flat coverage in introns, reflecting the coverage at the exon/intron boundaries, which is misleading. I tried to convert bam to BigWig format for coverage visualization, but it seems this problem carries on with the bedGraph and BigWig format. I wonder is there a way to work around this?

Thanks!

Ke

<image.png>

-- 
Ke Jiang, Ph.D. 
Post-doctoral Fellow
Delbruck Laboratory
Cold Spring Harbor Laboratory
Cold Spring Harbor, NY 11724

------------------------------------------------------------------------------
Meet PCI DSS 3.0 Compliance Requirements with EventLog Analyzer
Achieve PCI DSS 3.0 Compliant Status with Out-of-the-box PCI DSS Reports
Are you Audit-Ready for PCI DSS 3.0 Compliance? Download White paper
Comply to PCI DSS 3.0 Requirement 10 and 11.5 with EventLog Analyzer
http://pubads.g.doubleclick.net/gampad/clk?id=154622311&iu=/4140/ostg.clktrk_______________________________________________
Gmod-ajax mailing list
[hidden email]
https://lists.sourceforge.net/lists/listinfo/gmod-ajax


------------------------------------------------------------------------------
Meet PCI DSS 3.0 Compliance Requirements with EventLog Analyzer
Achieve PCI DSS 3.0 Compliant Status with Out-of-the-box PCI DSS Reports
Are you Audit-Ready for PCI DSS 3.0 Compliance? Download White paper
Comply to PCI DSS 3.0 Requirement 10 and 11.5 with EventLog Analyzer
http://pubads.g.doubleclick.net/gampad/clk?id=154622311&iu=/4140/ostg.clktrk
_______________________________________________
Gmod-ajax mailing list
[hidden email]
https://lists.sourceforge.net/lists/listinfo/gmod-ajax
Reply | Threaded
Open this post in threaded view
|

Re: tracks for bam coverage

Ke Jiang
Thanks, Vivek! '-split' did the trick! Actually, it was written upfront on the doc page for bedtools genomecov. I probably should have read the doc more carefully. However, the JBrowse 'FeatureCoverage' and 'SNPCoverage' still fill in the introns with coverage form partial reads spanning the junctions.

Ke

On Wed, Sep 24, 2014 at 3:20 PM, Krishnakumar, Vivek <[hidden email]> wrote:
Hi Ke,

Not sure why you’re seeing problems when visualizing RNA-seq coverage. On our Medicago JBrowse, we have several tracks which display read coverage, where the underlying data is in the form of BigWig files generated from BAM alignment files.

Here is an example region:
http://www.jcvi.org/medicago/jbrowse/?data=data%2Fjson%2Fmedicago&loc=chr1%3A26628677..26632655&tracks=gene_models%2CNodule%20RNA-seq%20Coverage&highlight=

Here is the protocol I use for the conversion from BAM to BigWig:
files=Mt*fastq
for file in $files; do
    outdir=`basename $file .fastq`
    if [ ! -d "$outdir" ]; then
        mkdir $outdir
    fi
    if [ ! -e "$outdir/accepted_hits.bam" ]; then
        qsub -N "$outdir" -P 0372 -cwd -pe threaded 16 -l medium -o $outdir/log -e $outdir/err \
            tophat2 -p 16 --output-dir ./$outdir \
                --min-intron-length 50 --max-intron-length 10000 \
                --library-type fr-unstranded \
                ./JCVI.Medtr.v4 $file
    else
        trackopts="name=\"$outdir RNA-seq coverage\""
        bedtools genomecov -ibam $outdir/accepted_hits.bam -g JCVI.Medtr.v4.fasta.sizes -bg -split > $outdir/$outdir.bedGraph
        bedGraphToBigWig $outdir/$outdir.bedGraph JCVI.Medtr.v4.fasta.sizes $outdir/$outdir.bw
    fi
done


I recall that it is important to specify the `-split` parameter for the bedtools genome coverage calculation step.

Hope this helps!
Vivek


On Sep 24, 2014, at 2:16 PM, Ke Jiang <[hidden email]> wrote:

Hi, 

I'm trying to visualize bam coverage on JBrowse. I noticed this potential problem (see the attached picture): when doing this for mRNA-Seq data, the reads that spanning more than one exons could be correctly visualized using 'Alignments2', the partial reads are connected by a thin black line. But if using 'SNPCoverage' or 'FeatureCoverage', the coverage between the partial reads (introns) are the same as the parital reads, instead of zero, which may lead the confusion of certain coverage in the intron. So on JBrowse with mRNA-Seq data, we don't see the 'pulse' of coverage on exons with zero coverage in introns, instead, we see flat coverage in introns, reflecting the coverage at the exon/intron boundaries, which is misleading. I tried to convert bam to BigWig format for coverage visualization, but it seems this problem carries on with the bedGraph and BigWig format. I wonder is there a way to work around this?

Thanks!

Ke

<image.png>

-- 
Ke Jiang, Ph.D. 
Post-doctoral Fellow
Delbruck Laboratory
Cold Spring Harbor Laboratory
Cold Spring Harbor, NY 11724

------------------------------------------------------------------------------
Meet PCI DSS 3.0 Compliance Requirements with EventLog Analyzer
Achieve PCI DSS 3.0 Compliant Status with Out-of-the-box PCI DSS Reports
Are you Audit-Ready for PCI DSS 3.0 Compliance? Download White paper
Comply to PCI DSS 3.0 Requirement 10 and 11.5 with EventLog Analyzer
http://pubads.g.doubleclick.net/gampad/clk?id=154622311&iu=/4140/ostg.clktrk_______________________________________________
Gmod-ajax mailing list
[hidden email]
https://lists.sourceforge.net/lists/listinfo/gmod-ajax




--
Ke Jiang, Ph.D. 
Post-doctoral Fellow
Delbruck Laboratory
Cold Spring Harbor Laboratory
Cold Spring Harbor, NY 11724


------------------------------------------------------------------------------
Meet PCI DSS 3.0 Compliance Requirements with EventLog Analyzer
Achieve PCI DSS 3.0 Compliant Status with Out-of-the-box PCI DSS Reports
Are you Audit-Ready for PCI DSS 3.0 Compliance? Download White paper
Comply to PCI DSS 3.0 Requirement 10 and 11.5 with EventLog Analyzer
http://pubads.g.doubleclick.net/gampad/clk?id=154622311&iu=/4140/ostg.clktrk
_______________________________________________
Gmod-ajax mailing list
[hidden email]
https://lists.sourceforge.net/lists/listinfo/gmod-ajax