ARTICLES
https://doi.org/10.1038/s41587-022-01311-4
nature
biotechnology
H Check for updates
OPEN
Scalable single-cell RNA sequencing from full transcripts with Smart-seq3xpress
Michael Hagemann-Jensen©1'2, Christoph Ziegenhain©12 and Rickard Sandberg©1H
Current single-cell RNA sequencing (scRNA-seq) methods with high cellular throughputs sacrifice full-transcript coverage and often sensitivity. Here we describe Smart-seq3xpress, which miniaturizes and streamlines the Smart-seq3 protocol to substantially reduce reagent use and increase cellular throughput. Smart-seq3xpress analysis of peripheral blood mononuclear cells resulted in a granular atlas complete with common and rare cell types. Compared with droplet-based single-cell RNA sequencing that sequences RNA ends, the additional full-transcript coverage revealed cell-type-associated isoform variation.
Many single-cell RNA sequencing (scRNA-seq) methods have been described and their respective strengths and weaknesses characterized1. However, researchers still face a compromise between methods with high cell throughput (that is, droplet or combinatorial indexing methods) but low transcript coverage and methods with high sensitivity and full-length transcript coverage (that is, plate-based methods). Based on Smart-seq3 (ref.2)—the method currently offering the highest information content per profiled cell—we systematically evaluated the feasibility of reducing volumes, reagents and experimental steps, without sacrificing data quality. This resulted in Smart-seq3xpress, a scalable nanoliter implementation of Smart-seq3 with throughput limited only by the available equipment (for example, polymerase chain reaction (PCR) instruments), and sequencing-ready libraries can be generated in a single workday.
We hypothesized that the Smart-seq3 chemistry would work in substantially lower volumes when covered in an inert hydrophobic substance ('overlay') (Fig. la). Using accurate non-contact nanoliter dispensers, we scaled the reaction volumes of the lysis, reverse transcription (RT) and pre-amplification PCR steps down to 1:2,1:5 and 1:10 of the established volumes and tested these conditions on K562 and HEK293FT cells (Fig. lb). For the lowest volumes, cells were sorted by fluorescence-activated cell sorting (FACS) into 300 nl of lysis buffer covered with 3 ul of Vapor-Lock, with subsequent additions of 100 nl for RT and 600 nl for PCR. After shallow sequencing, alignment and error correction of reads, we observed similar numbers of detected genes and molecules per cell at a sequencing depth of 100,000 reads per cell (Fig. lc and Extended Data Fig. la,b), confirming that reaction scaling is possible without compromising data quality or introducing unwanted variability. Further reduction of reaction volumes beyond 1:10 was also possible, although not further pursued as savings in reagents are diminishing and reactions may become vulnerable to variations in cell sorting fluid (~5nl). We hypothesized that the overlays would both protect the low reaction volumes from evaporation and provide a 'landing cushion for the FACS-sorted cells. Indeed, many overlays with varying chemical properties could be used with low-volume Smart-seq3 (Fig. Id), including silicone oils with high viscosities and hydrocarbons with higher freeze points. As expected, overlays did not interfere with the cDNA synthesis reaction when tested in larger volumes (Extended Data Fig. lc,d). Performing low-volume Smart-seq3 without an
overlay resulted in drastic losses (Fig. Id), explaining why earlier efforts to miniaturize plate-based scRNA-seq without overlay resulted in substantially decreased complexity3-5. Next, we investigated whether the time-consuming and plastics-consuming cDNA clean-up step could be omitted by instead diluting the cDNA now obtained in lower volumes. At equal sequencing depths, single-cell libraries generated with cDNA dilution and clean-up had no detectable differences (P = 0.89 and P = 0.71, respectively, f-test; Fig. le and Extended Data Fig. le,f). The lowered reaction volumes enabled better control of downstream reaction conditions, and we, therefore, investigated the relationship between RT and pre-amplification PCR volumes. Dilution of the RT reaction was found beneficial for optimal PCR performance, with the established 1.5X ratio of PCR to RT volume (Extended Data Fig. lg). However, PCR extension times could be reduced from 6 minutes to 4 minutes without complexity losses for longer transcripts (Extended Data Fig. lh-j).
A major cost in plate-based scRNA-seq is tagmentation, which needs to be performed individually on each cell, and, although reducing the amounts of commercial Tn5 has been suggested to cut reaction costs6'7, it is currently unclear to what degree Tn5 can be reduced without losing library complexity. We investigated how the cDNA input amounts influence library complexity, which revealed that 20-fold variations in cDNA amounts could be tolerated with minor effects on gene and unique fragment detection (Fig. lf,g). Library complexity was mainly a product of absolute Tn5 and cDNA amounts with little effect from varying reaction volumes (Extended Data Fig. 2a,b). The slightly reduced complexity at the highest cDNA amounts (Fig. If) likely resulted from insufficient tagmentation due to the limited Tn5 amounts used. Similar results were observed when varying the Tn5 amounts on a fixed amount of cDNA input (Extended Data Fig. 2c,d), and high-quality libraries were obtained using both commercial amplicon tagmentation mix (ATM) Tn5 and in-house8 Tn5 (Extended Data Fig. 2e). Altogether, these results demonstrated that tagmentation reactions were robust over large variations in cDNA input and Tn5 amounts and that considerable savings in resources can be made in reducing Tn5 amounts with only minimal effect on library complexity (Extended Data Fig. 2d).
Having demonstrated robust tagmentation over large ranges of cDNA, Tn5 and reaction volumes, we realized an opportunity to exclude several time-consuming and resource-intensive experimental steps, including excessive cDNA pre-amplification,
department of Cell and Molecular Biology, Karolínska Institutet, Stockholm, Sweden. 2These authors contributed equally: Michael Hagemann-Jensen, Christoph Ziegenhain. He-mail: rickard.sandberg@ki.se
1452
NATURE BIOTECHNOLOGY | VOL 40 | OCTOBER 2022 11452-1457 | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
concentration measurements, fragment length quality control (QC) traces and cDNA amount normalization (Fig. lh). Instead, cDNA products from fewer pre-amplification cycles could be directly tagmented without the above-mentioned steps—this strategy we termed Smart-seq3xpress. To explore this potential, we first generated libraries with low-volume RT (300 nl lysis volume) from HEK293FT cells and using a range of PCR pre-amplification cycles (10-20), which revealed very similar gene detection (Fig. li) without any need for additional enzymatic reaction clean-ups (Extended Data Fig. 2f,g). However, the resulting libraries were heavily biased toward 5' reads that contain the unique molecular identifier (UMI) at the expense of the internal reads important for full-transcript coverage scRNA-seq2. This resulted from inefficient tagmentation and the inability to modulate the ratio of UMTcontaining and internal reads by Tn5 amounts (Fig. lj). Specifically, the high salt concentration in the KAPA PCR buffer likely resulted in template blocking during tagmentation9, whereas other PCR polymerases were compatible with direct tagmentation (Extended Data Fig. 3a,b). Notably, tagmentation of SeqAmp and Platinum II amplified cDNA had considerably lowered fraction of 5' UMI reads, indicative of improved tagmentation compatibility (Extended Data Fig. 3c). Additional experiments identified SeqAmp as the preferred PCR polymerase, as it consistently performed well with improved gene and molecule detection over a range of template-switching oligonucleotide (TSO) and PCR primer amounts (Extended Data Fig. 3a-d).
To directly assess the accuracy of RNA counting in Smart-seq3xpress, we used the newly developed molecular spikes10 (UMIcountR R package). Whereas Platinum II had higher error rates, which resulted in inflated RNA counts (Extended Data Fig. 3e,f), SeqAmp with 2uJVI TSO and 1 uJVI of each PCR primer had high sensitivity and accuracy, and the fraction 5' UMI reads could be modulated as expected by varying the cDNA or Tn5 amounts (Fig. Ik).
Next, we realized that the Smart-seq3 TSO could mis-prime during RT to induce a strand invasion artifact1112 (Extended Data Fig. 4a,b). We screened TSO sequences of various designs and higher RT temperatures (Supplementary Note 1) for their effect on mis-priming, accuracy in UMTbased counting, gene and RNA molecule sensitivity and general quality (Fig. ll,m and Extended Data Fig. 4c,d). Although no TSO showed excellent performance
across all metrics, new TSO sequences with improved performance compared to Smart-seq3 and other TSOs12 were identified. Next, we validated candidate TSOs in primary cells to confirm sensitivity (Extended Data Fig. 5a) and fine-tuned oligonucleotide primer concentrations for optimal performance (Extended Data Fig. 5b). Thus, Smart-seq3xpress with the improved TSO has substantially reduced strand invasion in both peripheral blood mononuclear cells (PBMCs) (Extended Data Fig. 5c) and HEK293FT cells (Extended Data Fig. 5d,e) while maintaining accurate RNA counting (Fig. ll,m and Extended Data Fig. 4d).
Finally, we benchmarked both low-volume Smart-seq3 (lul, KAPA) and Smart-seq3xpress (12 PCR cycles, using either KAPA or SeqAmp with new and old TSO) against standard-volume Smart-seq3 (ref. 2). Notably, we observed an improved gene and molecule detection with SeqAmp-based Smart-seq3xpress (Fig. In and Extended Data Fig. 5f-h). Crucially, the material and resources needed to construct Smart-seq3xpress singe-cell libraries were ten-fold reduced, allowing researchers to substantially increase the cell numbers analyzed. Further streamlining and reduction of plastics consumables was achieved by collecting final libraries by centrifugation using a simple 3D-printed adapter (Extended Data Fig. 6a) and through contact-less combinatorial index dispensing or relying on tagmentation plates containing already dispensed desiccated index primers. Thus, the miniaturization and streamlining of Smart-seq3 was feasible at improved gene and molecular detection so that sequencing-ready Smart-seq3xpress libraries can be reached within a single workday (Extended Data Fig. 6b).
To showcase Smart-seq3xpress, we profiled 26,260 human peripheral blood mononuclear cells (hPBMCs) at an average depth of 258,000 read pairs per cell. Sequence data were demultiplexed, processed and quality controlled with zUMIs13, and, after stringent QC, Seurat14 was used for downstream analysis (Methods). The single-cell transcriptomes were visualized with uniform manifold approximation and projection (UMAP) and separated into 27 clusters (Fig. 2a) that were supported by all donors (Extended Data Fig. 7a) and protein staining from the index sorting (Extended Data Fig. 7b). Reconstruction of T cell receptor (TCR) sequences from the internal reads matched the identified T cell clusters (Fig. 2b), and the bulk of reconstructed TCR sequences corresponded to expected single-chain pairing (Fig. 2c). Consistently higher gene
->■
Fig. 1 | Scalable full-transcript coverage scRNA-seq with Smart-seq3xpress. a, Schematic of nanoliter cDNA synthesis reactions performed in wells of 384-well PCR plates with 3jil of hydrophobic overlay, b, Illustration of reduced-volume experiments with the lysis, RT and PCR volumes used, c, The number of genes detected per HEK293TF cell at each reaction volume, when sampling 100,000 sequencing reads (n = 100,19, 32 and 28 cells, respectively). P value represents a two-sided t-test between the 10-ul and 1-jil conditions, d, Influence of hydrophobic overlays on miniaturized cDNA synthesis (1 jil total volume). For each compound, boxes depict the number of genes detected per HEK293FT cell (n = 17, 34, 39, 28, 25, 24, 28, 38 and 70, respectively), subsampled at 200,000 sequencing reads per cell, e, Replacement of the bead-based cDNA cleanup by dilution in single HEK293FT (n = 58 and 52, respectively) cells. Box plots show the number of genes detected per cell and condition (at 100,000 reads) with P value for a two-sided t-test across conditions, f, Tagmentation complexity using 0.1 jil of ATM Tn5 enzyme per HEK293FT cell in relation to input cDNA. The median number of detected genes as a function of raw sequencing reads (n = 51, 53, 54, 53, 53 and 52 cells for 25, 50, 75,100, 200 and 500 pg, respectively), g, Tagmentation complexity for varying amounts of cDNA input. Complexity was summarized as unique aligned and gene-assigned UMI-containing read pairs per 400,000 raw reads and HEK293FT cell (n = 49, 51, 51, 50, 51 and 44). h, Schematic outline of the Smart-seq3 and Smartseq3xpress workflows, i, The number of genes detected with Smart-seq3xpress after variable amounts of pre-amplification PCR cycles. Median number of genes is reported as a function of raw sequencing reads in HEK293FT cells (n = 93, 98,108,113,102,114 and 118 cells for 10,12,13,14,15,16 or 20 cycles, respectively).], Fraction of UMI-containing reads to internal reads for HEK293FT cells prepared with Smartseq3xpress (KAPA HiFi; 12 PCR cycles), at a variable range of TDE1 Tn5 amounts (n = 64 cells each), k, Fraction of UMI-containing reads to internal reads for HEK293FT cells prepared with Smartseq3xpress (SeqAmp; 12 PCR cycles), at a variable range of TDE1 Tn5 amounts (n = 60 cells each). I,m, Optimization of RT and PCR conditions across 376 experimental conditions on HEK293FT cells. Colors indicate particular experimental conditions: Smart-seq3xpress with Smart-seq3 TSO (purple; n = 912), 52°C RT/alternate TSO implementation (yellow; n = 74), fixed spacer TSO variant (blue; n = 45), FLASH-seq TSO variant (green; n = 55), Smart-seq3xpress with improved TSO (pink; n = 63) and all other conditions (gray; n = 21,707). Scatter plots denote the level of artifactual TSO-UMI reads and RNA counting errors (I) as well as a percentage of ribosomal RNA (rRNA) mapped reads and number of detected genes in 100,000 reads after removal of strand invasion reads (m). n, Benchmarking of Smart-seq3 variants. Box plots show the number of genes detected per HEK293FT cell in full-volume Smart-seq3 (ref.2), low-volume Smart-seq3 and Smart-seq3xpress implementations, at the indicated read depths (n = 109-110,18-27, 9-170, 20-55 and 9-63 cells, depending on the cells available at the given sequencing depths). The box plots (in c, d, e, j, k and n) show the median and first and third quartiles as a box, and the whiskers indicate the most extreme data points within 1.5 lengths of the box. cSt, centistoke.
NATURE BIOTECHNOLOGY | VOL 40 | OCTOBER 2022 11452-1457 | www.nature.com/naturebiotechnology
1453
ARTICLES
NATURE BIOTECHNOLOGY
384-well plate
Hydrophobic . overlay (3 nl)
Nanoliter reaction volume
Reaction volumes (nl) Lysis: 3 1.5 RT: 1 PCR: 6 Total: 10
0.6
1.5 0.2
1.2
2
0.3 0.1 0.6 1
10,000 7,500 5,000
£ 8,000
o
o
S2 6,000
C
d) cn
£ 4,000
.Q
E
13
Z 2,000
HEK293FT P= 0.085
Tagmentation: input pre-amplified cDNA 25 pg 50 pg 75 pg
-•-100 pg-«- 200 pg -•- 500 pg
0 200,000 400,000 600,000 Number of reads
10 5 2 1 Total reaction volume (ni)
~ 320,000
x g
£ §> 280,000 E 2
o ~
£.§ 240,000
-Q 13
e3 No overlay b Alkane hydrocarbons ■ Silicone oils ■ Vapor-Lock
ä 12,000
S, 8,000
6,000
^ cí* ^ cí* íč^
, v vrtO" ri0" ri0" \>
' 12,500
0)
4^ I
25 50 75 100 200 500 Input cDNA (pg)
10,000
5,000
HEK293FT P =0.92
10,000
ä 6,000
2,000
Number of PCR cycles:
- 10 * 12* 13 * 14 * 15* 16
20
0 250,000 500,000 750,000 Number of reads
Smart-seq3
Standard volume
High PCR cycles ^_
cDNA pre-amplification bead clean-up
I
cDNA pre-amplification QC
cDNA normalization & dilution
Cell Isolation
Reverse transcription Template switching
[Tn5Motif][Tag] [U rGrGrG*.
C C C-
cDNA pre-amplification
Tagmentation
Sequencing
5' end reads (with UMI)
Smart-seq3xpress
Low volume + Overlay
*AAAAA
T[PCR]
Reduced PCR cycles
No clean-up, No QC, dilute cDNA & tagment directly
Spin out plates to pool
Internal reads
=> 0.75 .c
I
CO
1 0.50
0
KAPA 12 PCR cycles
0.7
Amount of TDE1 Tn5
=> 0.75
I
S 0.25
SeqAmp 12 PCR cycles
Amount of TDE1 Tn5
60
40
0 5 10 15 20 TSO artificial priming (%)
7,500
5,000
2,500
5 10 15 20 25 rRNA reads (%)
16,000 14,000 12,000 10,000 8,000 6,000 4,000 2,000
^ Smart-seq3 (low-vol, KAPA)
Ép Smart-seq3 (KAPA)
^ Smart-seq3xpress (KAPA)
^ Smart-seq3xpress (SeqAmp)
^ Smart-seq3xpress (SeqAmp, improved TSO)
• Smart-seq3xpress improved TSO • Fixed spacer TSO (ATAC) • All other conditions
• FLASH-seq TSO (CAGCA) 52 °C RT/alternate TSO • Smart-seq3 TSO
250k 500k 750k Number of reads
1454
NATURE BIOTECHNOLOGY | VOL 40 | OCTOBER 2022 11452-1457 | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
Monocytes (CD16+)
NK3/ILC
NK1
Monocytes, (CD14+)
Dendritic cells (conventional) 1
Cells with reconstructed T cell receptor
Dendritic cells (conventional) 2
CD8+ T cells
(EM)/ CD4+ Platelets CTL T cells (MAITJ \
4fr HSPCs CD8+T cells ?P(CM/EM) CD8+ T cells (naive)
CD4+T cells.4 (regulatory)
Dendritic cells (plasmacytoid)
Smar1-seq3xpress:KAPA Smar1-seq3xpress:SeqAmp TRG+TRD
T cells
(gamma-delta) 1 CD4+ T cells (CM/EM)
CD4+T cells (clonal)
CD4+ T cells (naive/CM) 2
fCD8+ T cells f (clonal)
' CD4/CD8 T cells (naive/CM) sT cells (gamma-delta) 2 s.CD4+ T cells CD4+ T cells (naive) (naive/CM) 1
-10
0
UMAP_1
10
15
* 0# , > -S> 0# a>
TCR chain pairing
B cells (memory) d
B cells (naive)
6,000 -
nes 5,000 -
cn 4,000 -
o 3,000 -
Ě E 2,000 -
z 1,000 -
0 -
B cell = 9.8 x 10~"
CD4 T cell P< 2.22 x 10~"
CDS T cell P< 2.22 x 10~"
NK cell P< 2.22 x 10""
■ Smart-seq2
■ Smart-seq3
■ Smart-seq3xpress
200 150 100 50 0
I TRAV8-2
TŘAV8-3
TRÁV8-4 TRAl/12-2
-0.5 0 0.5 1.0 1.5 Naive vs clonal CD4" T cells (fold change, log2)
CD4+ T cells (regulatory) T cells (prol.) Plalelels
Dendrilic cells (plasmacytoid) NK3/ILC NK2 NK1
CD4* T cells (naive/CM) 2 CD4* T cells (naive/CM) 1 CD4/CD8 T cells (naive/CM) CD8+ T cells (Naive) CD4+ T cells (Naive) B cells (Naive) B cells (Memory) T cells (MAIT) HSPCs
T-cells (gamma-della) 2 T-cells (gamma-della) 1 CD8+T cells (clonal) CD4+ T cells (clonal) Dendrilic cells (conventional) 2 Dendrilic cells (conventional) 1 CD8* T cells (EM)/CD4* CTL CD8*T cells (CM/EM) CD4*T cells (CM/EM) Monocytes (CD16+) Monocytes (CD14+)
o 15,000
o 10,000
B 5,000
o o
Q_
z
CO
Ť
I
t:
% \0 <> fď oř oř ,0 rÍJ S> Aj r> c% •Síl' % ?> ,■
Differential splicing
CD4 T cell
CD8 Tcell
60 40 20
&4
i j
i 10x Genomics v3.1 i Smart-seq3xpress
800
700
600
c 'CC 500
cn
O 400
m
_i 300
LU
200
100
0
-15-10-5 0 5 10 Effect size
ptprc
gusbp11
isg20
V.
1.0 0.8 0.6 0.4 0.2 0
NATURE BIOTECHNOLOGY | VOL 40 | OCTOBER 2022 11452-1457 | www.nature.com/naturebiotechnology
1455
ARTICLES
NATURE BIOTECHNOLOGY
<-
Fig. 2 | Application of Smart-seq3xpress to hPBMCs. a, Dimensional reduction (UMAP) of 26,260 hPBMC transcriptomes produced with Smartseq3xpress (KAPA, four donors; SeqAmp with improved TSO, three donors) colored and annotated by cell type. EM, effector memory; CM, central memory; NK, natural killer; ILC, innate lymphoid cell; HSPC, hematopoietic stem and progenitor cell; MAIT, mucosal-associated invariant T cell, b, Smartseq3xpress-based TCR reconstruction (TRaCeR) overlayed onto UMAP. c, QC of TCR reconstructions obtained with Scirpy, enumerating the number of T cells with certain types of TCR reconstructions, d, Benchmarking of Smart-seq2, Smart-seq3 and Smart-seq3xpress (SeqAmp, improved TSO) in primary hPBMCs. Each cell was downsampled to 100,000 reads, and the number of detected genes from exon-mapping reads is shown for representative cell types: B cells (n = 73, 366 and 859, respectively), CD4+T cells (n = 261,1,270 and 1,847, respectively), CD8+T cells (n = 76, 272 and 913, respectively) and NK cells (n = 73, 352 and 601, respectively). P values indicate results of two-sided t-tests between the Smart-seq3 and Smart-seq3xpress. e, Differential gene expression analysis (Wilcoxon test, Padj<0.01) between naive CD4~ T cell cluster (n = 2,476) and clonal CD4~ T cell cluster (n = 682). Indicated are the top five TCR genes driving the clonal CD4~ T cell cluster separation, f, Dot plot showing expression of selected marker genes for MAIT, gamma-delta and clonal CD4+/CD8+ T cells in all annotated clusters, with size of the dot denoting the detection of a gene within the cells of the cluster and color scale denoting the average expression level, g, Analysis of captured transcribed genetic variation in donor-matched Smart-seq3xpress and lOx Genomics 3' version 3.1 data. For each cell passing QC (n = 2,938 and 9,846, respectively), the number of SNPs with alternate allele coverage per cell are indicated (left) as well as the average SNP coverage normalized by the sequencing depth (right), h, Frequency of RNA-velocity-informative fully spliced reads in donor-matched Smart-seq3xpress and lOx Genomics 3' version 3.1 data. For each cell in representative cell types—B cells (n = 404 and 642, respectively), CD4+T cells (n = 1,320 and 1,317, respectively), CD8+T cells (n = 417 and 1,181, respectively) and NK cells (n = 441 and 498, respectively)—we summarized the percentage of reads spanning exon-exon junctions, with nominal Pvalues for a two-sided t-test. i, Differential splicing isoform analysis using BRIE2. Shown is a volcano plot of tested skipped-exon events color-coded by significant variation when testing for any cell type (ELBO gain >20). ELBO gain is a surrogate for the Bayes factor (that is, the likelihood ratio of two hypotheses). The x axis denotes the effect size on the distribution of PSI values in a cell type, j, Overlays of color-coded PSI values inferred for each sequenced cell (n = 26,260) in genes with significant cell type splicing variation (PTPRC, GUSBP11 and ISG20). The box plots (in d, g and h) show the median and first and third quartiles as a box, and the whiskers show the most extreme data points within 1.5 lengths of the box. HSPC, hematopoietic stem and progenitor cell.
detection was observed across cell types for Smart-seq3xpress in comparison to Smart-seq2 (ref.15) and Smart-seq3 (ref.2) (Fig. 2d). A recent study14 showed that scRNA-seq alone was not capable of separating certain cell types and states in hPBMCs without additional single-cell protein measurements. In a dataset with 200,000 hPBMC transcriptomes generated with lOx Genomics, 228 protein markers (CITE-seq) were used to distinguish unconventional T cell populations (mucosal-associated invariant T (MAIT) cells), gamma-delta T cells and effector memory T cells. All these cell types separated by the Smart-seq3xpress transcriptome data alone (Fig. 2a), exemplified with marker gene expression for MAIT cells, gamma-delta T cells and a CD4+ T cell population characterized only by their clonal expression of specific TCRs (Fig. 2e,f). Notably, this highly granular de novo cell type and state characterization was obtained from only 26,260 Smart-seq3xpress transcriptomes, and similar granularity remained at lower sequence depths (Extended Data Fig. 8).
Next, we performed a direct comparison of Smart-seq3xpress to droplet-based lOx Genomics (lOx Genomics 3' version 3.1; Methods) on a matched hPBMC donor. Tabulating the coverage over transcribed single-nucleotide polymorphism (SNP) positions—which is, for example, important for studies of allele-specific expression and cancer—showed that more SNPs were covered in the full-length Smart-seq3xpress data (9.2 million versus 2.1 million positions, over all cells). On a per-cell level, the number of SNPs with alternate allele coverage was approximately nine-fold higher with Smart-seq3xpress and with a three-fold-higher coverage per sequenced read (Fig. 2g). Similarly, full-length coverage with Smart-seq3xpress resulted in significantly (2-5-fold) increased read support over exon-exon and exon-intron splice junctions, which form the basis for RNA velocity inference16 (Fig. 2h). To perform a direct comparison of cluster granularity, we sampled an equal number of mean sequenced reads and cells for both lOx Genomics and Smart-seq3xpress and performed independent cluster analyses (Extended Data Fig. 9). Although both methods had approximately similar granularity at this number of cells and read numbers, cell type prediction based on Azimuth was slightly more consistent with clustering in Smart-seq3xpress (adjusted Rand index (ARI) of 0.59 and 0.49, respectively).
To study alternative splicing across the hPBMC Smart-seq3xpress data, we used BRIE2 (ref. 17) to identify 968 skipped exons with
significant variation in inclusion levels across cell types. (Fig. 2i). Overlaying cell-level percent spliced-in (PSI) estimates for alternative exons onto the UMAP revealed distinct splicing patterns across cell types. These included switching, graded changes, as well as heterogeneous cell types with variable inclusion levels (Fig. 2j), and, notably, they were independent of their gene-level expression levels between cell types (Extended Data Fig. 10). Thus, Smart-seq3xpress can generate unprecedented high-quality, full-transcript-coverage scRNA-seq data that enable scalable analyses of alternative splicing across complex human samples.
Large-scale efforts to enumerate cell types and states across human tissues and in model organisms are dominated by scRNA-seq methods that count RNAs by sequencing their ends. With the development of Smart-seq3xpress, we demonstrate a scalable solution for full-transcript-coverage scRNA-seq. Not only has Smart-seq3xpress overcome the main limitations of plate-based assays in terms of resources and material needed per cell, the data obtained also challenge the existing notions of the cell numbers required for efficient and finely resolved clustering of cells. Therefore, high-sensitivity scRNA-seq with isoform-specific and allele-specific resolution can, for the first time, be performed with Smart-seq3xpress at a scale suitable for large-scale cell atlas building.
Online content
Any methods, additional references, Nature Research reporting summaries, source data, extended data, supplementary information, acknowledgements, peer review information; details of author contributions and competing interests; and statements of data and code availability are available at https://doi.org/10.1038/ S41587-022-01311-4.
Received: 10 July 2021; Accepted: 8 April 2022; Published online: 30 May 2022
References
1. Mereu, E. et al. Benchmarking single-cell RNA-sequencing protocols for cell atlas projects. Nat. Biotechnol. 38, 747-755 (2020).
2. Hagemann-Jensen, M. et al. Single-cell RNA counting at allele and isoform resolution using Smart-seq3. Nat. Biotechnol. 38, 708-714 (2020).
3. Tabula Muris Consortium et al. Single-cell transcriptomics of 20 mouse organs creates a Tabula Muris. Nature 562, 367-372 (2018) .
1456
NATURE BIOTECHNOLOGY | VOL 40 | OCTOBER 2022 11452-1457 | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
4. Mayday, M. Y., Khan, L. M., Chow, E. D., Zinter, M. S. & DeRisi, J. L. Miniaturization and optimization of 384-well compatible RNA sequencing library preparation. PLoS ONE 14, e0206194 (2019).
5. Mamanova, L. et al. High-throughput full-length single-cell RNA-seq automation. Nat. Protoc. 16, 2886-2915 (2021).
6. Mora-Castilla, S. et al. Miniaturization technologies for efficient single-cell library preparation for next-generation sequencing. /. Lab. Autom. 21, 557-567 (2016).
7. Jaeger, B. N. et al. Miniaturization of Smart-seq2 for single-cell and single-nucleus RNA sequencing. STAR Protoc. 1, 100081 (2020).
8. Picelli, S. et al. Tn5 transposase and tagmentation procedures
for massively scaled sequencing projects. Genome Res. 24, 2033-2040
(2014) .
9. Nextera XT Library Prep: Tips and Troubleshooting. Illumina https://www.illumina.com/content/dam/illumina-marketing/documentsy products/technotes/nextera-xt-troubleshooting-technical-note.pdf
(2015) .
10. Ziegenhain, C, Hendriks, G.-J., Hagemann-Jensen, M. & Sandberg, R. Molecular spikes: a gold standard for single-cell RNA counting. Nat. Methods https://doi.org/10.1038/s41592-022-01446-x (2022).
11. Tang, D. T. P. et al. Suppression of artifacts and barcode bias in
high-throughput transcriptome analyses utilizing template switching. Nucleic Acids Res. 41, e44 (2013).
12. Hahaut, V, Pavlinic, D., Cowan, C. & Picelli, S. Lightning fast and highly sensitive full-length single-cell sequencing using FLASH-Seq. Preprint at https://www.biorxiv.org/content/10.1101/2021.07.14.452217vl (2021).
13. Parekh, S., Ziegenhain, C, Vieth, B., Enard, W. & Hellmann, I. zUMIs—a fast and flexible pipeline to process RNA sequencing data with UMIs. Gigascience 7, giy059 (2018).
14. Hao, Y. et al. Integrated analysis of multimodal single-cell data. Cell 184, 3573-3587 (2021).
15. Picelli, S. et al. Smart-seq2 for sensitive full-length transcriptome profiling in single cells. Nat. Methods 10, 1096-1098 (2013).
16. La Manno, G. et al. RNA velocity of single cells. Nature 560, 494-498 (2018).
17. Huang, Y. & Sanguinetti, G. BRIE2: computational identification of splicing phenotypes from single-cell transcriptomic experiments. Genome Biol. 22, 251 (2021).
Publisher's note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons. org/licenses/by/4.0/.
© The Author(s) 2022
NATURE BIOTECHNOLOGY | VOL 40 | OCTOBER 2022 11452-1457 | www.nature.com/naturebiotechnology
1457
ARTICLES
NATURE BIOTECHNOLOGY
Methods
Cell sources, culturing and sorting. HEK293FT cells (Invitrogen) were grown in DMEM medium (4.5 g L"1 of glucose and 6 mM L-glutamine, Gibco), supplemented with 10% FBS (Sigma-Aldrich), 0.1 mM MEM non-essential amino acids (Gibco), 1 mM sodium pyruvate (Gibco) and lOOugml"1 of pencillin-streptomycin (Gibco) at 37 °C. For sort of single cells, cells were harvested by incubation with TrypLE Express (Gibco). K562 (ATCC) cells were grown in RPMI medium, supplemented 10% FBS (Sigma-Aldrich) and 1% pencillin-streptomycin (Gibco). Frozen aliquots of 10 million hPBMCs from healthy individuals were purchased from Lonza, requiring healthy donors only. Written informed consent was obtained at sampling point from all donors by Lonza, and our analyses of hPBMCs were approved by the Regional Ethical Review Board in Stockholm, Sweden (2020-05070). hPBMCs were gently thawed and stained with PE mouse anti-human CCR7 (2-Ll-A,l:100), PE-Cy7 mouse anti-human CD4 (SK3,1:250), FITC mouse anti-human CD45RA (HI100,1:100), PerCP-Cy5.5/BB700 mouse anti-human CD8 (RPA-T8, 1:250) and PE-Cy5 mouse anti-human CD45RO (UCHL1, 1:250) (BD Biosciences) before sorting. For all cell types, dead cells were gated out after staining with propidium iodide (Thermo Fisher Scientific). Live, single cells were sorted into 384-well plates containing lysis buffer using a BD FACSMelody (BD FACSChorus version 1.3 software) equipped with a 100-um nozzle and plate cooling with index sorting on (BD Biosciences). After sorting, each plate was immediately spun down and stored at —80 °C.
Smart-seq3 library preparation. Full-volume Smart-seq3 library preparations were performed as previously described2. PCR was carried out using 20 cycles of amplification.
Low-volume Smart-seq3 and Smart-seq3xpress library preparation. For
experiments including overlays, including Vapor-Lock (Qiagen), silicon oils (Sigma-Aldrich) and tri-/te trade cane (Sigma-Aldrich), 3 ul of designated overlay was added to each well and stored at room temperature until use. Lysis buffer of various volumes (0.1-3 uL) was dispensed using either Formulatrix Mantis or Dispendix I.Dot One liquid dispenser to each well, all containing 0.1% Triton X-100, 5% PEG8000 adjusted to RT volume, 0.5 uM oligo(dT) adjusted to RT volume, 0.5 mM dNTPs each adjusted to RT volume and 0.5 U RNase Inhibitor (Takara, 40 U ul"1). After dispensing, lysis plates were briefly centrifuged down to ensure that lysis is properly collected and stored under the overlay. Stored plates of sorted cells were denatured at 72 °C for 10 minutes, followed by the addition of indicated volumes of RT mix; common for all is that the reagent concentrations are stable: 25 mM Tris-HCl ~pH 8.3 (Sigma-Aldrich), 30 mM NaCl (Ambion), 0.5mM GTP (Thermo Fisher Scientific), 2.5mM MgCl2 (Ambion), 8mM DTT (Thermo Fisher Scientific), 0.25 Uul"1 of RNase Inhibitor (Takara), 2uM TSO (5'-Biotin-AGAGACAGATTGCGCAATGNNNNNNNNrGrGrG-3'; IDT) and 2 Uul"1 of Maxima H Minus reverse transcriptase (Thermo Fisher Scientific). After RT mix dispensing, the plate was spun down to ensure merge of RT and lysis reactions. RT was performed at 42 °C for 90 minutes, followed by ten cycles of 50 °C for 2 minutes and 42 °C for 2 minutes. Indicated volumes of PCR master mix were dispensed, all containing constant reaction concentrations of lx KAPA HiFi PCR buffer (Roche), 0.3 mM dNTPs each (Roche), 0.5 mM MgCl2 (Ambion), 0.5 uM Smart-seq3 forward primer (5'-TCGTCGGCAGCGTCAGATGTGTATAAGAGACAGATTGCGCAATG-3'; IDT), 0.1 uM Smart-seq3 reverse primer (5'-ACGAGCATCAGCAGCATACGA-3'; IDT) and 0.02 Uul"1 of KAPA HiFi DNA polymerase (Roche). After dispensing, the plate was quickly spun down before being incubated in PCR as follows: 3 minutes at 98 °C for initial denaturation, 10-24 cycles of 20 seconds at 98 °C, 30 seconds at 65 °C and 2-6 min at 72 °C. Final elongation was performed for 5 minutes at 72°C. For conditions after cDNA pre-amplification clean-up: 100 nl of water, ExoSAP-IT express (Thermo Fisher Scientific) or 0.5UExoI (NEB) + 0.05 FastAP (Thermo Fisher Scientific) was dispensed per well and incubated at 37 °C for 15 minutes, followed by inactivation at 85 °C for 5 minutes.
For Smartseq3xpress with SeqAmp (Takara), lysis and RT was carried out with 0.125 uM or 0.5 uM oligodT30VN and 0.75 uM or 2 uM TSO unless otherwise indicated as described above. Original Smartseq3 TSO (5'-Biotin-AGAGACAGATTGCGCAATGNNNNNNNNrGrGrG-3'; IDT). Improved TSO (5'-Biotin-AGAGACAGATTGCGCAATGNNNNNNNN WWrGrGrG-3'; IDT). PCR mastermix was dispensed at 0.6 ul per cell containing lx SeqAmp PCR buffer, 0.025 Uul-1 of SeqAmp polymerase and 0.5uM/l uM Smartseq3 forward and reverse primer. After dispensing PCR mastermix, the plate was quickly spun down before being incubated as follows: 1 minute at 95 °C for initial denaturation, 6-18 cycles of 10 seconds at 98 °C, 30 seconds at 65 °C and 2-6 minutes at 68 °C. Final elongation was performed for 10 minutes at 72 °C. For Smartseq3xpress with NEBNext Ultra II Q5 Master Mix (NEB), PCR mastermix consisted of lx NEBNext Ultra II Q5 Master Mix and 0.5 uM/1 uM Smartseq3 forward and reverse primer and PCR was performed at 30 seconds at 98 °C for initial denaturation, 12 cycles of 10 seconds at 98 °C, 30 seconds at 65 °C and 6 minutes at 72 °C. Final elongation was performed for 5 minutes at 72 °C. For Smartseq3xpress with NEBNext Q5 Hot Start HiFi PCR Master Mix (NEB), PCR mastermix consisted of lx NEBNext Q5 Hot Start HiFi PCR Master Mix and 0.5 uM/1 uM Smartseq3 forward and reverse primer and PCR was performed
for 30 seconds at 98 °C for initial denaturation, 12 cycles of 10 seconds at 98 °C, 30 seconds at 65 °C and 1 minute at 65 °C. Final elongation was performed for 5 minutes at 65 °C. For Smartseq3xpress with Platinum SuperFi II DNA polymerase (Invitrogen), PCR mastermix consisted of lx SuperFi II Master Mix, 0.2 uM dNTPs and 0.5 uM/1 uM Smartseq3 forward and reverse primer and PCR was performed for 30 seconds at 98 °C for initial denaturation, 12 cycles of 10 seconds at 98 °C, 30 seconds at 60 °C and 6 minutes at 72 °C. Final elongation was performed for 5 minutes at 72 °C. For Smartseq3xpress with Platinum II Taq Hot Start DNA polymerase (Invitrogen), PCR mastermix consisted of lx Platinum II Taq Master Mix, 0.2uM dNTPs and 0.5uM/l uM Smartseq3 forward and reverse primer and PCR was performed for 2 minutes at 94 °C for initial denaturation, 12 cycles of 15 seconds at 94 °C, 30 seconds at 60 °C and 6 minutes at 68 °C. Final elongation was performed for 5 minutes at 68 °C.
A full and comprehensive protocol of Smart-seq3xpress has been deposited on protocols.io18.
After pre-amplification workflow. For regular Smart-seq3, pre-amplified cDNA libraries were purified with homemade 22% PEG beads at a ratio of 1:0.6. Library sizes were observed using Agilent Bioanalyzer High Sensitivity Chip, followed by concentration quantification using QuantiFlour dsDNA assay (Promega). cDNAwas subsequently diluted to 100 pgul"1 unless otherwise specified.
For low volume, pre-amplified cDNA libraries were diluted by the addition of 9 ul of water to each well, if not indicated otherwise, followed by a quick centrifugation. Library sizes were checked on an Agilent Bioanalyzer, using the high-sensitivity DNA chip; meanwhile, concentrations were quantified using QuantiFlour dsDNA assay (Promega). cDNA was normalized to 100pgul"1 if nothing else was specified.
For Smart-seq3xpress, pre-amplified cDNA libraries were diluted with the addition of 9 ul of water unless stated otherwise, before transferring 1 ul of diluted cDNA from each well into tagmentation.
Sequence library preparation for Smart-seq3xpress. Tagmentation was performed in 2 ul consisting of 1 ul of either diluted or normalized pre-amplified cDNA and 1 ul of lx tagmentation buffer (10 mM Tris pH 7.5, 5 mM MgCl2, 5% DMF), 0.025-0.5 ul of ATM (Illumina XT DNA sample preparation kit) or 0.0002-0.01 ul of tagmentation DNA enzyme 1 (TDElJllumina DNA sample preparation kit)). In the event of in-house Tn5, lx tagmentation buffer used consisted of 10 mM TAPS-NaOH pH 8.4, 5 mM MgCl2 and 8% PEG8000 and indicated amounts of 0.0005-0.01 uM in-house Tn5 enzyme. Samples were incubated at 55 °C for 10 minutes, followed by the addition of 0.5 ul of 0.2% SDS to each well. After addition of 1.5 ul/3.5 ul of custom Nextera index primers (0.5 uM) carrying 10-bp dual indexes, library amplification was started by the addition of 2/4ul of PCR mix (lx Phusion Buffer (Thermo Fisher Scientific), 0.01 Uul"1 of Phusion DNA polymerase (Thermo Fisher Scientific), 0.2 mM dNTP each) and incubated for 3 minutes at 72 °C; 30 seconds at 95 °C; 12-14 cycles of (10 seconds at 95 °C; 30 seconds at 55 °C; 30-60 seconds at 72 °C); and 5 minutes at 72 °C in a thermal cycler. Samples were pooled by spinning out each plate gently in a 300-ml robotic reservoir (Nalgene) fitted with a custom 3D-printed scaffold by pulsing to ~200g. The pooled library was purified with homemade 22% PEG beads at a ratio of 1:0.7.
lOx Genomics library preparation. After thawing the PBMC sample, we stained dead cells with propidium iodide (Thermo Fisher Scientific) and sorted 200,000 live cells into a 5-ml tube. After centrifugation, the cell suspension was resuspended and concentration of cells was determined using a Countess automated cell counter (Thermo Fisher Scientific). We loaded -13,000 cells for a target cell recovery of -8,000 cells and prepared libraries according to the lOx Genomics version 3.1 user guide. For both pre-amplification and post-fragmentation PCR, we applied 12 cycles of PCR.
Sequencing. Smartseq3 and Smartseq3xpress libraries were sequenced on a Illumina NextSeq 500 (Illumina NextSeq Control Software 2.2.0) or MGI DNBSEQ G400RS platform (version 1.1.0.108 software). For NextSeq runs with Smart-seq3, denatured libraries were loaded on HighOutput version 2.5 cartridges at 2.1-2.3 pM. For G400RS runs, libraries were created using phosphorylated index primers or subjected to five cycles of adapter conversion PCR using the MGIEasy Universal Library Conversion Kit (MGI) and subsequently circularized from 1 pmol of dsDNA according to the manufacturers protocol. Next, 60 fmol of circular ssDNA library pools were used for DNA nanoball (DNB) making using a custom rolling-circle amplification primer (5'-TCGCCGTATCATTCAAGCAGAAGACG-3'). DNBs were loaded on FCL flow cells (MGI) and sequenced using SE100, PE100 or PE150 cartridges using custom sequencing primers (Read 1: 5'-TCGTCGGCAGCGTCAGATGTGTA TAAGAGACAG-3'; MDA: 5'-CGTATGCCGTCTTCTGCTTGAATGATA CGGCGAC-3', Read 2: 5'-GTCTCGTGGGCTCGGAGATGTGTATAAGAGA CAG-3'; i7 index: 5'-CCGTATCATTCAAGCAGAAGACGGCATACGA GAT-3'; i5 index: 5'-CTGTCTCTTATACACATCTGACGCTGCCGACGA-3'). lOx Genomics version 3.1 libraries were sequenced on a NextSeq 500 according to
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
manufacturer specifications (1.8pM loading concentration; HighOutput version 2.5 150-cycle kits, 28-8-92 cycles for readl, indexl and read2).
Primary data processing. zUMIs13 version 2.8.2 or newer was used to process raw FASTQ files. Reads were filtered for low-quality barcodes and UMIs (4 bases < phred 20, 3 bases < phred 20, respectively) and UMI-containing reads parsed by detection of the pattern (ATTGCGCAATG) while allowing up to two mismatches. Reads were mapped to the human genome (hg38) using STAR version 2.7.3, and error-corrected UMI counts were calculated from Ensembl gene annotations (GRCh38.95). zUMIs was also used to downsample cells to equal raw sequencing depth to facilitate method benchmarking.
Analysis of Smartseq3xpress hPBMC data. Cells were filtered for low-quality libraries, requiring (1) more than 50% of read pairs mapped to exons+introns, (2) more than 20,000 read pairs sequenced, (3) more than 500 genes (exon+intron quantification) detected per cell and (4) less than 15% of read pairs mapped to mitochondrial genes. Furthermore, a gene was required to be expressed in at least ten cells. Analysis was done using Seurat v4.0.1 (ref.14). Data were normalized ('LogNormalize'), scaled to 10,000 and total number of counts, and mitochondrial fraction was regressed out. Using the Seurat integration function, the donor effect from the seven different donors in the dataset was removed. The top 10,000 variable genes were considered and 35 principal components for shared nearest neighbor (SNN) neighborhood construction and UMAP dimensionality reduction. Cell clusters were produced using Louvain algorithm at a resolution of 0.8. Cell types were identified by using the R package Presto (Wilcoxon & AUC, version 1.0.0). For the Azimuth predictions, a QC-filtered count matrix was uploaded to the Azimuth web-based application and processed according to the Azimuth app. For the direct donor comparison with lOx Genomics version 3.1 data, read counts from only donor 7 were downsampled to similar median sequencing depth as the comparable lOx dataset and quality filtered as follows: at least 10,000 read pairs, more than 50% of read pairs mapped to exons+introns and less than 15% read pairs mapped to mitochondrial genes. A gene was required to be expressed in at least ten cells. Data were subset to 3,000 cells for analysis in Seurat. Data were LogNormalized, scaled to 10,000 and mitochondrial genes regressed out. Default Seurat settings were used for neighborhood construction and dimensionality reduction. Cell clusters were assigned using the Louvain algorithm at a resolution of 0.8. Cell type identification was performed as above.
Analysis of lOx Genomics version 3.1 donor 7 hPBMC data. Raw sequencing data in FASTQ format were processed using zUMIs version 2.9.3 with automatic barcode detection based on the lOx Genomics version 3.1 allow-list. After completion, we exported full count tables including empty droplets and assigned ambient RNA droplets and real cells using the CellBender (version 0.2.0) re move-background function19. To filter for doublets, the CellBender output.h5 file was used with Solo20 (version 0.6). Additional doublets were discarded by manually inspecting the distribution of total UMI counts per droplet and discarding those greater than 45,000. For downstream analysis in Seurat, a low-quality filter was applied based on requiring at least 10,000 read pairs and less than 10% read pairs mapped to mitochondrial genes. A gene should be expressed in at least ten cells to be included. A subset of 3,000 cells out of 6,483 passing QC was used for direct comparison to Smart-seq3xpress. Seurat was run at default settings using SCTransform, and cell clusters were assigned at resolution 0.8 using the Louvain algorithm. Cell types were identified using Presto (Wilcoxon & AUC) together with reference-based approach performed by the Azimuth app.
TCR reconstruction. TCR sequences were reconstructed using TraCeR version 0.6.0 (ref.21) run in the teichlab/tracer Docker environment and using the --loci A B G D --species Hsap flags. Scirpy22 (version 0.8.0) was used to summarize and QC the output from TraCer.
Molecular spike data processing and analysis. Molecular spike data were extracted from aligned zUMIs BAM files and analyzed using the UMIcountR package10 (https://github.com/cziegenhain/UMIcountR, version 0.1.1). After loading the data using the extract_spike_daf function, overrepresented spikes were discarded with a read cutoff of 25 and higher. We next used molecular spike observations across all cells and conditions with at least five reads per molecule to sample 26 ground truth mean expression levels from 1 to 316 molecules per cell using the subsample_recompute' function. We then plotted the mean counting difference shaded by the standard deviation.
Identification of TSO strand invasion artifact. To identify UMI reads with the TSO mis-primed artifact, we loaded sequencing reads into R using Rsamtools (version 2.6.0). In the case of paired-end sequencing, only first-in-pair reads were selected using the appropriate SAM flags. The strand orientation of the mapped reads was also determined from SAM flags. Then, we extracted a 20-bp window of genome sequence upstream of the read start position on the positive strand (+stranded mappings) or downstream of the read start position+read length on the negative strand (—stranded mappings) using the BSgenome package (version 1.62.0, human hg38). Afterwards, we checked for presence of the UMI sequence
(with or without addition of GGG overhang) in the genomic window using Rs fuzzy string matching function (allowing 0,1 or 2 mismatches). This identification procedure of artifactual UMI reads was also implemented in Python3 to process aligned BAM files and remove all artifactual reads/read pairs (available on GitHub: https://github.com/cziegenhain/pyTSOfilter).
Isoform-based analysis. For analysis of skipped-exon (SE) isoform differences, we retrieved annotations from GenCode (Human v39) and produced the SE annotation in GFF file format using BRIEkit-event (version 0.2.2). We filtered SE events using BRIEkit-event-filter with the following criteria: (1) retain SE events on autosomes and X/Y chromosomes; (2) SE events not overlapped by any other AS-exon; (3) surrounding introns are no shorter than a fixed length (100 bp);
(4) presence of specific splice sites (that is, surrounded by AG-GT); and
(5) SE events have a minimum distance (lObp) from transcription start site or transcription termination site. Next, we summarized the coverage over the filtered SE events for each cell using the brie-count command from BRIE2 (version 2.0.6) using per-cell demultiplexed, aligned and TSO-artifact-filtered (see above) BAM files as input. The resulting count files in h5ad format were used as input for the Bayesian regression-based inference of PSI values and variable splicing detection over cell types. We applied the aggregated imputation mode introduced by BRIE2 to fit the gene-wise prior distribution through aggregation of data over all cells for each gene. Default settings for Monte Carlo EM were applied. Genes were filtered by requiring at least 50 counts, ten unique counts and at least 30 cells with unique counts. The minimum required minor isoform frequency was 0.001 (default settings). For variable splicing detection, we annotated each cell with a binarized dummy factor of cell type identity (Seurat clustering; Louvain resolution 2.0), removing the most common cell type to avoid collinearity of the design matrix. We loaded the resulting h5ad file into scanpy23 (version 1.8.2)
for visualization of PSI values. For selection of SE events with significant cell type difference, we selected the highest evidence lower bound (ELBO) value per SE for each of the cell type LRT indices. Gene model plots to visualize significant cell-type-variable SE events were generated using the Gviz (version 1.38.1L, https://link.springer.com/protocol/10.1007%2F978-l-4939-3578-9_16) and rtracklayer (version 1.54.0, https://academic.oup.com/bioinformatics/ article/25/14/1841/225816) R packages.
SNP and junction coverage analysis. Coverage over transcribed SNPs was analyzed per cell using the cellsnp-lite24 package (version 1.0.0) over the most common human polymorphisms (1000 Genomes Project minor allele frequency >0.005, 36 million positions). We applied default settings of minimum aggregated count over cells 20 and minimum MAPQ for read filtering of 20 (essentially discarding multimapping reads due to the mapping quality encoding of the STAR aligner). Coverage on RNA velocity informative positions was tabulated from zUMIs output BAM files. Fully spliced exon-exon reads were identified by the presence of splicing in their CIGAR value and exclusive assignment to exon regions, whereas nascent (that is, unspliced or partially spliced) exon-intron spanning reads were identified by the overlap with both exonic and intronic regions of the same gene.
Reporting summary. Further information on research design is available in the Nature Research Reporting Summary linked to this article.
Data availability
Sequencing data have been deposited at ArrayExpress under the following accession numbers: E-MTAB-11488, E-MTAB-11452 and E-MTAB-11467.
Code availability
Processing of Smart-seq3xpress data was implemented in zUMIs (https://github. com/sdparekh/zUMIs). Code to filter reads with the TSO strand invasion artifact was implemented in a stand-alone script pyTSOfilter (https://github.com/ cziegenhain/pyTSOfilter).
References
18. Hagemann-Jensen, M., Ziegenhain, C. & Sandberg, R. Smart-seq3xpress. protocols.io https://www.protocols.io/view/smart-seq3xpress-bwh4pb8w (2022).
19. Fleming, S. J., Marioni, J. C. & Babadi, M. CellBender remove-background: a deep generative model for unsupervised removal of background noise from scRNA-seq datasets. Preprint at https://www.biorxiv.org/ content/10.1101/791699vl (2019).
20. Bernstein, N. J. et al. Solo: doublet identification in single-cell RNA-seq via semi-supervised deep learning. Cell Syst. 11, 95-101 (2020).
21. Stubbington, M. J. T. et al. T cell fate and clonality inference from single-cell transcriptomes. Nat. Methods 13, 329-332 (2016).
22. Sturm, G. et al. Scirpy: a Scanpy extension for analyzing single-cell T-cell receptor-sequencing data. Bioinformatics 36, 4817-4818 (2020).
23. Wolf, F. A., Angerer, P. & Theis, F. J. SCANPY: large-scale single-cell gene expression data analysis. Genome Biol. 19, 15 (2018).
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
ARTICLES
NATURE BIOTECHNOLOGY
24. Huang, X. & Huang, Y. Cellsnp-lite: an efficient tool for genotyping single cells. Bioinformatics 37, 4569-4571 (2021).
Acknowledgements
We are grateful to G.-J. Hendriks for 3D printing the adapter for library pooling, I. Dumitru for help with lOx Genomics library preparations and M. Seifert and L. Engstrand for support with sequencing. We thank the Sandberg Laboratory for fruitful discussions. This work was supported by grants to R.S. from the Swedish Research Council, the Knut and Alice Wallenberg Foundation and the Goran Gustafsson Foundation. C.Z. was supported by an EMBO long-term fellowship (ALTF 673-2017).
Author contributions
Developed and optimized the method: M.H.-J. Designed the experimental infrastructure for miniaturization and automation: C.Z. Performed experiments and analyzed data: C.Z. and M.H.-J. Conceived the study and wrote the manuscript: C.Z., M.H.-J. and R.S. Supervised the work: R.S.
Funding
Open access funding provided by Karolinska Institute.
Competing interests
M.H-J. and R.S. are inventors on the patent relating to Smart-seq3 that is licensed to Takara Bio USA. C.Z. declares no competing interests.
Additional information
Extended data is available for this paper at https://doi.org/10.1038/s41587-022-01311-4.
Supplementary information The online version contains supplementary material available at https://doi.org/10.1038/s41587-022-01311-4.
Correspondence and requests for materials should be addressed to Rickard Sandberg.
Peer review information Nature Biotechnology thanks Mario Suva and the other, anonymous, reviewer(s) for their contribution to the peer review of this work.
Reprints and permissions information is available at www.nature.com/reprints.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
a
^10,000
T3 03
cd
>- 8,000
cd -q
E 2,000
K562
10 5 2 1
Total reaction volume (pL)
k562
^ <á^ ^ ^ i£ ^
co ■d
ra
cu
cd 6,000
O 4,000
i—
cd
= 2,000
K562 P = 0.75
hek293t
TT
k562
^3
10000-
8000- t
? 6000-
4000-
2000-
10 pL 10 mL no overlay with overlay
10000
8000
2000
K562 p-0.1
J—J
3
6000 4000
10 (JL 10 ML no overlay with overlay
g
■a 12,500
I :§ io,ooo-
^ S 7,500.
§ — 5,000-O
"= « 10,000
I £ 8,000
'1 6,000
5 -1 4 000
90,000
I 30,000 0
^ 0.25X ^ 0.75X ^ 1.25X ^ 0.50X ^ 1.0X ^ 1.50x
pcr:rt volume
pi? '
in
jtF=*f=f
200,000
Number of reads
T T
—
PCR extension time
2min ■ 3min 4min ■ 6min
> 4000 bp
3500-4000 bp
3000-3500 bp
£. 2500-3000 bp
c 2000-2500 bp £
< 1750-2000 bp
or
e 1500-1750 bp >-
XI
■o 1250-1500 bp c
J 1000-1250 bp I 750-1000 bp 500-750 bp 250-500 bp 0-250 bp
ř
500 1000 1500
Number of detected genes (per bin and cell
> 4000 bp 3500-4000 bp 3000-3500 bp 2500-3000 bp 2000-2500 bp 1750-2000 bp 1500-1750 bp 1250-1500 bp 1000-1250 bp 750-1000 bp 500-750 bp 250-500 bp 0-250 bp
pcr extension time 0 3min ^ 4min 0 6min
i
f
Í
500 1000 1500 2000
Number of detected genes (per bin and cell)
Extended Data Fig. 11 See next page for caption.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
ARTICLES
NATURE BIOTECHNOLOGY
Extended Data Fig. 11 Optimization of low-volume cDNA synthesis, (a) Reduction of reaction volumes in single K562 cells. Shown are the number of genes detected per cell at each reaction volume when sampling 100,000 sequencing reads (n = 63, 39, 55, 53 cells, respectively). P-value represents a t test between the 10 j-iL and 1 jiL conditions, (b) Coefficient of variance of cells at scaled volumes in both HEK293FT (n = 100,19, 32, 28 cells, respectively) and K562 cell lines (n = 63, 39, 55, 53 cells, respectively), (c-d) For both (c) HEK293FT and (d) K562 cells the influence of creating standard volume Smart-seq3 libraries without (HEK293FT n = 37; K562 n = 26) and with an overlay (VaporLock; HEK293FT n = 15; K562 n = 31) was compared. Boxplots show genes detected and p-values show result of a two-sided t-test. (e) Replacement of the bead-based cDNA cleanup by dilution in single K562 (n = 57, 38, respectively) cells. Shown are the number of genes detected per cell for each condition at 100,000 reads with a p-value for a two-sided t-test within cell types, (f) Coefficient of variance between cells having received cDNA clean-up or dilution in HEK293FT (n = 58, 52, respectively) and K562 (n = 57, 38, respectively), (g) For indicated ratios of PCR volume to RT volume, boxplots show genes detected, genes with UMIs detected and UMIs captured, downsampled by sequenced reads. For PCR ratios O.lx (n = 43), 0.2x (n = 48), 0.3x (n = 48), 0.4x (n = 47), 0.5x (n = 47), 0.6x (n = 45) HEK293FT cells were analyzed, (h) Boxplots for each PCR extension time (2, 3, 4, 6 min) using KAPA HiFi Hot Start polymerase are shown as detected genes binned by their transcript length at 250,000 reads per cell, (n = 29, 24, 29,10; respectively) (i) Number of genes detected binned by transcript length for each extension time 3 min (n = 64), 4 min (n = 64), 6 min (n = 60) using SeqAmp polymerase at 250,000 reads per cell. The boxplots (in a-i) show the median, first and third quartiles as a box, and the whiskers show the most extreme data points within 1.5 lengths of the box.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
a 320,000
0.1|jLATMTn5
number of cells =
50 75 100 200 500
cDNA input (pg) Input scaled (reaction volume: 2uL)
E o
co
number of cells =
22 24 39 43 44 47
E±3
number of cells =
0.025 0.05 0.075 0.1 0.2 0.5
ATM Tn5 amount (uL)
— O)
a. cu £ *=
81;
£•=> ra a)
I Uii#
0.0250.050.075 0.1 0.2 0.5
ATM Tn5 (uL)
g
Volume scaled (100 pg cDNA input)
Density:
CD -0 275,000 Q. Si E J
8 "250,000
2 i
225,000
0.4 0.6 1 2
0.02 0.025 0.05 0.075
Reaction volume (uL), ATM Tn5 amount O.luLATM Tn5, 100pg cDNA
8 u 260,000 5 5 240,000 220,000
10,000
It
e a
I in-house Tn5
100 pg cDNA B ATM Tn5
Tn5
exosap19ul
noc!ean19ul
exofastap19ul
exosap19ul
noclean19ul
exofastap19ul
exosap9fjl
noclean9ul
Extended Data Fig. 2 | See next page for caption.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
ARTICLES
NATURE BIOTECHNOLOGY
Extended Data Fig. 2 | Investigating tagmentation complexity, (a) Systematic investigation of tagmentation complexity was performed by varying cDNA input with constant Tn5 amount, varying input and Tn5 amounts, varying Tn5 amounts with constant cDNA input, scaling of reaction volumes and Tn5 amounts with constant cDNA input, scaling of reaction volumes with constant cDNA and Tn5 amounts. For each boxplot, shown is the library complexity in terms of unique gene-assigned UMI-anchored read pairs (unique per-molecule cut-sites) from 400,000 raw sequencing reads. Each condition contains between 22 and 73 HEK293FT cells as annotated above each box. (b) Concordance of gene expression levels between HEK293FT cells tagmented using 0.05(.iL ATM Tn5 (n = 23 cells) and 0.1 jiL ATM Tn5 (n = 11 cells) (mean UMI counts over 15,915 genes), (c) Tagmentation complexity using 100 pg cDNA per single HEK293FT cell in relation to the enzyme amount (ATM Tn5). For each dot, the median number of detected genes is calculated from the indicated number of raw sequencing reads and plotted from n = 55, 54, 56, 58, 52, 51 cells (0.025 (il_, 0.05 mJ_, 0.075jiL, 0.1 \iL, 0.2 nl, 0.5 \iL). (d) For varying amounts of Tn5 enzyme (see (c)), tagmentation complexity was summarized as unique aligned and gene-assigned UMI-containing read-pairs per 400,000 raw reads per HEK293FT cell (n = 53, 46, 56, 57, 52, 50, respectively), (e) The use of in-house Tn5 relative to performance of ATM Tn5 was compared in HEK293FT cells. Each boxplot shows the number of genes detected at 500,000 raw reads (n = 182, 52, 24, 226, 41, 53, 50, respectively), (f) Influence of post PCR clean-up for Smartseq3xpress. Genes detected after treating preamplified cDNA libraries to reduce "contaminants", such as dNTPs and oligonucleotides, before going into tagmentation. All conditions (dilution alone (n = 152,129), ExoSAP IT-express ("ExoSAP"; n = 254,178), or Exol + Fast-AP (n = 259,167)) were either diluted in 9 or 19 jiL of water, (g) Influence of post PCR clean-up for UMI captures. Pairwise comparisons of mean expression estimates (mean UMIs per gene) for evaluated clean-up conditions. The boxplots shown in a, d, e and f indicate the median, first and third quartiles as a box, and the whiskers show the most extreme data point within 1.5 lengths of the box.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
KAPA .1 . ......J . ..i.J . ..I...J NEBNexl J . ..I...J . ..i.J . ..I...J Platinumll J . ..i.J . ..I...J . ,.i...J
05U«rall f' * t*.......\ ..'fiHw- SeqAmp \ ......■! ......J ......A 1 ......4 .......1 ......i
X110000' (D
"Ö 5000 O
"S o
T3 CO CD
KAPA NEBNext Platinumll
✓ l..líf«rf..*J / 1 ř.H...d M 1 ..,J .£ ..,J
Q5Ultrall SeqAmp SuperFJII
/ LiftiJ, *íf?.,iJ / / IM s&J
tjp s«f <ď tjftjř a? (J? áfttf # ° ^
Sequenced read pairs
Sequenced read pairs
o c 40 S, 30 KAPA NEBNext Platinumll
Ul Ju S 20 I 10 3 ° - ........ **-
Ü 40 Q. 05Ultral[ SeqAmp SuperFill
Q. 30 1 20 ^ 0 rÖ rö rö •'3*
75
T3 50 0>
Q.25
O.
(0
E c
=) 75
c£ 50
25
..Jt".Ji
rS$> JBM^
Sequenced read pairs
KAPA NEBNext Platinumll
i ...ij ...iJ ...iJ i .....jTffli i ...ü ...,jS1
Q5Ullrall SeqAmp SuperFill
:/&
I ->-i-J .....~l .....-I I ...ij .TEfSii I ,„|J ,„1,4 „ 40 ~" 20
Q
100 1 10 100 1 10 Mean expression level (molecules)
i
„100000 = 10000
5 1000
- 100
=100000 I 10000 - 1000 100
KAPA NEBNext Platinumll
F 4* F- ^ % J> F
F...U ...d F....J ...ü ...Ü F....J ...u
Q5Ullrall SeqAmp SuperFill
F« F I'
oř
..,|J n.lJ
(*■ cÖ1 cÖ' sff# fÄtfctf é" & O*
Sequenced read pairs
d
r - 1
........................ _____________________________ -......T-A-.f*g«j-3g^-: - - - - -
PCR polymerase —
KAPA — Platinumll SeqAmp NEBNext Q5Ultrall - SuperFill
Template switching oligo concentration
0.5 nM 10 uM 2.0 uM
£ 5000
0 0.2 0.4 0.6 0 0.2
Sequence depth (M reads)
Template switching oligo concentration
0 5 M 1.0 |jM 2.0 pM
s
or
Sequence depth (M reads)
, Ö KAPA E5 Platinumll Ö SeqAmp
PCR polymerase ^ NEBNext Ö QSUItrall E5 SuperFMI
Transitions Transversions
Extended Data Fig. 3 | See next page for caption.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
ARTICLES
NATURE BIOTECHNOLOGY
Extended Data Fig. 3 | Performance of pre-amplification polymerases for Smartseq3xpress. (a) Mapping statistics of the six different tested commercial polymerases: KAPA HiFi Hot Start (n = 384), NEBNext Q5 ("NEBNext", n = 384), NEBNext Ultra II Q5 ("Q5Ultrall", n = 384), Platinum II (n = 384), SuperFi II (n = 384), SeqAmp (n = 360). Dotplots show percentage of mapped reads to exons, introns, intergenic regions and unmapped with the sequenced read depth, (b) Number of Genes and UMIs detected per cell relative to sequencing depth for all six polymerases tested, (c) Fraction of UMI-containing reads versus internal reads for each of the tested polymerases, (d) Lines show median number of genes detected and median number of unique tagmentation sites after downsampling of sequenced reads using 0.5 jiM, 1.0 jiM and 2.0 jiM TSO concentration in RT, combined with either 0.5 jiM or 1 jiM of each forward and reverse PCR primer. The replicate numbers for each of the panels are: TSO = 0.5 jiM, PCR = 0.5 jiM: n = 58, 57, 62, 63, 52, 55 cells; TSO = 0.5nJVi, PCR = 1.0nJVi: n = 51, 60, 61, 45, 51, 55 cells; TSO = 1.0nJVi, PCR = 0.5nJVi: n = 57, 58, 61, 52, 51, 57 cells; TSO = 1.0nJVi, PCR = 1.0nJVi: n = 51, 56, 63, 49, 52, 55 cells; TSO = 2.0(iM, PCR = 0.5(iM: n = 62, 59, 61,53,54, 51 cells; TSO = 2.0 (iM, PCR = 1.0 (iM: n = 49, 53, 59, 50, 55, 49 cells (KAPA, Platinumll, SeqAmp, NEBNext, Q5Ultrall, SuperFill, respectively), (e) Molecular spike-ins were used to assess the accuracy of each polymerase in mRNA molecule counting, based on the counting difference between internal molecular spikes counts and Smart-seq3xpress UMIs at indicated amounts of TSO and PCR primer concentrations. Colored lines indicate the mean counting difference between the unique spike identifiers and quantified UMIs when sampling the spike at the indicated mean expression levels for each of the polymerases, shaded by the standard deviation. The counting differences are expressed in absolute deviance or relative to the mean molecule number, (f) Rate of base conversions in aligned reads relative to the reference genome. For every polymerase, we compute the average fraction of transitions and transversions, shown as boxplots over all cells, n = 328, 367, 315, 343, 302, 322 (for KAPA, Platinumll, SeqAmp, NEBNext, Q5Ultrall, SuperFill). Boxplots indicate the median, first and third quartiles as a box, and the whiskers show the most extreme data point within 1.5 lengths of the box.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
NNNNNNNNGGG I I I I I I I I I I I I I I I I I I I I I I
Alignment start position
Reference sequence
screen forTSO-UMI sequence in 20 bp upstream window (=< 1 error)
1000 1500 Distance (bp)
o 20 E
£ 15
E 0 E 20
Colums: OligodT, Rows: Fwd/Rev PCR primer
TSO concentration •0.75[iM 6 1uM 0.0625pM
0.125|jM
number of cells
number of cells
96 96 120 96
128 96 96 120 96 96 120
96 96 120 193 192 256 192 191 256
■---*
number of cells
4^//////^/////////,
0'
nnnnnnn ntwwwc nnnnnnnnatnacn
FLASH-seq: nnnnnnnncatca
NNNNN NNN ATAC NNNNNNNNATACN NNNNNNNNNATN NNNNYRNN NWWA NNNNNNNNWWWC
NNNNYflNNNWWC
NNNNNNNNATWW NNNNNNWWWWC NNNNWNNNWWC
FLASH-seq: nnnnnnnncagca
NNNNWNNNWWA NNNNNNNNNATW NNNNN NN WWWC NNNNNNN N WWC NNNNNNNWWWH NNNNWNNNWWW NNNNNNWHATN NNNNNNWWWH
NNNNWNNNANC NNNNWNNNWNA NNNNATNNNNC NNNNNNNNWWW
NNNNNNWWWW NNNNNNNWATN NNNNYRNN N WW NNNNNNN N AWW NNNNN NN WWI I NNNNN]NNWWM NNNNNNN WWC
NN NM NN.NN I:
S m a rt-seq 3x p re ss I m proved:
NNNNWNNNWW NNNNNNHHHH (0.75 uM) NNNNNNNNAC NNNNNNNANC
NNNNNNTNNC NNNNNNANTN NNNNATNNNN NNNNYRNNNW NNNNNNANNC NNNNNNATNN
TSO artificial priming (%)
Number of genes {100,000 reads)
TSO RNA counting errors Cells
UMI Coding capacity (k) 250 500 750 1000
Extended Data Fig. 4 | See next page for caption.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
ARTICLES
NATURE BIOTECHNOLOGY
Extended Data Fig. 4 | TSO strand-invasion artifact and improved TSO designs, (a) Schematic representation of the search procedure for artifactual TSO-UMI reads. For every aligned 5' read, a 20 bp window of reference genome sequence upstream of the alignment start position is considered. Within this sequence, we search for the UMI sequence allowing up to 1 mismatch, (b) Strand invasion leads to shortened captured molecules. We grouped reads by the presence or absence of TSO-UMI match in the 20 bp upstream window and retrieved the closest annotated transcription start site (TSS) of the assigned gene. Shown is the distance to TSS as the cumulative percentage of reads analyzed, (c) Influence of TSO design (x-axis), oligo-dT primer amount (columns) and forward/reverse PCR primer amounts (rows) on the occurrence of strand-invasion artifacts (y-axis). TSO concentration is indicated in color (red: 0.75 uM, teal: 1 pM). Replication was performed over many cells per condition as annotated above each respective box. (d) Shown are all evaluated UMI sequences incorporated into the Smart-seq3 TSO with their base composition in terms of random or stable bases. For each TSO, we show sequencing metrics in HEK293FT cells (numbers of cells per condition annotated in the Figure), in terms of the frequency of artificial TSO priming, Number of genes detected after discarding TSO primed molecules (100,000 raw reads) and the accuracy of the UMI counting as assayed by molecular spikes. Every box is colored by the coding capacity of the associated random bases in the TSO. Boxplots in c and d indicate the median, first and third quartiles as a box, and the whiskers show the most extreme data point within 1.5 lengths of the box.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
0 0.25 0.5 0.75 1 Number of reads (M)
S—
Z§ 2.000
0 0.25 0.5 0.75 1 Number of reads (M)
o u 20000-
lo
- NNNNNNNNC TSO (0.125MM OligodT)
Smartseq3xpress lmproved_TSO(0.125|jM OligodT) NNNNNNNNWW
- NNNNNNWWWW_TSO (0.125pM OligodT) NNNNNNNWWM TSO (0.125pM OligodT)
- NNNNNNNWM TSO (0.125pM OligodT)
- Smartseq3 Original TSO (0.5uM OligodT) NNNNNNNN
0 0.25 0.5 0.75 1 Number of reads (M)
PBMCs PBMCs
Smartseq3xpress Improved TSO (0.75 l>M) Smartseq3 Original TSO {1 pM)
Number of UMI-reads
J?
Number of UMI-reads
PBMCs
4^
4P
HEK293FT
HEK293FT
HEK293FT
12000
tß 10000
125,000
0
S 100,000 .a
1 75,000 Z
50,000 25,000 0
M Low-volume Smart-seq3 (KAPA)
O Smart-seq3 (KAPA)
£3 Smart-seq3xpress (KAPA)
■ Smart-seq3xpress (SeqAmp)
E^3 Smart-seq3xpress (SeqAmp, Improved TSO)
0 100,000 200,000 300,000 400,000 500,000 Number of sequenced read pairs
g
E 0.£
m 2 1e+01 CT Q.
« iE i. cl ÍĎ E
nneighbors I
1e-01 1e+01
Smart-seq3
Extended Data Fig. 5 | See next page for caption.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
ARTICLES
NATURE BIOTECHNOLOGY
Extended Data Fig. 5 | Evaluation of candidate reaction conditions in human PBMCs and HEK293FT cells, (a) Candidate TSO sequences and previous TSO ("Smartseq3 Original") were evaluated in human PBMC samples with a change to 0.125jiM oligo-dT primer. At indicated sequencing depths, we show the median number of detected genes (left), median number of detected genes after discarding TSO-priming artifact reads (middle), and median number of detected UMIs after discarding TSO-priming artifact reads (right), (b) Investigation of the optimal amount of oligo-dT primer (colors) in the context of new TSO ("Smartseq3xpress improved"; left; n = 74, 71, 74, 83 cells for 0.0625 uM, 0.125 uM, 0.25 uM, 0.5 uM, respectively) and previous TSO ("Smartseq3 Original"; right; n = 293, 286, 313, 319 cells for 0.0625 uM, 0.125 uM, 0.25 jiM, 0.5 uM, respectively) in PBMCs. Shown are the median number of detected genes per cell after discarding TSO-priming artifact reads, (c) Frequency of TSO artifact reads in PBMC cells for new TSO ("Smartseq3xpress Improved"; left; n = 64, 63, 53, 58 cells, respectively) and previous TSO ("Smartseq3 Original"; right; n = 266, 256, 262, 251 cells, respectively), with colors denoting the amount of oligo-dT primer used, (d) Benchmarking of new Smart-seq3xpress Improved TSO in HEK293FT cells. At the indicated sequencing depth, we show the number of genes in internal+UMI reads (left) and TSO-artifact filtered UMI reads (middle), (e) The reduction in the occurrence of TSO primed strand-invasion artifacts is shown as a boxplot for n = 94, 330 HEK293FT cells (improved TSO, original TSO). (f) Benchmarking of Smartseq3 variants and Smart-seq3xpress iterations. Shown are the number of UMIs captured in HEK293FT cells in the full-volume Smart-seq3 (n = 110 cells), low-volume Smart-seq3 (n = 27 cells) and Smart-seq3xpress iterations (n = 170, 55, 63 cells for KAPA, SeqAmp and SeqAmp improved TSO, respectively) at the indicated read depths, (g) Reproducibility over cells visualized by cell-to-cell correlation for Smartseq3 and Smartseq3xpress (n = 107, 62, respectively). Two-sided t-test p-value < 2.2e-16. (h) Representative correlation in captured molecules between a Smartseq3 cells and a Smartseq3xpress cell. Boxplots in c, e, f and g show the median, first and third quartiles as a box, and the whiskers show the most extreme data point within 1.5 lengths of the box.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY ARTICLES
a
b
-i-1-1-1-1-1-1-1-1-1-1-1—
01 23456789 10 11
Time (hours)
Extended Data Fig. 6 | Experimental improvements to the Smart-seq3xpress workflow, (a) Introduction of a 3D-printed adapter to facilitate pooling of final libraries by centrifugation. Pictures showing reservoir and 3D printed frame/holder and assembly to facilitate pooling of 384 well plates quickly via gentle centrifugation. (b) For Smart-seq3 and Smart-seq3xpress (42 °C and 52 °C reverse transcription implementations), we show the workflow in the library preparation process with the associated timings. Any steps that require hands-on work are shaded with pattern.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
ARTICLES
NATURE BIOTECHNOLOGY
donor
• donorl donor3 • donor6
Naive CD4+ T_1 CD4+ TCMrTEM_1 Naive CD8+ T_1 CD4+ TCMrTEM_2 Inter/Memory B_1 NK_3 NK_1 NK_2 Naive B
Clonal CD4+ T (TEM)_2 Clonal CD4+ T (TEM)_1 Inter/Memory B_2 MAIT CD8+ TEM_2 CD4/CD8 T Naive CD16+ Mono_1 CD8+ TEM_1 CD16+ Mono_2 Clonal CD4+ T (TEM)_3 NMLC Memory B_1 Clonal CD8+ T (TEM)_1 Tregs cDC_1 CD4+ TCM/TEM_3 CD14+ Mono_2 CD14+ Mono_1 gdT_1 Naive CD8+ T_2 Naive CD8+ T_3 NK/CTL pDC Memory B_2 Proliferating cells Platelets gdT_2 Plasma Cells Naive CD4+ T_2 cDC_2 HSPCs
I--< ■ ~|-«* •
■ ,-f<
•-• f- - ■• -
4.---
••-('■•
-4 -4-
•-•^ - » — -•
+—■ + »■»•
h-—--
; -
:
+-
4 —
:
5 10 15
% of cells in cluster
b
Annotated Subset2 (CD4,CD8) Subset3 (CD45RA, CCR7)
-10 0 10 -10 0 10 -10 0 10
UMAP 1 UMAP 1 UMAP 1
Extended Data Fig. 7 | Analysis of PBMC samples using Smart-seq3xpress. (a) Distribution of cell type abundances as a percentage of all cells from each of the 7 donors, (b) FACS-based (index sorting) data overlayed onto UMAP embeddings. "Top Left" shows annotated UMAP based on gene expression. "Top middle", shows classification of all sorted live cells based on their expression of CD4 and CD8. (DnT = double negative CD4- and CD8-, DpT = Double positive CD4 + and CD8 + , NA = cells without staining). "Top Left" all recorded, and index sorted live cells are categorized based on their expression of CD45RA and CCR7 overlayed over UMAP embeddings (CM = CD45RA- CCR7 + , EM = CD45RA- CCR7-, TEMRA = CD45RA + CCR7-, TN = CD45RA + CCR7 + , NA = cells without stainings). "Lower panels" show protein levels for antibody-staining against CD45RA, CD45RO and CCR7 overlayed over UMAP embeddings of sequenced transcriptomes as log-scaled mean fluorescence intensity (mFI).
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
Extended Data Fig. 8 | See next page for caption.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
ARTICLES
NATURE BIOTECHNOLOGY
Extended Data Fig. 8 | Downsamplingof Smart-seq3xpress PBMC data. Sequencing data generated for 26,260 human PBMCs was downsampled from a median depth of 258,000 read pairs per cell to steps of 75% (-193,000 reads), 50% (-129,000 reads), 25% (-64,000 reads), 10% (-26,000) of reads while retaining the relative abundances of per-cell coverage, (a) At each of the downsampling depths, we repeated the Seurat workflow of normalization, clustering and dimensionality reduction using UMAP. (b) We tracked the assigned cluster identity of every cell in the dataset over the downsampling depths and visualized the flow of clustered cells by connecting cluster labels (nodes) with line widths scaled to the number of cells.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
3 Smart-seq3xpress (downsampled - sequence depth matched}
Clusters
-5 0 5 10 15 UMAP_1
Annotation subset
-5 0 5 10 15 UMAP_1
Azimuth
-5 0 5 10 15 UMAP_1
|^ 10x Genomics (downsampled - sequence depth matched)
a. 0 <
1 -5
-10 -5 0 5 10
UMAPJ
10 5
CM
CL 0
Annotation subset ™..Sk
TEU Naive CO'
8* TEM TCM/TEM
10
5
eg
i ° 1 -5
Azimuth
CDS TEH
—10 -5 0 5 10
UMAP_1
-10 -5 0 UMAPJ
C Smart-seq3xpress (all)
Clusters
10H
0 10 UMAP 1
Annotation
atmcDa+ tew :D4+TCMITEM_2
0 10 UMAP 1
Azimuth
0 10 UMAP 1
d 10x Genomics (all)
Clusters
10
-15 -10 -5 0 5
UMAP_1
COG* TCM/TEM Naiva CD4+ TÍCOJ+ TCM/TEM
Uneyj T Naive
-5 0 5
UMAP_1
y gdT NK_CD56tnighl
-15 -10 -5 0 5
UMAP_1
predicted.celltype.l2
• ASDC
B intermediate
• B memory
• B naive CDUMono
• CD16Mono CD4 CTL
• CD4 Naive CD4 Proliferating
• CD4TCM CD4TEM CD8 Naive
• CD8 Proliferating CD8TCM
• CD8TEM
UMAP_1
Extended Data Fig. 9 | See next page for caption.
• CDC1
• CDC2
• dnT
• Erytri gdT HSPC
• ILC
MAIT NK
• NK Proliferating
• NK_CD56brigrit pDC
Plasmablast
• Platelet Treg
predicted.celltype.l3
• ASDC_pDC CDS Naive • gdT 2
• B intermediate kappa • CDS Naive_2 gdT_3
• B intermediate lambda CD8 Proliferating gdT 4
B memory kappa • CD8TCM 1 • HSPC
• B memory lambda • CD8TCM 2 • ILC
• B naive kappa CD8TCM 3 MAIT
• B naive lambda CDS TEM_1 NK Proliferating
• CD14 Mono • CDSTEM 2 • NK 1
• CD16 Mono CD8 TEM_3 • NK_2
CD4 CTL CD8 TEM_4 • NK 3
CD4 Naive • CD8TEM 5 • NK 4
• CD4 Proliferating • CDSTEM 6 NK CD56bright
CD4TCM 1 • cDC1 • pDC
• CD4TCM 2 • cDC2 1 • Plasma
CD4TCM_3 cDC2_2 • Plasmablast
• CD4TEMJ dnT_2 • Platelet
CD4TEM 2 • Eryth • Treg Memory
• CD4TEM_3 • gdT 1 • Treg Naive
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
ARTICLES
NATURE BIOTECHNOLOGY
Extended Data Fig. 9 | Direct comparison of Smart-seq3xpress and lOx Genomics 3' v3.1 data. Human PBMCs from the same donor were processed and sequenced using both protocols. We sampled the same number of sequenced reads (approx. median 140,000 per cell) and cells (3,000) and applied an identical data analysis workflow using Seurat (see Methods), (a) Shown are the UMAP embeddings for Smartseq3xpress identifying the number of clusters (left), annotated cluster names (middle) and cluster identities by reference-based identification via Azimuth cell type predictions (level 2 annotations) (b) 10xv3.1 UMAP embeddings showing clusters found (left), annotated cluster names (middle) and reference based predicted annotations by Azimuth (prediction level 2). (c-d) Shows the full the donor matched data without downsampling or cell subsetting after QC for both Smartseq3xpress (n = 3187) and 10xv3.1 (n = 6483). (left) clusters detected, (middle) annotated cell types based on clusters, (right) reference-based azimuth cell type predictions (prediction level 2). (e-f) Reference-based annotation of all human PBMC generated in this study using Smartseq3xpress by Azimuth (Hao et al., 2020) using the human PBMC reference dataset available at httpsyyazimuth.hubmapconsortium.org. (e) Smartseq3xpress UMAP is colored and annotated by Azimuth cell type annotation level 2. (f) Smartseq3xpress UMAP is colored and annotated by Azimuth cell type annotation level 3.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
NATURE BIOTECHNOLOGY
ARTICLES
CD27-AS1
ETHE1
POP5
._ -0.8 i
"O
s_
~ P-0.6 q- I
w
H-'
c cu Ü
1—
03 CL
PTPRC
15 10
Q.
<
I 5
-5 -10
o- C
15
10-
C\J
5
<
0
-5-
-10-
C\J Q.
<
15 10
5
0--5 -10
-10
-10
0 10 UMAP_1
ISG20
0 10 UMAP_1
GUSBP11
-10
0 10 UMAP_1
g '(/> in
0J
±-
Cl X (D (D C CD O
■
CM
Q.
<
15 10 5 0 -5 -10
Q.
<
15 10-
5
0 -5 -10
c\j
Q.
<
15 10 5 0 -5 -10
POP5
Q. X (D
(D C CD O
-10
0
UMAP_1 CD27-AS1
10
-10
0
UMAP_1 ETHE1
10
9
-10
0
UMAP_1
10
■
Q. X CD
CD
CD O
g en
CD
Q. X CD
CD != CD O
1
Extended Data Fig. 10 | PBMC splicing analysis using Smart-seq3xpress. (a) Percent-spliced-in (PSI) for skipped exon events in the CD27-AS1, ETHE1 and POP5 genes are overlayed as color scale on the UMAP embedding of sequenced PBMC cells (n = 26,260). (b) For genes with significant cell-type specific splicing patterns (shown in Figures. 3j, Supplementary Fig. 18a), scaled expression levels (log-normalized read counts) are annotated as colors.
NATURE BIOTECHNOLOGY | www.nature.com/naturebiotechnology
nature research
Reporting Summary
Corresponding author(s): Rickard Sandberg
Last updated by author(s): Mar 23, 2022
Nature Research wishes to improve the reproducibility of the work that we publish. This form provides structure for consistency and transparency in reporting. For further information on Nature Research policies, see our Editorial Policies and the Editorial Policy Checklist.
Statistics
For all statistical analyses, confirm that the following items are present in the figure legend, table legend, main text, or Methods section.
n/a
□
□
□
□
□
□
□
□
□
Confirmed
^1 The exact sample size (n) for each experimental group/condition, given as a discrete number and unit of measurement
^1 A statement on whether measurements were taken from distinct samples or whether the same sample was measured repeatedly
K-pi The statistical test(s) used AND whether they are one-or two-sided
^ Only common tests should be described solely by name; describe more complex techniques in the Methods section.
[5^1 A description of all covariates tested
[5^1 A description of any assumptions or corrections, such as tests of normality and adjustment for multiple comparisons
prpi A full description of the statistical parameters including central tendency (e.g. means) or other basic estimates {e.g. regression coefficient)
^ AND variation (e.g. standard deviation) or associated estimates of uncertainty (e.g. confidence intervals)
prpi For null hypothesis testing, the test statistic (e.g. F, t, r) with confidence intervals, effect sizes, degrees of freedom and P value noted
^ Give P values as exact values whenever suitable.
^ For Bayesian analysis, information on the choice of priors and Markov chain Monte Carlo settings
~] For hierarchical and complex designs, identification of the appropriate level for tests and full reporting of outcomes
[5^1 Estimates of effect sizes (e.g. Cohen's d, Pearson's r), indicating how they were calculated
Our web collection on statistics for biologists contains articles on many of the points above.
Software and code
Policy information about availability of computer code
Data collection The following software was used in data collection: FACS sorter software (BD FACSChorus 1.3), sequencer software (lllumina NextSeq Control Software 2.2.0 and MGI DNBSEQG400RS 1.1.0.108).
Data analysis
[Analysis was performed using zUMIsv2.8.2 or zUMIs v2.9.3f, STAR v2.7.3a, Seurat v4.0.1, presto (vl.0.0), Cellbender, Solo v0.6, TRaCeR vO.6.0, Rsamtoolsv2.6.0, UMIcountR vO.1.1, cellSNP-lite vl.0.0, BRIEkit vO.2.2, BRIE2 V2.0.6., CellBenderv0.2.0,BSgenome vl.62.0, Scanpy vl.8.2, Scirpy vO.8.0, Gviz vl.38.lL, rtracklayer vl.54.0, Further details are listed in the Methods section.
zUMIs is available at https://github.com/sdparekh/zUMIs. Code to filter reads with the TSO strand invasion artifact is implemented in a standalone script pyTSOfilter (https://github.com/cziegenhain/pyTSOfilter).
For manuscripts utilizing custom algorithms or software that are central to the research but not yet described in published literature, software must be made available to editors and reviewers. We strongly encourage code deposition in a community repository (e.g. GitHub). See the Nature Research guidelines for submitting code & software for further information.
Data
Policy information about availability of data
All manuscripts must include a data availability statement. This statement should provide the following information, where applicable: -Accession codes, unique identifiers, or web links for publicly available datasets -A list of figures that have associated raw data -A description of any restrictions on data availability
Sequencing data have been deposited at ArrayExpress, European Bioinformatics Institute, under the following accession numbers; E-MTAB-11488, E-MTAB-11452, I E-MTAB-11467. Human genome build hg38 fasta files and gene annotation in GTF format (Grch38.95) were obtained from Ensembl.
Field-specific reporting
Please select the one below that is the best fit for your research. If you are not sure, read the appropriate sections before making your selection. Life sciences ] Behavioural & social sciences ] Ecological, evolutionary & environmental sciences
For a reference copy of the document with all sections, see nature.com/documents/nr-reporting-summary-flat.pdf
Life sciences study design
All studies must disclose on these points even when the disclosure is negative,
Sample size Sample sizes were not predetermined using statistical analysis. They were determined based on allowable size within a reasonable budget cost for preparing libraries and sequencing
(
Data exclusions Single-cell RNA-seq data were filtered according to established criteria. Cutoffs are listed where appropriate. Further data exclusions were not performed.
Replication
For each experimental condition, a large number of single cells was sequenced to ensure reprodocibility. Sample sizes are clearly indicated throughout.
Randomization Not relevant because FACS sorting of individual cells into random wells of microplates. Blinding (investigators were not blinded to groups of samples as it was not practically feasible.
Reporting for specific materials, systems and methods_
We require information from authors about some types of materials, experimental systems and methods used in many studies. Here, indicate whether each material, system or method listed is relevant to your study. If you are not sure if a list item applies to your research, read the appropriate section before selecting a response.
Materials & experimental systems Methods
n/a Involved in the study n/a Involved in the study
□ £3 Antibodies □ ChlP-seq
□ E^l Eukaryoticcelllines ~\ Flow cytometry
] Palaeontology and archaeology ] MRI-based neuroimaging
~~\ Animals and other organisms
□ E^l Human research participants
~\ Clinical data
] Dual use research of concern
Antibodies
Antibodies used
Validation
. All antibodies used are from BD Bioscience (Pharmingen & Horizon). PE Mouse Anti-Human CCR7 (Cat.no: 566742, Clone: 2-L1-A, Lot:1006133), PE-Cy7 Mouse Anti-Human CD4 (Cat. no: 557852 Clone:SK3, Lot: 3060697), FITC Mouse Anti-Human CD45RA (Cat. no: 561882, Clone:HllOO, Lot: 9301282), PerCP-Cy5.5 / BB700 Mouse Anti-Human CD8 (Cat. no: 566451, Clone:RPA-T8, Lot: 1011278), PE-Cy5 Mouse Anti-Human CD45RO (Cat.no: 561888, Clone:UCHLl, Lot: 9101742)
All antibodies are validated for reactivity against Human, and applicable to FACS sorting according to manufacturer BD Bioscience.
Eukaryotic cell lines
Policy information about cell lines Cell line source(s)
Authentication
Mycoplasma contamination
HEK293FT: Thermo Fisher, K562: DSMZ (Braunschweig, Germany)
HEK293FTand K562 were authenticated by PCR-single-locus-technology (Eurofins Forensik)
HEKF293FT and K562 were confirmed free of mycoplasma contamination (Eurofins)
Commonly misidentified lines HEK293FTand K562 were used but authenticity was confirmed, both are not part of common misidentified cell lines (see
(See ICLAC register)
I above).
2
Human research participants
Policy information about studies involving human research participants
Population characteristics
Recruitment
Ethics oversight
No specific population selections or criterias were used, other than overall healthy donors. The 7 healthy donors, include 4 I male, 2 female, 1 N/A. Age range 21-47
Frozen Aliquots of cryopreserved human PBMCs were bought from Lonza.
Etikprovningsmyndigheten (Sweden) Dnr 2020-05070
Note that full information on the approval of the study protocol must also be provided in the manuscript.