Main

It has long been recognized that in biological systems ranging from the Drosophila embryo to the hepatic lobule, many genes need to be properly regulated in space for the system to function. To study the spatial patterns of gene expression, many different spatial transcriptomics methods, which produce spatially localized quantification of messenger RNA (mRNA) transcripts as proxies for gene expression, have been developed. Thanks to growing interest in the field, several reviews have been written in the past 5 years, providing overviews of experimental techniques for data collection1,2, and describing how such techniques can be applied to specific biological systems, for example tumors3, brain4, and liver5. These reviews typically begin with either laser capture microdissection (LCM) or single-molecule fluorescent in situ hybridization (smFISH) in the late 1990s, although the quest to profile the transcriptome in space is much older.

Unlike the previous reviews, this paper presents a database of literature dating back to 1987 comprehensively documenting the historical evolution and current development in data collection and analysis in spatial transcriptomics. In addition, we have analyzed the literature metadata from the database to show trends in the field. Key highlights from the database and analyses are presented in this paper, and more details are presented in our book-length supplement: https://pachterlab.github.io/LP_2021/. Section and figure numbers of the supplement in this paper refer to those in the DOI PDF version, while those in the online HTML version are subject to change, as it is continuously updated to reflect changes in the field. This database was curated by searching keywords such as “spatial transcriptomics” and “Visium” on PubMed and bioRxiv, and manually screening literature citing influential papers in the field. Literature metadata collected include the date published or posted and the institution of the first author. In addition, metadata for publications concerning new datasets include the species and tissue from which the data were collected, the experimental techniques used to collect the data, and the programming languages used to analyze the data. Metadata for publications concerning new data-analysis methods include the programming languages used in the implementation, the code repository of the implementation, and whether the code is packaged and documented. The database is continuously updated by manually screening RSS feeds from PubMed and bioRxiv for relevant keywords, or by submission via a Google Form.

Prequel era

By “spatial transcriptomics”, we mean attempts to quantify mRNA expression of large numbers of genes within the spatial context of tissues and cells. Some important technologies enabling spatial transcriptomics date back to the 1970s (Chapter 2 in Supplementary Information). Various forms of in situ hybridization (ISH) have been used for a long time to visualize gene expression in space. Radioactive ISH was first introduced in 1969, visualizing ribosomal RNA6 and DNA7 in Xenopus laevis oocytes, and was first used to visualize transcripts of specific genes (globin) in 1973 (ref. 8) (Fig. 1a). Non-radioactive fluorescent or colorimetric ISH was developed in the 1970s and the early 1980s, improving spatial resolution, enabling three-dimensional (3D) staining, and shortening required exposure times9,10 (Fig. 1a). Early ISH was performed in tissue sections, making it challenging to apply to blastulas and to reconstruct 3D tissue structures; whole-mount ISH (WM ISH) was first introduced in Drosophila in 1989 (ref. 11) and was soon adapted to other species, such as mice, in the early 1990s (ref. 12).

Fig. 1: Timelines of major events.
figure 1

a, Development of prequel era technologies. References: 1969 radioactive ISH of ribosomal RNA (rRNA)6,7, 1973 radioactive ISH of goblin mRNAs8, 1977 FISH of rRNA10, 1982 immunological FISH with biotin-labeled probe9, 1982 FISH of actin mRNA135, 1987 Drosophila enhancer trap14, 1989 WM ISH in Drosophila11, 1989 ES cell enhancer and gene trap in mice15, 1991 in situ reporter in Caenorhabditis elegans136. b, Major WM ISH atlases and gene expression pattern databases. References: 1994 scaling up WM ISH in C. elegans137, 1995 first mouse WM ISH138, 1998 AXelDB139, 1999 mouse: GXD140, 2000 Maboya Gene Expression patterns and Expression Sequence Tags (MAGEST)141, 2001 the Nematode Expression Pattern Database (NEXTDB)142, 2001 Ghost143, 2002 GenePaint144, 2002 D. melanogaster: Berkeley Drosophila Genome Project (BDGP)24, 2003 Medaka Expression Pattern Database (MEPD)145, 2003 Zebrafish Information Network (ZFIN)31, 2004 Gallus Expression In Situ Hybridization Analysis (GEISHA)25, 2005 miRNA atlas29, 2006 Allen26, 2006 BDTNP146, 2007 Fly-FISH147, 2007 Xenbase148, 2011 mouse Genitourinary Development Molecular Anatomy Project (GUDMAP)27, 2017 LungMAP28, 2020 ZEBrA149. c, Development of current-era technologies and their notable precursors, colored by type of technology. References: 1976 LCM17, 1988 ligase-mediated single -nucleotide variant (SNV) detection150, 1989 single-cell cDNA amplification151,152, 1989 FISH with combinatorial barcoding21, 1995 cDNA microarray153, 1996 commercial LCM18,19, 1998 smFISH23, 1999 LCM + microarray154, 2002 combinatorial FISH22, 2008 RNA-seq155, 2012 Tomo-array156 (the cDNA microarray predecessor of Tomo-seq), 2013 high-throughput RCA + ISS59, 2014 seqFISH20, 2015 MERFISH50, 2016 ST87, and 2019 GeoMX DSP44.

Another strand of development in early spatial transcriptomics was the enhancer and gene trap screen, which was developed in the 1980s when DNA sequencing throughput was increasing13 and metazoan genomes were newly opened frontiers. The first screens in Drosophila14 and mice15 were performed in the late 1980s in order to visualize expression of untargeted, and often previously unknown, genes. With increasing throughput, enhancer and gene traps became the technology of choice for spatial transcriptomics in the 1990s, until the rise of WM ISH in the late 1990s, which leveraged automation. WM ISH also avoided the need for transgenic lines, and was facilitated by the availability of reference genomes in the early 2000s for computational probe design. Although now eclipsed by newer methods, enhancer trap, gene trap, and in situ reporter methods have been used to build reference databases of gene expression and enhancer usage patterns in transgenic lines throughout the 2000s and 2010s16.

The foundation for many current-era technologies was built in the decades between the 1970s and the 2000s (Fig. 1c). For example, ultraviolet (UV) laser was first used to cut tissue in 1976 (ref. 17). Popular infrared (IR) and UV LCM systems were first reported in 1996 (refs. 18,19) and were soon commercialized. Some highly multiplexed smFISH technologies, such as sequential FISH (seqFISH)20, rely on combinatorial barcoding, that is encoding each gene with a combination of colors so transcripts of more genes with easily discernible colors (up to 5) can be quantified simultaneously. Combinatorial barcoding was first reported in immunological DNA FISH in 1989 (ref. 21) and was first used for transcripts in 2002 (ref. 22). The first unequivocal demonstration of smFISH showing each mRNA molecule as a spot was reported in 1998 (ref. 23). Highly multiplexed smFISH would not have been possible without the development of these technologies.

WM ISH was the technology of choice in the late 1990s and the 2000s, before the rise of highly multiplexed, high-resolution, and more quantitative technologies, and has been used to create gene expression atlases in embryos of several species such as Drosophila melanogaster24, Mus musculus, and Gallus gallus25; in various mouse organs such as the brain26, genitourinary tract27, and lung28; and for specific types of genes, such as microRNAs (miRNAs)29 (Fig. 1b). For miRNAs and many species other than mice and humans, the only spatial transcriptomics resources currently available are, for the most part, WM ISH atlases. Model-organism databases collecting proliferating gene expression patterns from various sources were also established in this period, such as Gene Expression Database (GXD)30 and Zebrafish Information Network31 (Fig. 1b). The golden age of WM ISH seems to have ended in the 2010s (Fig. 1b), perhaps due to some of the technology’s disadvantages, such as requiring stereotypical tissue structure, the need for thousands of animals to generate an atlas, and the largely qualitative nature of results.

Early motivating applications for spatial transcriptomics included identification of genes with restricted patterns that indicated function in development, identification of novel cell-type markers, and identification of novel cell types not evident from tissue morphology14,15. In the 1980s and 1990s, analyses were typically done manually, although more recently automated methods have been developed (Chapter 3 in Supplementary Information). Convergence of strands of technologies, including more powerful computing infrastructure, decreasing cost of sequencing, and the generation of more quantitative data, have mainstreamed and revolutionized spatial transcriptomics and opened up new possibilities. However, the legacy of the prequel era still lives on in the usage of prequel resources, such as referencing the Allen Brain Atlas (ABA)32 and the Allen Mouse Common Coordinate Framework33, and in institutions such as the Allen Brain Institute and the Jackson Laboratory, which are contributing to the current era of research34,35.

Data collection

Current-era technologies broadly fall into five categories in terms of how spatial information is acquired: region of interest (ROI) selection (Section 5.1 in Supplementary Information), smFISH (Section 5.2 in Supplementary Information), in situ sequencing (ISS) (Section 5.3 in Supplementary Information), next-generation sequencing (NGS) with spatial barcoding (Section 5.4 in Supplementary Information), and methods not requiring a priori spatial locations (Section 5.6 in Supplementary Information). Developers of such technologies often seek to enable a trifecta of transcriptome-wide profiling, single-cell resolution, and high gene-detection efficiency. Although this achievement appears to be increasingly within reach, current-era technologies are characterized by trade-offs between these goals.

ROI selection

Spatial locations can be obtained by selection and isolation of ROIs of known locations and shapes, which can be performed by physical (Section 5.1.3 in Supplementary Information) and optical marking of ROIs for isolation (Section 5.1.4 in Supplementary Information). The isolated ROIs can then be analyzed with complementary DNA (cDNA) microarray or RNA sequencing (RNA-seq), or dissociated into single cells for single-cell RNA-seq (scRNA-seq).

Physical microdissection includes LCM, 2000s voxelation36, and Tomo-seq37, which sections a tissue with a cryotome along an axis of interest, followed by RNA-seq on each section. Since 1999, by far the most widely used microdissection technology is LCM, which has been used in various biological fields, such as oncology, neuroscience, immunology, developmental biology, and botany (see Chapter 6 in the Supplementary Information for topic modeling of PubMed and bioRxiv LCM literature). In LCM, ROIs in the tissue section are dissected by either UV laser cutting (lasers manufactured by Zeiss and Leica) or fusion of tissue with a membrane by IR laser (manufactured by Arcturus, Fig. 2a); the two are combined in recent versions of Arcturus, in which IR fusion removes the ROI cut using UV. Combining LCM and Tomo-seq, the spatial transcriptome in 3D can be profiled as in geographical position sequencing (Geo-seq)38, albeit with limited spatial resolution. An innovative physical microdissection method is STRP-seq39, which slices adjacent tissue sections into stripes at different angles and reconstructs gene expression patterns in 3D with an algorithm inspired by ray-based computerized tomography. Manual dissection is commonly used to profile gene expression along one spatial axis of interest in plants40.

Fig. 2: Schematics of common current-era technologies.
figure 2

a, IR LCM. b, GeoMX DSP. The purple circle in step 2 is the UV-illuminated ROI. c, seqFISH barcoding and error correction scheme: if signal from one round of hybridization is missing, the remaining rounds can still uniquely identify the gene barcoded. d, MERFISH Hamming distance 4 barcoding and error correction scheme: from the design of the barcodes, if signal from one round of hybridization is missing, the correct barcode can be recovered. If two rounds are missing, the remaining signals are equidistant to two different barcodes so the original barcode cannot be recovered. e, Cartana ISS with cPAL sequencing: many copies of the gene barcode are made with RCA for signal amplification, which are then sequenced in situ with cPAL. The orange line stands for the RCA amplicon. Short blue lines stand for the gene barcode. Brown stands for the probe; bases not labeled are degenerate. Gray stands for primer matching constant region. f, NGS barcoding techniques. In Visium, the spots are arranged in a hexagonal grid, 100 µm apart center to center and 55 µm in diameter. In DBiT-seq, positional barcodes are deposited in microfluidic channels and spatial resolution is determined by the width (down to 10 µm) and spacing of the channels. In Slide-seq, barcoded beads 10 µm in diameter are spread in a single layer on a slide. In XYZeq, spatial barcodes are conferred on multiple cells in wells 500 µm in diameter, which are then dissociated for scRNA-seq. In Seq-Scope, the tissue is mounted on a repurposed Illumina flow cell with barcoded polony spots 0.6 µm apart on average. For Visium and Slide-seq, the lines represent oligonucleotides attached to the slide or bead. For DBiT-seq, red and green lines represent the flow in microfluidic channels carrying barcoding oligonucleotides. For Seq-Scope, the tissue (pink block) is mounted on repurposed Illumina flow cell with bridge-amplified polonies each with its own spatial barcode represented by different colors. For XYZeq, different colors of the cells represent different spatial barcodes in the microwells, and the cells are dissociated for scRNA-seq. t-SNE, t-distributed stochastic neighbor embedding. g, Data-analysis workflow: upstream analysis is technology specific, and includes image processing for smFISH and ISS-based technologies, and FASTQ file processing, quality control of the gene count matrix, and data normalization for NGS-based technologies. Non-spatial scRNA-seq data can be integrated by mapping cells to locations with landmark genes in the smFISH or ISS data or deconvolving cell types in Visium spots. Downstream analyses tend to be technology agnostic, and include finding spatially variable genes, transcriptionally defined spatial regions, and cell–cell interactions. Created with BioRender.com.

Optical marking of ROIs includes NICHE-seq41, which uses two-photon irradiation to mark ROIs in tissue from transgenic mice expressing photoactivatable green fluorescent protein (PA-GFP), and then uses fluorescence-activated cell sorting (FACS) to isolate cells with activated PA-GFP for scRNA-seq. Similar to NICHE-seq but without transgenic mice is spatially photoactivatable color encoded cellular address tags (SPACECAT)42, which stains cultured live cells or organoids with photocaged fluorophores and photoactivates ROIs for FACS and scRNA-seq. Also using photocaging, ZipSeq43 attaches anchor oligonucleotides with photocaged overhangs to tissue with antibodies or lipid insertion, and adds spatial ‘zipcodes’ to photoactivated ROIs hybridizing to the overhangs. A more popular commercial optical ROI-selection technique is the GeoMX Digital Spatial Profiler (DSP)44 and whole-transcriptome atlas (WTA)45 of Nanostring (Fig. 2b), which shines UV light on ROIs to release photo-cleavable gene barcodes for quantification with either nCounter or NGS. As GeoMX uses predefined gene panels rather than poly-A capture, Nanostring provides the Cancer Transcriptome Atlas (CTA) gene panel with over 1,800 genes, as well as human and mouse whole-transcriptome panels with over 18,000 genes.

Single-molecule FISH

Chronologically, the next technology developed in the current era is highly multiplexed single-molecule FISH (smFISH), which began with a 2012 prototype (seqFISH) that relied on super-resolution microscopy (SRM) to simultaneously profile 32 genes in yeast by hybridizing probes with different colors to transcripts, and then deducing the relative locations of the colors present46. SRM is no longer needed; in 2014, seqFISH20 was published, in which one color per gene is visualized per round of hybridization, and the probes are stripped before the next round for the next color in the barcode. All transcripts of the same gene have the same barcode. Four colors and 8 rounds of hybridization (48 = 65,536) are more than enough to encode all genes in the human or mouse genome. In practice, an error-correcting round of hybridization is performed, so that genes can still be distinguished if signal from one round of hybridization is missing47 (Fig. 2c). More recently, in a version of seqFISH based on RNA sequential probing of targets (RNA SPOTs)48, the ‘colors’ themselves are one-hot encoded by a sequence of hybridizations, expanding the palette to 20 ‘colors’ per channel and enabling the profiling of 10,000 genes49.

Another smFISH technique is multiplexed error-robust FISH (MERFISH)50, which uses a different barcoding strategy, in which each gene is encoded by a binary code. The color codes in each experiment must be separated by a Hamming distance (HD) of four to allow for correction of missing signal in one round, and by two to identify error without the facility for correcting it (Fig. 2d). The length of barcodes can be increased to encode 10,000 genes51. As only the fluorophores are removed but the probes are not stripped, numerous rounds of hybridization in MERFISH are less time consuming than those in seqFISH. Most other smFISH-based techniques, such as hybridization-based ISS (HybISS)52 and split-FISH53, use either seqFISH-like or MERFISH-like barcoding.

smFISH faces a number of challenges, which have been addressed by various methods. Signal-to-noise ratio can be improved with rolling circle amplification (RCA)52, branched DNA (bDNA)54, hybridization chain reaction (HCR)47, primer exchange reaction55, and tissue clearing56. With an increasing number of genes profiled, the transcript spots are increasingly likely to overlap, causing optical crowding. This can be mitigated by expansion microscopy (ExM)57, only imaging a subset of probes at a time and using computational super-resolution49, imaging highly expressed genes without combinatorial barcoding50, and computationally resolving overlapping spots58.

In situ sequencing

ISS methods yield spatial transcriptome information by sequencing, typically by ligation (SBL), gene barcodes (targeted), or short fragments of cDNAs (untargeted) in situ. Such methods rely on ligase joining only two pieces of DNA—a primer with known sequence and a probe—if they match the template, and non-matching probes are washed away. The probes used are degenerate except for one or two query bases encoded by a color. RCA is commonly used for signal amplification. The 2013 ISS59, later commercialized by Cartana, and barcoded oligonucleotides ligated on RNA amplified for multiplexed and parallel in situ analyses (BOLORAMIS)60 use one query base per probe, as in combinatorial probe anchor ligation (cPAL)61, to sequence gene barcodes (Fig. 2e). In cPAL, each probe queries one base in the gene barcode. Fluorescence ISS (FISSEQ)62 and a later adaptation with ExM, called ExSeq63, use SOLiD, which uses two query bases per probe to sequence circularized and RCA-amplified cDNAs. In spatially resolved transcript amplicon readout mapping (STARmap)56, gene barcodes are sequenced by sequencing with error-reduction by dynamic annealing and ligation (SEDAL), in which SOLiD-like two query bases are used to reject error, but one-base encoding can also be used. Barcode analysis by sequencing (BAR-seq) also RCA amplifies probes with gene barcodes, but uses sequencing by synthesis (SBS) instead of SBL to sequence the barcodes64.

NGS with spatial barcoding

Spatial locations of transcripts can also be preserved by capturing the transcripts from tissue sections on in situ arrays. Such arrays can be manufactured by printing spot barcodes, unique molecular identifiers (UMIs), and poly-T oligonucleotides on commercial microarray slides to capture polyadenylated transcripts, as in the spatial transcriptomics (ST) and Visium technologies (Fig. 2f). They can also be Drop-seq-like beads65 with split pool barcodes, UMIs, and poly-T oligonucleotides spread on slides in a single layer (for example, Slide-seq66) or confined in wells etched on the slides (for example, high-definition spatial transcriptomics (HDST)67), with bead barcodes subsequently located using in situ SBL. Alternatively, in deterministic barcoding in tissue for spatial omics sequencing (DBiT-seq)68, an array is generated by microfluidic channels, which are used to deposit one type of barcode in one direction and then another in a perpendicular direction, with the orthogonal barcodes ligated so each spot can be identified with a unique pairwise combination. While NGS barcoding techniques are typically designed for 3′-end Illumina sequencing, Visium has been adapted to Nanopore long-read sequencing69.

NGS barcoding techniques have been applied to large areas of tissue33, and their use is increasing (Fig. 4b). Nevertheless, they do not have single-cell spatial resolution. The commonly used Visium has spots in a hexagonal array 100 µm center to center and 55 μm in diameter (Fig. 2f). Bead diameter is 10 μm in Slide-seq and 2 μm in HDST (Fig. 2f). Slide-seq and HDST use bead sizes smaller than single cells, but they may not always provide single-cell resolution because one bead can span two or more cells. Resolution of DBiT-seq is determined by channel width (either 50, 25, or 10 μm, Fig. 2f). More recently, the spot size can be reduced to below 1 μm, with RCA-amplified DNA nanoballs as small as 0.22 μm across, with spot barcodes deposited in wells that are 0.5 or 0.715 μm apart in Stereo-seq70, and in Seq-Scope polymerase colonies (polonies) with spatial barcodes ~0.6 μm center to center on an Illumina flow cell that has been repurposed to capture transcripts from tissue sections71 (Fig. 2f). Another polony-based method, PIXEL-seq, achieves a spot diameter of about 1.22 μm, but unlike in the flow cell, polony (or DNA cluster)-indexed library-sequencing (PIXEL-seq) does not have much spacing around each polony72. Techniques such as XYZeq73 and sci-Space74 have been developed to dissociate the single cell or nuclei in spatially barcoded spots for scRNA-seq, so the data have single-cell transcriptomic, but not spatial, resolution (Fig. 2f).

De novo reconstruction of spatial information

Some technologies have been developed to preserve information necessary to computationally reconstruct spatial gene expression patterns without knowing or collecting spatial locations. One such technology is DNA microscopy75,76, which records proximity between cDNAs. This information can be used to reconstruct relative locations of transcripts. At the cellular level, gene expression in rare cell types can be reconstructed by deliberately assaying multiplets and then mapping them to locations in a spatial reference on the basis of gene expression of cells from common cell types attached to cells from the rare cell types77. Variants of the term “spatial transcriptomics” have also been used to describe techniques localizing transcripts to organelles (for example, APEX-seq78), although no spatial coordinates are recorded.

Multi-omics

The transcriptome is only one aspect of cell function. Other aspects, such as the proteome, neuronal connectome, and 3D chromatin conformation are also important to cell function, and some methods have been developed to profile them along with the transcriptome in the same cells (Section 5.8 of Supplementary Information). For the proteome, oligonucleotide-tagged antibodies are used to detect proteins of interest, and the oligonucleotide signifying the protein species can be detected with smFISH-based methods. Such antibody panels have been combined with transcriptomics, such as in DBiT-seq68, SM-Omics79, GeoMX DSP44, and MERFISH80. With the oligonucleotide barcode, over 100 antibodies can be used, such as when using all available antibody panels for GeoMX DSP. For 3D chromatin conformation, MERFISH and seqFISH+ have been adapted to visualize chromatin structure, by targeting DNA genomic loci81 or introns of nascent transcripts81,82. For the neuronal connectome, multiplexed transcript quantification can also be combined with neuron projection tracing. For instance, cholera toxin subunit b (CTb) retrograde tracing has been used in conjunction with MERFISH to visualize axons83. Also, BAR-seq was originally designed to use ISS for axon tracing by sequencing neuron-specific barcodes introduced by a virus injected into the brain, but was later adapted to sequence gene barcodes64 as well. In addition, while not an -ome per se, electrophysiology has been recorded prior to transcriptome profiling in the same cells, such as with a patch–clamp in explanted human neurons, followed by HCR–smFISH84, and with extracellular electrodes in cultured cardiomyocytes, followed by STARmap in electro-seq85.

Comparison across categories

In this section, we discuss trade-offs, among high detection efficiency, transcriptome-wide profiling, high spatial resolution, and sometimes larger tissue area, made by different types of technologies, as well as practical factors relevant to selection of technology, such as FFPE compatibility and cost/usability.

Detection efficiency

Detection efficiency is commonly estimated by performing non-barcoded smFISH with near 100% sensitivity for select marker genes on the same cell type and comparing the average number of transcripts detected for each gene per cell for techniques in which cells can be segmented, or per unit tissue area for techniques without single-cell resolution. For NGS-based techniques with UMI, sometimes the number of UMIs and genes detected per cell or unit area is compared with that of other techniques with UMI. Note that comparisons of efficiencies are confounded by different tissues and methods used to estimate efficiencies in different studies and by different sequencing depths in NGS.

Highly multiplexed smFISH techniques tend to excel in this area, with ~95% for Hamming distance 4 MERFISH86 compared with non-barcoded smFISH; multiple rounds of hybridization tend to decrease the efficiency, in part because barcodes with incorrigible errors are discarded. NGS barcoding techniques tend to have lower efficiency. For select genes in the same tissue type, ST detected around 6.9% as many UMIs as transcript spots detected by non-multiplexed smFISH per unit area analyzed87, comparable to the detection efficiency of scRNA-seq per cell analyzed. Visium’s efficiency seems to be moderately higher than that of ST, and DBiT-seq’s is even higher, at ~15.5% per area compared with smFISH68. Efficiencies of the submicrometer techniques, in the number of UMIs per unit area in the same tissue, might be comparable to that of Visium72. ISS tends to be less efficient, in part because of inefficiency of reverse transcription (RT) and SBL. Whereas the detection efficiency of scRNA-seq techniques is between 3% and 25% (refs. 65,88,89,90,91), the detection efficiencies of Cartana ISS and FISSEQ92 are ~5% and ~0.005% respectively, with STARmap being only marginally better than scRNA-seq. However, compared with smFISH, ExSeq claims up to 62% efficiency per cell for genes tested63. Newer technologies tend to skip RT and make ligation of the padlock probe on an RNA template more efficient, such as in BOLORAMIS and hybridization-based RNA ISS (HybRISS)93, or to substitute SBL with seqFISH-like barcoding, as in HybISS, to improve detection efficiency.

Transcriptome-wide profiling

Techniques not targeting specific genes with a panel of known probes are transcriptome wide, such as ROI selection followed by NGS, and NGS barcoding, where NGS is performed on poly-A captured transcripts, as well as untargeted ISS, such as FISSEQ and untargeted ExSeq. However, these transcriptome-wide techniques tend to have lower detection efficiency. It is possible to use certain techniques that require gene probe panels to quantify transcripts of over 10,000 genes, such as seqFISH+, MERFISH, and GeoMX WTA, although unlike in NGS, novel transcripts not targeted by the probes cannot be detected. While GeoMX WTA has been used in some studies outside Nanostring, where GeoMX originated94, the number of overall genes profiled with smFISH-based techniques per dataset has not increased over time (Fig. 3g). Instead, in studies using smFISH- and ISS-based techniques, a smaller number of genes is profiled, and the smFISH or ISS dataset is complementary to a transcriptome-wide scRNA-seq dataset95.

Fig. 3: Current-era metadata.
figure 3

a, Number of institutions that have published papers or preprints with each technique, excluding LCM literature too vast to be manually curated. Only techniques used by at least three institutions are shown. b, Number of publications for each healthy organ in humans (male shown here, as there is no study on healthy female-specific organs in humans at present). c, Number of publications for pathological organs in humans (female shown here, but there are two studies on prostate cancer). d, Number of publications per species. e, Number of publications per healthy organ in mice. f, Number of publications for pathological organs in mice. g, Number of genes per dataset over time. Gray ribbon in g and h stands for 95% confidence interval. The slope is not significantly different from 0 in g (t test). In g and h, the y axis is log-transformed. h, Total number of cells per study profiled by smFISH-based techniques over time among studies that reported the number of cells. IceFISH, intron chromosomal expression FISH; C-FISH, consecutive FISH; MOSAICA, multi-omic single-scan assay with integrated combinatorial analysis; SGA, spatial genomic analysis; corrFISH, correlation FISH; EASI-FISH, expansion-assisted iterative FISH; par-seqFISH, parallel seqFISH; CISI, composite in situ imaging; SCRINSHOT, single-cell-resolution in situ hybridization on tissue; coppaFISH, combinatorial padlock-probe-amplified FISH. i, Number of publications for data collection using each of the fice most popular programming languages for downstream data analysis. j, Number of publications for data analysis using each of the five most popular programming languages for package development. In both i and j, each icon stands for 20 publications. Note that multiple programming languages can be used in one publication.

The number of genes that can be detected by highly multiplexed smFISH is limited by optical crowding, and expansion microscopy was used to address this issue in MERFISH and ExSeq. However, expansion reduces the amount of tissue covered per field of view, thus limiting imaging throughput.

Spatial resolution

smFISH- and ISS-based techniques have single-cell and single-molecule resolution, although cell segmentation can be challenging. In addition, smFISH- and ISS-based techniques can be applied to cleared thick tissue sections80, although the number of genes profiled in this case is much smaller than in most two-dimension (2D) highly multiplexed smFISH studies. All other types of techniques require tissue sections and are thus limited to 2D, or 3D with z resolution limited to section thickness, which is usually at least 10 μm for frozen sections. Although there are submicrometer-resolution NGS barcoding techniques, and the ROIs of LCM and GeoMX can in principle be single-cell or smaller, these types of techniques, as they are most commonly used, tend to have lower spatial resolution, such as 55 μm for Visium and several hundred micrometers across for GeoMX (for example 700 × 800 μm in ref. 94), owing to insufficient sensitivity of transcript detection at single-cell or subcellular resolution96.

Tissue area

Overall, techniques with lower detection efficiencies tend to be better at profiling larger tissue areas, and for smFISH, there seems to be a trade-off between the number of cells and the number of genes. In current-era spatial transcriptomics, a tissue section several millimeters across, such as a substantial portion of a mouse brain coronal section, which can fit into a Visium or ST tissue capture area, is considered large, and increasing tissue area and sequencing depth for sensitivity would increase sequencing cost. Cartana ISS and HybISS have also been used to profile large areas of tissue several millimeters across, but only around 100 genes97. An advantage of HybISS here is strong RCA signal and less optical crowding, thanks to lower detection efficiency facilitating lower magnification (×20; MERFISH uses ×60) and thus faster imaging. While most highly multiplexed smFISH datasets remain at hundreds of genes (Fig. 3g), among studies that reported the number of cells, the total number of cells per study has increased (Fig. 3h, P < 0.001, two-sided t-test). ROI-selection techniques are generally used for small numbers of ROIs, as it’s labor intensive to select very large numbers of ROIs and process them separately without spatial barcoding. However, when high spatial resolution is not as crucial or practical, ROIs with very low resolution can be selected to cover more tissue, as in the LCM dataset in the Allen Human Brain Atlas98.

Usability

While most techniques were originally developed for frozen sections, some are compatible with FFPE, which, as this is a common tissue archive, may at times be the only type of tissue available. Among smFISH-based techniques, ACD’s RNAscope99 is FFPE compatible but can profile only 12 genes at a time in FFPE, compared with 48 in frozen sections. Among NGS barcoding techniques, Visium100 and DBiT-seq101 are FFPE compatible, but owing to crosslinking and RNA fragmentation in archival storage, detection efficiency as number of UMIs and genes detected per spot in FFPE tissues is about 5 to 10 times lower than in their frozen counterparts. LCM has long been applied to FFPE tissues, even at single-cell resolution with the sensitive SMART-3Seq102. GeoMX is not only FFPE compatible, but also predominantly used on pathological human FFPE tissues (Figure 5.8 in Supplementary Information).

While many new techniques have been developed, most never spread beyond their institutions of origin (Figure 4.9 in Supplementary Information). Among those that have spread far, the most popular tend to have commercial platforms, such as LCM, 10X Visium and its precursor ST, Cartana ISS (acquired by 10X), and Nanostring GeoMX (Fig. 3a). In addition, many major institutions have core facilities for NGS, if not LCM, Visium, and GeoMX (for example, the TPCL at the University of California, Los Angeles, and the Advanced Genomics Core at University of Michigan, Ann Arbor), reducing the cost of purchasing new equipment and training personnel in individual laboratories. Tomo-seq has also spread, perhaps because of its ease of implementation with standard equipment. In contrast, smFISH-based techniques have not spread as widely thus far, perhaps due to the complicated home-built fluidic system, long imaging time, terabytes of images, and expensive probes. However, some smFISH techniques are being commercialized with automated imaging and fluidic platforms, such as MERFISH, commercialized as MERSCOPE by Vizgen, and another smFISH-based technique, Resolve Biosciences’s Molecular Cartography platform. In addition, Rebus Esper can be programmed to automate different smFISH technologies and can process images online as in Illumina sequencing, and has been used to automate cyclic-ouroboros smFISH (osmFISH)103. With the new automated commercial platforms, the popularity of smFISH-based techniques might rise, especially if such platforms are adopted by core facilities.

Data analysis

The processing and analysis of high-throughput spatial transcriptomics data requires new methods and tools, especially for problems such as image preprocessing, spatial reconstruction of scRNA-seq data, cell-type deconvolution of NGS barcoding data, identification of spatially variable genes, and inference of cell–cell interactions (Fig. 2g).

Upstream

Upstream data analysis converts raw data into forms more amenable to biological interpretation and is dependent on the data-collection technology.

For smFISH- and ISS-based methods, the raw data consist of images of fluorescent spots, which must be processed to identify transcript spots, match spots to genes, and assign spots to cells (Section 7.1 of Section Information). smFISH and ISS studies often use classical image-processing tools, such as top-hat filtering, to remove background, translation to align images from different rounds of hybridization, and watershed for cell segmentation47,56,86. Machine learning in Ilastik, deep learning packages like DeepCell104, and alternative tools incorporating scRNA-seq data105 can also be used for cell segmentation. However, without visualizing the plasma membrane, accuracy of cell segmentation is limited. Some analyses, such as identification of tissue regions, can be performed without cell segmentation105. Until 2019, image processing was typically performed with poorly documented and technique-specific code written in the proprietary language MATLAB, but more recently, such code is increasingly written in the open-source language Python. The package starfish106 was developed as an attempt to provide a unified and well-documented user interface to process images from different techniques, such as seqFISH, MERFISH, and ISS, but it has not been widely adopted.

Improvements in scRNA-seq technology have inspired new methods for leveraging the complementary nature of high-resolution transcriptome quantification with spatial transcriptomics data. For smFISH and ISS data that are not transcriptome wide, expression patterns of genes not profiled in the spatial data can be imputed with scRNA-seq data, either by mapping dissociated scRNA-seq cells to the spatial reference or by directly imputing gene expression in space using expression profiles from scRNA-seq (Section 7.3 in Supplementary Information). Cells can be mapped to spatial locations on an existing spatial dataset with genes shared by the two datasets, with an ad hoc score favoring similarity between cell and location107 or via optimal transport modeling108. While ad hoc scoring is simple to implement, the results tend to be qualitative. Gene expression in space can also be imputed from scRNA-seq without explicitly mapping scRNA-seq cells to locations. A common approach is to project the spatial and scRNA-seq data into a shared low-dimensional and batch-free latent space, and to subsequently estimate gene expression by projecting the spatial cells into the latent space. Examples of this approach include Seurat3 (ref. 32) and gimVI109. These methods may also be used to add spatial context to single-cell multi-omics data when spatial techniques for some of the multi-omics data are not available.

In spatial data that are not single-cell resolution, such as those derived from ST and Visium, scRNA-seq data can inform cell-type composition of the spots or voxels (Section 7.4 of Supplementary Information). Negative binomial models and non-negative least squares (NNLS) are common principles underlying cell-type deconvolution methods. Negative binomial models are typically parameterized with rate and dispersion, and the rate is modeled as a weighted sum of cell-type signatures from scRNA-seq, with scaling factors for library size and technology sensitivity; the non-negative weights may be normalized to sum up to 1 as cell-type proportions per spot. Negative-binomial-based methods include stereoscope110 and cell2location111. Simpler than negative binomial, gene expression is modeled as Poisson instead in RCTD112. Cell-type deconvolution can also be performed by modeling gene expression at each spot as a weighted sum of cell-type signatures outside the rate parameter of negative binomial distributions, and the weights are inferred with NNLS. For example, AdRoit113 uses the means of negative binomial distributions fitted to spot gene expression and to scRNA-seq cell-type signatures. The cell-type signatures can be non-negative matrix factorization (NMF) cell factors from scRNA-seq assigned to cell types, as in NMFreg66 and SPOTlight114. The cell-type weights can be regularized or thresholded to limit the number of cell types assigned to each spot. Parallels can also be drawn between cell-type deconvolution and topic modeling in text mining; cell types are analogous to topics, and genes are analogous to words. Latent Dirichlet allocation (LDA) from topic modeling has been adapted to cell-type deconvolution, such as in spatial transcriptomics deconvolution by topic modeling (STRIDE)115 and STdeconvolve116; the latter is unsupervised and does not require a scRNA-seq reference.

Downstream

Downstream analyses most often apply to the gene count matrix and cell or spot locations, and are thus largely independent of data-collection technologies.

Given the relevance of scRNA-seq to spatial data, and how spatial data are often analyzed like scRNA-seq data in exploratory data analysis (EDA), popular scRNA-seq EDA ecosystems, such as Seurat32, SCANPY (which spatical single-cell analysis in Python (Squidpy) is built on)117, and SingleCellExperiment (extended by SpatialExperiment)118, have added functionalities for spatial data, such as updates to data containers and functions to facilitate visualization of gene expression and cell or spot metadata at spatial locations (Section 7.2 of Supplementary Information). EDA packages dedicated to spatial data with beautiful graphics and good documentation have also been written, such as Giotto119 and STUtility120. Seurat and Giotto also implement basic methods to identify spatially variable genes. In addition, Giotto implements methods to identify cell-type enrichment in ST and Visium spots, to identify gene coexpression and association between gene expression and cell-type colocalization, and to identify spatial regions121.

Spatially variable genes are genes whose expression is associated with spatial location (Section 7.5 of Supplementary Information). Three approaches are commonly used for these genes: Gaussian process regression (GPR)122 and its generalization to Poisson123 and NB124, Laplacian score125, and Moran’s I. GPR-based methods model normalize gene expression or the rate parameter of Poisson or NB gene expression as a GPR and find whether the model better describes the data with the spatial term than without. Laplacian-score-based methods identify genes whose expression better reflects the structure of a spatial neighborhood graph. The locations of cells can also be modeled as a spatial point process with gene expression as marks; spatially variable genes can be identified as marks associated with locations126. Fitting GPR models to numerous genes can be time consuming, especially when a Bayesian approach with Markov chain Monte Carlo is used. Permutation testing used in Laplacian-score-based methods can also be time consuming. As both GPR- and Laplacian-score-based methods seek to identify spatial autocorrelation, sometimes the classic spatial autocorrelation metric Moran’s I is directly used to identify spatially variable genes, as in Seurat v3 and above. MERINGUE127 uses a local version of Moran’s I. Moran’s I and its significance testing are implemented in established geospatial packages and are easy and fast to run, but may have less statistical power than model-based methods123.

Spatial information also enables identification of potential cell–cell interactions (Supplementary Section 7.8). This is commonly done with knowledge of ligand–receptor (L–R) pairs, and can test which L–R pairs are more likely to be expressed in neighboring cells or spots128 or whether two cell types each expressing the ligand and the receptor are more likely to colocalize127. The cross-type L function from a spatial point process can be used to find cell types that colocalize129. Expression of genes of interest can also be modeled, including a term for cell–cell colocalization; a gene is considered associated with cell–cell colocalization if the model better describes the data with this term than without130.

There are many other types of downstream analysis that are useful for spatial transcriptomics analysis, including identification of archetypal gene patterns (Section 7.6 of Supplementary Information), spatial regions defined by the transcriptome (Section 7.7 of Supplementary Information), inferring gene–gene interactions (Section 7.9 of Supplementary Information), subcellular transcript localization (Section 7.10 of Supplementary Information), and gene expression imputation from H&E images (Section 7.11 of Supplementary Information).

Trends in the spatial transcriptomics field

The quality versus quantity trade-off inherent in existing technologies means that there is no single “best” solution currently available, and the difficulty in implementing methods has resulted in many technologies never spreading beyond their institutions of origin. LCM, Visium, ST, GeoMX DSP, and Tomo-seq have been the most widely adopted (Fig. 3a), and in almost all cases in the United States and western Europe (Figures 4.12, 5.27, and 5.33 in Supplementary Information). In terms of tissues analyzed, multiplexed current-era techniques have been used widely to characterize human tissues131, tumors87 (especially breast tumors), and pathological tissues that don’t necessarily have a stereotypical structure132 (Fig. 3b,c). In the SARS-CoV-2 pandemic, GeoMX DSP has been used for spatial transcriptomic profiling in lung autopsies of people who died due to COVID-19 (ref. 94).

Some of the processed data, and associated spatially variable genes, can be downloaded and visualized from SpatialDB133. Excluding LCM literature too vast to manually curate, the vast majority of current-era studies were performed in either humans or mice (Fig. 3d), and the brain is the most studied healthy organ while the lungs (particularly due to COVID-19) and breast tumors are also often studied in humans (Fig. 3b,e,f). In particular, the international project Brain Research through Advancing Innovative Neurotechnologies (BRAIN) Initiative—Cell Census Network (BICCN) is constructing a multi-modal atlas for human, mouse, and non-human primate brains, including spatial data such as those from MERFISH and seqFISH34.

All packages mentioned in the “Data analysis” section are open source and written in languages such as R, Python, and Julia. Downstream analyses in studies primarily concerning new data anddata-analysis packages predominantly use open-source programming languages, such as R, Python, and C++ (Fig. 3i,j). While MATLAB is still popular, its use has not risen, as was the case for R and Python (Figure 7.12 in Supplementary Information). While R is more popular for downstream analyses and EDA, Python, and C++ are more popular for package development (Fig. 3i,j). Most of the packages are not hosted on standard repositories, such as the Comprehensive R Archive Network (CRAN), Bioconductor, pip, or conda (Figure 7.13 in Supplementary Information). While most packages using R, Python, and C++ are well-documented, many MATLAB packages are not (Figure 7.12 in Supplementary Information). The standard repositories and documentation make packages more usable; this is discussed in more detail in Section 7.12 of the Supplementary Information.

Future perspective

While technologies of the past are rapidly depreciating, the ideas and methods that underlie them are fundamental to current-era spatial transcriptomics. The field has dramatically expanded over the past 5 years (Fig. 4a), with a plethora of new techniques and the popularization of Visium driving growth (Fig. 4b and Figures 4.9, 5.38, and 8.1 of Supplementary Information).

Fig. 4: Growth of the current era.
figure 4

a, Number of publications over time for current-era data collection and data analysis. Bin width is 120 days; the curves drop because the plot was made at the beginning of a new bin. Non-curated LCM literature is excluded. b, The data collection curve in a, broken down by category of techniques. The colors are stacked and sorted in descending order of total number of publications using techniques in that category.

What lies ahead of the rising curves (Fig. 4)? First, more can be done to improve data-collection techniques. For example, most current-era techniques require tissue sections. Highly multiplexed whole-mount smFISH and tissue clearing protocols, and more efficient computational tools that will align multiple sections that may come from multiple individuals or even developmental stages, should be developed to extend current-era techniques to 3D and to spatiotemporal analysis. Future techniques may also extend the current era from the scale of millimeters to centimeters and across other modalities, such as epigenomics and metabolomics, to give a fuller picture of cellular function. Furthermore, smFISH and ISS techniques, with signal amplification to reduce the number of probes per transcript, can be adapted to target isoform-specific exons or untranslated regions, rather than all transcripts of a gene.

Second, current-era data have not yet been integrated into comprehensive databases. Prequel databases, such as GXD and e-Mouse Atlas and Gene Expression (EMAGE)134, include data from multiple sources and can be queried by gene symbol and developmental and spatial ontologies. In addition, ABA26 and EMAGE aligned ISH images to common coordinates and can be queried with expression patterns. While some current-era authors provide online interactive visualization of datasets from their studies33, comprehensive databases integrating, querying, and visualizing data from multiple sources, as in the prequel era, have not yet been developed. Furthermore, while prequel ontologies are still used in current-era studies, such ontologies may be improved with the transcriptome-wide quantitative data from the current era.

Third, outside of LCM, the current era is highly focused on humans and mice, while potential spatial transcriptomics investigations of other species, such as plants and invertebrates, lag behind. Technological modernization of prequel consortia for organisms other than humans and mice holds much promise for the development of useful spatial transcriptomics atlases.

Fourth, an open-source, well-documented, interoperable, and scalable workflow with an integrated, easy-to-use interface would greatly simplify spatial transcriptomics data collection and analysis. At present, for tasks beyond EDA, users still often need to learn new syntax, convert object types, and even learn new languages to use some data-analysis tools. Finally, our survey of methods shows that spatial transcriptomics methods need to be more open and accessible so that they become adopted around the world and are not restricted to elite Western institutions.