call peak for all replicates in all conditions using macs2

Question

[QX]

Hi all,

I am doing ATAC-seq analysis, where I want to make a matrix of all consensus peaks from all replicates in all conditions. Can I just call all of them in one call using macs2?

macs2 callpeak --verbose 3 -t *.bed \
    -f BEDPE \
    -g ${macs_genomesize} -B \
    -q ${macs_q} \
    --extsize 200 --nomodel --shift -100 --nolambda \
    --keep-dup all --outdir ${peak_all_BEDPE_dir} -n "peak_all" --call-summits

or do I need to perform intersect for replicates and then merge for all conditions? what is a more properly to deal with multiple replicates and conditions

[rfran010]

I would call peaks for individual replicates and conditions rather than pool them together. This will also help you check replicate agreement.

[QX]

after calling individually, how can you call the consensus peaks for all conditions? can you share the command for that?

[zau saa]

Fully agree to rfran010's reply. And, compared to intersect peaks for all replicates, I think you may set consensus peaks for each condition to intervals which are more than 500bp and identified as peaks by at least 50% replicates.

[QX]

do you know how to set this 500 and 50%? is that a macs2 option?

[zau saa]

They're two independent steps. You need to call peak for each replicates and then use bedtools genomecov to compute coverage of each interval.

https://bedtools.readthedocs.io/en/latest/content/tools/genomecov.html