Hi,
I am trying peak calling for H3K27ac ChIP-seq for enhancer and super enhancer analysis. Initially, I was planning to use default narrowPeak mode in macs2. However, I looked that some publications use --broad parameter and others use --no-model (they don't specify if they are using narrow/broad, but as no flag is specified I am assuming narrowPeak). I also found this analysis on github by crazyhottommy, and wanted to give a try and evaluate performance. Paired end fastq files were aligned using bowtie2, duplicates and no/multi mapper were removed using sambamba.
I used following commands for narrow, broad and broad no-model runs. Narrow:
macs2 callpeak -t ${ID}.sorted.dedup.bam \
-c INPUT.sorted.dedup.bam \
-f BAMPE -q 0.01 -g hs -n ${ID} --tempdir $tmpd 2> ${ID}_macs2.log
Broad:
macs2 callpeak -t ${ID}.sorted.dedup.bam \
-c INPUT.sorted.dedup.bam \
--broad -g hs -f BAMPE --broad-cutoff 0.1 \
-n ${ID}_b --tempdir $tmpd 2> ${ID}_macs2b.log
Broad No model
macs2 callpeak -t ${ID}.sorted.dedup.bam \
-c INPUT.sorted.dedup.bam \
--broad -g hs -f BAMPE --broad-cutoff 0.1 --nomodel --extsize 146 \
-n ${ID}_bnm --tempdir $tmpd 2> ${ID}_macs2bnm.log
Upon analysis, I see more peaks in broad and broad no-model, but I get exactly the same peaks for both. Am I not supposed to get at least some difference? Just wanted to see if I am doing something wrong?
I checked that the same thing happened with default call (narrow peak). Is it because I am using -f BAMPE (I have paired end data)? Does MACS2 not take into account --nomodel flag with BAMPE?
Thanks!