Single-Cell RNA Sequencing Reveals Multiple Pathways and the Tumor Microenvironment Could Lead to Chemotherapy Resistance in Cervical Cancer

Background Cervical cancer is one of the most common gynecological cancers worldwide. The tumor microenvironment significantly influences the therapeutic response and clinical outcome. However, the complex tumor microenvironment of cervical cancer and the molecular mechanisms underlying chemotherapy resistance are not well studied. This study aimed to comprehensively analyze cells from pretreated and chemoresistant cervical cancer tissues to generate a molecular census of cell populations. Methods Biopsy tissues collected from patients with cervical squamous cell carcinoma, cervical adenocarcinoma, and chronic cervicitis were subjected to single-cell RNA sequencing using the 10× Genomics platform. Unsupervised clustering analysis of cells was performed to identify the main cell types, and important cell clusters were reclustered into subpopulations. Gene expression profiles and functional enrichment analysis were used to explore gene expression and functional differences between cell subpopulations in cervicitis and cervical cancer samples and between chemoresistant and chemosensitive samples. Results A total of 24,371 cells were clustered into nine separate cell types, including immune and non-immune cells. Differentially expressed genes between chemoresistant and chemosensitive patients enriched in the phosphoinositide 3-kinase (PI3K)/AKT pathway were involved in tumor development, progression, and apoptosis, which might lead to chemotherapy resistance. Conclusions Our study provides a comprehensive overview of the cancer microenvironment landscape and characterizes its gene expression and functional difference in chemotherapy resistance. Consequently, our study deepens the insights into cervical cancer biology through the identification of gene markers for diagnosis, prognosis, and therapy.


INTRODUCTION
As the fourth most common gynecological malignant tumor, cervical cancer is a leading cause of cancer-related deaths among women and poses a serious threat to the health of women worldwide (1). In 2018, approximately 570,000 new cases of cervical cancer and 311,000 deaths from this cancer were reported (1). Paclitaxel, cisplatin, carboplatin, or a combination of these agents is the front-line treatment for cervical cancer (2). However, the efficacy of current chemotherapeutic agents is limited, with relatively low response rates of 29%-63% because of chemotherapy resistance (3). The combination of paclitaxel and cisplatin is one of the most commonly utilized regimens in the metastatic disease setting (4). However, in actual clinical treatment, tumor cells often develop resistance.
The tumor microenvironment (TME) comprises various cell types [such as fibroblasts, endothelial cells (ECs), and immune cells] and extracellular components (such as cytokines, hormones, extracellular matrix, and growth factors), which surround tumor cells as a vascular network. The TME not only plays a pivotal role during tumor initiation, progression, and metastasis, but it also has profound effects on therapeutic efficacy. TME-mediated chemotherapy resistance is a result of complex crosstalk between tumor cells and their surrounding environment (5,6).
For example, the TME and therapeutic response can be induced by soluble factors secreted by tumors. The adhesion of tumor cells to stromal fibroblasts can also affect responses to chemotherapy (7), and immune cells also play an important role in improving and obstructing therapeutic efficacy (7). The interaction between chemotherapy sensitivity and TME is a complex phenomenon. Cancers can develop remarkable resistance to various treatments that target different molecular pathways (8).
Research studies have shown that the cellular and molecular mechanisms underlying the development of resistance are multifactorial and include genetic and epigenetic alterations, cell detoxification, and abnormal drug efflux and accumulation (9,10). However, the molecular mechanisms underlying the occurrence and development of resistance are poorly understood. Thus, there is an urgent need to identify the basic factors that determine chemotherapy resistance in cancer.
Previous studies on molecular mechanisms and chemotherapy resistance in cervical cancer patients have focused mainly on bulk genomic or transcriptome profiling methods and in-situ hybridization techniques (11)(12)(13). Consequently, studies of chemoresistance mechanisms based on the signatures of distinct cell populations are obscure (10)(11)(12). Single-cell RNA sequencing (scRNA-seq) techniques have emerged as a powerful tool for analyzing the spectrum of cell populations in tissues. Furthermore, these techniques have been widely used for elucidating the complex subpopulations in tissues of organs such as the lung (14), heart (15), and brain (16), as well as in various cancers including melanoma (17), ovarian (18), and colon (19) cancers (20,21).
Although studies have used scRNA-seq on the cervix uteri (22) and drug-resistant cell lines (12), to the best of our knowledge, scRNA-seq profiling of human cervical cancer tissues has not been reported to date. Although scRNA-seq is increasingly being adopted, its application to tumors has been limited to several types, but not cervical cancer. Because of this limited elucidation of human tumors and the lack of TME profiling, the intratumoral transcriptomic heterogeneity of the most common cancer in women is largely unknown.
Exploring the molecular mechanism of chemotherapy resistance is important in the development of strategies to overcome tumor resistance and provides a theoretical basis for reversing tumor resistance and improving cancer chemotherapy efficacy. The development of chemotherapy combinations that could simultaneously target tumor cells and the microenvironment may represent a solution to overcome chemotherapy resistance. This study aimed to analyze cells from pretreated and chemoresistant cervical cancer tissues at a much higher scale, to generate a molecular census of cell populations.
Furthermore, we sought to uncover the cell heterogeneity using unbiased scRNA-seq techniques. Consequently, we performed an scRNA-seq survey of 23,444 cells from five tissues from pretreated cervical cancer patients and constructed a single-cell transcriptome atlas for early malignancy of cervical cancer ( Figure 1). Our study provides novel insights into the heterogeneity of cervical cancer at the single-cell level and will serve as a valuable resource for understanding chemotherapy resistance mechanisms in tumor progression.

Patients and Tumor Specimens
Four female patients with a pathologic diagnosis of cervical cancer and one female patient diagnosed with cervicitis were enrolled at Renmin Hospital of Wuhan University, Wuhan ( Table 1). All enrolled patients signed the written consent, and this study was approved by the Institutional Review Board (IRB) of Renmin Hospital of Wuhan University (IRB no. WDRY2021-K014). Fresh tumor samples (at least 1.5 cm 3 ) were surgically resected from all enrolled patients. None of the patients were treated with chemotherapy prior to tissue sample collection. After sample resection, three of the four cervical cancer patients were treated with chemotherapy, and one patient was a responder and the two were resistant to chemotherapy.

Human Cervical Cancer Tissue Cell Dissociation
All fresh cervical cancer tissues were transferred from the operating room to the dissociation laboratory in cold Hank's balanced salt solution (HBSS) medium supplemented with 1% penicillin-streptomycin within 30 min. Samples were gently washed in phosphate-buffered saline (PBS) after removing the adipose tissue and minced into pieces of approximately 1 mm 3 using an Iris scissors. Enzymatic digestion was performed using the MACS tumor dissociation kit (Miltenyi Biotec) according to FIGURE 1 | Single-cell transcriptome regulatory network of the tumor microenvironment (TME) and chemoresistance in cervical cancer. In total, 24,371cells were clustered into nine separate cell types, including immune cells (T cells, B cells, and myeloid) and non-immune cells [fibroblast cells and endothelial cells (ECs), epithelial cells, and smooth muscle cells]. Among these, fibroblast cells formed five distinct subtypes, and cluster 4 (C4) contained proliferative fibroblast cells enriched in cancer samples. ECs comprised two subclusters: C0 corresponded to tumor-associated ECs. T and B cells formed five and three subclusters, respectively, where C1 of T cells expressed higher levels of immune checkpoint molecules (PDCD1) than the other clusters did and C4 highly expressed proliferating cancer marker genes, MKI167 and BIRC5. Subpopulations of B cells were not strongly affected by drug resistance. Differentially expressed genes (DEGs) between chemoresistant and chemosensitive patients were enriched in phosphoinositide 3-kinase (PI3K)/AKT pathway involved in tumor development, progression, and apoptosis, which might lead to chemoresistance. the instructions of the manufacturer. The cell suspension was further filtered through a 40-mm cell-strainer (BD) and centrifuged at 400g for 10 min to remove the cell aggregates. After removing the supernatant, the pelleted cells were resuspended in 2 ml red blood cell lysis buffer (SolarBio), and then resuspended in PBS with 10% fetal bovine serum (FBS). The viability of the obtained single-cell suspension was detected using a hemocytometer with trypan blue (0.4%, 420301, Gibco).
Single-Cell RNA Capturing, Library Preparation, and Sequencing All single-cell capturing and downstream library construction were performed using the 10× Chromium Single Cell platform (Single Cell 3 library and Gel Bead kit v.3). Briefly, the cell suspensions at a concentration of 1,000 cells/ml were loaded into a 10x Genomics microfluidics chip and encapsulated with barcoded gel beads to generate gel beads in emulsion (GEM). Reverse transcription of generated droplets was performed at 53°C for 45 min. cDNA was amplified for 12 cycles total on T100 Thermal Cycler (Bio-Rad). Then, RT-cDNA was recovered using Recovery Agent provided by 10× Genomics and purified with DynaBeads MyOne Silane Beads (Thermo Fisher Scientific) as outlined in the user guide. Subsequently, purified cDNAs were amplified and cleaned up with the SPRIselect Reagent Kit (Beckman Coulter, USA). Quantification of constructed libraries of each sample was detected using a Qubit2.0 Fluorometer (Invitrogen) before pooling. Pooled libraries were sequenced to a depth of an average of 50,000 reads per cell on a NovaSeq 6000 sequencer (Illumina, San Diego) at 2 × 150 bp sequencing model (23)(24)(25).

scRNA-Seq Data Analysis
scRNA-seq raw data were demultiplexed to FASTQ files (observed average read depth per cell was~1.6 million reads) and aligned to an indexed ensembl_92-GRCh38.92 RefSeq genome to generate gene expression matrices using 10× Genomics pipeline CellRanger v.2.1.0 (25). The number of unique molecular identifiers (UMIs), the number of genes, and the percentage of mitochondrial genes were examined for quality control. Cells expressing <500 or >4,000 genes (potential cell duplets) and gene expression not detected in fewer than three cells were trimmed from the library.
Cells containing >10% mitochondrial genes were also discarded because of their poor cell viability (Table S1). We detected 92,675 genes in 21,433 cells from five samples. After data normalization, variably expressed genes were normalized and scaled, where single-cell gene counts were normalized to the total gene counts presented in that cell at a normalized expression between a low cutoff of 0.0125 and a high cutoff of 3 and a quantile-normalized variance >0.5, using the Seurat R package. The resulting gene expression values were transformed into a log space.

Major Cell-Type Clustering
Principal component analysis (PCA) was used to reduce the dimensionality of the results of variably expressed genes based on the JackStraw function. Then, the first 10 principal components were selected as a statistically significant input for further twodimensional visualization using t-distributed stochastic neighbor embedding (t-SNE) plots (RunTSNE function, the default setting). Cell clusters were annotated and identified to known cell types using specific marker genes identified using the Seurat "FindAllMarkers" function with the default setting (26).

Marker Gene and Cell-Type Identification
Cluster-specific genes were acquired using the Seurat native FindMarkers function with a log-fold change threshold of 0.25. Receiver operating characteristic (ROC) analysis was used to identify cell markers. To further characterize these cell types in each cluster, we used the automated annotation tool SingleR (27) and manually checked using known cell surface markers based on related references.

Identification of Differential Gene Expression and Cell Function in Different Samples
Differential expression analysis was performed after normalization and removal of the batch effects of total genes from the specific cluster from different samples. The function of FindAllMarkers wrapped in the R package of Seurat was used to identify differentially expressed genes for each of cell clusters compared to others. Then, the functions of FindIntegrationAnchors and IntegrateData wrapped in the R package of Seurat were employed to remove the batch effect referred to the standard procedure. We detected differentially expressed genes (DEGs) with a false discovery rate (FDR) <10e-6 and abs (log2 ratio) >2. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses were performed to investigate the cell functional status. DEGs and molecular regulators in the clusters were investigated using GO and KEGG pathway analyses, respectively. GO and KEGG terms with an FDR <0.05 were considered significantly enriched (28).

Single-Cell Profiling of Cervical Cancer and Chronic Cervicitis Tissues
The single-cell atlas of the cervical tissues was characterized using five biopsy samples comprising three cervical squamous cell carcinomas (CQ), one cervical adenocarcinoma (CA), and one chronic cervicitis (CE). Two of the CQ patients were chemoresistant, one was chemosensitive, and three were infected with human papillomavirus. Each sample was processed to isolate single cells without prior selection of cell types, and then we performed scRNA-seq using a 10× Genomics Chromium platform to generate RNA-seq data. After quality filtering, 24,371 high-quality cells from five cervical biopsy samples with a median of 1,303-2,214 genes per cell were analyzed ( Figure 2). Subsequently, the cells were further identified to be nine separate cell types, including fibroblast cells (34.17%), epithelial cells (24.07%), smooth muscle cells (4.3%), ECs (4.26%), and immune cells (33.2%).
The main cell types were identified based on gene expression patterns obtained using dimensionality reduction and unsupervised cell clustering with the described Seurat pipeline. Based on the marker genes of each cell cluster, nine distinct cell clusters were assigned to known cell lineages, which mainly comprised immune and non-immune cells ( Figure 3A). The t-SNE plot also showed distinct clustering according to the tumor origin ( Figure 3B). The heatmap showed differentially expressed marker genes in nine clusters ( Figure 3C).

Reclustering and Differential Gene Profiles of Fibroblast and ECs in Cervicitis and Cervical Cancer Samples
We detected 8,329 fibroblast cells among the five samples and most were found in the CA and the CE samples, which accounted for 95.1% of the total fibroblast cells of the five samples ( Figure 4A). To gain a better understanding of these cell types, we performed a reclustering of 8,329 fibroblasts and 1,037 ECs and assigned each of these subclusters based on known cell markers.
In our study, five distinct subtypes of fibroblast cells were identified. Cluster 0 (C0), C1, and C3 contained cells from all samples, but C2 and C4 were strongly enriched in the CA sample ( Figure 4B). C0 and C3 were similar, corresponding to cells expressing extracellular matrix (ECM) molecules, such as collagen type III alpha 1 chain (COL3A1), COL6A3, and COL1A1, and these cells represent a population of collagengenerating fibroblasts. C1 demonstrated differential activation of FOS, heat shock protein 90 (HSP90s), and ETS-related gene-1 (ERG1).
Cells in C2 exhibited differentially elevated expression of genes involved in translation initiation [ribosomal protein L10 (RPL10) and RPS3] (37) and iron metabolism regulation [ferritin light chains (FTLs)] (38). The final subpopulation C4 was speculated to be a group of proliferative fibroblast cells, which has not been previously recognized, based on the expression of a relatively high level of proliferative genes such as the S100 (s100B) and CDH19 genes (18, 39) ( Figure 4C).
To better understand the differential gene expression profiles between cervicitis and cervical cancer samples, bulk differential gene expression analysis was performed for cervical cancer and cervicitis samples. The two-tailed Wilcoxon rank-sum test, which is implemented in R, was used to conduct the bulk differential gene expression analysis. The volcano map ( Figure 4D) and the heat map ( Figure 4E) showed significant diversity in the expressed genes across the two samples.. For example, wnt2, a member of the WNT gene family, is highly expressed in cervical cancer samples, and WNT signaling is normally involved in the development and progression of various cancers. We also found that the extracellular sulfatase, SULT2, was highly expressed in the cancer samples, which is consistent with the findings of other studies (40). Upregulation of SULF2 gene expression is associated with proliferation, invasion, and migration of cervical cancer cells by its regulation of the extracellular signal-related kinase (ERK)/AKT signaling pathway (41). The ECM is important in tumor genesis and progression and fibroblast-activating protein (FAP), which plays a crucial role in ECM production and remodeling genes and is highly upregulated in tumor fibroblast cells (42).
C4 mainly contained proliferative fibroblast cells, which were mostly from the CA sample. Furthermore, we studied the gene profiles and performed functional enrichment, which demonstrated that several crucial genes were upregulated or downregulated in this cluster. The functional enrichment analysis results indicated that C4 was enriched in cell adhesion molecule binding and collagen-binding growth factor binding in molecular functions ( Figure 4F).
As shown in Figures 5A-C, the 1,037 ECs were clustered in two separate subclusters. C0 corresponded to tumor-associated [heparan sulfate proteoglycan 2 (HPSG2) and plasmalemmal vesicle-associated protein (PLVAP)] and C1 ECs were blood ECs (FLT). Cells from sample 2 were spread in C0 but not in C1. PLVAP is located in the fenestral diaphragm and is considered to play a role in the passage of proteins through the fenestrae (43,44). The generation of Plvap-deficient mice has highlighted the  structural role of PLVAP in the maintenance of size-selective permeability in fenestrated endothelia (45,46). Moreover, Plvapdeficient mice that survived postnatally showed growth retardation (46). The bulk differential gene expression analysis between cervicitis and cervical cancer samples showed that several cancer marker genes, which have been reported in various malignant tumors responsible for angiogenesis, metastasis, and invasion, were significantly elevated in cancer ECs ( Figure 5D). For example, we detected elevated levels of phospholipase D1 (PLD1), a key enzyme involved in lipid metabolism, indicating that abnormal lipid metabolism might be involved in the tumorigenesis and progression of cervical cancer (47,48).

Acquired Resistance Was Associated With Differential Immune Cell Subpopulation Distribution and mRNA Expression
Both adaptive immune cells (T and B lymphocytes) and innate immune cells (such as macrophages, mast cells, neutrophils, dendritic cells, and natural killer cells) have critical roles in the TME and are considered to interact with tumor cells by direct contact or through different chemokine and cytokine signaling pathways that regulate the response of tumors to therapy. To explore the diversity of immune cells in cervical cancer, we extracted 4,961 T cells, 1,480 B cells, and 1,476 myeloid cells from three individual CQ patients and reclustered them (Figures 6, 7). Furthermore, the main T-cell clusters could be  Figure 5A] (26,29). Cells from all clusters were detected in all patients. Cells from all five clusters were detected in both chemoresistant and chemosensitive patients, except for C1, which was predominantly derived from chemoresistant patients. It is worth noting that there were high levels of immune checkpoint molecules, including approved target programmed cell death 1 (PDCD1) and other targets that are currently undergoing clinical   (49), indicating that the cytotoxic activities of CD8+ are significantly curtailed by high-level expression of the checkpoint gene. The robust CD8+ response was confirmed to play a critical role in a mammary tumor model treated with a HER2 inhibitor (50). Additionally, the average expression of the proliferating cancer marker genes, MKI167 and BIRC5, was significantly enriched in C4 ( Figure 6D). These two cancer marker genes, MKI167 and BIRC5, have been reported to be highly expressed in cervical cancer in many studies, where they are associated with the cell cycle pathway and cell apoptosis inhibition, respectively, and are both involved in breast cancer pathogenesis (51)(52)(53). B cells (1,480) were detected and reclustering from C6 and C7 revealed three subclusters, and cell type assignments were performed based on known marker genes. C0 corresponded to follicular B cells based on marker MS4A1, whereas C1 corresponded to plasma B cells based on MZB1, CD38, and IGHG1 ( Figures 7A, C) (29,54). C2 contained lower quality B cells and showed no B-cell markers; thus, it was not further analyzed. One sample was mostly from patient 4, and two others were mostly from patients 1 and 5. The t-SNE profile failed to show the separation of subclusters between chemoresistant and chemosensitive patients, indicating that the B-cell subpopulation was not strongly affected by drug resistance.

DEG Pathways Between Chemoresistant and Chemosensitive Patients
To investigate the enriched functions of DEGs between chemoresistant and chemosensitive samples, signaling pathway analysis was performed among subpopulations of epithelial, T, and B cells. In total, 3,157, 1,399, and 1,141 DEGs were annotated into 327, 312, and 312 pathways for epithelial, T, and B cells, respectively, in the KEGG database (Tables S2-S4). The top 25 annotation results, which were classified according to the pathway types, are shown in Figure 8. Several signaling pathways that are closely correlated with chemotherapy were enriched in all three subpopulations, including the PI3K/AKT and mitogen-activated protein kinase (MAPK) signaling pathways. The PI3K/AKT signaling pathway is involved in tumor development, progression, cellular survival, and apoptosis, and its correlation with chemoresistance has been presented in numerous studies (55)(56)(57) and cervical cancer (12). Another chemotherapy resistance-associated pathway, the MAPK signaling pathway, was significantly enriched and upregulated in chemoresistant patients and is responsible for the regulation of cell migration, survival, proliferation, and progression (58). Activation of survival pathways involving PI3K/Akt and MAPK caused by integrins and soluble factors secreted in the TME results in elevated expression of antiapoptotic proteins, leading to cell viability and drug resistance (7).

DISCUSSION
Advanced scRNA-seq technology has allowed the comprehensive investigation of tumor heterogeneity gene expression differences with single-cell resolution. Although scRNA-seq profiles of other gynecological cancers, including breast and ovarian cancers, have been extensively studied (18,34,59), profiles of the microenvironment of cervical cancer, particularly with respect to chemotherapy resistance, have not been demonstrated.
To the best of our knowledge, this is one of the first few studies to comprehensively characterize a single-cell atlas of cervical cancer patients in China and worldwide (60).
Here, we present a comprehensive catalog of cell types and subtypes, including fibroblast, endothelial, and immune cells in CQ, adenocarcinoma, and non-tumorous CE biopsy samples at single-cell resolution. In this study, we isolated 24,371 cells from five cervical biopsy samples and identified nine separate cell types, including fibroblast, epithelial, smooth muscle, and immune cells as well as ECs. Then, we subclustered the main cell types into different subpopulations based on the expression of the marker gene. We then studied pathway signatures and activities of distinct subpopulations that represent different biological and molecular entities between cervical cancer and cervicitis samples, and assessed the differential gene expression and signaling pathways between chemoresistant and chemosensitive cervical cancer patients.
Fibroblasts are considered a heterogeneous population in many cancers, but the extent of heterogeneity remains unclear, especially in cervical cancer because they are difficult to culture and highly dependent on context (26). Furthermore, fibroblasts are a large source of growth factors and cytokines, including stromal cell-derived factor 1 (SDF1), hepatocyte growth factor (HGF), and vascular endothelial growth factor (VEGF), which all promote tumor growth and contribute to chemotherapy resistance (61)(62)(63). Although most clusters were composed of cells originating from different samples, the fibroblast cells and ECs were observably enriched in CA and CE patients.
The significant shift might have occurred because the condition of the patients who provided the cervical squamous cell carcinoma samples was likely more serious, which would result in the recruitment of more immune cells to exert their antitumor function. In addition, upregulated expression of WNT2 in tumor fibroblast cells, which was observed in this study, could also result in tumor growth and promotion of invasion by activating the canonical WNT/b-catenin signaling pathway (64,65). A previous study also revealed that downregulating WNT2 significantly suppressed cell motility and invasion and reversed epithelial-mesenchymal transition (EMT) progression in cervical cancer (66). By comparing the differential gene expression in fibroblasts and ECs between cervical cancer and cervicitis patients, we found that several cancer marker genes, which have been reported in various malignant tumors responsible for angiogenesis, metastasis, and invasion, via various metabolic pathways related to cervical cancer, were significantly elevated in cancer ECs. For instance, genes involved in translation initiation such as RPL10 and RPS3 (37) and iron metabolism regulation such as FTLs (38), which have been reported to be involved in increased cell viability, migration, and invasion in different cancers, were elevated in the identified C2. This observation indicates that these cells are highly active and associated with tumorigenesis. In addition, the differential activation of FOS, HSP90s, and ERG1 exhibited by the cells in C1 indicated a strong correlation with myogenesis and fibrosis (37,67). Tumor cells have been reported to trigger the immunosuppressive activities of tumor-associated ECs that influence therapeutic responses (68). ECs selectively permit transmigration of immunosuppressive myeloid cells from the blood to the tumor, thereby impairing antitumor immune responses (69). Moreover, ECs involved in the TME may also suppress T-cell function through the expression of different inhibitory molecules, such as programmed cell death ligand 1 (PDL1) and PDL2 (70).
Checkpoint immunotherapy is rarely used in cervical cancer, and there is a poor understanding of the outcome of immunotherapy. The results of this study were consistent with this phenomenon, and we found that the immune checkpoint molecules, PDCD1, LAG3, and HAVCR2, were highly expressed in CD8+ cells. This observation indicates that cytotoxicity was inhibited by high checkpoint expression (54) and that immunotherapy could also be useful in cervical cancer by targeting these pathways. This result is in agreement with the good clinical results observed in immunotherapy of cervical cancer (71).
The results of the analysis of differentially enriched signaling pathways did not show significant differences in chemoresistant patients. However, the PI3K/AKT pathway enriched most DEGs, which suggests that the activation state of this common pathway, which is involved in tumor development, progression, and apoptosis, might contribute to the development of resistance to chemotherapy in cervical cancer patients. Inhibition of certain components of the PI3K/AKT and MAPK pathways not only enhances chemotherapy efficacy in cervical cancer, but also has the potential to overcome resistance (72,73).
There are some limitations to our study that need to be addressed in future studies. First, because sample collection was challenging, the number of samples in this study was relatively small, especially the chemoresistant and chemosensitive samples, which indicates that our observations cannot reflect the comprehensive profiles of the TME and generalize the gene expression differences in cervical cancer and its chemoresistant features.

CONCLUSION
In conclusion, we comprehensively described the profiles of immune and non-immune cells in cervical cancer and cervicitis samples at single-cell resolution. Furthermore, we compared the distinctive features of signaling pathways among subpopulations from chemoresistant and chemosensitive samples in this study. Our results of the analysis of extensive bioinformatics data demonstrated the scope and potential impact of heterogeneity and suggest that single-cell profiling could identify and characterize clinically important subpopulations to develop successful targeted treatments. Our findings also indicated the need for comprehensive single-cell gene profiling and characterization of heterogeneous tumors such as chemoresistant cervical cancer. Finally, by assessing the cell subpopulations, differential gene expression between chemoresistant and chemosensitive samples, and distinct signaling pathways, we expect that our findings will provide novel and deeper insights into human cervical cancer and facilitate the advancement of its diagnosis and treatment.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding authors.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by the Institutional Review Board of Renmin Hospital of Wuhan University (Institutional Review Board no. WDRY2021-K014). The patients/participants provided their written informed consent to participate in this study.

ACKNOWLEDGMENTS
We would like to thank Professor Congyi Zheng for his helpful suggestions regarding sample sources. We also thank all the individuals who participated in this study and donated samples.