# Bio Applied Scatac Chromatin

> TF-IDF normalize scATAC-seq peak matrices and run LSI/SVD, dropping depth-correlated component 1, via SnapATAC2 or Signac/Seurat. Use when clustering 10x fragments.tsv.gz, running LSI/UMAP on ATAC data, or linking peaks to genes.

- Skill: `pavel-kravchenko/bio-applied-scatac-chromatin` (Agent Skill)
- Install (CLI): `npx skillmds@latest add pavel-kravchenko/bio-applied-scatac-chromatin`
- Raw SKILL.md: https://api.skillmd.com/api/skills/pavel-kravchenko/bio-applied-scatac-chromatin/raw
- Safety review: pending
- Works with: Claude Code, Claude.ai, OpenAI Codex
- Category: Web & Frontend
- Author: pavel-kravchenko (https://skillmd.com/u/pavel-kravchenko)
- Updated: 2026-09-17
- Page: https://skillmd.com/skills/pavel-kravchenko/bio-applied-scatac-chromatin

---


# Single-Cell ATAC-seq: Chromatin Accessibility

**References:** [Signac](https://stuartlab.org/signac/) | [SnapATAC2](https://kzhang.org/SnapATAC2/) | [ArchR](https://www.archrproject.com/) | [chromVAR](https://www.nature.com/articles/nmeth.4401)

scATAC-seq measures chromatin accessibility at single-cell resolution — a direct readout of where transcription factors bind and regulatory decisions are made, upstream of expression changes.

## When to Use

- Processing a 10x Genomics scATAC-seq or multiome `fragments.tsv.gz` file into a cell-by-peak matrix.
- Normalizing sparse binary accessibility data (library-size normalization doesn't apply) and clustering cells.
- Diagnosing why PC1/LSI1 correlates with sequencing depth instead of biology.
- Linking distal peaks to genes (co-accessibility, peak-gene links) or scoring TF motif activity (chromVAR-style).
- QC'ing an ATAC library with TSS enrichment score or nucleosome fragment-length banding.

## Version Compatibility

Python ≥3.10, SnapATAC2 ≥2.6, anndata ≥0.10, scikit-learn ≥1.3, umap-learn ≥0.5. R ≥4.3, Signac ≥1.13, Seurat ≥5.0. Cicero (R/monocle3) for co-accessibility. chromVAR ≥1.24 (Bioconductor) for motif deviations.

## Prerequisites

`pip install snapatac2 anndata scikit-learn umap-learn scipy pandas matplotlib` or in R `install.packages("Signac"); BiocManager::install("chromVAR")`. Assumes familiarity with `bio-single-cell-preprocessing` (AnnData basics) and `bio-genome-intervals-bed-file-basics` (BED coordinates). A pseudo-bulk peak set (e.g. from MACS3) is a prerequisite input for the peak-by-cell matrix.

## Technology

**Tn5 tagmentation:** Hyperactive Tn5 inserts into nucleosome-free (accessible) regions, fragmenting DNA and ligating sequencing adapters. For 10x: isolate nuclei first (removes mtDNA), then add Tn5 inside droplets.

**Fragment file format** (TSV.gz, Tabix-indexed):
```text
chr1    10000    10200    ACGTCAGTACGT-1    1
```
Columns: chrom, start, end, cell barcode, read count. Compact and cell-indexed — faster than BAM for per-cell operations.

**QC metrics:**
1. **TSS enrichment score** — ratio of coverage at TSS vs. flanking background. Score > 4 acceptable; > 8 high quality.
2. **Nucleosome banding** — fragment length histogram shows mono- (147 bp), di- (294 bp), tri-nucleosome peaks + sub-nucleosomal (<147 bp) free DNA.
3. **Fragments per cell** — 10,000–50,000 unique fragments typical for 10x.

| Feature | Bulk ATAC-seq | scATAC-seq |
|---------|--------------|------------|
| Input | 50,000–500,000 cells | Single nuclei |
| Sparsity | Dense (~100% at peaks) | Sparse (binary: 0 or 1) |
| Peak calling | Directly on BAM | Requires pseudo-bulk |

## Peak Calling and TF-IDF Normalization

**Why not call peaks per cell:** each cell has only 0–1 reads per position — insufficient signal. Solution: pseudo-bulk peak calling (cluster cells → merge fragments per cluster → run MACS3 → build a master peak set).

**TF-IDF** handles sparse binary ATAC data where library-size normalization is meaningless:
- TF = fragments in peak / total fragments in cell (depth normalization)
- IDF = log(1 + n_cells / cells_with_peak_accessible) (downweights ubiquitous peaks)
- Final = TF × IDF — emphasizes cell-type-specific peaks

SnapATAC2 applies TF-IDF to raw fragment counts; ArchR binarizes (0/1) by default. Both work; binary is more robust to PCR duplicates.

**Goal:** build a cell-by-peak matrix and cluster cells by chromatin accessibility.
**Approach:** load fragments, bin the genome into tiles, select variable features, run TF-IDF + spectral (LSI) embedding, then UMAP/Leiden.

```python
import snapatac2 as snap

def cluster_scatac(fragment_file: str, bin_size: int = 500, n_comps: int = 30):
    """Build a tile matrix from a fragments.tsv.gz, normalize, and cluster.

    Returns an AnnData with obsm['X_spectral'] (LSI) and obs['leiden'].
    """
    data = snap.read(fragment_file, backed="r")
    snap.pp.make_tile_matrix(data, bin_size=bin_size)
    snap.pp.select_features(data, n_features=50000)
    snap.tl.spectral(data, n_comps=n_comps)  # TF-IDF + SVD (LSI) internally
    snap.pp.knn(data)
    snap.tl.leiden(data)
    snap.tl.umap(data)
    return data
```

### Signac (R)
```r
library(Signac); library(Seurat)

chrom_assay <- CreateChromatinAssay(
  counts = peak_counts, fragments = "fragments.tsv.gz",
  genome = "hg38", min.cells = 10
)
seurat_obj <- CreateSeuratObject(counts = chrom_assay, assay = "ATAC")
seurat_obj <- RunTFIDF(seurat_obj)
seurat_obj <- FindTopFeatures(seurat_obj, min.cutoff = "q75")
seurat_obj <- RunSVD(seurat_obj)
DepthCor(seurat_obj)  # LSI_1 usually shows high depth correlation -> drop it
seurat_obj <- RunUMAP(seurat_obj, reduction = "lsi", dims = 2:30)
```

## LSI Dimensionality Reduction

**Why PCA fails:** scATAC matrices are binary/sparse; PC1 correlates with read depth (technical artifact), not biology.

**LSI** (SVD on the TF-IDF matrix): compute TF-IDF, apply truncated SVD, **discard component 1** (almost always depth-correlated), then use components 2–30 for UMAP/clustering.

**Goal:** confirm LSI1 is a depth artifact and remove it before downstream analysis.
**Approach:** correlate each LSI component against total fragment count per cell; drop any component with |r| > 0.5 (normally just component 1).

```python
import numpy as np
from sklearn.decomposition import TruncatedSVD
from scipy.stats import pearsonr


def tfidf_lsi(counts: np.ndarray, n_components: int = 30):
    """TF-IDF normalize a binary cell-by-peak matrix and run LSI (SVD).

    counts: (n_cells, n_peaks) binary/count matrix.
    Returns (X_lsi, depth_corr) — depth_corr[i] is Pearson r of LSI comp i
    with per-cell fragment depth; drop components with |r| > 0.5.
    """
    n_cells = counts.shape[0]
    cell_totals = counts.sum(axis=1, keepdims=True) + 1e-6
    tf = counts / cell_totals
    cells_with_peak = (counts > 0).sum(axis=0) + 1
    idf = np.log1p(n_cells / cells_with_peak)
    tfidf = tf * idf[np.newaxis, :]

    svd = TruncatedSVD(n_components=n_components, random_state=42)
    x_lsi = svd.fit_transform(tfidf)

    depth = counts.sum(axis=1)
    depth_corr = np.array([pearsonr(x_lsi[:, i], depth)[0] for i in range(n_components)])
    return x_lsi, depth_corr


# Example: drop components correlated with depth, keep the rest
# x_lsi, corr = tfidf_lsi(counts_matrix)
# keep = np.where(np.abs(corr) <= 0.5)[0]
# x_lsi_clean = x_lsi[:, keep]
```

## Co-accessibility and Peak-Gene Links

**Cicero** (Pliner et al. 2018): trains graphical LASSO on pseudo-cells to compute co-accessibility scores (0–1) within 500 kb. Score > 0.25 suggests a regulatory connection.

### Peak-gene linking (Signac)
```r
seurat_obj <- LinkPeaks(
  object = seurat_obj,
  peak.assay = "ATAC",
  expression.assay = "RNA",
  distance = 5e5
)
```

**Validation:** Hi-C/HiChIP (3D proximity), eQTL overlap, ChromHMM enhancer states, ENCODE cCREs.

**Limitations:** co-accessibility requires >500 cells; correlation ≠ causation; cell-type mixing creates spurious co-accessibility.

## Pitfalls

- **Coordinate systems**: BED uses 0-based half-open; VCF/GFF use 1-based inclusive — mixing causes off-by-one errors.
- **Depth confound**: never use raw LSI1/PC1 for clustering without checking `DepthCor`/depth correlation first.
- **Batch effects**: always check for batch confounding (e.g. by sample/lane) before interpreting biological signal.
- **Multiple testing**: apply FDR correction (Benjamini-Hochberg) when testing thousands of peaks or peak-gene links simultaneously.

## See Also

- `bio-single-cell-scatac-analysis` — broader single-cell ATAC preprocessing workflow.
- `bio-atac-seq-atac-peak-calling` — MACS3 pseudo-bulk peak calling.
- `bio-atac-seq-motif-deviation` — chromVAR TF motif deviation scores.
- `bio-workflows-multiome-pipeline` — joint RNA+ATAC multiome pipeline.