Skip to main content

ATAC-seq Accessibility Heatmap

Author

[Editor] ZCode (GLM) โ€” AI-generated, independently AI-reviewed;

Modified

2026-11-05

๐Ÿค– AI Skill โ€” Download the Bizard unified skill for your AI assistant โฌ‡๏ธ Download Skill ZIP
๐Ÿ”Œ MCP Server โ€” Let AI agents (Claude Code, Cursor, โ€ฆ) query Bizard live: claude mcp add bizard -- npx @openbiox/bizard-mcp โ€” or the free hosted instance (no install, faster from mainland China): claude mcp add --transport http bizard http://biotree.top:38123/bizard-mcp/mcp ๐Ÿ“– MCP Docs
๐Ÿ“– Cite Bizard โ€” if Bizard helps your research, please cite: Li, K., Zheng, H., Huang, K., et al. (2026). Bizard: A Community-Driven Platform for Accelerating and Enhancing Biomedical Data Visualization. iMetaMed, e70038. ๐Ÿ”— DOI / BibTeX

Chromatin-accessibility experiments (ATAC-seq, DNase-seq) are reported as aggregate heatmaps over genomic regions: rows are peaks, columns are binned signal per condition, sorted so induced and unresponsive regions separate [1]. This tutorial draws that layout on simulated bigWig-style binned matrices โ€” the plotting half of the deeptools workflow (align โ†’ MACS2 โ†’ multiBigwigSummary โ†’ this figure), with the locus-level difference view the reference paper highlights.

Example

AtacSeqPlot DEMO

Setup

  • System Requirements: Cross-platform (Linux/MacOS/Windows)

  • Programming Language: R

  • Dependent packages: ComplexHeatmap; ggplot2; dplyr; tidyr; circlize

# Installing necessary packages
if (!requireNamespace("ComplexHeatmap", quietly = TRUE)) {
  BiocManager::install("ComplexHeatmap")
}
if (!requireNamespace("ggplot2", quietly = TRUE)) {
  install.packages("ggplot2")
}
if (!requireNamespace("dplyr", quietly = TRUE)) {
  install.packages("dplyr")
}
if (!requireNamespace("tidyr", quietly = TRUE)) {
  install.packages("tidyr")
}
if (!requireNamespace("circlize", quietly = TRUE)) {
  install.packages("circlize")
}

# Load packages
library(ComplexHeatmap)
library(ggplot2)
library(dplyr)
library(tidyr)
library(circlize)
sessioninfo::session_info("attached")
โ”€ Session info โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€
 setting  value
 version  R version 4.6.2 (2026-10-29)
 os       Ubuntu 26.04.1 LTS
 system   x86_64, linux-gnu
 ui       X11
 language (EN)
 collate  C.UTF-8
 ctype    C.UTF-8
 tz       UTC
 date     2026-11-05
 pandoc   3.7.0.2 @ /usr/bin/ (via rmarkdown)
 quarto   1.11.9 @ /usr/local/bin/quarto

โ”€ Packages โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€
 ! package        * version    date (UTC) lib source
 P circlize       * 0.4.18     2026-04-04 [?] RSPM
 P ComplexHeatmap * 2.28.0     2026-04-28 [?] Bioconductor 3.23 (R 4.6.2)
 P dplyr          * 1.2.1      2026-04-03 [?] RSPM
 P ggplot2        * 4.0.3.9000 2026-11-03 [?] Github (tidyverse/ggplot2@b638d65)
 P tidyr          * 1.3.2      2025-12-19 [?] RSPM

 [1] /home/runner/.cache/R/renv/library/Bizard-5140bdc6/linux-ubuntu-resolute/R-4.6/x86_64-pc-linux-gnu
 [2] /home/runner/work/_temp/Library
 [3] /opt/R/4.6.2/lib/R/library

 * โ”€โ”€ Packages attached to the search path.
 P โ”€โ”€ Loaded and on-disk path mismatch.

โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€

Data Preparation

In practice the input is the matrix multiBigwigSummary BED-file writes (one column per bigWig, one row per peak). We simulate 120 peaks measured in four conditions โ€” WT and KO macrophages, each with and without LXR agonist treatment โ€” with a treatment-responsive block and an unresponsive block.

set.seed(701)
n_induced <- 50; n_flat <- 70
n_bin <- 30
conditions <- c("WT", "WT_treated", "KO", "KO_treated")

## binned ChIP/ATAC signal: rows = peaks, columns = 30 bins x 4 conditions
shape_vals <- dnorm(seq(-3, 3, length.out = n_bin))              # peaked profile
mk_profile <- function(peak_effect, n) {
  do.call(cbind, lapply(conditions, function(cn) {
    eff <- peak_effect *
      ifelse(cn == "WT_treated", 1, ifelse(cn == "KO_treated", 0.25, 0))
    do.call(cbind, lapply(seq_len(n_bin), function(b)
      0.6 + eff * shape_vals[b] + rnorm(n, sd = 0.12)))          # n peaks x 1 bin
  }))
}
peaks_induced <- mk_profile(2.2, n_induced)
peaks_flat    <- mk_profile(0.0, n_flat)
colnames(peaks_induced) <- colnames(peaks_flat) <-
  paste0(rep(conditions, each = n_bin), "_b", rep(seq_len(n_bin), times = 4))

peak_mat <- rbind(peaks_induced, peaks_flat)
rownames(peak_mat) <- c(paste0("Induced_peak_", seq_len(n_induced)),
                        paste0("Flat_peak_", seq_len(n_flat)))
peak_group <- factor(c(rep("GW3965-induced", n_induced), rep("Unresponsive", n_flat)),
                     levels = c("GW3965-induced", "Unresponsive"))
dim(peak_mat)
[1] 120 120

Visualization

1. Group-Mean Accessibility Profiles

prof <- as.data.frame(peak_mat) |>
  mutate(group = peak_group) |>
  group_by(group) |>
  summarise(across(where(is.numeric), mean), .groups = "drop") |>
  pivot_longer(-group, names_to = "key", values_to = "signal") |>
  mutate(condition = sub("_b[0-9]+$", "", key),
         bin = as.integer(sub(".*_b", "", key)))

ggplot(prof, aes(bin, signal, colour = condition, fill = condition)) +
  stat_summary(geom = "line", linewidth = 0.9) +
  facet_wrap(~group, nrow = 1) +
  scale_colour_manual(values = c(WT = "#2874C5", WT_treated = "#2874C5",
                                 KO = "#B2182B", KO_treated = "#B2182B"),
                      labels = c(WT = "WT", WT_treated = "WT + GW3965",
                                 KO = "KO", KO_treated = "KO + GW3965"),
                      name = NULL) +
  scale_fill_manual(values = c(WT = "#2874C5", WT_treated = "#2874C5",
                               KO = "#B2182B", KO_treated = "#B2182B")) +
  labs(x = "Position over peak (binned)", y = "Mean accessibility") +
  guides(fill = "none") +
  theme_classic(base_size = 11) +
  theme(legend.position = "top", legend.key.size = unit(3, "mm"))
Figureย 1: Mean binned accessibility per condition for the induced region set: treatment raises WT strongly, KO only weakly.

The WT + GW3965 curve separates from KO + GW3965 in the induced set but not in the unresponsive set โ€” the treatment effect the reference figure reports.

2. Accessibility-Difference Heatmap

col_fun <- circlize::colorRamp2(c(0, 1.2, 2.8), c("#F7F7F7", "#FDAE61", "#B2182B"))
row_ord <- order(peak_group,                                  # group first
                 rowMeans(peak_mat[, grep("KO_treated_b", colnames(peak_mat))] -
                            peak_mat[, grep("KO_b", colnames(peak_mat))]),
                 decreasing = TRUE)
mat_sorted <- peak_mat[row_ord, ]
grp_sorted <- peak_group[row_ord]

col_names <- rep(seq_len(n_bin), times = 4)
ha_top <- HeatmapAnnotation(
  Condition = rep(conditions, each = n_bin),
  col = list(Condition = c(WT = "#BDD7E7", WT_treated = "#2171B5",
                           KO = "#FCAE91", KO_treated = "#B2182B")),
  show_legend = FALSE)

Heatmap(mat_sorted, name = "Accessibility",
        top_annotation = ha_top,
        col = col_fun,
        cluster_rows = FALSE, cluster_columns = FALSE,
        row_split = grp_sorted, column_split = factor(rep(conditions, each = n_bin), conditions),
        column_gap = unit(1, "mm"),
        show_column_names = FALSE, show_row_names = FALSE,
        border = TRUE)
Figureย 2: Log2 accessibility heatmap (rows = peaks sorted within each group, columns = 30 bins x 4 conditions), split by region group.
Tip

Practice tips

  • Sort rows by the biological contrast (here KO-treated minus KO) inside each group โ€” clustering would scramble the induced/unresponsive split the figure exists to show.
  • column_split + column_gap keeps the four condition blocks visually separate; the top annotation maps each block to its condition colour.
  • Upstream, the real matrices come from multiBigwigSummary BED-file; normalise samples first (deeptools bigwigCompare or RPGC), or the condition colours will reflect depth, not accessibility.

Application

Aggregate accessibility heatmaps are the standard figure for reporting condition-specific chromatin opening (drug treatment, knockout, lineage transitions). The same layout serves ChIP-seq signal over peak sets.

Reference

[1] Ramรญrez F, Dรผndar F, Diehl S, Grรผning BA, Manke T. deepTools: a flexible web platform for exploring high-throughput sequencing data. Nucleic Acids Res. 2014;42(W1):W187-191. doi:10.1093/nar/gku365

[2] Xiaofan Lu, et al.ย FigureYa: A Standardized Visualization Framework for Enhancing Biomedical Data Interpretation and Research Efficiency. iMetaMed. 2025. doi:10.1002/imm3.70005