What is the input format to peak calling?

[adamklie]

I've never felt comfortable with what combination of -f/--format parameter and file type to feed to MACS callpeak. I've seen multiple GitHub issues (macs3-project/MACS#145, ENCODE-DCC/atac-seq-pipeline#265, stuart-lab/signac#682), Twitter threads (https://twitter.com/XiChenUoM/status/1336658454866325506) and even google group discussions (https://groups.google.com/g/macs-announcement/c/4OCE59gkpKY/m/v9Tnh9jWriUJ?pli=1) that touch on this (non-exhaustive list). From my understanding, the ENCODE "best practice" is to call peaks on tagAlign files with -f BED parameter (and some other important parameter choices).

I've been rummaging through the SnapATAC2 code base to try to figure out the format of the files on which peaks are called, but my rust is pretty shaky. It looks like the export_tags function is called here:

SnapATAC2/snapatac2-python/python/snapatac2/tools/_call_peaks.py

Line 129 in 17537d6

fragments = _snapatac2.export_tags(adata, tmpdirname, groupby, replicate, max_frag_size, selections)

I also can't figure out what format MACS3 PeakDetect is expecting.

@kaizhang can you provide some details on this.

As an aside. it would be great to have the best practices for this formalized somewhere (e.g. a spec for tagAlign would be a good start!), though I know peak calling has it's issues. Folks with expertise, feel free to chime in!

[kaizhang]

In MACS, there are different parsers for different format that are used to parse the read information and to create the pileup signals (FWTrack object). The pileup signals are then given to PeakDetect to call peaks.

SnapATAC2 exports fragments (paired reads), converts fragments to single-ended reads (to be consistent with most literature), and finally creates FWTrack object directly from single-ended reads and call peaks (bypassing MACS's built-in parser for efficiency reason).

[adamklie]

When you say single-ended reads, do you mean insertion sites? Otherwise this makes sense.

[kaizhang]

Yes, I meant insertion sites

[kaizhang]

In MACS, there are different parsers for different format that are used to parse the read information and to create the pileup signals (FWTrack object). The pileup signals are then given to PeakDetect to call peaks.

SnapATAC2 exports fragments (paired reads), converts fragments to single-ended reads (to be consistent with most literature), and finally creates FWTrack object directly from single-ended reads and call peaks (bypassing MACS's built-in parser for efficiency reason).