Visualizing Coverage Of Ngs Reads Across Reference Genome

Question

[nicolazilio]

I am interested in finding out at how well my genome is covered (depth of coverage) by the reads I got out on an Illumina run to determine if there are loci that are under-represented, or missing (i.e. large deletion) or over-represented. Is there any script that can do that? I already aligned the reads using BWA.

Thank you.

Answer 1

[KCC]

You probably have a file in SAM format. Convert to BAM (samtools):

samtools view -bS -F 4 data.sam > data.bam

Convert from BAM to BED (bedtools command):

bamToBed -i data.bam > data.bed

BED to Sorted BED (unix system command):

sort -k1,1 data.bed > data.sorted

Sorted BED to BedGraph (bedtools command):

genomeCoverageBed -bga -i data.sorted -g chrom.sizes > data.bedGraph

Convert to bigwig for fast viewing (UCSC utilities):

bedGraphToBigWig data.bedGraph chrom.sizes data.bw

Finally, use a program like IGV to view it. On a Mac, it's as easy as selecting the right genome, then drag and drop the bigwig file into the window. You can now explore your underrepresented features.

Note that the chrom.sizes file is a tab delimited file with the format: <chromosome name><tab><number of bases in chromosome>

[matted]

You can skip all the bed stuff and just load a sorted and indexed bam file directly into IGV. It shows a coverage track along with the reads, so it should give you the same result with fewer steps.

[Geparada]

Here are other ways to do it! http://onetipperday.blogspot.com/2012/07/three-ways-to-convert-bambed-file-to.html

Answer 2

[swbarnes2]

If you want to eyeball coverage, one easy way is with IGV, from the Broad institute.

[Charles Warden]

IGV works well, but make sure you have a coverage file. This can be created using igvtools (should be under "Run --> Run igvtools...")

[swbarnes2]

Like I said, you can eyeball coverage. Just throw the indexed .bam in there.

Answer 3

[nicolazilio]

Thank you, that worked!