Single-Cell RNA Sequencing Analysis With AI

The standard single-cell pipeline is six steps long and two of them silently decide how many cell types you will report. Knowing which two is most of what separates a defensible analysis from a picture.

The matrix you start with

After demultiplexing and alignment, the data is a matrix of counts: rows are cells, columns are genes, and the entry is how many unique transcripts of that gene were captured in that cell. A typical droplet-based experiment gives something like 10,000 cells by 20,000 genes.

The matrix is overwhelmingly zero — commonly 90 to 95 percent — and the reason for the zeros is the thing to be clear about, because it decides how you treat them. Some are biological: the gene is not expressed in that cell. Some are technical: the gene is expressed, but capture efficiency for any single transcript is low, so a lowly expressed gene is often missed entirely. The two are not distinguishable from the matrix, and pipelines that treat all zeros as biological absence and pipelines that impute all of them as missing data both introduce artefacts. The mainstream position is to leave the zeros alone and use methods that tolerate them.

Total counts per cell also vary by an order of magnitude, driven by cell size, RNA content and capture efficiency rather than by biology. Every normalisation step below exists to stop that variation dominating.

Quality control and normalisation

Three filters are near-universal, and each targets a specific failure of the droplet chemistry.

• Mitochondrial fraction. A cell whose membrane broke before encapsulation loses cytoplasmic RNA while mitochondrial RNA stays put, so a high fraction of mitochondrial reads is the signature of a dying or damaged cell. The threshold is tissue-dependent — cardiomyocytes are legitimately mitochondria-rich — so a single cut-off copied from a tutorial will delete a real population.
• Detected genes per cell. Very low counts indicate an empty droplet containing only ambient RNA; very high counts, roughly double the mode, indicate two cells in one droplet. Doublet detection tools do better than a threshold by simulating artificial doublets from the observed data and scoring each real cell against them.
• Ambient RNA. Free-floating transcripts from lysed cells contaminate every droplet, so a marker gene from one cell type appears at low level everywhere. This is estimated from empty droplets and subtracted, and if you skip it, every cluster looks slightly like every other one.

Normalisation then divides each cell’s counts by its total, scales to a fixed size factor — 10,000 is the conventional choice — and takes a natural log of one plus the result. The plus-one keeps zeros at zero. The log is doing real work: expression is roughly multiplicative, so a gene going from 10 to 20 counts and one going from 100 to 200 are the same biological event, and only on a log scale are they the same distance.

Two reductions doing different jobs

Next, select highly variable genes — typically 2,000 — by fitting the relationship between a gene’s mean and its variance across cells and keeping the genes whose variance most exceeds what their mean predicts. Genes that are uniformly high or uniformly absent carry no information about which cell is which, and including them adds noise dimensions.

Then principal component analysis on the scaled, variable-gene matrix, keeping something like 50 components. This is the first reduction and it is linear and reversible; each component has loadings you can inspect to see which genes drive it, which is how you catch a component that is really cell-cycle phase or dissociation stress. See what PCA is actually computing for the mechanism.

The second reduction is UMAP or t-SNE down to two dimensions, and its only job is the picture. This is the most misread object in the field. Distances on a UMAP plot are not meaningful: the algorithm preserves local neighbourhood structure and explicitly does not preserve global distances, so two clusters being far apart on the plot says nothing about how different they are, and the apparent size of a cluster says nothing about its variance. Nothing downstream should be computed from the two-dimensional coordinates. Clustering happens on the principal components, not on the embedding you are looking at.

Graph clustering, and the resolution knob

Clustering is done on a graph rather than on the coordinates. Build a k-nearest-neighbour graph in the principal-component space — k around 15 by convention — then optionally weight edges by shared neighbours, then partition the graph by community detection. The Leiden algorithm is the current default, having replaced Louvain after it was shown that Louvain can produce internally disconnected communities.

These methods take a resolution parameter, and it is the knob that decides your answer. Higher resolution yields more, smaller communities; lower yields fewer, larger ones. There is no value of it that is correct, because there is no threshold in the data at which a cell state becomes a cell type. Doubling resolution will split your T cells into subtypes, and doubling it again will split those. Anyone reporting “we identified 14 cell types” has reported a choice of resolution as much as a property of the tissue, and the honest version of that sentence says which resolution and shows stability across a range of them.

Turning clusters into cell types

A cluster is a set of cell indices; making it a cell type requires outside knowledge. Two routes are used, and they fail differently.

Marker-based annotation ranks genes by differential expression between one cluster and the rest, then matches the top genes against a curated marker list. It is interpretable and it depends entirely on the marker list being right for your tissue and species. Reference-based mapping instead projects your cells into the space of an annotated atlas and transfers labels by nearest neighbours. It is fast and consistent and it cannot discover anything the reference does not contain — a genuinely novel population gets confidently labelled as the nearest thing in the atlas, with no signal that anything went wrong. Running both and examining where they disagree is more informative than either.

Three traps in the standard pipeline

• Double dipping. Clustering on expression and then testing for differential expression between those clusters, on the same data, uses the data twice. The clusters were constructed to maximise separation, so the p-values from that test are anticonservative by a wide margin — this is a known selective-inference problem, and the corrections for it involve sample splitting or explicitly conditioning on the clustering. Marker gene p-values from a naive pipeline are descriptive, not inferential, and should be reported as such.
• Batch effects that look like biology. Cells processed on different days or from different donors separate cleanly, and that separation is often larger than the biological signal. Integration methods correct for it, but every one of them can also remove real differences if the biological variable is confounded with the batch — which it is, if each donor is one condition. The design decision precedes the analysis and cannot be fixed afterwards.
• Dissociation artefacts. Getting cells out of a tissue requires enzymatic digestion, which itself induces stress-response genes. A cluster defined by immediate-early genes may be an artefact of the protocol rather than a state in the tissue, and the check is whether it appears in nuclei preparations of the same tissue.