I'm trying to count paired end ATAC-seq alignments that fall within peaks specified by a BED file for differential accessibility analyses across two conditions. Peaks were called using MACS2 with --shift -75 --extsize 150 --nomodel --keep-dup all --SPMR, to center reads around sites of transposase insertion. My question is about applying featureCounts to identify the number of transposition events in each BED-specified interval: I'm not supplying the -p arg to featureCounts because I don't care where the fragments are, but rather where the transposition sites fall. Assuming this is a correct way to think about this, subread v1.6.2 provides a --read2pos 5 option for featureCounts. Should I be using this to specifically count where the transposition events are occurring, rather than the overlap of library fragments over intervals? Alternatively, to keep things consistent with the MACS2 peak calling, is it better to call featureCounts with the same --shiftsize and --extsize parameters I used with MACS2?

I don't know what are the transposition events you are looking for. But if you are looking for the cleavage sites in the open chromatin regions, you can use the start position of reads to search such sites. The --read2pos 5 option in featureCounts can help you to achieve this.

The shifting and extending parameters in featureCounts and MACS2 have different meanings. I don't think values you provided for these parameters in one program can guide you to set parameter values in the other program.