Q: MACS3 Broad Peak Cutoff Analysis

[infurnoheat]

Hi,

We are currently working on paired-end data generated from CUT&Tag across all 4 histones. Our sample call for our broad peaks was:

macs3 callpeak -f BAMPE -g hs --keep-dup all --broad --nolambda --nomodel -q 0.01 --broad-cutoff 0.1

As recommended in the documentation, we employed the elbow method cutoff analysis for our data, resulting in the above cutoff values. The average of H3K27me3, H3K36me3, and H3K9me3 was about 2. As such, I used 0.01 for the q cutoff. However, I am unsure if that value should be applied to the q-cutoff or the broad-cutoff (which is currently just the default). Looking at our data, we do have a higher number of peaks than expected for H3K36me3 and H3K27me3.

It is important to note that in our data, many of the avelpeak elbow points differed significantly from the rest of the data or were undeterminable. If you have any further advice or recommendations, that would be greatly appreciated.

[taoliu]

@infurnoheat Your usage is right that 0.01 should be applied to -q. However for general 'broad peak' calling, I suggest you avoid using '--broad' since it can be confusing and the resulting nested peak may not be what you want (also confusing). Instead, you can just use macs3 callpeak -f BAMPE -g hs --keep-dup all --nolambda --nomodel -q 0.01 in your case. When you pick the cutoff, the avelpeak (average length of peak) is one factor, but you may also need to add some other factors to help make decision. For example, I know H3K36me3 normally can be found at the exon region of expressed genes, so the npeaks (number of peaks) should be around 10^6 (average 10 exons per gene and about 10k genes expressed in the cell). I don't have recommendation for the repressive mark H3K27me3 and it's always a challenge to get good results from H3K27me3 ChIP-seq experimentally and computationally. Not sure if CUT&Tag H3k27me3 data can be much cleaner than ChIP-seq.