Single‐cell transcriptomics provides insights into the origin and immune microenvironment of cervical precancerous lesions
Authors: Chunbo Li, Keqin Hua
Journal: Cancer Communications (2023)

Dear editors, Despite efforts to implement vaccination and screening programs, cervical cancer (CC) remains a public health problem [1]. Over 90% of cases are associated with persistent human papillomavirus (HPV) infections. HPV integrates its DNA into the basal cells of the transformation zone, the region where the epithelium transitions from the columnar epithelial cells of the endocervix to the stratified squamous epithelium cells of the ectocervix, leading to the production of proteins (E6 and E7) that eventually cause dysplasia [2, 3]. Squamous cell carcinoma (SCC) and cervical adenocarcinoma (AD) are the two most common types of HPV-related CC [4]. SCC and AD have a long-term stage of preinvasive (PRE) lesion, which are named high-grade squamous intraepithelial lesion (HSIL) or adenocarcinoma in situ (AIS), respectively [5]. However, no studies were conducted to determine why some individuals with HPV infections develop HSIL while others develop AIS. To reveal the different characteristics of normal cervical progression to HSIL or AIS, we performed single-cell RNA-sequencing (scRNA-seq) on 8 samples, including 2 HPV-positive (HPV-P) normal cervical, 2 HPV-negative (HPV-N) normal cervical, 2 HSIL (PRE-HSIL) and 2 AIS (PRE-AIS) samples (Supplementary Table S1). The study protocols can be found in the Supplementary Materials. A total of 80,238 cells were obtained, with an average gene number of 2,260. Then, we acquired 33 clusters (Supplementary Figure S1A), which could be divided into 10 major cell types through the expression of typical marker genes (Supplementary Figure S1B-C). Based on the comparison of cell numbers, we found that, compared to HPV-N, HPV-P led to an increase in immune cells, such as B cells, natural killer (NK) cells/T cells, and neutrophils (Supplementary Figure S1D). When the lesion progressed to a preinvasive lesion, the prominent characteristic was the increase of epithelial cells in both AIS and HSIL. To determine the different cellular origins of AIS and HSIL, we reclustered the extracted 28,364 epithelial cells and acquired 13 clusters. The split Uniform Manifold Approximation and Projection (UMAP) map revealed that most epithelial cells in cluster 3 were from AIS, whereas those in clusters 4 and 5 were from HSIL (Figure 1A-B). Surprisingly, epithelial cells in clusters 3, 4 and 5 were found in lower numbers in the HPV-N or HPV-P samples, indicating that these cells represented the specific subpopulations of AIS and HSIL, respectively. Gene Ontology (GO) analysis revealed that cells in cluster 3 expressed high levels of cell cycle-related signaling pathways, such as the mitotic cell cycle and chromosome condensation (Supplementary Figure S2). Cells in cluster 4 were related to epithelial cell differentiation and actin cytoskeleton organization, whereas cells in cluster 5 exhibited humoral immune response, negative regulation of endopeptidase activity, inflammatory response and regulation of apoptotic signaling pathway (Supplementary Figure S2). Single-cell transcriptomics identifies the origin and immune microenvironment of cervical precancerous lesions. (A) UMAP showing that the number of epithelial cells in cluster 3 was higher in AIS and epithelial cells in cluster 4 and 5 were higher in HSIL compared to HPV-P and HPV-N group. Each dot represents a single cell, color-coded by cell clusters. X- (UMAP 1) and Y-axis (UMAP 2) represented the location of epithelial cells in the four groups after dimension reduction and unsupervised clustering (B) Bar graph showing the proportion of the 13 epithelial cell clusters quantified in four groups. (C) 2D graph of the pseudotime-ordered epithelial cells from HPV-N, HPV-P, PRE-AIS and PRE-HSIL samples, respectively. In the top panel, X-axis indicates the pseudotime, and Y-axis indicates the cell density distribution along the trajectory. In the bottom panel, X- (Component 1) and Y-axis (Component 2) label represented the location of epithelial cells on the pseudotime trajectories in the four groups. In the right panel, the branch trajectory plot with dark blue to light blue colors indicated the pseudotime of single cells. (D) Heatmap showing the dynamic changes in gene expression along the pseudotime. (E) Bubble plot showing the expression of markers genes in epithelial cells in the five cell types. The fraction of cells expressing genes was indicated by the size of the circle, and their scaled expression levels were indicated by the color of the circle. (F) UMAP of five types of epithelial cells. X- (UMAP 1) and Y-axis (UMAP 2) represented the location of five epithelial cell types. (G) Heatmap showing the enriched hallmark gene signatures of epithelial cells in five epithelial cell types. (H) The proportion of the five epithelial cell types in the four groups. The green arrow indicated the proportion of undefined Epi (I) UMAP of 13 NK/T cell clusters. X- (UMAP 1) and Y-axis (UMAP 2) represented the location of NK/T cells. (J) The proportion of 13 NK/T cell types in the four groups. (K) Boxplots comparing cytotoxicity, exhausted, naïve, proliferative, regulatory scores for NK/T cells between HPV-N, HPV-P, PRE-AIS and PRE-HSIL. (***P < 0.001; n = 10, 512). (L) Bubble plot showing the expression of the representative genes of NK/T cells in the four groups. The fraction of cells expressing genes was indicated by the size of the circle, and their scaled expression levels were indicated by the color of the circle. (M) Heatmap depicting the enriched hallmark gene signatures of NK/T cells in the four groups. (N) Schematic illustration of the development of human cervical AIS and HSIL. Abbreviations: UMAP, uniform manifold approximation and projection; AIS, adenocarcinoma in situ; HSIL, high-grade squamous intraepithelial lesion; HPV-P, human papillomavirus-positive normal cervix; HPV-N, human papillomavirus-negative normal cervix; 2D, 2-ddemension; PRE, preinvasive; NK, natural killer; To reveal the functional changes during lesion progression, we projected all epithelial cells onto pseudotime trajectories (Figure 1C). This unsupervised approach identified continuous cell states and formed two distinct trajectories, starting from state 1 and gradually progressing toward states 2 and 3, revealing a common origin with divergent fates (Figure 1C and Supplementary Figure S3). State 1 mainly comprised epithelial cells from HPV-N and HPV-P normal cervix, which were then separated into two branches. One branch described the developmental pathway of cells derived mainly from HSIL (cell fate 1, state 2), whereas the other branch was occupied by cells from AIS (cell fate 2, state 3). The branched heatmap showed that significant genes could be divided into four different clusters according to the gene expression dynamics (Figure 1D). GO analyses revealed that the highly expressed genes (gene set 2) in cell fate 1 were related to keratinocyte differentiation, epidermis development, and epidermal cell differentiation (Pathway 2). Meanwhile, immune regulation-related pathways, such as leukocyte-mediated immunity, complement activation and phagocytosis, were up-regulated from the normal cervix to the HSIL (Pathway 3). The highly expressed genes (gene set 1) in cell fate 2 had high enrichment of response to ketone, cellular response to metal ion and RNA splicing (Pathway 1). Together, these results demonstrated the different functional characteristics of epithelial cells between HSIL and AIS. Based on top differentially expressed genes (DEGs) and functional characteristics, we classified all epithelial cells (Epi) into five types: proliferative Epi (C3), differentiated Epi (C4, 6, and 10), inflammatory Epi (C2 and 7), metabolic Epi (C0, 1, 5, 8, 9, and 11) and undefined Epi (C12) (Figure 1E-F). Gene Set Variation Analysis (GSVA) confirmed their specific functional characteristics (Figure 1G). We compared the percentages of the five epithelial types in different groups and found that AIS had the highest percentage of proliferative Epi cells, while HSIL had the highest percentage of differentiated Epi cells (Figure 1H). Many studies have reported that epithelial-mesenchymal transition (EMT) is closely related to SCC progression and metastasis [6]. Herein, we found that differentiated Epi cells had a higher EMT score than other cell types (Supplementary Figure S4A-B). Similarly, cells in cluster 4 had the highest EMT score compared with other epithelial cell clusters (Supplementary Figure S4C), indicating that EMT played a key role in regulating the progression of lesions from HPV infection to HSIL. In our study, we identified 13 clusters (Supplementary Figure S5A-B), including 6 CD8+ T cell clusters (C0, 2, 3, 4, 10, and 12), 3 CD4+ T cell clusters (C1, 6, and 8), and 4 NK cell clusters (C5, 7, 9, and 11) (Figure 1I and Supplementary Figure S5C-E). Different groups had different percentages of 13 clusters (Figure 1J). Our results showed that, compared to the HPV-N cervix, the HPV-P cervix had a higher infiltration of cytotoxicity CD8+ T cells and effector memory CD8+ T cells but a lower infiltration of naïve CD4+ T cells. Similarly, HPV-P had higher cytotoxic and exhausted scores but lower naïve scores than HPV-N (Figure 1K), indicating that HPV infection could trigger an immune response. Compared with HSIL, AIS was associated with higher percentages of cytotoxic CD8+ T cells and exhausted CD8+ T cells, but a lower percentage of NK cells, indicating that both types were dominated by different effector cells. HSIL had the highest cytotoxic and exhausted scores compared to other groups (Figure 1K), indicating that HSIL had the best response rate to immunotherapy. More importantly, AIS and HSIL strongly expressed various immune-related genes and proliferative signaling pathways (Figure 1L-M), indicating that abnormal proliferation of epithelial cells contributed to the complexity of the microenvironment. Then, we performed a pseudotime analysis to understand the underlying evolution of the cellular status of CD8+ T, CD4+ T and NK cells (Supplementary Figure S6A-C). Interestingly, by combining the clustering and pseudotime analyses findings, we observed a gradual transition of CD8+ T and NK cells towards the subpopulations with immune enhancement from HPV infection to HSIL or AIS. Such an activation status for immune cells was indicated by the upregulation of GNLY, GZMB, FCER1G, and NKG7 and the downregulation of IL7R, CD69, GPR183 and NR4A1 in HSIL and AIS (Supplementary Figure S6D-E). More importantly, HSIL and AIS had similar immune activation states, which were characterized by the high infiltration of proliferative NK cells, regulatory CD4+ T cells, cytotoxic CD8+ T cells and exhausted CD8+ T cells. However, compared with AIS, HSIL was associated with higher infiltration of GNLY+ NK cells, naïve NK cells, and effector memory CD8+ T cells (Supplementary Figure S7). Our findings suggested that AIS might originate from proliferative epithelial cells, whereas HSIL might originate from differentiated epithelial cells. Further, HPV infection could affect the immune microenvironment, and HSIL and AIS exhibited different characteristics of immune activation (Figure 1N). Collectively, the results of this study provide new insights into the progression of normal cervix to the AIS and HSIL stages after HPV infection. Chunbo Li and Keqin Hua designed this study; Chunbo Li collected the samples, performed the experiments, analyzed the data, and drafted the manuscript. All authors participated in the writing and have read and approved the final manuscript. The authors would like to thank the patients, data managers, and all study investigators for their contributions to the study. The authors declare that they have no competing interests. This study was supported by the National Natural Science Foundation of China (Grant No. 82173188) to Keqin Hua, the Clinical Research Plan of Shanghai Hospital Development Central (SHDC) (SHDC2020CR1045B, SHDC2020CR6009) to Keqin Hua, the Shanghai Municipal Health Commission (20194Y0085) to Chunbo Li, and the Shanghai “Rising Stars of Medical Talent” Youth Development Program (SHWSRS2020087) to Chunbo Li. Not applicable This study was approved by the Ethical Committee of our hospital (2022-143). The patient provided written informed consent. All data in the study are available from the corresponding author upon reasonable request. Please note: The publisher is not responsible for the content or functionality of any supporting information supplied by the authors. Any queries (other than missing content) should be directed to the corresponding author for the article.