Nature Biotechnology
\n (2026) Cite this article<\/p>\n
Single-cell omics methods provide sparse, noisy measurements of high-dimensional cell states, whose underlying distributions remain poorly understood, raising an urgent need for exploratory analysis and visualization methods. However, current methods are ad hoc and uninterpretable, and they distort the structure in the data. We overcome these challenges by representing data on trees and present Bonsai, a method that reconstructs the most likely tree relating any set of high-dimensional objects with arbitrary heterogeneous measurement noise. Bonsai automatically regularizes noise, accurately recovers differentiation trajectories, preserves high-dimensional distances and improves nearest-neighbor identification. When applied to blood cell data, Bonsai not only accurately recovers known lineage relationships but also discovers a subtype of natural killer (NK) cells deriving from the myeloid lineage, pinpointing genes distinguishing myeloid NK from lymphoid NK cells. Bonsai has no tunable parameters, integrates downstream exploratory analyses methods with Bonsai-scout and scales to large datasets, making it applicable to visualizing the structure in any set of high-dimensional objects.<\/p>\n
Single-cell omics methods promise to transform our understanding of the gene regulatory processes underlying cellular dynamics. Major efforts have measured the states of single cells using scRNA-seq for mRNA expression levels1 and scATAC-seq for chromatin accessibility2,3, including atlases of model organisms4,5,6.<\/p>\n
However, a consensus on how to analyze such data is yet to emerge; more than 1,750 scRNA-seq analysis tools were published in the last 8\u2009years7, almost as many as the total number of scRNA-seq papers (Supplementary Fig. 1)7,8. Analysis of these data not only faces the major technical challenges that single-cell omics data are very high-dimensional and sparse, in addition to having highly heterogeneous noise levels that vary over several orders of magnitude across measurements; the more fundamental problem is that we do not know what structure to expect in the data. For example, we do not know to what extent cells occur in discrete \u2018cell types\u2019 versus being distributed on a continuous manifold or what topology and dimensionality this manifold of cell states may have. Distributions of points in high-dimensional spaces can have extremely complex structures for which we lack intuition.<\/p>\n
There is, thus, an urgent need for methods that can explore the complex high-dimensional structure of these data. This requires methods that can visualize local and global relationships in the data without distortion and without hallucinating structure that does not exist. Unfortunately, currently popular methods such as t-distributed stochastic neighbor embedding (t-SNE)9 and uniform manifold approximation and projection (UMAP)10 are well known to fail in this regard11,12,13. In the absence of methods that reliably represent the structure in the data without distortion, the current practice in the field is to implement visualization methods in an almost trial-and-error manner14,15,16, tweaking their tunable parameters until the visualizations match prior biological knowledge or preconceived expectations. However, this practice hinders making truly novel observations or falsifying strongly held beliefs.<\/p>\n
Here we show that all these challenges can be overcome by representing single-cell omics data using tree structures. First, there are many examples where relationships between high-dimensional objects were efficiently captured by hierarchical representations (for example, the early taxonomies of organisms17). Trees are used ubiquitously to describe relationships between DNA sequences18 and hierarchical structures have also proven effective for representing high-dimensional objects in machine learning19,20.<\/p>\n
Second, as cells from a single organism are related through a lineage tree of cell divisions, gene expression patterns of single cells have in fact diverged along the branches of a tree structure. In addition, we show that, because of a feature of high-dimensional spaces that we call the \u2018blessing of dimensionality\u2019, distances between objects in high-dimensional spaces can generically be accurately represented along the branches of a tree. As trees can always be displayed in two dimensions, this allows for distortion-free visualizations of the relationships between cells.<\/p>\n
Although conceptually similar to hierarchical clustering21 and the phylogeny reconstruction problem22, developing a method for relating a set of single-cell gene expression states into a most likely tree structure faced several challenges. These included developing a probability model for movement through the high-dimensional gene expression space, using this to derive a likelihood of a given tree given the expression states of cells at its leaves and developing methods for rigorously accounting for the complex heterogeneous noise properties of scRNA-seq data. In addition, we had to develop algorithms that efficiently find the maximum-likelihood tree for datasets with many thousands of cells.<\/p>\n
We here present Bonsai, a Bayesian method that takes any set of objects with estimated coordinates in a high-dimensional continuous space, together with individual error bars on each estimated coordinate of each object, and reconstructs the most likely tree structure relating the objects at its leaves. Bonsai is derived from first principles with minimal assumptions and without any tunable parameters. Using extensive tests on real and simulated scRNA-seq data, we show that, in contrast to existing methods, Bonsai accurately represents the structure in the data on all scales. Rather than just clustering cells into different \u2018cell types\u2019, Bonsai not only infers the trajectories along which gene expression states have differentiated but also accurately preserves pairwise distances between all cells at all scales. Moreover, applying Bonsai to scRNA-seq data of cord-blood cells, it not only recovers the known differentiation hierarchy of blood cell types but also uncovers lineage relationships that, as far as we know, are previously undescribed.<\/p>\n
In addition to providing Bonsai as a standalone tool23, we implemented an automated pipeline for scRNA-seq analysis as a webserver (https:\/\/bonsai.unibas.ch) that starts from an mRNA count matrix, normalizes the data using Sanity24, identifies clusters of statistically indistinguishable cells using Cellstates25 and then runs Bonsai to reconstruct a tree relating the cells. Lastly, to facilitate exploring Bonsai\u2019s results, we also provide Bonsai-scout, an interactive app that can be used to view the tree using different layouts, to zoom in on specific parts, to define clusters of cells by finding maximally separated clades, to overlay gene expression on the tree and to find marker genes that most distinguish different clades of cells. Example results of our integrated scRNA-seq pipeline on cell atlases of human and mouse4,5 are provided as community resources.<\/p>\n
As it is generally impossible to represent the distances between objects in a high-dimensional space by distances in a two-dimensional image, currently popular visualization methods focus solely on preserving nearest-neighbor relationships at a particular chosen scale, for example, such that the k-nearest-neighbor cells in the high-dimensional space are also near in the two-dimensional image26. In contrast, we here aim to provide visualizations that capture the full structure in the high-dimensional data by reconstructing a tree with cells at its leaves such that, for each pair of cells, the path along the branches of the tree corresponds to the most likely trajectory through gene expression space connecting the cells. Equivalently, going out from the root of the tree to the cells at the leaves, the trajectories along the branches capture the most likely trajectories along which the gene expression states have diverged into the states at the leaves.<\/p>\n
To define the likelihood of a tree, we need to assign probabilities to the gene expression changes along its branches and this may seem like a poorly defined problem. It is often imagined that gene expression trajectories are constrained to some manifold embedded in the high-dimensional gene expression space but we know virtually nothing about the structure or dimensionality of this manifold. However, as any smooth manifold is locally isomorphic to flat Euclidean space, at least for nearby cells, their Euclidean distance in gene expression space may be a good approximation to their distance along the manifold. As we outline in the next section, if we assume that cells have been sampled sufficiently densely that, for each individual branch of the tree, the Euclidean distance in gene expression space is a good approximation of the distance along the unknown underlying manifold, this suffices to uniquely define the likelihood of any tree connecting the gene expression states of the cells at the leaves.<\/p>\n
Instead of defining the structure in the data in terms of the most likely trajectories through gene expression space connecting the cells, it is more common to define the quality of a visualization by how well it conserves the true distances between the cells. This raises the following question: To what extent does the tree of most likely trajectories also preserve cell-to-cell distances? Strikingly, we see below that, for high-dimensional data, the reconstructed trees in fact also accurately conserve cell-to-cell distances. We hypothesize that this is a consequence of the fact that, in high-dimensional spaces, virtually all directions are mutually orthogonal. If the underlying manifold along which gene expression states diverge is sufficiently high-dimensional, then the movements along the branches of the tree are all mutually orthogonal so that pairwise squared Euclidean distances sum along the tree\u2019s branches. Note that this is reminiscent of what happens in the evolution of DNA sequences. Even though the evolution of protein-coding sequences is constrained to the complex subspace of sequences that code for a functional protein, this space is sufficiently high-dimensional that Hamming distances between pairs of sequences generally accurately reflect the time since the sequences diverged from a common ancestor. More surprisingly, we see below that this preservation of distances even extends to high-dimensional datasets that were not in fact generated by a tree process.<\/p>\n
As we focus on the analysis of scRNA-seq data, we will in the following presentation assume that the objects of our dataset are cells and their positions in high-dimensional space are their gene expression states, even though our methods apply much more generally.<\/p>\n
For a given dataset D, consisting of measurements of the gene expression states of a collection of cells, we need to define a model that assigns a likelihood P(D\u2223T, t) to any tree, as defined by its topology T and branch lengths t. This likelihood should reflect the probability for the gene expression states of the cells at the leaves to have diverged along the branches of this tree. The optimal representation of the dataset then corresponds to the tree that maximizes this likelihood. The likelihood naturally decomposes into factors corresponding to probabilities for the movements along the branches of the tree and factors for the probabilities of the measurements given the positions of the cells at the leaves:<\/p>\n