Bio Single Cell Markers Annotation

Version Compatibility

Reference examples tested with: scanpy 1.10+, Seurat 5.0+, anndata 0.10+

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 to verify parameters

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

Marker Gene Detection and Manual Annotation

"Find marker genes for my clusters" -> Rank genes that separate each cluster from the rest, then map clusters to cell types using canonical markers.

• Python: sc.tl.rank_genes_groups() -> filter by effect size + fraction expressing -> adata.obs[...].map(labels)
• R: Seurat::FindAllMarkers() -> filter by avg_log2FC + pct.1/pct.2 -> RenameIdents()

Governing principle

Marker detection is descriptive ranking, NOT inference. Two distinct questions get sloppily called "DE" and must never be conflated: (1) marker detection - "which genes are higher in cluster X vs the rest?" is a ranking/annotation task where the cell is the unit and Wilcoxon is an acceptable heuristic; (2) condition DE - "which genes change in cell type X between treated and control?" is a population claim that requires biological replicates, where the unit of replication is the sample/donor, not the cell. Question 2 must use pseudobulk (aggregate raw counts per sample x cell type, then DESeq2/edgeR/limma-voom); treating cells as replicates is pseudoreplication and inflates false positives by orders of magnitude (Squair 2021). Post-clustering marker p-values are double-dipping. Clusters were defined to maximize between-group separation, so testing those same clusters for markers tests a hypothesis built from the data used to test it. The Wilcoxon/t-test null assumes fixed a-priori labels; under a single homogeneous population the statistic does not follow its nominal null, type-I error approaches 1 as resolution rises, and BH correction does nothing because the p-values are invalid before correction. Cluster-marker p-values are descriptive labels, never evidence that a cluster is a real cell type. Rank and filter markers by effect size and fraction-expressing, not by p-value; a gene can be "significant" at p=1e-40 (n is thousands) yet useless as a marker (60% in-group vs 55% out-group).

This skill covers marker discovery for clusters plus manual labeling. Automated reference-based label transfer (SingleR, CellTypist, Azimuth, scANVI) lives in single-cell/cell-annotation. Cross-condition compositional change lives in single-cell/differential-abundance.

Choosing a marker / DE method

Marker tools (rank_genes_groups, FindMarkers) will technically run a treatment-vs-control contrast cell-by-cell and return tidy tiny p-values. That is statistically invalid for a population claim. The tool not stopping the user is why this error is so common. When methods compete, verify current defaults against installed docs.

Defaults that bite (verify before trusting tutorials)

Scanpy marker detection

Goal: Rank cluster-specific markers and filter them by specificity, not p-value alone.

Approach: Run Wilcoxon explicitly (scanpy's default is t-test), pull results to a DataFrame with pts=True for in/out fraction, then keep genes with a large positive log fold change and a high in-group / low out-group fraction.

import scanpy as sc

adata = sc.read_h5ad('clustered.h5ad')

sc.tl.rank_genes_groups(adata, groupby='leiden', method='wilcoxon', pts=True, corr_method='benjamini-hochberg')
markers = sc.get.rank_genes_groups_df(adata, group=None)

specific = markers[(markers['logfoldchanges'] > 1) & (markers['pct_nz_group'] > 0.5) & (markers['pct_nz_reference'] < 0.25)]
print(specific.groupby('group').head(10)[['group', 'names', 'logfoldchanges', 'pct_nz_group', 'pct_nz_reference']])

Seurat marker detection

Goal: Rank markers per cluster and keep specific ones for labeling.

Approach: Run FindAllMarkers with only.pos=TRUE, install presto so Wilcoxon is fast, then rank within cluster by avg_log2FC and require a pct.1-pct.2 gap.

library(Seurat)
library(dplyr)

all_markers <- FindAllMarkers(seurat_obj, only.pos = TRUE, logfc.threshold = 0.25, min.pct = 0.1)

specific <- all_markers %>%
    filter(p_val_adj < 0.05, avg_log2FC > 1, (pct.1 - pct.2) > 0.2) %>%
    group_by(cluster) %>%
    slice_max(n = 10, order_by = avg_log2FC)
print(specific)

Seurat v5 lowers thresholds to 0.1/0.01, so explicit logfc.threshold=0.25 and a pct.1-pct.2 filter restore a marker-grade (specific) shortlist from a permissive run.

Gene signature scoring

Goal: Score each cell for a curated panel without library-size confounding.

Approach: Both tools subtract an expression-binned control set; thresholds are dataset-relative and must never be ported as absolute cutoffs.

t_cell_panel = ['CD3D', 'CD3E', 'CD4', 'CD8A', 'CD8B']
sc.tl.score_genes(adata, gene_list=t_cell_panel, ctrl_size=50, n_bins=25, score_name='T_cell_score')

seurat_obj <- AddModuleScore(seurat_obj, features = list(c('CD3D', 'CD3E', 'CD4', 'CD8A', 'CD8B')), ctrl = 100, name = 'T_cell_score')

scanpy uses 25 control bins, Seurat uses 24 by default (both follow Tirosh 2016) - a real cross-ecosystem non-reproducibility source for small panels.

Cell-cycle scoring

Goal: Assign each cell an S and G2/M score and a phase, to diagnose (and optionally regress) cell-cycle-driven structure.

Approach: Score the Tirosh S and G2/M gene panels; both tools ship the lists. Regression is optional and confounded with biology (cycling is a real state in proliferating populations) - diagnose first and regress only when the cycle is a confound, not reflexively.

sc.tl.score_genes_cell_cycle(adata, s_genes=s_genes, g2m_genes=g2m_genes)

seurat_obj <- CellCycleScoring(seurat_obj, s.features = cc.genes.updated.2019$s.genes, g2m.features = cc.genes.updated.2019$g2m.genes)

Provide s_genes/g2m_genes as the Tirosh 2016 panels (Seurat's cc.genes.updated.2019 exposes both lists directly); scanpy ships no built-in list, so load the panels from the reference or a regev-lab gene file.