Do I need an input control for ATAC-seq peak calling?

No, there's no equivalent to a ChIP-seq input or IgG control in ATAC-seq because there's no antibody pull-down to background-correct against. Instead, the defense against artifacts is ENCODE blacklist filtering, FRiP and TSS enrichment thresholds, and a manual browser check of your top peaks against expected biology.

What FRiP score counts as good for ATAC-seq?

ENCODE's ATAC-seq standards call FRiP above 0.3 ideal for bulk ATAC-seq, with 0.2-0.3 acceptable but indicating elevated background. Below 0.2, most reads fall outside any called peak, and you should treat the peak set as unreliable until you've checked library complexity and fragment size.

Why are my strongest ATAC-seq peaks sitting in centromeres or repeat regions?

This is the classic signature of mappability and repeat-content artifacts: reads from repetitive or high-copy-number DNA pile up and confidently mismap to a single reference locus, which MACS's background model can't distinguish from real accessibility. Intersecting your peaks against the ENCODE blacklist for your exact genome build removes the worst of these before they reach downstream analysis.

Should I use narrow or broad peaks for ATAC-seq?

Default narrowPeak mode (no --broad) is what almost every ATAC-seq pipeline uses, because accessible chromatin regions are focal, nucleosome-free windows, closer in shape to a sharp TF footprint than to a broad histone domain. Reach for broad-peak mode only if you're specifically studying large accessible domains, and check the signal shape in a browser first rather than assuming.

MACS2 or MACS3 for ATAC-seq peak calling?

MACS3 is the actively maintained successor with Python 3 support and adds an ATAC-specific HMM-based caller (hmmratac) alongside standard callpeak. If you don't specifically need hmmratac, the standard BAMPE narrowPeak workflow works the same way in both versions; the practical reason to move is maintenance and long-term support rather than a documented accuracy jump.

Sanity check · ATAC-seq

How to Call Peaks You Can Trust in ATAC-seq

No input control exists for ATAC-seq, so ENCODE blacklists, FRiP, and a five-minute browser check of your top peaks are what stand between you and peaks called in centromeres.

By Ming "Tommy" Tang, Director of Bioinformatics in Big Pharma · Reviewed September 2026 · 5 min read

You ran macs2 callpeak on your ATAC-seq BAM the way you'd run it for ChIP-seq, except there was no -c control.bam to give it, because ATAC-seq has no input sample. The caller ran without complaint, handed back tens of thousands of peaks passing your q-value cutoff, and the first thing you notice when you sort by score is that your strongest "peak" sits inside a centromeric satellite repeat, not anywhere near a gene you'd expect for your sorted cell population.

Miss this before running differential accessibility or motif enrichment and you'll spend days chasing signal that isn't there: "differential" peaks driven by mappability differences between conditions, motif hits explained by Tn5 sequence bias rather than TF occupancy, and eventually a collaborator asking why your top hit is annotated to nothing. None of this shows up in a read-count QC report. A library can pass alignment rate and depth and still be unusable.

This page walks through the checks that catch a bad ATAC-seq peak set before it reaches your differential test: fragment size and TSS enrichment, which tell you if the library itself is any good; FRiP, which tells you if your peak calls actually capture where the signal is; and blacklist filtering plus a manual browser check, which catch what statistics alone won't. Work through them in order and you'll know within the hour whether your peak set is trustworthy or needs to be recalled.

What it looks like when it's happening

macs2 or macs3 callpeak run without a -c control finishes with tens of thousands of peaks, but the top peaks by score sit in centromeric satellite repeats or on chrM
FRiP score comes back below 0.2
Fragment size histogram is flat or missing the nucleosome-periodic ~147-294 bp bump
TSS enrichment score is below 5 (hg38) or 10 (mm10)
The same handful of peak coordinates show up in every sample and condition regardless of biology
IGV shows your highest-scoring peaks sitting over segmental duplications or low-mappability regions
Motif enrichment on the peak set returns no expected transcription factor motifs despite thousands of 'significant' peaks
Differential accessibility results are dominated by peaks near assembly gaps or unannotated regions

Why it happens

There is no equivalent of a ChIP-seq input or IgG control in ATAC-seq, because there's no antibody pull-down to background-correct against. Everything the assay sees becomes a candidate peak. Tn5 cuts nucleosome-free regions, but it also cuts nonspecifically wherever DNA looks accessible for non-biological reasons: high-copy-number regions, mitochondrial DNA fragments integrated into the nuclear genome (NUMTs), satellite repeats, and segmental duplications where reads confidently mismap to the wrong locus. Without a background sample to subtract, MACS's local background model calibrates against the library's own genome-wide read distribution, so if the library itself carries chrM contamination or repeat artifacts, that noise gets baked into what statistically looks "significant."

The ENCODE blacklist is the closest thing ATAC-seq has to a control. It was built by scanning many independent experiments across cell types and assays for regions with anomalous, unstructured, or consistently high signal regardless of biology. If a region is a hotspot everywhere, it isn't real accessibility, it's a technical byproduct of repeat content, low mappability, or assembly gaps. That's an empirically defined danger list, not a matched biological negative, but it's what stands between you and peaks called on artifacts.

Peak-calling parameters add a second failure mode. MACS2/3 build a model by default assuming ChIP-seq fragment structure, paired peaks flanking a true binding site. That model doesn't fit ATAC, where insertions happen directly at accessible DNA rather than around it. --nomodel with --shift/--extsize forces pileups centered on the Tn5 cut site instead of an inferred fragment center. But that recipe is written for single-end, cut-site input. Feed MACS true paired-end fragments with -f BAMPE and the tool already has the real insert size for every fragment; by its own documentation it won't let you override --shift away from 0. Mixing the two conventions, common in copy-pasted tutorials, produces a command MACS partially ignores.

FRiP failures have a simpler mechanism: a library where too many fragments come from PCR duplicates or open, non-regulatory chromatin spreads reads genome-wide instead of concentrating them in narrow accessible windows. The peak caller still finds loci where local density edges out the genome-wide average, and those pass a q-value threshold, but the threshold is calibrated against a noisy background, not a clean one.

The checks

Run them in order. Each one tells you what healthy looks like and what the problem looks like.

0/8 checked · saved in this browser

01Alignment rate, chrM fraction, and duplication rate
Run samtools flagstat on the aligned BAM, then count reads mapping to the chrM contig specifically before you do anything else. Pull NRF, PBC1, and PBC2 library-complexity metrics from the filtered BAM.
bash
```
samtools flagstat aligned.bam
samtools view -c aligned.bam chrM
samtools view -c aligned.bam
```
Healthy
>95% of reads aligned, chrM reads a minor fraction removed before peak calling, NRF >0.9, PBC1 >0.9, PBC2 >3.
Red flag
Alignment rate below 80%, a large chrM fraction inflating apparent depth, or NRF/PBC1 below 0.9 meaning much of your coverage is PCR duplicates collapsing onto a handful of positions.
02Fragment size distribution shows nucleosome periodicity
Generate an insert-size histogram from the paired-end BAM and look for the nucleosome-periodic pattern: a nucleosome-free peak below ~100 bp, a mononucleosome peak in the 147-294 bp range, then fainter di- and tri-nucleosome bumps.
bash
```
picard CollectInsertSizeMetrics I=aligned.bam O=insert_size_metrics.txt H=insert_size_histogram.pdf
```
Healthy
A histogram with a clear nucleosome-free peak, a mononucleosome bump around 147-294 bp, and successively smaller multi-nucleosome bumps.
Red flag
A flat or unimodal distribution with no periodicity, meaning Tn5 tagmentation didn't respect nucleosome structure, usually a failed, over-, or under-tagmented library.
03TSS enrichment score
Aggregate ATAC signal around a reference set of TSS positions and compute the ratio of signal at the center to signal in the flanking background, using the ENCODE TSS enrichment approach or an equivalent metagene profile centered on a TSS BED file.
Healthy
hg38: >7 ideal, 5-7 acceptable. mm10: >15 ideal, 10-15 acceptable.
Red flag
hg38 score below 5 or mm10 below 10, meaning signal isn't concentrating at active promoters, usually pointing to dead or degraded nuclei going into the tagmentation reaction.
04Confirm which MACS mode you actually ran
Check whether your peak-calling command used true paired-end fragments (-f BAMPE) or single-end cut-site mode (-f BED/BAM with --nomodel --shift -100 --extsize 200), and don't mix them: MACS documentation states --shift cannot be set to a nonzero value when the format is BAMPE or BEDPE.
bash
```
# paired-end, true fragments (recommended default)
macs3 callpeak -f BAMPE -t ATAC.bam -g hs -n test -B -q 0.01

# single-end / cut-site-centered mode
macs2 callpeak -f BAM --nomodel --shift -100 --extsize 200 -t file.bam -g hs -n peaks -q 0.05
```
Healthy
The command matches the input format: BAMPE with no --shift override, or single-end BAM/BED explicitly paired with --nomodel --shift -100 --extsize 200.
Red flag
A copy-pasted command combining -f BAMPE with --shift -100 --extsize 200, producing peaks that silently aren't shifted the way you assumed.
05FRiP score
Count reads overlapping your called peaks and divide by total usable reads in the library, post-dedup and post-chrM removal.
bash
```
reads_in_peaks=$(samtools view -c -L peaks.narrowPeak aligned.bam)
total_reads=$(samtools view -c aligned.bam)
echo "scale=3; $reads_in_peaks / $total_reads" | bc
```
Healthy
FRiP above 0.3 for a solid bulk ATAC-seq library; 0.2-0.3 acceptable but flags elevated background.
Red flag
FRiP below 0.2: most sequencing lands outside any peak, meaning your 'significant' peaks are riding on top of genome-wide noise.
06Blacklist overlap on your peak set
Download the ENCODE blacklist BED matching your exact genome build from the Boyle-Lab Blacklist repository, then intersect it against your peaks and check both how many peaks it removes and whether any top-scoring peaks are in the list.
bash
```
bedtools intersect -v -a peaks.narrowPeak -b blacklist.bed > peaks.filtered.narrowPeak
wc -l peaks.narrowPeak peaks.filtered.narrowPeak
```
Healthy
Blacklisted regions account for a small fraction of total peaks (well under a few percent), and none of the top 20 peaks by score.
Red flag
A meaningful chunk of peaks disappear after filtering, or your single highest-scoring peak overlaps a blacklisted region.
07Browser check of top-scoring peaks against expected biology
Sort your peak file by score (or -log10(q) in narrowPeak format) and load the top 20-50 alongside your BAM/bigWig tracks in IGV. Ask whether they sit over promoters, enhancers, or genes plausible for your sorted cell population.
Healthy
Top peaks correspond to accessible regions near genes with a defensible biological story for your cell type or treatment.
Red flag
Top peaks sit in centromeric satellite arrays, segmental duplications, or gene deserts with no annotation nearby.
08Replicate concordance via IDR
Call peaks separately per replicate and on the pooled data, run IDR between replicate pairs, and compare the rescue ratio and self-consistency ratio.
Healthy
Both ratios under 2, meaning replicates agree on which regions are accessible.
Red flag
Ratios above 2, or replicate peak sets barely overlapping, meaning your 'differential accessibility' calls may just be replicate noise.

What to do about it

Filter every peak set against the ENCODE blacklist for your exact genome build

When: Before any peaks feed into counting, differential testing, or motif enrichment, for any ATAC-seq or ChIP-seq experiment.

Download the matching blacklist BED from the Boyle-Lab Blacklist repository (dm3, dm6, ce10, ce11, mm10, hg19, hg38 available) and run bedtools intersect -v -a peaks.bed -b blacklist.bed > peaks.filtered.bed.

Caveat: The blacklist is build-specific. Running the hg19 blacklist against hg38 peaks (or vice versa) silently filters almost nothing, and you'll believe the peak set is clean when it isn't.

Remove chrM reads and PCR duplicates before peak calling, not after

When: At the alignment/filtering stage of the pipeline, before Tn5 shift correction and peak calling.

Filter chrM reads during BAM filtering (exclude the mitochondrial contig with samtools view) and mark/remove duplicates with Picard MarkDuplicates, then recompute NRF/PBC1/PBC2 on the filtered BAM before calling peaks.

Caveat: Doing this after peak calling means your peak scores, FRiP, and TSS enrichment were all computed on a contaminated read pool, so you have to rerun the whole QC chain anyway. Do it first.

Match MACS shift/extsize settings to your actual input format

When: Calling peaks with MACS2 or MACS3 on paired-end ATAC-seq BAMs, especially when copying a command from a tutorial.

For true paired-end fragments, use -f BAMPE with no --shift override, e.g. macs3 callpeak -f BAMPE -t ATAC.bam -g hs -n test -B -q 0.01. For the classic Tn5 cut-site-centered pileup, convert to single-end BED/BAM input and use --nomodel --shift -100 --extsize 200, since MACS documentation states shift can't be overridden away from 0 in BAMPE/BEDPE mode.

Caveat: The two modes produce differently shaped and differently wide peaks. Don't call some samples in one mode and others in the other within the same comparison; the widths won't be comparable.

Build a consensus peak set from replicates with IDR

When: Any cross-condition comparison where you have, or can get, at least two biological replicates per condition.

Call peaks per replicate and on pooled reads, run IDR between replicate pairs, and keep only peaks meeting the rescue and self-consistency ratio cutoffs (<2) as your consensus set for counting and differential testing.

Caveat: IDR needs at least two real biological replicates; with n=1 per condition you have no way to distinguish reproducible signal from library-specific artifact.

Eyeball your top-scoring peaks in a browser before trusting the list

When: Right after peak calling, before building a count matrix or running motif enrichment, on any new dataset or cell type.

Sort peaks by score, load the top 20-50 with bigWig tracks in IGV, and check them against genes expected to be accessible in that cell population, the same way you'd check an ER ChIP for signal at TFF1 or GREB1.

Caveat: Doesn't scale past a handful of samples and won't catch subtle genome-wide biases; it's a fast sanity check for gross artifacts, not a replacement for FRiP, TSS enrichment, and blacklist filtering.

When not to "fix" it

Not every peak in an unusual location is an artifact worth discarding. Segmental duplications and structurally complex loci, like parts of the HLA region or amplified loci in cancer and immortalized cell lines, can carry genuine biological accessibility differences that happen to overlap regions flagged by mappability-based blacklists. If your question is specifically about copy-number-altered or structurally variable loci, blanket blacklist removal can throw away the signal you're studying; check whether the region's problem is mappability-driven misassignment or a true amplification in your sample, ideally with orthogonal WGS/CNV evidence, before you drop it. Likewise, a muted or asymmetric fragment-size periodicity isn't automatically a failed library in very low-input or FACS-sorted samples with only a few thousand cells. Before rejecting the library, check FRiP and TSS enrichment together; those two agreeing the signal is real matters more than fragment-size shape alone in a low-input regime.

Five things experienced analysts do here

Pin the blacklist file version to your genome build in the pipeline config, not just by species name, so a hg38-vs-patch mismatch can't slip through unnoticed.
Compute FRiP and TSS enrichment before you even look at total peak count; MACS will happily hand you hundreds of thousands of 'significant' peaks from a low-complexity library.
Log peak counts before and after blacklist filtering and chrM/duplicate removal for every sample, so a sudden jump signals a build mismatch or a broken filtering step immediately.
Rerun the blacklist intersect any time you change genome build or patch version; this is one of the most common silent reproducibility breaks in ATAC pipelines.
Never skip the manual browser look at your top 20-50 peaks, even when every automated metric passes; it catches per-sample weirdness that aggregate QC numbers average away.

Questions people ask

Do I need an input control for ATAC-seq peak calling?: No, there's no equivalent to a ChIP-seq input or IgG control in ATAC-seq because there's no antibody pull-down to background-correct against. Instead, the defense against artifacts is ENCODE blacklist filtering, FRiP and TSS enrichment thresholds, and a manual browser check of your top peaks against expected biology.
What FRiP score counts as good for ATAC-seq?: ENCODE's ATAC-seq standards call FRiP above 0.3 ideal for bulk ATAC-seq, with 0.2-0.3 acceptable but indicating elevated background. Below 0.2, most reads fall outside any called peak, and you should treat the peak set as unreliable until you've checked library complexity and fragment size.
Why are my strongest ATAC-seq peaks sitting in centromeres or repeat regions?: This is the classic signature of mappability and repeat-content artifacts: reads from repetitive or high-copy-number DNA pile up and confidently mismap to a single reference locus, which MACS's background model can't distinguish from real accessibility. Intersecting your peaks against the ENCODE blacklist for your exact genome build removes the worst of these before they reach downstream analysis.
Should I use narrow or broad peaks for ATAC-seq?: Default narrowPeak mode (no --broad) is what almost every ATAC-seq pipeline uses, because accessible chromatin regions are focal, nucleosome-free windows, closer in shape to a sharp TF footprint than to a broad histone domain. Reach for broad-peak mode only if you're specifically studying large accessible domains, and check the signal shape in a browser first rather than assuming.
MACS2 or MACS3 for ATAC-seq peak calling?: MACS3 is the actively maintained successor with Python 3 support and adds an ATAC-specific HMM-based caller (hmmratac) alongside standard callpeak. If you don't specifically need hmmratac, the standard BAMPE narrowPeak workflow works the same way in both versions; the practical reason to move is maintenance and long-term support rather than a documented accuracy jump.

Guide · How to Call Peaks You Can Trust in ChIP-seq
Guide · How to Call Peaks You Can Trust in Single-Cell ATAC-seq
Glossary · Count matrix
Glossary · Peak calling
Glossary · Q-value

Sources

ATAC-seq Data Standards and Processing Pipeline (ENCODE4) — FRiP, TSS enrichment, NRF/PBC, and peak count thresholds cited in this guide
MACS3 callpeak Documentation — Definitions of --shift, --extsize, --nomodel and the BAMPE/BEDPE shift restriction
The ENCODE Blacklist: Identification of Problematic Regions of the Genome — Blacklist region sources: NUMTs, repeats, low-mappability regions, assembly gaps
ATAC-seq Data Analysis, NBIS Epigenomics Workshop 2025 — Practical peak-calling parameter walkthrough referenced for shift/extsize usage

Part of the Peak calling controls and blacklists series.