Visualizing RNA-seq Data with pheatmap: Advanced Features for Annotations, Gaps, and Custom Colors

By Dr. Zubair Khalid, DVM, MS, PhD · 2026-08-20

Key Takeaways

• pheatmap enhances RNA-seq visualization by integrating experimental metadata directly into heatmaps via annotation_col and annotation_row parameters, enabling simultaneous display of expression patterns and sample/gene characteristics. This is crucial for interpreting complex experimental designs involving multiple treatment groups, time points, or biological replicates.
• Hierarchical clustering in pheatmap benefits from domain-specific choices; correlation-based distances (e.g., 1-Pearson correlation) are often superior to Euclidean distance for RNA-seq data, as they focus on expression profile shapes rather than absolute magnitudes, revealing biologically relevant co-expression patterns. Ward's method or average linkage are recommended for merging clusters to minimize variance or provide a balance between cluster compactness and separation.
• Visual separation of experimental groups or gene clusters is achieved through gaps_col and gaps_row parameters, which insert white space based on column/row order. This is particularly effective when combined with manual column ordering to delineate treatment groups or when applied to row clusters identified by dendrogram cutting (cutree_rows) to highlight functional modules.
• Custom color palettes and breaks are essential for accurate representation of RNA-seq expression data, especially for log2-transformed values or z-scores. Diverging palettes centered at zero (e.g., blue-white-red) are recommended for log2 fold changes to clearly distinguish up- and down-regulation, while careful use of the breaks parameter can optimize color mapping for skewed data distributions.
• Reproducibility in pheatmap generation necessitates meticulous documentation of all parameters, including input matrix preprocessing (normalization methods like DESeq2's regularized log), gene selection criteria (e.g., top variable genes, significant DE genes), clustering methods, annotation color assignments, and color scale breaks. This ensures that figures can be accurately recreated and validated.

RNA sequencing produces high-dimensional expression matrices that require careful visualization for biological interpretation. The pheatmap package in R provides a practical solution for creating publication-ready heatmaps with sample annotations, gene groupings, and custom color schemes. This article addresses the specific challenge of moving beyond default heatmap output to build figures that convey experimental design, sample metadata, and expression patterns in a single coherent graphic. The content applies to biology students, researchers, laboratory professionals, and life-science practitioners who have completed differential expression analysis and need to communicate results effectively.

The Role of Heatmaps in RNA-seq Analysis Workflows

Heatmap visualization sits at the intersection of quality assessment and biological interpretation in RNA-seq projects. Before differential expression testing, heatmaps of sample-to-sample distances or top variable genes reveal batch effects, outlier samples, and overall data structure. After differential expression analysis, heatmaps of significant genes across conditions provide an intuitive summary of expression patterns that tables of log2 fold changes cannot match.

The practical value of advanced pheatmap features becomes clear when you consider the typical RNA-seq experiment. A study comparing treated and control samples might include biological replicates, multiple time points, and several treatment groups. The expression matrix itself contains thousands of genes, but the genes of interest after filtering may number in the hundreds. A default heatmap of this matrix would produce a figure with no sample grouping information, no gene clustering structure, and a default color scale that may obscure biologically meaningful differences.

The pheatmap package addresses these limitations through parameters that control annotation display, gap placement, and color mapping. These features allow researchers to embed experimental metadata directly into the figure, separate gene clusters visually, and choose color scales that match the perceptual needs of the data. The result is a figure that answers questions about sample relationships and gene expression patterns simultaneously.

RNA-seq analysis requires proficiency with computational approaches to manage technical issues and large data sizes, as noted in a practical guide for molecular biologists encountering RNA-seq analysis for the first time [<a href="#ref-1">1</a>]. Visualization tools like pheatmap form part of this computational skill set, enabling researchers to move from raw count matrices to interpretable biological conclusions. The broader bioinformatics training ecosystem, including resources from the European Bioinformatics Institute, emphasizes practical analysis education that includes data visualization competencies [<a href="#ref-2">2</a>].

Understanding the pheatmap Input Structure

The pheatmap function accepts a numeric matrix where rows represent genes or features and columns represent samples. For RNA-seq data, this matrix typically contains normalized expression values such as log2-transformed counts per million, variance-stabilized counts, or regularized log counts from DESeq2. The choice of input values affects the visual output substantially, and researchers should document this choice in their analysis records.

Matrix Preparation and Normalization Considerations

Raw count data should not be passed directly to pheatmap. The variance structure of raw counts is dominated by highly expressed genes, and the range of values spans several orders of magnitude. Normalization methods address library-specific compositional differences, as demonstrated in studies using the Trimmed Mean of M-values method for RNA-seq data [<a href="#ref-3">3</a>]. After normalization, a log transformation compresses the dynamic range and makes expression differences visually interpretable.

For heatmap purposes, the most common input is a matrix of log2-transformed normalized counts. Some researchers prefer z-scores calculated across rows, which centers each gene's expression at zero and scales by the standard deviation. Row scaling is particularly useful when the goal is to compare expression patterns across genes with very different baseline expression levels. The choice between raw log2 values and row-scaled values depends on whether absolute expression differences or relative patterns matter more for the biological question.

The normalization step carries consequences for downstream interpretation. Studies that apply robust normalization methods before visualization produce heatmaps where the visual patterns reflect biological variation instead of technical artifacts [<a href="#ref-3">3</a>]. Researchers should record which normalization method was applied and how the transformed values were calculated, because this information affects how readers interpret the color scale.

Selecting Genes for Visualization

The full expression matrix contains too many rows for meaningful heatmap visualization. A heatmap of 20,000 genes produces a figure where individual gene patterns cannot be distinguished. Practical selection strategies include using the top variable genes across samples, the significant genes from differential expression analysis, or a curated set of genes related to the biological pathway of interest.

For differential expression results, the selection typically follows significance thresholds and fold-change cutoffs. The statistical approach matters here, as methods like Generalized Linear Models with Quasi-Likelihood handle the discrete and overdispersed nature of RNA-seq data more appropriately than simple t-tests [<a href="#ref-3">3</a>]. After identifying significant genes, researchers often further filter by the magnitude of expression change to focus on biologically meaningful differences.

The number of genes selected affects both the visual clarity and the computational performance of pheatmap. A matrix of 50 to 200 genes produces a readable figure where individual gene rows can be examined. Larger selections require larger figure dimensions and may benefit from clustering to reveal structure within the gene set.

At a Glance: Key pheatmap Parameters for RNA-seq Visualization

The following table summarizes the primary parameters covered in this article, their purpose, and typical use cases in RNA-seq analysis.

Building the Base Heatmap with Clustering Controls

The default pheatmap call clusters both rows and columns using hierarchical clustering with Euclidean distance and complete linkage. These defaults work for many datasets, but RNA-seq data often benefits from explicit choices about distance metrics and clustering methods.

Distance Metrics and Linkage Methods

Euclidean distance treats each gene's expression profile as a point in sample space and measures the straight-line distance between points. Correlation-based distances, such as 1 minus Pearson correlation, focus on the shape of expression profiles instead of their absolute magnitude. For RNA-seq data where genes have very different baseline expression levels, correlation-based distances often produce more biologically meaningful clusters.

The linkage method determines how clusters are merged during hierarchical clustering. Complete linkage uses the maximum distance between points in two clusters, which tends to produce compact clusters. Ward's method minimizes the total within-cluster variance and often produces visually appealing results for expression data. Average linkage represents a compromise between the extremes of single and complete linkage.

The clustering parameters should be recorded in the analysis documentation. Different choices can produce different cluster assignments, and reproducibility requires that these choices be explicit. The pheatmap function allows users to supply precomputed distance matrices or clustering objects, which provides additional control for advanced users.

Row and Column Ordering

By default, pheatmap orders rows and columns according to the hierarchical clustering result. This ordering places similar expression profiles adjacent to each other, revealing blocks of co-expressed genes and groups of samples with similar expression patterns.

For column ordering, researchers may prefer to specify the order manually to match the experimental design. This approach places all control samples together, followed by treatment samples, regardless of the clustering result. Manual ordering is useful when the experimental design has a natural structure that should be reflected in the figure.

The gaps_col parameter works with column ordering to insert visual breaks between groups of columns. When combined with manual column ordering, gaps can separate treatment groups, time points, or any other experimental factor. The gaps_row parameter provides the same functionality for rows, allowing separation of gene clusters or functional categories.

Adding Sample Annotations with annotation_col

The annotation_col parameter transforms a basic heatmap into a figure that communicates experimental design. This parameter accepts a data frame where each row corresponds to a sample and each column corresponds to an annotation variable such as treatment group, time point, or patient identifier.

Building the Annotation Data Frame

The annotation data frame must have row names that match the column names of the expression matrix. Each column of the annotation data frame becomes a colored annotation bar displayed above the heatmap. The pheatmap function automatically assigns colors to factor levels, and users can customize these colors with the annotation_colors parameter.

For RNA-seq experiments, common annotation variables include treatment condition, replicate number, batch identifier, and sample quality metrics. Including batch information in the annotation is particularly valuable because it allows readers to assess whether sample clustering correlates with technical factors instead of biological ones.

The annotation data frame should be constructed carefully to ensure that factor levels are ordered correctly. The order of levels determines the order of the legend, and for ordinal variables such as time points, the levels should be specified in the correct sequence.

Customizing Annotation Colors

The annotation_colors parameter accepts a named list where each element corresponds to an annotation variable. Each element is itself a named vector mapping factor levels to colors. This parameter provides control over the color scheme of the annotation bars, allowing researchers to use consistent colors across multiple figures in a publication.

Color choices for annotations should consider colorblind accessibility. Red-green color schemes are problematic for readers with common forms of color vision deficiency. Alternatives such as blue-orange or viridis-based palettes provide better accessibility while maintaining visual distinction between categories.

The annotation legend appears alongside the heatmap and shows the mapping between colors and factor levels. This legend is essential for interpreting the figure, and its presence should be verified in the final output.

Creating Visual Gaps with gaps_col and gaps_row

The gaps_col and gaps_row parameters insert white space between groups of columns or rows in the heatmap. These gaps serve a visual organization function, separating experimental groups or gene clusters so that the figure communicates structure more effectively.

Using gaps_col for Sample Group Separation

The gaps_col parameter accepts a numeric vector specifying the positions where gaps should be inserted. The positions refer to the column indices after any clustering or manual ordering has been applied. For example, if the columns are ordered with three control samples followed by three treatment samples, gaps_col = 3 inserts a gap after the third column.

When columns are clustered, the gap positions must correspond to the clustered order. This requirement makes manual column ordering more practical when gaps are needed, because the researcher controls the exact column sequence. The combination of manual ordering and gaps produces a figure where each experimental group appears as a distinct block.

Using gaps_row for Gene Cluster Separation

The gaps_row parameter works analogously for rows. After row clustering, the researcher can inspect the cluster assignments and insert gaps between clusters. This approach separates groups of co-expressed genes visually, making it easier to identify functional modules within the heatmap.

Determining gap positions for rows requires examining the clustering result. The pheatmap function returns the clustering information when the return_value parameter is set to "cluster". This output includes the row order and cluster assignments, which can be used to calculate gap positions.

A practical approach is to cut the dendrogram at a specified height to obtain a fixed number of clusters, then calculate the cumulative sizes of these clusters to determine gap positions. This method produces consistent gaps that correspond to actual gene clusters instead of arbitrary positions.

Custom Color Palettes for Expression Data

The color parameter controls the gradient used to represent expression values in the heatmap. The default palette in pheatmap uses a purple-to-orange gradient, but RNA-seq data often benefits from custom palettes that match the data characteristics and publication requirements.

Choosing Appropriate Color Scales

For log2-transformed expression data, a diverging color scale with a neutral color at the midpoint works well. The midpoint should correspond to the center of the data distribution, which may be near zero for row-scaled data or near the median expression level for unscaled data.

Common choices include blue-white-red gradients, where blue represents low expression and red represents high expression. The colorRampPalette function in R allows users to create custom gradients with any number of intermediate colors. For example, a gradient from dark blue through white to dark red provides strong visual contrast for extreme values while maintaining a neutral center.

For data with a natural zero point, such as log2 fold changes, the color scale should be centered at zero. This centering ensures that genes with no expression change appear in the neutral color, while up-regulated and down-regulated genes appear in contrasting colors.

Breaking the Color Scale

The breaks parameter allows researchers to specify exact break points for the color gradient. This parameter accepts a numeric vector that defines the boundaries between color intervals. Using breaks provides control over the mapping between data values and colors, which is particularly useful when the data distribution is skewed.

For example, if most expression values fall between -2 and 2 but a few extreme values reach -6 and 6, the default color scaling would compress the common values into a narrow color range. Specifying breaks that concentrate color transitions in the -2 to 2 range would improve visual discrimination for the majority of genes.

The breaks vector must have a length that matches the number of colors in the color vector. Specifically, if there are n colors, there must be n+1 breaks. This requirement ensures that each color interval has defined boundaries.

Advanced Annotation Features for Complex Experimental Designs

RNA-seq experiments often involve multiple grouping variables that cannot be captured by a single annotation bar. The pheatmap package supports multiple annotation columns, allowing researchers to display treatment, time point, batch, and other variables simultaneously.

Multiple Annotation Columns

The annotation_col data frame can contain multiple columns, each representing a different annotation variable. The pheatmap function displays each column as a separate annotation bar above the heatmap. The order of bars follows the column order in the annotation data frame.

When multiple annotations are displayed, the color assignments for each variable are independent. This independence allows researchers to use consistent colors for the same factor levels across different annotation variables. For example, if treatment appears in two annotation columns, the same color can be used for the same treatment level in both columns.

The annotation_colors list must include entries for each annotation variable. If a variable is missing from the list, pheatmap assigns default colors. For publication figures, explicit color specification ensures consistency and accessibility.

Annotation for Rows

The annotation_row parameter provides the same functionality for rows as annotation_col provides for columns. This parameter accepts a data frame where each row corresponds to a gene and each column corresponds to an annotation variable such as gene family, pathway membership, or functional category.

Row annotations are displayed to the left of the heatmap. This placement allows researchers to see the functional category of each gene alongside its expression pattern. The combination of row annotations and row clustering can reveal whether genes in the same functional category share expression patterns.

For RNA-seq data, common row annotations include gene type, chromosome location, and pathway membership. These annotations add biological context to the heatmap without requiring additional figure panels.

Practical Workflow for Creating Publication-Ready Heatmaps

The process of creating an advanced pheatmap figure follows a systematic workflow that moves from data preparation through figure customization to output verification.

Step 1: Prepare the Expression Matrix

Start with a normalized expression matrix where rows are genes and columns are samples. Apply any filtering to select genes of interest, such as significant differential expression results or highly variable genes. Transform the data as needed, either through log2 transformation or row scaling.

Verify that the matrix contains no missing values, as pheatmap cannot handle NA entries. If missing values exist, address them through imputation or gene filtering before proceeding.

Step 2: Construct Annotation Data Frames

Create the annotation_col data frame with row names matching the column names of the expression matrix. Include all experimental variables that should appear in the figure. Create the annotation_row data frame if row annotations are needed.

Check that factor levels are ordered correctly and that all levels have defined colors in the annotation_colors list.

Step 3: Determine Clustering and Ordering

Decide whether to use default clustering or manual ordering. For column ordering, consider whether the experimental design has a natural structure that should be reflected in the figure. For row ordering, decide whether clustering or manual ordering based on gene categories is more appropriate.

If gaps are needed, determine the gap positions based on the final ordering. For clustered rows, cut the dendrogram to obtain clusters and calculate cumulative sizes.

Step 4: Configure Colors and Breaks

Select a color palette appropriate for the data type and distribution. Define breaks if the data distribution requires custom color mapping. Verify that the number of colors and breaks are compatible.

Step 5: Generate and Inspect the Figure

Run the pheatmap function with all parameters specified. Inspect the output figure for readability, checking that annotations are legible, gaps appear in the correct positions, and the color scale represents the data appropriately.

Save the figure in a publication-quality format such as PDF or PNG with sufficient resolution. Record all parameters used in the analysis documentation for reproducibility.

Parameter Comparison for Common RNA-seq Visualization Scenarios

The following table compares parameter configurations across common RNA-seq visualization scenarios, helping researchers select appropriate settings for their specific goals.

Records and Measurements for Reproducible Heatmap Generation

Reproducibility in bioinformatics requires documentation of all analysis parameters, including those used for visualization. The specific parameters passed to pheatmap should be recorded in the analysis script or notebook, along with the version of R and the pheatmap package.

Documenting Parameter Choices

The analysis documentation should include the input matrix and its preprocessing steps, the gene selection criteria, the distance metric and clustering method, the annotation variables and their colors, the color palette and breaks, and the gap positions. This documentation allows another researcher to reproduce the exact figure from the same input data.

Version control for analysis scripts is essential for reproducibility. Tools like Git provide a record of changes to analysis code over time, and platforms like GitHub enable sharing of analysis code with collaborators and reviewers [<a href="#ref-4">4</a>]. The Carpentries lessons provide foundational training in these computing practices [<a href="#ref-4">4</a>].

Workflow management systems offer another layer of reproducibility for complete analysis pipelines. Community-driven pipeline frameworks such as nf-core provide standardized workflow documentation and configuration practices that support reproducible analysis [<a href="#ref-5">5</a>]. While these frameworks focus on the full analysis pipeline instead of individual visualization steps, their documentation standards illustrate the level of detail expected for reproducible computational work [<a href="#ref-5">5</a>].

Verifying Figure Accuracy

Before including a heatmap in a publication or report, verify that the figure accurately represents the underlying data. Check that sample labels match the correct columns, annotation colors correspond to the correct factor levels, and the color scale spans the appropriate range of values.

For figures that will undergo peer review, consider whether the raw data and analysis code should be made available. Public repositories for sequencing data, such as those maintained by the National Center for Biotechnology Information, provide a venue for sharing the underlying data [<a href="#ref-6">6</a>]. The NCBI maintains databases and search systems that support data sharing and reuse [<a href="#ref-6">6</a>].

Common Failure Patterns in pheatmap Visualization

Several recurring problems appear when researchers create heatmaps with pheatmap. Recognizing these failure patterns helps avoid wasted effort and produces better figures.

Mismatched Row Names in Annotation Data

A frequent error occurs when the row names of the annotation data frame do not match the column names of the expression matrix. This mismatch produces an error or a figure with missing annotations. The solution is to verify that the row names of annotation_col match the column names of the expression matrix exactly, including any naming conventions.

Color Scale Obscuring Data Patterns

The default color scale may not suit the data distribution, resulting in a figure where most cells appear in similar colors. This problem occurs when the data has outliers that compress the range of common values. The solution is to examine the distribution of expression values and adjust the color breaks accordingly.

Unreadable Row and Column Labels

When the heatmap contains many genes or samples, the default label size may produce overlapping or unreadable text. The fontsize_row and fontsize_col parameters control label sizes independently. Increasing figure dimensions or reducing font sizes can resolve readability issues.

Clustering That Separates Known Groups

Sometimes the clustering algorithm separates samples that should group together based on the experimental design. This result may indicate a batch effect or other technical variation that should be investigated. The annotation display helps identify whether clustering correlates with technical factors instead of biological ones.

Gap Positions That Do Not Match Visual Groups

When gaps are specified without verifying the final column or row order, the gaps may appear in unexpected positions. This problem occurs when clustering reorders samples or genes after the gap positions were calculated. The solution is to determine gap positions after the final ordering is established, either by using manual ordering or by extracting the clustering result first.

Integration with RNA-seq Quality Control

Heatmap visualization plays a role in RNA-seq quality control, complementing other assessment methods. Sample-level heatmaps can reveal outliers, batch effects, and unexpected sample relationships before downstream analysis proceeds.

Sample-to-Sample Distance Heatmaps

A heatmap of sample-to-sample distances provides an overview of sample relationships. Samples from the same experimental group should cluster together, while samples from different groups should separate. Unexpected patterns may indicate sample mislabeling, contamination, or technical variation.

The distance matrix for this visualization is typically calculated from normalized expression data using Euclidean or correlation-based distances. The resulting heatmap uses a color scale where low distances appear in one color and high distances in another.

Top Variable Gene Heatmaps

A heatmap of the most variable genes across samples reveals the major sources of expression variation. If the top variable genes separate samples by experimental group, the biological signal is strong. If the separation correlates with batch or processing date, technical variation may dominate.

This visualization is often performed before differential expression analysis to assess data quality. The results inform decisions about whether additional normalization or batch correction is needed.

Quality Control in Training and Pipeline Contexts

Quality control practices for RNA-seq data are emphasized across bioinformatics training resources. The Galaxy Training Network provides accessible workflow training that includes quality assessment steps for sequencing data [<a href="#ref-7">7</a>]. These training materials demonstrate how visualization tools fit into the broader quality control workflow, from raw read assessment through expression quantification [<a href="#ref-7">7</a>].

Similarly, the Bioconductor project maintains packages and workflows for genomic analysis that include visualization components [<a href="#ref-8">8</a>]. The official documentation for Bioconductor packages provides installation guidance and reproducible analysis examples that researchers can adapt for their own data [<a href="#ref-8">8</a>].

Limitations of Heatmap Visualization for RNA-seq Data

Heatmaps provide a valuable visualization tool, but they have inherent limitations that researchers should understand when interpreting results.

Loss of Quantitative Precision

A heatmap displays relative expression levels through color intensity, but the human eye cannot accurately recover exact values from color alone. For precise quantitative comparisons, researchers should consult the underlying expression values or supplementary tables. The heatmap serves as a pattern-discovery tool instead of a quantitative reporting format.

Sensitivity to Gene Selection

The genes included in a heatmap determine the patterns visible in the figure. Different gene selection criteria can produce very different visual impressions of the same dataset. Researchers should document their selection criteria and consider whether the selected genes represent the biological question of interest.

Clustering Instability

Hierarchical clustering results can be sensitive to the choice of distance metric, linkage method, and the specific samples included in the analysis. Small changes in the input data can produce different cluster assignments. This instability means that cluster boundaries should not be overinterpreted as definitive biological groupings.

Computational Scaling

Very large expression matrices can be slow to process and produce figures that are difficult to interpret. The practical limit for readable heatmaps is on the order of a few hundred genes. For larger gene sets, alternative visualization approaches such as pathway-level summaries may be more appropriate.

Interpretation Context

Heatmap patterns require biological context for meaningful interpretation. A cluster of co-expressed genes may indicate shared regulatory mechanisms, but confirming this interpretation requires additional analysis such as pathway enrichment or transcription factor binding analysis. The heatmap generates hypotheses instead of confirming biological mechanisms.

Professional Escalation Criteria for Visualization Problems

Some visualization problems indicate underlying data issues that require attention beyond adjusting pheatmap parameters. Recognizing these situations helps researchers address root causes instead of symptoms.

Escalate When Clustering Reveals Unexpected Sample Relationships

If sample clustering consistently separates samples in ways that do not match the experimental design, investigate potential causes. Sample mislabeling, contamination, and batch effects can all produce this pattern. Consult with colleagues or bioinformatics support to determine whether the issue reflects a technical problem or a genuine biological signal.

Escalate When Annotation Reveals Batch Effects

If the annotation display shows that samples cluster by batch instead of by treatment, the data may require batch correction or additional normalization. This situation warrants discussion with a bioinformatician or statistician before proceeding with downstream analysis.

Escalate When Color Scaling Cannot Reveal Patterns

If no color scale or break configuration produces a figure with interpretable patterns, the data may have quality issues. Investigate the distribution of expression values, check for outliers, and verify that normalization was performed correctly.

Escalate When Heatmap Patterns Contradict Statistical Results

If the heatmap shows clear expression differences between groups but the differential expression analysis reports no significant genes, the discrepancy may indicate a problem with the statistical model or the normalization approach. This situation requires consultation with a statistician familiar with RNA-seq analysis.

Frequently Asked Questions

What is the difference between annotation_col and annotation_row in pheatmap?

The annotation_col parameter adds annotation bars above the heatmap for each sample column, while annotation_row adds annotation bars to the left of the heatmap for each gene row. Both parameters accept data frames where row names match the column names or row names of the expression matrix respectively. The annotation_col data frame has one row per sample and one column per annotation variable, while annotation_row has one row per gene and one column per annotation variable.

How do I choose the right color palette for RNA-seq expression data?

The choice of color palette depends on the data type and the message the figure should convey. For log2-transformed expression values, a diverging palette with a neutral center works well. For row-scaled z-scores, a blue-white-red gradient centered at zero provides clear visual distinction between low and high expression. Consider colorblind accessibility when selecting colors, and verify that the palette renders correctly in the final output format.

Why do my sample annotations not appear in the heatmap?

Missing annotations usually result from mismatched row names between the annotation data frame and the expression matrix. Verify that the row names of annotation_col match the column names of the expression matrix exactly. Also confirm that the annotation data frame contains no missing values and that all factor levels have defined colors in the annotation_colors list.

How do I insert gaps between treatment groups in my heatmap?

Use the gaps_col parameter with a numeric vector specifying the column positions where gaps should appear. The positions refer to the column order after any clustering or manual ordering. For reliable gap placement, order the columns manually to match the experimental design, then specify gap positions at the boundaries between groups.

What is the best way to cluster genes in a heatmap?

The best clustering approach depends on the data and the biological question. Correlation-based distances often work well for expression data because they focus on pattern similarity instead of absolute magnitude. Ward's linkage method tends to produce compact, visually interpretable clusters. Experiment with different combinations and inspect the results to determine what works best for your data.

How many genes should I include in a heatmap?

The number of genes should balance biological relevance with visual clarity. A heatmap of 50 to 200 genes produces a readable figure where individual gene patterns can be distinguished. Larger selections require larger figure dimensions and may obscure individual gene patterns. Select genes based on differential expression significance, variability, or pathway membership instead of including all genes in the dataset.

Can I use pheatmap for single-cell RNA-seq data visualization?

Pheatmap can visualize single-cell RNA-seq data, but the scale of such data often requires aggregation before visualization. Common approaches include averaging expression across cell clusters or selecting marker genes for visualization. The same annotation and gap features apply, with cell cluster identity serving as a natural annotation variable. Single-cell studies frequently integrate multiple datasets and require careful quality control before visualization, as demonstrated in curated single-cell resources for immune checkpoint blockade research [<a href="#ref-9">9</a>].

How do I save a pheatmap figure for publication?

Use the pdf or png functions to open a graphics device, run the pheatmap function, then close the device with dev.off. Specify appropriate dimensions to ensure readability, and consider the target journal requirements for figure resolution and format. For large heatmaps, vector formats like PDF provide the best quality.

Related Bioinformatics Analysis Context

Heatmap visualization does not operate in isolation within the RNA-seq analysis workflow. The interpretation of heatmap patterns gains meaning when connected to the broader analytical context, including differential expression testing, pathway analysis, and integration with other data types.

Connecting Heatmap Patterns to Differential Expression Results

The genes displayed in a post-analysis heatmap typically come from differential expression testing. The statistical methods used to identify these genes influence which genes appear and therefore shape the visual patterns. Methods that account for the discrete and overdispersed nature of RNA-seq counts, such as Generalized Linear Models with Quasi-Likelihood, provide more reliable gene rankings than approaches assuming normal distributions [<a href="#ref-3">3</a>]. Researchers should ensure that the gene selection for heatmap visualization follows sound statistical practice.

Integrating Heatmap Visualization with Multi-omics Data

RNA-seq analysis increasingly integrates with other molecular data types. Studies combining RNA-seq with chromatin accessibility data, for example, use integrated analysis platforms to jointly interpret transcriptomic and regulatory information [<a href="#ref-10">10</a>]. Heatmap visualization of expression data can be complemented by similar visualizations of chromatin accessibility or methylation data, allowing researchers to compare patterns across molecular layers.

The interpretation of expression patterns in heatmaps also benefits from knowledge of the underlying biological system. Studies of specific diseases or organisms often identify gene signatures that appear as coordinated expression patterns in heatmaps. For example, transcriptomic analyses of cancer samples have identified gene expression signatures associated with metastasis and treatment response [<a href="#ref-3">3</a>][<a href="#ref-11">11</a>]. These signatures often appear as distinct blocks in heatmaps when the genes are clustered appropriately.

Reproducibility Standards for Visualization Code

The computational reproducibility of heatmap figures depends on the same standards that apply to other analysis steps. Training resources from The Carpentries emphasize version control and reproducible computing practices that apply directly to visualization code [<a href="#ref-4">4</a>]. Researchers should maintain their pheatmap scripts under version control and document the package versions used.

Workflow management frameworks provide structured approaches to reproducible analysis. The nf-core community maintains documentation for pipeline standards that emphasize reproducibility and configuration management [<a href="#ref-5">5</a>]. While these pipelines focus on primary analysis, the documentation principles apply to visualization steps as well.

Related Bioinformatics Guides

• RNA-Seq Databases: Accessing and Using Public RNA-Seq Data
• RNA-Seq vs ChIP-Seq: Complementary Approaches for Gene Regulation
• RNA-Seq Data Analysis in Galaxy: A User-Friendly Platform
• RNA-Seq Data Analysis Workflow: From Raw Reads to Insights
• RNA-Seq Visualization: Volcano Plots, Heatmaps, and PCA

Related Clinical & Scientific Guides

• A Practical Guide to Detecting Antimicrobial Resistance Genes in Shotgun Metagenomic Data
• Computational Immunology: Modeling the Immune System
• How to Set Hard Filters for Germline Variant Calling: A Practical Guide to GATK Best Practices

References and Further Reading

This article is educational and does not replace validated analysis plans, institutional policy, clinical interpretation, or specialist review.