bio-data-visualization-oncoprint-mutation-matrices

What this skill does

## Version Compatibility

Reference examples tested with: ComplexHeatmap 2.18+, maftools 2.18+, comut 0.0.3+, MAFtools requires R 4.0+; comut.py requires pandas 2.0+, matplotlib 3.8+.

Before using code patterns, verify installed versions match. If versions differ:
- Python: `pip show <package>` then `help(module.function)` to check signatures
- R: `packageVersion('<pkg>')` then `?function_name`

If code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.

# OncoPrint and Mutation Matrix Plots

**"Plot mutations across a cohort"** -> Render a gene-by-sample matrix where each cell stacks colored rectangles encoding alteration class (missense, truncating, splice, copy-gain, copy-loss, fusion). Sort samples by burden, optionally split by clinical group, and overlay co-mutation / mutual-exclusivity annotations. OncoPrint (Cerami 2012 *Cancer Discov* 2:401; canonical at cBioPortal) is the genre-defining visualization.

- R: `ComplexHeatmap::oncoPrint`, `maftools::oncoplot`
- Python: `comut.CoMut`, `cbioportal`-style implementations

## The Single Most Important Modern Insight -- Cell Stacking Encodes Multiple Alterations Per Cell

OncoPrint differs from a generic heatmap because each cell can encode multiple alterations simultaneously through *stacked* rectangles. A patient with both a missense and a copy-gain in TP53 shows one cell with two overlapping colored rectangles (e.g., green diamond inside red square). This stacking is the whole point — it preserves the multi-modal alteration landscape that flattening to a single category destroys.

In ComplexHeatmap's `oncoPrint`, the `alter_fun` argument is the rendering specification: a named list of functions, one per alteration class, each drawing its rectangle inside the cell. Get this right and the figure works; get it wrong and overlapping alterations are invisible.

## Decision Tree by Cohort and Question

| Question | Sort by | Display |
|----------|---------|---------|
| Which genes are most altered? | Gene frequency (default) | Bar above samples (sample TMB); bar right of genes (gene frequency) |
| Per-patient burden patterns | Sample burden | TMB bar on top; sample-name labels |
| Subtype-driver enrichment | Clinical group then burden | `column_split` by group; per-group frequency right bar |
| Mutual exclusivity (BRAF vs NRAS) | Custom (alphabetic-by-mutation pattern) | Memo sort; overlay log10(OR) heatmap |
| Co-occurrence (TP53 + MYC) | Custom | Same pattern; positive OR coloring |
| Driver vs passenger comparison | Two panels | Concatenate two oncoPrints horizontally |

## ComplexHeatmap::oncoPrint -- Canonical Implementation

**Goal:** Render a cohort mutation matrix with stacked alteration-class encoding, sample annotations, and a sample-sorted, gene-frequency-ranked layout.

**Approach:** Convert the MAF/variant table to a gene-by-sample matrix of `;`-delimited alteration strings; define `alter_fun` rendering one rectangle per class; pass to `oncoPrint()` with column annotations.

```r
library(ComplexHeatmap)
library(circlize)

# Input: matrix where each cell is a string like 'Missense;Amp' or '' for no alteration
# Rows = genes; columns = samples

# Color per alteration class
col <- c('Missense'   = '#56B4E9',
         'Truncating' = '#000000',
         'Splice'     = '#CC79A7',
         'Amp'        = '#D55E00',
         'HomDel'     = '#0072B2',
         'Fusion'     = '#009E73')

# alter_fun -- one function per class, each drawing inside the cell
alter_fun <- list(
    background = function(x, y, w, h)
        grid.rect(x, y, w - unit(0.5, 'mm'), h - unit(0.5, 'mm'),
                  gp = gpar(fill = '#EEEEEE', col = NA)),
    Amp = function(x, y, w, h)
        grid.rect(x, y, w - unit(0.5, 'mm'), h - unit(0.5, 'mm'),
                  gp = gpar(fill = col['Amp'], col = NA)),
    HomDel = function(x, y, w, h)
        grid.rect(x, y, w - unit(0.5, 'mm'), h - unit(0.5, 'mm'),
                  gp = gpar(fill = col['HomDel'], col = NA)),
    Missense = function(x, y, w, h)
        grid.rect(x, y, w - unit(0.5, 'mm'), h * 0.5,
                  gp = gpar(fill = col['Missense'], col = NA)),
    Truncating = function(x, y, w, h)
        grid.rect(x, y, w - unit(0.5, 'mm'), h * 0.33,
                  gp = gpar(fill = col['Truncating'], col = NA)),
    Splice = function(x, y, w, h)
        grid.rect(x, y, w - unit(0.5, 'mm'), h * 0.25,
                  gp = gpar(fill = col['Splice'], col = NA)),
    Fusion = function(x, y, w, h)
        grid.points(x, y, pch = 17, size = unit(2, 'mm'),
                    gp = gpar(col = col['Fusion'])))

# Clinical column annotation
ha_clin <- HeatmapAnnotation(
    Subtype = clinical$subtype,
    Stage   = clinical$stage,
    col = list(Subtype = c(Luminal='#0072B2', Basal='#D55E00', HER2='#009E73'),
               Stage   = c(I='#FFFFCC', II='#FED976', III='#FD8D3C', IV='#BD0026')))

oncoPrint(mat,
          alter_fun = alter_fun,
          col = col,
          top_annotation = ha_clin,
          column_title = 'TCGA-BRCA mutation landscape',
          row_names_gp = gpar(fontsize = 8),
          pct_gp = gpar(fontsize = 7),
          show_pct = TRUE,
          remove_empty_columns = FALSE,
          remove_empty_rows = FALSE)
```

## maftools::oncoplot -- Faster Onboarding

For TCGA-style MAF files, `maftools::oncoplot` is the lower-friction option:

```r
library(maftools)
maf <- read.maf(maf = 'tcga.maf', clinicalData = clinical)
oncoplot(maf = maf,
         top = 20,                            # top 20 mutated genes
         clinicalFeatures = c('Subtype', 'Stage'),
         annotationColor = list(Subtype = c(Luminal='#0072B2', Basal='#D55E00'),
                                 Stage = c(I='#FFFFCC', IV='#BD0026')),
         sortByAnnotation = TRUE,
         removeNonMutated = FALSE)
```

maftools defaults handle alteration-class colors, sample sorting, and percentage bars automatically. Customization is more limited than ComplexHeatmap.

## Mutual Exclusivity and Co-Occurrence

```r
# maftools provides somaticInteractions
si <- somaticInteractions(maf = maf, top = 20,
                          pvalue = c(0.05, 0.01),
                          fontSize = 0.7)
# Plot returns a matrix of -log10(p) with sign by direction (+ co-occur, - mutex)
```

Mutual-exclusivity testing on small cohorts (N < 50) is underpowered; reported "significant" mutex on n=20 with 2 mutations each is uninterpretable. Aggregate to larger cohorts (TCGA + ICGC pan-cancer) or report effect size with CI rather than p-value.

**Fisher exact vs DISCOVER:** standard 2x2 Fisher tests sample-mutation pairs, ignoring per-gene mutation rate background. DISCOVER (Canisius 2016 *Genome Biol* 17:261) models per-tumor mutation probability and is preferred for pan-cancer analyses where mutation rate varies 100× across samples.

## comut.py -- Python Equivalent

```python
import comut
import pandas as pd

# Long-format: columns = sample, category (gene), value (alteration class)
toy_comut = comut.CoMut()
toy_comut.add_categorical_data(
    data=mutation_long_df,
    name='Mutations',
    category_order=top_genes,
    value_order=['Truncating', 'Missense', 'Splice', 'Amp', 'HomDel'],
    mapping={'Truncating': '#000000', 'Missense': '#56B4E9',
             'Splice': '#CC79A7', 'Amp': '#D55E00', 'HomDel': '#0072B2'})

toy_comut.add_categorical_data(
    data=clinical_long_df,
    name='Subtype',
    mapping={'Luminal': '#0072B2', 'Basal': '#D55E00'})

toy_comut.add_continuous_data(
    data=tmb_long_df,
    name='TMB',
    mapping='viridis',
    value_range=(0, 30))

toy_comut.plot_comut(figsize=(12, 8))
toy_comut.figure.savefig('comut.pdf', dpi=300, bbox_inches='tight')
```

## Per-Method Failure Modes

### Alterations flattened to a single class

**Trigger:** Reducing each cell to a single most-severe alteration, losing the stack.

**Mechanism:** Loses the multi-alteration biology (e.g., MYC amp +