I have called peaks using MACS2 and now I'm trying to run ChIPQC bioconductor package. Can someone please tell me should I use sorted BAM files or unsorted BAM files? I have created sorted BAM files and from that indexes for BAM files. When I was trying run ChIPQC earlier, I got an error. Then I figured out that is because I don't have index files. Then I created indexes after sorting BAM files.
You need to provide sorted bam files, ideally with an index. If you run
ChIPQCsample, it calls
ChIPQC:::sampleQC which iterates through the chromosomes or contigs, using the index, and the absence of the index will crash it.
We can use an example, first in linux:
wget ftp://hgdownload.cse.ucsc.edu/goldenPath/hg19/encodeDCC/wgEncodeHaibTfbs/wgEncodeHaibTfbsA549RxlchV0416101Etoh02AlnRep1.bam -O Etoh02AlnRep1.bam samtools index ./Etoh02AlnRep1.bam # s samtools view -uh Etoh02AlnRep1.bam chr1 > sorted.bam samtools index sorted.bam # quick way to unsort samtools sort -n -o unsorted.bam sorted.bam # it fails to index samtools index unsorted.bam [E::hts_idx_push] Unsorted positions on sequence #1: 246814554 followed by 42922635 samtools index: failed to create index for "unsorted.bam"
Then in R:
library(ChIPQC) ChIPQCsample("unsorted.bam") Bam file has 25 contigs Error in value[[3L]](cond) : failed to open BamFile: failed to load BAM index file: unsorted.bam ChIPQCsample("sorted.bam",annotation=NULL) Bam file has 25 contigs Calculating coverage histogram for chr1 Calculating SSD for chr1 Calculating unique positions per strand for chr1 Calculating shift for chr1 300 / 300 Counting reads in features for chr1 ....