Lesson 1 — Introduction to Scanpy & AnnData

What is Scanpy?

Scanpy (Single-Cell ANalysis in Python) is the primary Python toolkit for analysing single-cell RNA-sequencing (scRNA-seq) data. It was developed at the Helmholtz Centre Munich and published in 2018 by Wolf et al. in Genome Biology. Today it is the most-used scRNA-seq framework in the Python ecosystem, with a vibrant extension package landscape including scvi-tools, squidpy, cellrank, and more.

Why Scanpy instead of Seurat?
Seurat (R) is the dominant framework in R-centric labs. Scanpy is its Python equivalent. Both implement the same biological workflow: QC → normalise → HVGs → PCA → UMAP → cluster → annotate. Choosing Scanpy means you stay inside the Python data-science ecosystem (NumPy, Pandas, scikit-learn, matplotlib), making integration with machine-learning pipelines effortless. For a plant genomics researcher, Scanpy is also straightforward to combine with GWAS summary statistics and custom Python scripts.

The Scanpy workflow follows eight canonical steps identical to Seurat's — only the syntax changes:

Step	Seurat function	Scanpy equivalent
Load data	`Read10X()`	`sc.read_10x_mtx()`
QC filtering	`subset()`	`sc.pp.filter_cells()`
Normalise	`NormalizeData()`	`sc.pp.normalize_total()`
HVGs	`FindVariableFeatures()`	`sc.pp.highly_variable_genes()`
Scale	`ScaleData()`	`sc.pp.scale()`
PCA	`RunPCA()`	`sc.tl.pca()`
UMAP	`RunUMAP()`	`sc.tl.umap()`
Cluster	`FindClusters()`	`sc.tl.leiden()`

The AnnData Object

The heart of Scanpy is the AnnData (Annotated Data) object. It is a single Python object that holds everything about your experiment: raw counts, QC metrics, dimensionality reductions, cluster labels, and metadata — all linked together.

Why one object for everything?
In classical bulk RNA-seq (DESeq2), you typically keep a counts matrix, a colData table, and results tables as separate R objects. With single-cell data — where you may have 10,000–50,000 cells — managing separate objects becomes error-prone. AnnData enforces a linked structure: the row (cell) axis of every annotation is guaranteed to match the row axis of the count matrix. You cannot accidentally apply a cluster label from sample A to the cells of sample B.

AnnData object with n_obs × n_vars = 2700 × 32738 adata.X — count matrix (2700 cells × 32738 genes) # sparse CSR matrix adata.obs — cell metadata (DataFrame, 1 row per cell) # n_genes, pct_mt … adata.var — gene metadata (DataFrame, 1 row per gene) # n_cells, highly_variable … adata.obsm — cell embeddings (dict of arrays) # PCA, UMAP coords adata.obsp — cell graphs (sparse matrix) # kNN adjacency adata.varm — gene embeddings (dict of arrays) # PCA loadings adata.uns — unstructured (dict) # colours, params adata.layers — extra matrices (dict of matrices) # raw counts backup

Think of AnnData like a spreadsheet where rows are cells and columns are genes, but with extra filing cabinets attached for dimensionality reductions (obsm), graphs (obsp), and free-form notes (uns).

💡

obs vs var: obs = observations = cells. var = variables = genes. This follows the statistics convention where "observations" are samples and "variables" are features. In scRNA-seq, each cell is one observation and each gene is one variable.

Installing Scanpy

Why a dedicated conda environment?
Scanpy depends on specific versions of NumPy, anndata, leidenalg, igraph, and umap-learn. Mixing these with your existing system Python or other conda environments can cause version conflicts that are very hard to debug. Creating a clean scanpy-env is the professional standard — it matches exactly what you would do on an HPC cluster or in a Docker container for a published analysis.

📁 Run from: any terminal (home directory is fine)

Bash — Create environment

# Create a dedicated conda environment for Scanpy
conda create -n scanpy-env python=3.11 -y

# Activate it
conda activate scanpy-env

# Install Scanpy and core dependencies
pip install scanpy

# Install Leiden algorithm for clustering (required for sc.tl.leiden)
pip install leidenalg igraph

# Install Jupyter for interactive analysis
pip install jupyter notebook

# Verify installation
python -c "import scanpy as sc; print(sc.__version__)"
# Expected: 1.10.x (or similar)

⚠️

leidenalg is not installed automatically by Scanpy. If you run sc.tl.leiden() without it, you get an ImportError. Always install leidenalg and igraph alongside Scanpy. Alternatively use sc.tl.louvain() (fewer dependencies), but Leiden is the modern standard.

Prefer a single conda command? Use this all-in-one approach:

Bash — All-in-one conda install

conda create -n scanpy-env -c conda-forge -c bioconda \
    scanpy leidenalg python-igraph jupyter -y

First Import & Settings

Why set global settings at the top?
sc.settings controls how Scanpy saves figures, how verbose it is, and what figure format it uses. Setting these at the very top of your script means every plot generated later automatically goes to the right folder with the right resolution. This is essential for reproducibility — when you re-run the script six months later, you get identical output files in the same location.

📁 Run from: ~/scrna-scanpy-pbmc/

Python

# ── Standard Scanpy imports ─────────────────────────────────
import scanpy as sc
import anndata as ad
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# ── Global settings ─────────────────────────────────────────
sc.settings.verbosity = 3             # 0=errors only, 3=hints (recommended)
sc.settings.set_figure_params(
    dpi=100,                             # screen resolution
    dpi_save=300,                        # publication resolution
    figsize=(6, 5),
    facecolor='white'
)
sc.settings.figdir = './figures/'     # all sc.pl plots saved here

# ── Reproducibility: fix the random seed ────────────────────
np.random.seed(42)

# ── Print version info ──────────────────────────────────────
sc.logging.print_header()
# scanpy==1.10.x  anndata==0.10.x  umap==0.5.x  numpy==1.26.x ...

💡

sc.logging.print_header() prints the exact version of every key dependency. Always include this at the top of published analysis scripts — reviewers can reproduce your environment precisely by reading these version numbers.

Loading PBMC Data

We use the PBMC 3k dataset — 2,700 peripheral blood mononuclear cells from a healthy donor, sequenced with 10x Genomics Chromium. This is the "Hello, World" of scRNA-seq: small enough to process in minutes on a laptop, well-documented, with known cell types (T cells, B cells, NK cells, monocytes).

What is the 10x MTX format?
The 10x Chromium pipeline (Cell Ranger) outputs three files into filtered_feature_bc_matrix/:
• matrix.mtx.gz — the sparse count matrix in Market Exchange Format
• barcodes.tsv.gz — one cell barcode per row (these become row names = cell IDs)
• features.tsv.gz — one gene per row with Ensembl ID and gene name
sc.read_10x_mtx() reads all three files and assembles them into an AnnData object automatically.

📁 Run from: ~/scrna-scanpy-pbmc/

Python — Load PBMC 3k (built-in)

# ── Option A: Scanpy's built-in PBMC 3k loader ──────────────
# Downloads ~22 MB to ~/.cache/scanpy/ on first run
adata = sc.datasets.pbmc3k()

# AnnData object with n_obs x n_vars = 2700 x 32738

To load your own Cell Ranger output (what you would do in a real project):

Python — Load from Cell Ranger folder

# ── Option B: Read from Cell Ranger output folder ───────────
adata = sc.read_10x_mtx(
    './data/raw/pbmc3k/filtered_feature_bc_matrix/',
    var_names='gene_symbols',   # use gene names, not Ensembl IDs
    cache=True                # cache to disk for faster re-loading
)

# Make gene names unique (some symbols appear at two loci)
adata.var_names_make_unique()

💡

var_names='gene_symbols' uses human-readable gene names (CD3D, MS4A1) instead of Ensembl IDs (ENSG00000167286). For sorghum data, gene names follow the Sorghum bicolor annotation convention (e.g. Sobic.001G000100) but the parameter works identically.

Inspecting the Object

Why inspect before analysing?
Your first action after loading any dataset should always be to examine its shape and structure. In bulk RNA-seq you would run dim(dds) or head(colData). In Scanpy the equivalent is printing the AnnData object and exploring its slots. This confirms the data loaded correctly, tells you how many cells and genes you have, and reveals whether metadata is already present.

📁 Run from: ~/scrna-scanpy-pbmc/

Python — Inspect AnnData

# ── Basic summary ────────────────────────────────────────────
print(adata)
# AnnData object with n_obs x n_vars = 2700 x 32738

# ── Shape: (cells, genes) ────────────────────────────────────
print(f"Cells: {adata.n_obs}, Genes: {adata.n_vars}")
# Cells: 2700, Genes: 32738

# ── Cell metadata (obs) — empty at this stage ────────────────
print(adata.obs.head())
# Empty DataFrame (QC metrics added in Lesson 2)

# ── Gene metadata (var) ──────────────────────────────────────
print(adata.var.head())
#                  gene_ids        feature_types
# MIR1302-2HG    ENSG00000243485   Gene Expression
# FAM138A        ENSG00000237613   Gene Expression

# ── First few cell barcodes ──────────────────────────────────
print(adata.obs_names[:5])
# Index(['AAACATACAACCAC-1', 'AAACATTGAGCTAC-1', ...])

# ── Count matrix type ────────────────────────────────────────
print(type(adata.X))
# <class 'scipy.sparse.csr_matrix'>

# ── Sparsity ─────────────────────────────────────────────────
sparsity = 1.0 - adata.X.nnz / (adata.n_obs * adata.n_vars)
print(f"Sparsity: {sparsity:.1%}")
# Sparsity: ~94.4%

# ── Memory comparison ────────────────────────────────────────
dense_mb = (adata.n_obs * adata.n_vars * 4) / 1e6
print(f"Dense would use: {dense_mb:.0f} MB — sparse uses ~8 MB")

The ~94 % sparsity is normal and expected for scRNA-seq data. Each cell only expresses a fraction of all genes. This is why Scanpy stores adata.X as a sparse matrix (scipy CSR format) — a dense float32 matrix would use ~354 MB; the sparse version uses ~8 MB.

⚠️

Always call adata.var_names_make_unique() right after loading. Some genomes have duplicate gene symbols (two loci annotated as TBCE, for example). Scanpy will raise cryptic errors later in the pipeline if gene names are not unique.

Quick Reference

Slot	What it stores	Example access
`adata.X`	Count matrix (cells × genes), sparse	`adata.X.toarray()`
`adata.obs`	Cell metadata DataFrame	`adata.obs['n_genes']`
`adata.var`	Gene metadata DataFrame	`adata.var['highly_variable']`
`adata.obsm`	Cell embeddings (PCA, UMAP)	`adata.obsm['X_pca']`
`adata.obsp`	Cell-cell adjacency (kNN graph)	`adata.obsp['connectivities']`
`adata.uns`	Unstructured metadata dict	`adata.uns['leiden_colors']`
`adata.layers`	Extra count matrices	`adata.layers['raw_counts']`
`adata.n_obs`	Number of cells (int)	`print(adata.n_obs)`
`adata.n_vars`	Number of genes (int)	`print(adata.n_vars)`
`adata.obs_names`	Cell barcode index	`adata.obs_names[:5]`
`adata.var_names`	Gene name index	`adata.var_names[:5]`

Exercises

Print the AnnData summary

After loading the PBMC 3k dataset with sc.datasets.pbmc3k(), print the full AnnData summary. How many cells and genes does it contain? Which slots are populated at this stage?

Show answer

Python

adata = sc.datasets.pbmc3k()
print(adata)
# n_obs x n_vars = 2700 x 32738
# Only adata.X and adata.var (gene_ids, feature_types) are populated.
# adata.obs is empty — QC metrics are added in Lesson 2.

Explore sparsity and memory

Calculate the sparsity of the PBMC count matrix. Then compare the hypothetical dense size (float32) against the actual sparse memory usage.

Show answer

Python

sparsity = 1.0 - adata.X.nnz / (adata.n_obs * adata.n_vars)
print(f"Sparsity: {sparsity:.1%}")  # ~94.4%

dense_mb = (adata.n_obs * adata.n_vars * 4) / 1e6
print(f"Dense: {dense_mb:.0f} MB")    # ~354 MB

sparse_mb = (adata.X.data.nbytes + adata.X.indices.nbytes + adata.X.indptr.nbytes) / 1e6
print(f"Sparse: {sparse_mb:.1f} MB")  # ~8-12 MB

Access a single gene's counts

Extract the counts for CD3D (T-cell marker) across all 2,700 cells. Print the first 10 values and the mean expression.

Show answer

Python

cd3d_counts = adata[:, 'CD3D'].X.toarray().flatten()
print(cd3d_counts[:10])
print(f"Mean CD3D: {cd3d_counts.mean():.4f}")
# Low (~0.3) because most cells are not T cells

Set up your project folder

Create the standard directory structure for this Scanpy module, then verify your installation produces the expected version header.

Show answer

Bash

mkdir -p ~/scrna-scanpy-pbmc/{data/raw,figures,results,scripts}
conda activate scanpy-env
python -c "import scanpy as sc; sc.logging.print_header()"