🧬 CORE OMICS DISCIPLINE 05 DNA Methylation (WGBS/EM-Seq) • ATAC-Seq • CUT&Tag • Chromatin ⏱️ Comprehensive Production Blueprint • MACS3 • Bismark • TOBIAS Footprinting

Epigenomics: DNA Methylation, Chromatin Accessibility & Histone Modification Dynamics

Epigenomics investigates heritable biochemical modifications to DNA and chromatin proteins that regulate gene expression without altering the underlying nucleotide sequence. By mapping 5-methylcytosine (5mC), open euchromatin (ATAC-seq), and post-translational histone marks (CUT&Tag), epigenomics decodes how a single genome gives rise to hundreds of distinct cellular lineages.

1. Biological Foundations 2. Experimental Technologies 3. Visual Architecture (SVG) 4. Computational Pipeline & Code 5. Real-World Case Study 6. Software Directory & Links 7. Pitfalls & QC Standards 8. Exercises & Datasets

🧬 Module 1: Biological Foundations of Epigenetic Regulation

Chromatin Architecture

Nuclear DNA does not exist as naked strands; it is wrapped around octamers of core histone proteins (two copies each of H2A, H2B, H3, and H4) to form nucleosomes. The structural packaging of chromatin regulates the physical access of transcription factors and RNA Polymerase II to genetic loci.

DNA Methylation (5mC & 5hmC)

DNA methyltransferases (DNMT1, DNMT3A/B) transfer methyl groups onto cytosine C5 in CpG dinucleotides. While ~70–80% of human CpGs are methylated, unmethylated CpG islands at gene promoters enable transcription. Ten-eleven translocation (TET) enzymes oxidize 5mC to 5hmC to initiate active DNA demethylation.

The Histone Code

Post-translational modifications on histone N-terminal tails act as regulatory beacons:
• H3K4me3: Canonical mark of active transcription initiation promoters.
• H3K27ac: Active promoters and distal enhancers (distinguishes active from primed enhancers).
• H3K27me3: Polycomb repressive complex 2 (PRC2) facultative silencing.
• H3K9me3: Constitutive heterochromatin (centromeres, satellite DNA).

Chromatin Accessibility & Footprints

Open euchromatin regions (~2–3% of the genome) represent nucleosome-depleted regions (NDRs) bound by pioneer transcription factors. The precise physical obstruction caused by a bound TF protects the underlying base pairs from enzymatic cleavage, producing a digital genomic footprint.

🔬 Module 2: Epigenomic Profiling Assays & Chemistries

Assay Protocols

Epigenomic assays employ distinct enzymatic, chemical, and antibody-directed tagmentation strategies:

Assay	Target Mechanism & Reaction	Input Requirement	Resolution	Key Computational Tools
WGBS / EM-Seq	Sodium bisulfite or enzymatic deamination (APOBEC3A + TET2) converts unmethylated C → U (reads as T), while 5mC remains C.	10–100 ng DNA	Single-base CpG resolution	Bismark, bwa-meth, methylKit, DSS
ATAC-Seq	Hyperactive Tn5 transposase inserts Illumina sequencing adapters directly into accessible, open chromatin.	50,000 cells (down to single-cell)	Nucleosome & TF footprint scale	Bowtie2, MACS3, TOBIAS, HMMRATAC
CUT&Tag / CUT&RUN	Protein A-Tn5 fusion protein guided by primary antibody directly tagments chromatin surrounding specific histone marks in permeabilized nuclei.	1,000–50,000 cells (Ultra-low input)	High signal-to-noise peaks	Bowtie2, SEACR, MACS3, deepTools

🖼️ Module 3: Chromatin Landscape & Epigenomic Signal Integration

Scientific Diagram

Figure 5.1: Epigenomic Regulatory Profiling — (1) Tn5 transposase insertion into open chromatin; (2) Chemical/enzymatic conversion of unmethylated cytosines to uracil; (3) Multi-track alignment demonstrating concordance between ATAC accessibility, promoter H3K4me3, and enhancer H3K27ac peaks.

💻 Module 4: Computational Pipeline & Production Code

ATAC-Seq & Bismark

Part A: ATAC-Seq / CUT&Tag Alignment, Tn5 Shifting & Peak Calling (Bash)

#!/usr/bin/env bash
set -euo pipefail

# 1. Trimming Nextera Transposase Adapters
fastp -i atac_R1.fastq.gz -I atac_R2.fastq.gz \
  -o atac_trimmed_R1.fastq.gz -O atac_trimmed_R2.fastq.gz \
  --adapter_sequence=CTGTCTCTTATACACATCT --thread 16

# 2. Bowtie2 Alignment with fragment length constraint up to 2000 bp
bowtie2 -p 16 --very-sensitive --no-discordant --no-mixed -X 2000 \
  -x /genomes/GRCh38/bowtie2_index/GRCh38 \
  -1 atac_trimmed_R1.fastq.gz -2 atac_trimmed_R2.fastq.gz | \
  samtools view -b -q 30 -F 1804 - | \
  samtools sort -@ 8 -o atac_aligned.bam -

# 3. Purge Mitochondrial DNA (chrM) & Mark Duplicates
samtools idxstats atac_aligned.bam | cut -f 1 | grep -v chrM | xargs samtools view -b atac_aligned.bam > atac_no_chrM.bam
picard MarkDuplicates I=atac_no_chrM.bam O=atac_dedup.bam M=dup_qc.txt REMOVE_DUPLICATES=true
samtools index atac_dedup.bam

# 4. MACS3 Peak Calling with Paired-End Mode (-f BAMPE)
macs3 callpeak \
  -t atac_dedup.bam \
  -f BAMPE \
  -g hs \
  -n sample_atac \
  -q 0.01 \
  --keep-dup all \
  --outdir ./macs3_peaks/

# 5. Generate bigWig normalized coverage track for IGV
bamCoverage -b atac_dedup.bam -o atac_norm_coverage.bw --normalizeUsing RPKM -p 16

Part B: Bismark Bisulfite Alignment & CpG Methylation Extraction (Bash + R)

# 1. Bisulfite 3-letter reference genome preparation
bismark_genome_preparation --bowtie2 /genomes/GRCh38/

# 2. Align trimmed paired-end reads to in-silico bisulfite-converted genome
bismark --bowtie2 -p 8 --multicore 4 /genomes/GRCh38/ \
  -1 wgbs_val_1.fq.gz -2 wgbs_val_2.fq.gz \
  -o ./bismark_aligned/

# 3. Mark PCR duplicates in bisulfite alignments
deduplicate_bismark -p --bam ./bismark_aligned/wgbs_val_1_bismark_bt2_pe.bam

# 4. Extract single-base CpG methylation ratios & coverage bedGraph
bismark_methylation_extractor -p --comprehensive --no_overlap --bedGraph \
  --counts --buffer_size 10G --multicore 8 \
  wgbs_val_1_bismark_bt2_pe.deduplicated.bam

# R Analysis with methylKit: Call Differentially Methylated Regions (DMRs)
# library(methylKit)
# myobj <- methRead(list("sample1.cov", "sample2.cov"), sample.id=list("S1","S2"), assembly="hg38", context="CpG")
# myDiff <- calculateDiffMeth(myobj)
# dmrs <- getMethylDiff(myDiff, difference=25, qvalue=0.01)

🏥 Module 5: Real-World Case Study

Enhancer Reprogramming

Case: De Novo Super-Enhancer Assembly Driving Oncogenic TAL1 in T-ALL

Biological Problem: In T-cell acute lymphoblastic leukemia (T-ALL), the *TAL1* master proto-oncogene is aberrantly overexpressed in patients who lack chromosomal translocations or gene amplifications. Paired ATAC-seq, H3K27ac CUT&Tag, and whole-genome sequencing were deployed to resolve non-coding regulatory architecture.

Genetic Lesion: 12-bp Non-Coding Indel Insertion (chr1:47,713,290) creating a de novo MYB binding motif.

Epigenomic Signal: Emergence of a massive 18-kb H3K27ac super-enhancer cluster >8 kb upstream of TAL1 TSS.

ATAC-Seq Footprint: High-occupancy TF footprinting confirmed continuous physical recruitment of MYB, CBP, and RUNX1.

Transcriptional Outcome: 45-fold transcriptional induction of TAL1 causing total differentiation arrest.

Translational Resolution: Targeted disruption of the non-coding insertion via CRISPR-Cas9 or treatment with CDK7/THZ1 transcriptional kinase inhibitors abolished the super-enhancer and extinguished leukemic cell viability.

📦 Module 6: Essential Software Directory & Official Links

Verified Tool Stack

MACS3

Peak Calling

Model-based Analysis of ChIP-Seq and ATAC-Seq for identifying enriched binding sites and open chromatin.

Official MACS3 GitHub →

Bismark

Bisulfite Aligner

Flexible tool to map bisulfite-treated sequencing reads and determine cytosine methylation states.

Babraham Bismark Portal →

TOBIAS Footprinting

TF Footprints

Transcription factor footprinting framework correcting for Tn5 transposase sequence insertion bias in ATAC-seq.

Looso Lab TOBIAS →

deepTools

ChIP Matrix & Heatmaps

User-friendly suite for normalizing BAM coverage tracks, plotting metaplots around TSSs, and heatmap clustering.

deepTools Documentation →

methylKit

DMR Calling

R package for DNA methylation analysis, filtering CpG coverage, clustering, and calling differentially methylated regions.

Bioconductor methylKit →

SEACR (CUT&Tag Peak Calling)

Ultra-Low Background

Sparse enrichment analysis for CUT&RUN and CUT&Tag data tailored to extremely low background non-specific signals.

Fred Hutch SEACR →

⚠️ Module 7: Bioinformatic Pitfalls & Quality Control Standards

Best Practices

1. Mitochondrial Read Contamination in ATAC

Because mitochondrial DNA is completely circular and naked (lacking nucleosomes), Tn5 transposase preferentially tagments chrM, frequently consuming 30–70% of sequencing reads. Use mild detergent lysis or bioinformatic chrM purging.

2. Tn5 Insertion Shifting Requirement

Tn5 binds as a homodimer with a 9-bp steric duplication between cleavage events. To identify true single-base transcription factor cut centers, forward strand reads MUST be shifted +4 bp and reverse strand reads shifted -5 bp.

3. Incomplete Bisulfite Conversion

If bisulfite conversion does not reach >99% efficiency, unmethylated cytosines will fail to convert to uracil and be erroneously called as methylated, globally corrupting methylation profiles. Always spike unmethylated Lambda phage DNA as a conversion control.

📝 Module 8: Self-Assessment Exercises & Practice Datasets

Hands-on Mastery

📦 Curated Benchmark: Chr1 Bismark Coverage & ATAC Peak Matrix

Size: 8.4 MB • Bismark .cov.gz and MACS3 .narrowPeak benchmarks for epigenomic practice.

Access Practice Dataset →

Q1: What is the biological difference between H3K4me1, H3K4me3, and H3K27ac?

Answer: H3K4me3 localizes sharply at actively transcribed transcription start sites (promoters). In contrast, H3K4me1 marks enhancer regions (both primed and active). H3K27ac distinguishes active enhancers from inactive or primed enhancers (which are H3K4me1+ but H3K27ac-). When a region exhibits both H3K4me1 and H3K27ac, it represents a bonafide active regulatory enhancer.

Q2: Why does ATAC-seq fragment size distribution exhibit periodic 200 bp banding?

Answer: Nucleosomes are spaced at approximately ~147 bp of DNA wrapped around the core octamer plus ~50 bp of linker DNA, creating a ~200 bp periodic repeat. Tn5 can only cut within the linker DNA between nucleosomes. Consequently, fragment size analysis reveals a dominant nucleosome-free peak (<100 bp), followed by mononucleosome (~200 bp), dinucleosome (~400 bp), and trinucleosome (~600 bp) periodic distributions.

Q3: Why can't standard BWA or Bowtie2 align bisulfite-treated DNA without specialized index preparation?

Answer: Bisulfite treatment chemically converts unmethylated cytosines into uracil (which amplifies as Thymine). As a result, the read alphabet is effectively reduced to a 3-letter code (A, G, T) with very few C's. Aligning these reads against a standard 4-letter genome leads to massive mismatch penalties across every unmethylated cytosine. Bismark solves this by performing in-silico C-to-T conversion of both the reference genome and sequencing reads before three-letter alignment.