Single‐cell RNA sequencing of cells from fresh or frozen tissue reveals a signature of freezing marked by heightened stress and activation

Thawing of viably frozen human tissue T cells, ILCs, and NK cells and subsequent single-cell RNA sequencing reveals that recovery of cellular subclusters is variably impacted. While freeze-thawing does not alter the transcriptional profiles of cells, it upregulates genes and gene pathways associated with stress and activation.


Technical Report
Single-cell RNA sequencing of cells from fresh or frozen tissue reveals a signature of freezing marked by heightened stress and activation Single-cell RNA sequencing (scRNAseq) is increasingly used to decipher immune cell functions within tissues.Sequencing human tissue lymphocytes is often performed using frozen samples because of logistical constraints or to allow for the batch processing of multiple samples.While published work has shown freezethawing to impact bulk RNA sequencing, scRNAseq seems less impacted by it [1][2][3][4][5].However, an outstanding methodological question is if there is an impact on scRNAseq of tissue-derived immune cells.We investigated this using immune cells isolated from tonsils, which contain both mucosal and lymphoid structures, and found the data between the two are largely comparable but cells are in a heightened state of activation and stress after freeze-thawing.We generated singlecell suspensions from tonsils and immediately sorted and sequenced T cells, natural killer (NK) cells, and innate lymphoid cells (ILCs) while also freezing half of the cell suspension in FBS/10% DMSO.The frozen cells were subsequently thawed, sorted, and sequenced (Supporting information Fig. S1A).After sequencing, low-quality cells were filtered using standard quality control (QC) metrics and remaining cells were randomly downsampled to match the tonsil and condition with the lowest cell number (Supporting information Fig. S1B).For both tonsils, more cells were initially recovered from freeze-thawed samples, perhaps because fragile cells will die during thawing, leaving a more robust population for subsequent sorting and sequencing.Quality control revealed that the largest difference after freeze-thawing was in the number of unique genes detected per cell, which was higher in cells from fresh tissues both pre-and post-QC filtering (Fig. 1A-D, Supporting information Fig. S1C-F).The data were integrated using Harmony to reduce donor-but not condition-dependent effects.Clusters were annotated using differentially expressed (DE) genes (Fig. 1E-F) [6][7][8][9][10].The top DE genes per cluster were largely unaffected by freeze-thawing, showing only minor differences in magnitude and not in patterns of expression (Fig. 1E).However, we found the overall composition of cell clusters to be variably affected by freeze-thawing (Fig. 1G).We confirmed that neither tonsil had an outsized impact on clustering and by projecting integrated clusters onto the non-integrated UMAP (Uniform Manifold Approximation and Projection), we confirmed that integration had no major impacts on DE gene identification between conditions (Supporting information Fig. S2).To further investigate compositional differences, we performed subcluster analysis of ILC1/NK cells, T cells, and ILC2/ILC3/ILCp and found γδ T, T C , innate-like T cells, and NKp44 + ILC3 to be the most reduced during the freeze-thawing process (Supporting information Fig. S3).
To investigate how cell states differ between conditions in a cluster agnostic fashion, we used differential abundance (DA) testing [11] which revealed 7 DA populations (Fig. 1H), of which some derived disproportionally from the two tonsils (Supporting information Fig. S4).However, all of them consisted of only a small fraction of total cells in each cluster (Fig. 1I).We next investigated the transcriptional profiles of the DA populations using DE analysis comparing each DA population to the rest of the cells in the clusters from which they derive.DA populations enriched after freeze-thawing were largely marked by increased expression of genes associated with activation and stress (Table S1).
To further explore differences in gene expression after freeze-thawing, we investigated a curated list of transcription factors, cell surface markers, lineage defining genes, and a randomly selected subset of the lowliest expressed genes as they are likely the most susceptible to differences caused by the degradation and loss of mRNA.We found only minor differences in the magnitude of expression of some genes (Fig. 2A and Supporting information Fig. S5A).To better understand the observed increase in activation markers in differential abundance testing and which genes were most impacted by freeze-thawing, we performed DE analysis comparing the two conditions on a per cluster basis to account for the skewed compositions.Note that 215 genes across eight of the 10 cell clusters were differentially expressed excluding ribosomal, mitochondrial, and non-protein coding genes (Fig. 2B and Supporting information Table S2).All differences were smaller than onefold (Supporting information Table S2 and Fig. S5B-C) and freeze-thawing downregulated more genes (n = 145) than it upregulated (n = 70) (Fig. 2B).Several DE genes were shared across clusters (Fig. 2C and Supporting information Table S2).These shared genes showed no overall trends of enrichment toward ILC/NK cells or T cells (Fig. 2C).Pathway and gene set enrichment analysis of DE genes revealed an increase of pathways associated with activation and cellular stress (Fig 2D and E).Using a previously published gene set associated with stress [12], we assigned each cell a stress score and confirmed the freeze-thawed  cells score higher, confirming that freezethawing before scRNAseq places the cells in a heightened state of stress and activation (Fig. 2F and Supporting information Fig. S5D).
Single-cell sequencing technologies are increasingly improving while becoming more readily available, simpler, and cheaper to perform.Therefore, we must also be fully cognizant of how they can be applied and their limitations to avoid erroneous conclusions.An outstanding technical question is the impact of freezethawing on cells derived from human tissue.We found that freeze-thawing had a small impact on cellular composition and sequencing library quality and that gene expression differences were minor.Further, the changes were limited to the magnitude of expression while overall patterns of cell-defining genes remain unchanged.Many of the genes upregulated by freezethawing are markers of stress and activation, similar to reported increases after enzymatic digestion [12,13], suggesting that this transcriptional signature can also be derived from or amplified by freezethawing.Therefore, while there are no compelling arguments to avoid the use of freeze-thawed ILCs, NK cells, or T cells for scRNAseq studies, there are important aspects to keep in mind, including the altered recovery of cellular populations and a heightened state of transcriptional activation and stress.Validation of these findings in lymphocytes from tissues other than tonsils is now eagerly awaited.

Figure 1 .
Figure 1.Quality control, cluster composition, and differential abundance testing.(A) Violin plot showing amount of RNA per cell in fresh or freeze-thawed cells after QC. (B) Violin plot showing unique genes per cell in fresh or freeze-thawed cells after QC. (C) Violin plot showing percent mitochondrial genes per cell in fresh or freezethawed cells after QC. (D) Violin plot showing percent ribosomal genes per cell in fresh or freeze-thawed cells after QC. (E) Heatmap showing the average expression in freeze-thawed cells, left half of the columns, and fresh cells, right half of the columns, of the top 10 DE genes per cluster.(F) UMAP with cell types annotated by DE genes.(G) Bar plot shows the cluster composition of the clusters from F in fresh or freeze-thawed cells, and numbers show the number of cells recovered from fresh or freeze-thawed tonsils.(H) UMAP with cell types annotated by differentially abundant (DAseq) cell populations between conditions.Gray cells are similar between conditions, cells colored in shades of blue are more abundant in the freeze-thawed set, and cells colored in shades of red are more abundant in the fresh data set.(I) Bar plot showing the cluster composition of the clusters from (F) by differentially abundant cell population.Gray represents the fraction of each cluster that are similar between conditions, shades of blue represent cell states enriched in the freeze-thawed, and shades of red cell states enriched in the fresh, data set.The total percentage of differentially abundant cell populations is marked on the bars.Statistics: Black bars in (A-D) show the geometric mean.

Figure 2 .
Figure 2. Gene expression differences by cluster and pathway and enrichment analysis.(A) Heatmap showing the average expression in freeze-thawed cells, left half of the columns, and fresh cells, right half of the columns, of curated lists of genes made up of lineage markers, transcription factors, or surface markers.(B) Bar plot showing the number of genes up-or downregulated by freeze-thawing in each cluster.(C) Schematic plot and table showing genes commonly up-or downregulated by freeze-thawing in each cluster.(D) KEGG pathway analysis heatmap generated using all genes that are differentially expressed between preparation condition.Pathways are shown on the left-hand side and which genes are included in the pathways are along the bottom.(E) Lollipop plot showing the top five gene sets found enriched in genes DE between conditions in gene set enrichment analysis using the MSigDB Hallmark gene sets.Gene sets which are enriched are listed on the left.X axis is the normalized enrichment score.Circle size indicates the number of genes from the gene set which are found differentially expressed between conditions and the color indicates if the enrichment is significant.(F) Violin plots split by condition, fresh or freeze-thawed, showing the stress scores for all cells in each cluster.Statistics: Black bars in (F) show the geometric mean.