visualization-workflow

Visualization Workflow for ENCODE Data

When to Use

• User wants to create genome browser visualizations, heatmaps, or signal track plots from ENCODE data
• User asks about "visualization", "genome browser", "deeptools", "heatmap", "signal track", or "IGV"
• User needs to generate publication-ready figures from ChIP-seq, ATAC-seq, or other genomic data
• User wants to compare signal profiles across conditions, tissues, or histone marks
• Example queries: "visualize H3K27ac signal at promoters", "create a heatmap of ChIP-seq signal", "set up a UCSC track hub for my data"

Help the user create informative, publication-quality visualizations of ENCODE genomic data. This skill covers four major visualization approaches: deepTools heatmaps and profiles, IGV genome browser views, UCSC track hubs for sharing, and publication-quality static plots using R and Python. Visualization is not decorative -- it is an essential analytical step that reveals patterns invisible in summary statistics and validates computational findings.

Literature Foundation

| Reference | Journal | Key Contribution | DOI | Citations | |-----------|---------|-----------------|-----|-----------| | Ramirez et al. (2016) | Nucleic Acids Research | deepTools2: next-generation server for deep-sequencing data analysis; heatmaps, profiles, correlation, PCA | 10.1093/nar/gkw257 | ~3,800 | | Robinson et al. (2011) | Nature Biotechnology | Integrative Genomics Viewer (IGV): interactive exploration of large genomic datasets | 10.1038/nbt.1754 | ~10,000 | | Kent et al. (2002) | Genome Research | The Human Genome Browser at UCSC: foundation for track-based genomic visualization | 10.1101/gr.229102 | ~8,000 | | Ramirez et al. (2014) | Nucleic Acids Research | deepTools: flexible platform for exploring deep-sequencing data; original computeMatrix/plotHeatmap framework | 10.1093/nar/gku365 | ~2,500 | | Amemiya et al. (2019) | Scientific Reports | ENCODE Blacklist: comprehensive identification of artifact regions to exclude from visualization | 10.1038/s41598-019-45839-z | ~1,372 | | Wickham (2016) | Springer | ggplot2: Elegant Graphics for Data Analysis; grammar of graphics for genomic visualization | ISBN: 978-3-319-24277-4 | ~30,000+ |

Part 1: deepTools Heatmaps and Profiles

deepTools (Ramirez et al. 2014, 2016) is the standard toolkit for visualizing ChIP-seq and ATAC-seq signal across genomic regions. The core workflow is: compute a signal matrix, then render it as a heatmap or profile plot.

1a. computeMatrix: Building the Signal Matrix

computeMatrix extracts signal values from bigWig files across a set of genomic regions. Two modes are available:

reference-point mode -- centers the signal on a single anchor point (e.g., TSS, peak summit):

# Signal centered on peak summits, +/- 3kb
computeMatrix reference-point \
    -S H3K27ac_fc.bigWig H3K4me3_fc.bigWig ATAC_fc.bigWig \
    -R peaks.bed \
    --referencePoint center \
    -b 3000 -a 3000 \
    --binSize 50 \
    --missingDataAsZero \
    --sortRegions descend \
    --sortUsing mean \
    -o matrix_refpoint.gz \
    -p 8

scale-regions mode -- scales all regions to uniform length (e.g., gene bodies):

# Signal across scaled gene bodies with 2kb flanks
computeMatrix scale-regions \
    -S H3K36me3_fc.bigWig RNA_signal.bigWig \
    -R genes.bed \
    --regionBodyLength 5000 \
    -b 2000 -a 2000 \
    --binSize 50 \
    --missingDataAsZero \
    -o matrix_scaled.gz \
    -p 8

When to use which mode:

• reference-point: TF ChIP-seq peaks, ATAC-seq summits, TSSs, enhancer centers -- any feature defined by a point
• scale-regions: gene bodies, broad histone domains (H3K27me3, H3K36me3), TADs -- features with variable length

1b. plotHeatmap: Rendering the Matrix

plotHeatmap -m matrix_refpoint.gz \
    -o heatmap.png \
    --colorMap RdYlBu_r \
    --whatToShow "heatmap and colorbar" \
    --sortRegions descend \
    --sortUsing mean \
    --heatmapHeight 15 \
    --heatmapWidth 4 \
    --zMin 0 --zMax 10 \
    --samplesLabel "H3K27ac" "H3K4me3" "ATAC" \
    --regionsLabel "Peaks" \
    --dpi 300

Clustering: To reveal sub-patterns within peak sets:

plotHeatmap -m matrix_refpoint.gz \
    -o heatmap_clustered.png \
    --kmeans 4 \
    --colorMap viridis \
    --zMin 0 --zMax 10 \
    --outFileSortedRegions clusters.bed \
    --dpi 300

The --outFileSortedRegions flag exports the cluster assignments as a BED file, enabling downstream analysis of each cluster separately.

Recommended color maps by mark type: | Mark Type | Recommended colorMap | Rationale | |-----------|---------------------|-----------| | Active marks (H3K27ac, H3K4me3) | Reds, YlOrRd | Warm colors for activation | | Repressive marks (H3K27me3, H3K9me3) | Blues, PuBu | Cool colors for repression | | Accessibility (ATAC, DNase) | Greens, YlGn | Distinct from histone colors | | Multi-mark comparison | viridis, inferno | Perceptually uniform, colorblind-safe |

1c. plotProfile: Average Signal Plots

Profile plots show the average signal across all regions, useful for comparing samples:

plotProfile -m matrix_refpoint.gz \
    -o profile.png \
    --perGroup \
    --plotTitle "Signal at H3K27ac peaks" \
    --yAxisLabel "Fold change over input" \
    --samplesLabel "H3K27ac" "H3K4me3" "ATAC" \
    --dpi 300

Use --perGroup when you have multiple region sets (e.g., active vs poised enhancers) and want separate profile lines for each group.

1d. Signal Correlation and PCA

Before making complex visualizations, verify that replicates correlate and conditions separate:

# Build correlation matrix
multiBigwigSummary bins \
    -b sample1.bw sample2.bw sample3.bw sample4.bw \
    --labels Rep1 Rep2 Rep3 Rep4 \
    --binSize 10000 \
    -o results.npz \
    -p 8

# Correlation heatmap
plotCorrelation -in results.npz \
    --corMethod pearson \
    --whatToPlot heatmap \
    --plotFile correlation.pdf \
    --skipZeros

# PCA plot
plotPCA -in results.npz \
    --plotFile pca.pdf \
    --labels Rep1 Rep2 Rep3 Rep4

Part 2: IGV Visualization

The Integrative Genomics Viewer (Robinson et al. 2011) provides interactive, locus-level inspection of ENCODE data. IGV is essential for validating computational findings at individual loci.

2a. Loading ENCODE Files in IGV

ENCODE data can be loaded directly from URLs without downloading:

1. Open IGV and select the correct genome (hg38 for GRCh38, mm10 for mouse)
2. File > Load from URL > paste the ENCODE file download URL
3. For bigWig files, IGV streams data on-the-fly (no full download needed)

Recommended file types for IGV: | File Type | IGV Display | Best For | |-----------|------------|----------| | bigWig (fold change over control) | Continuous signal track | Viewing signal intensity | | bigBed (IDR thresholded peaks) | Discrete interval track | Viewing peak locations | | BAM (alignments) | Read pileup + coverage | Inspecting read-level evidence |

2b. Batch Screenshots with IGV

For systematic locus-level visualization across many genes, use IGV batch scripting:

new
genome hg38
load https://www.encodeproject.org/files/ENCFF.../@@download/ENCFF....bigWig
load https://www.encodeproject.org/files/ENCFF.../@@download/ENCFF....bigBed
snapshotDirectory /path/to/output/
goto chr11:2,159,779-2,161,209
snapshot INS_locus.png
goto chr7:44,182,955-44,184,393
snapshot GCK_locus.png
goto chr17:40,927,190-40,928,775
snapshot HNF1B_locus.png

Run with: igv.sh -b batch_script.txt

2c. IGV.js for Web-Based Viewing

For sharing interactive browser views without requiring local IGV installation:

<div id="igv-div"></div>
<script src="https://cdn.jsdelivr.net/npm/igv@2.15.0/dist/igv.min.js"></script>
<scr