If you don't remember your password, you can reset it by entering your email address and clicking the Reset Password button. You will then receive an email that contains a secure link for resetting your password
If the address matches a valid account an email will be sent to __email__ with instructions for resetting your password
Division of Gastroenterology and Hepatology, Key Laboratory of Gastroenterology and Hepatology, Ministry of Health, Shanghai Institute of Digestive Disease, Renji Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai, China
State Key Laboratory of Natural Medicines, Research Center of Biostatistics and Computational Pharmacy, China Pharmaceutical University, Nanjing, China
Division of Gastroenterology and Hepatology, Key Laboratory of Gastroenterology and Hepatology, Ministry of Health, Shanghai Institute of Digestive Disease, Renji Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai, China
Division of Gastroenterology and Hepatology, Key Laboratory of Gastroenterology and Hepatology, Ministry of Health, Shanghai Institute of Digestive Disease, Renji Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai, China
State Key Laboratory of Natural Medicines, Research Center of Biostatistics and Computational Pharmacy, China Pharmaceutical University, Nanjing, China
Division of Gastroenterology and Hepatology, Key Laboratory of Gastroenterology and Hepatology, Ministry of Health, Shanghai Institute of Digestive Disease, Renji Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai, China
Division of Gastroenterology and Hepatology, Key Laboratory of Gastroenterology and Hepatology, Ministry of Health, Shanghai Institute of Digestive Disease, Renji Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai, China
Division of Gastroenterology and Hepatology, Key Laboratory of Gastroenterology and Hepatology, Ministry of Health, Shanghai Institute of Digestive Disease, Renji Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai, China
Division of Gastroenterology and Hepatology, Key Laboratory of Gastroenterology and Hepatology, Ministry of Health, Shanghai Institute of Digestive Disease, Renji Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai, China
Fang-Rong Yan, State Key Laboratory of Natural Medicines, Research Center of Biostatistics and Computational Pharmacy, China Pharmaceutical University, Nanjing 211198, China.
State Key Laboratory of Natural Medicines, Research Center of Biostatistics and Computational Pharmacy, China Pharmaceutical University, Nanjing, China
Haoyan Chen, State Key Laboratory for Oncogenes and Related Genes, Shanghai Cancer Institute, Renji Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai 200032, China.
State Key Laboratory for Oncogenes and Related Genes, Shanghai Cancer Institute, Renji Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai, China
Division of Gastroenterology and Hepatology, Key Laboratory of Gastroenterology and Hepatology, Ministry of Health, Shanghai Institute of Digestive Disease, Renji Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai, China
Approximately one-third of colorectal cancers develop from serrated lesions (SLs), including hyperplastic polyp (HP), sessile serrated lesion (SSL), traditional serrated adenoma (TSA), and SSL with dysplasia (SSLD) through the serrated neoplasia pathway, which progresses faster than the conventional adenoma-carcinoma pathway. We sought to depict the currently unclarified molecular and immune alterations by the single-cell landscape in SLs.
Methods
We performed single-cell RNA sequencing of 16 SLs (including 4 proximal HPs, 5 SSLs, 2 SSLDs, and 5 TSAs) vs 3 normal colonic tissues.
Results
A total of 60,568 high-quality cells were obtained. Two distinct epithelial clusters with redox imbalance in SLs were observed, along with upregulation of tumor-promoting SerpinB6 that regulated ROS level. Epithelial clusters of SSL and TSA showed distinct molecular features: SSL-specific epithelium manifested overexpressed proliferative markers with Notch pathway activation, whereas TSA-specific epithelium showed Paneth cell metaplasia with aberrant lysozyme expression. As for immune contexture, enhanced cytotoxic activity of CD8+ T cells was observed in SLs; it was mainly attributable to increased proportion of CD103+CD8+ tissue-resident memory T cells, which might be regulated by retinoic acid metabolism. Microenvironment of SLs was generally immune-activated, whereas some immunosuppressive cells (regulatory T cells, anti-inflammatory macrophages, MDK+IgA+ plasma cells, MMP11-secreting PDGFRA+ fibroblasts) also emerged at early stage and further accumulated in SSLD.
Conclusion
Epithelial, immune, and stromal components in the serrated pathway undergo fundamental alterations. Future molecular subtypes of SLs and potential immune therapy might be developed.
The study identified epithelial cells specific to serrated lesions (SLs) with dysregulated expression of genes related to redox balance, and changes specific to SLs in nonparenchymal compartment (eg, accumulation of tissue-resident memory T cells and PDGFRA+ fibroblasts) by depicting the single-cell landscape of 16 SLs.
Colorectal cancer (CRC) ranks as the third most common malignant tumor and the third leading cause of cancer deaths around the globe.
CRC originates from premalignant lesions through 2 major pathways: conventional adenoma-carcinoma sequence and alternative serrated tumorigenesis pathway. Characterized by saw-toothed appearance of crypts, serrated lesions (SLs) account for precancerous lesions of 15% to 35% sporadic CRC.
According to the latest World Health Organziation classification, SLs consist of the following categories based on pathological features: hyperplastic polyp (HP), sessile serrated lesion (SSL), SSL with dysplasia (SSLD), traditional serrated adenoma (TSA), and unclassified.
Proximal (right-sided) microvesicular HPs are precursors of SSLs. It is now recognized that SSLs and TSAs have the malignant potential to develop into serrated CRC.
SSL locates in proximal colon with predominant BRAFV600E mutation, which can progress into BRAF-mutant right-sided CRC with either microsatellite instability-high (MSI-H) or microsatellite stable (MSS).
Traditional serrated pathway-associated colorectal carcinoma: morphologic reappraisal of serrated morphology, tumor budding, and identification of frequent PTEN alterations.
Of note, some case studies revealed that lesions arising from serrated pathway showed accelerated carcinogenesis, with reports of rapid transformation from SSL to invasive CRC within several months,
Furthermore, SSL is easily missed during colonoscopy examination due to its flat appearance. As such, SSLs have become a major cause of interval CRCs after a negative colonoscopy.
its major molecular alterations, and its interaction with the immunological environment are not well-characterized. Understanding the pathogenesis of serrated neoplasia is vital for precision treatment and early intervention for CRC. However, a key challenge in further elucidating the serrated pathway is that it is difficult to clearly distinguish serrated pathway-derived invasive CRC from those originated from conventional adenomas, owing to the loss of serrated morphology at the late stage.
To fulfill this gap, in this study, we focused on well-characterized early lesions of the serrated subtypes and performed single-cell RNA sequencing (scRNA-seq) of 16 SL samples, for more in-depth understanding of pathogenesis of serrated colorectal neoplasia and the serrated pathway to CRC.
Results
Single-cell Transcriptomic Atlas of Human SLs
To explore transcriptional alterations in serrated tumorigenesis, we generated an scRNA-seq profile from 3 normal contrasts (NCs), 4 right-sided microvesicular HPs, 5 SSLs, 2 SSLDs, and 5 TSAs samples using 10× Genomics sequencing (Figure 1, A). The clinical characteristics and endoscopic imaging of these patients are shown in Supplementary Table 1 and Figure 1, B–E. All serrated lesions used for analysis used in our study were sporadic.
Figure 1General clinical and single-cell transcriptional profile of SL samples included in this study.A, Schematic diagram of scRNA-seq workflow in 16 SLs. B, Representative endoscopic image and the hematoxylin and eosin (H&E) staining of HP samples in this study. C, Representative endoscopic image and the H&E staining of SSL samples in this study. D, Endoscopic image and the H&E staining of SSL with dysplasia in this study. E, Representative endoscopic image and the H&E staining of TSA samples in this study. F, UMAP plot showed 6 color-coded global cell types of 60,568 cells analyzed. G, Fractions of epithelial, stromal, myeloid, T, B, and mast cells among normal, HP, SSL, SSLD, and TSA samples. H, Overview of sequencing cell number, cell proportion, and anatomical region of each sample.
After human genome mapping and quality filtering, we obtained a total of 60,568 high-quality cells, of which 11,661 cells were from NC, 12,465 cells were from HP, 16,713 cells were from SSL, 3003 cells were from SSLD, and 16,726 cells were from TSA (Figure 1, F–H; Figure 2, A).
Figure 2Dissection of the various cell types in human SLs with scRNA-seq.A, Violin plot showed the expressed genes (nGene) and read counts (transcript, nUMI). B, UMAP plot showed Louvain clusters of 60,568 cells analyzed. C, UMAP plot showed distribution of cells from 19 samples. D, Feature plots illustrated canonical markers of global cell types. E, Heatmap illustrated the top 10 DEGs (adjusted P < .05; |log2FC| > 0.25) in each global cell type.
We identified 16 clusters in total through Louvain clustering and uniform manifold approximation and projection (UMAP) for dimensionality reduction and visualization (Figure 2, B). We then reclassified these Louvain clusters into global cellular compartments using graph-based clustering and known classic markers (epithelial cells: EPCAM; T cells: CD3D; B cells: MS4A1, CD79A; myeloid cells: CD14, CD163; stromal cells: COL1A1, PECAM1). Mast cells were removed from analysis as fewer than 1000 cells were detected. Each cell type was from multiple patients, which indicated correction of the batch effect (Figure 2, C–D). The top 10 differentially expressed genes (DEGs; adjusted P < .05; |log2FC| > 0.25) between major cell types further confirmed the clustering results (Supplementary Table 2 and Figure 2, E).
Identification of an Epithelial Cell Population Featuring Redox Imbalance Specifically Existing in SLs
In this study, we identified 20 main epithelial cell types in serrated polyps, after the reclustering of 7294 epithelial cells (Figure 3, A–C and Supplementary Table 3). Seven types of normal cells in large intestine were designated: BEST4+ colonocytes (high expression of BEST4, OPOT2, LYPD8, CA4, CA7), mature colonocytes (BEST4lo GUCA2Bhi, and SLC26A3hi), goblet cells (high expression of TFF3, SPINK1, REG4, AGR2, MUC2, and WFDC2), transit amplifying cells (high expression of mitochondrial ATP synthase component ATP5MC1, OLFM4+, and PCNA+), intermediate cells (expressed both absorptive and secretory markers and MKI67), tuff cells (HPGDS+ and PTGS1+), and enteroendocrine cells (SCGN+ and PCSK1N+).
Figure 3Epithelial cell subclusters in SLs.A, UMAP representation of dimensionally reduced data following graph-based clustering with marker-based cell type assignments (left) and proportion of epithelial cell types in NC and SLs (right). B, Dot plot depicting expression levels of canonical colonocyte marker genes together with the percentage of cells expressing the marker. C, UMAP plot showing Seurat clusters of epithelial cells and association of Seurat clusters with defined epithelial cell types. D, Expression of GUCA2A and GUCA2B in epithelia of different groups. E, UMAP plot showing the sample origin of epithelial cells. F, Volcano plot depicting the differentially expressed genes between Epi-SL cluster and normal colonocytes (BEST4+ colonocytes, mature colonocytes, and goblet cells). G, Expression of LYZ, RNF43, and MYC in the epithelia of different groups.
In SLs, we observed a decreased proportion of these normal cell components, with lower expression of endogenous hormones guanylin (GUCA2A) and uroguanylin (GUCA2B), which were involved in intestinal fluid homeostasis and normal secreting function (Figure 3, D). BEST4/OTOP2 cells were a cluster of newly identified pH-sensing absorptive colonocyte expressing calcium-sensitive chloride channel BEST4,
which was almost absent in SLs, suggesting impaired ion secretion and intestinal barrier function.
Interestingly, a distinct epithelial subpopulation was found predominantly in SLs but was rarely present in the normal tissue (Figure 3, E). Thus, we named these SL-specific epithelia as Epi-SL (cluster 1, 2, 5, 6, 9, 12, 15, 16, 17). We compared DEGs between Epi-SL and normal colonocytes (BEST4+ colonocytes, mature colonocytes, and goblet cells). As shown in Figure 3, F–G, top overexpressed genes of Epi-SL included OLFM4 (intestinal stem cell marker), GSTP1 (an antioxidant enzyme), PKM (a key rate-limiting enzyme in glycolysis), and lysozyme (LYZ). Interestingly, LYZ was a known marker of Paneth cell, which predominantly existed in the small intestine and was nearly absent in normal large intestine epithelium.
Of note, several significantly upregulated genes in Epi-SL, such as SOX9, MYC, GOLM1, LCN2, S100A11, FOXQ1, ANXA2, ENO1, and PKM (Supplementary Table 4), have been established to play oncogenic roles in colon tumorigenesis,
indicating the carcinogenic nature of Epi-SL cluster.
To further characterize Epi-SL cells, which represent altered epithelium in premalignant colon SLs, we focused on genes significantly upregulated in this subpopulation compared with normal absorptive and secretory colonocytes (adjusted P < .05; log2FC > 0.25). Gene Ontology (GO) enrichment analysis showed that pathways related to response to oxidative stress and oxidoreductase activity were significantly enriched in Epi-SL cluster (all adjusted P < .05) (Figure 4, A). These findings suggested that redox imbalance might play a role in serrated tumorigenesis. In Epi-SL, upregulation of both oxidant (DUOX1, DUOX2, DUOXA1, DUOXA2, GPX2, AQP5) and antioxidant genes (GSTP1, NQO1, NQO2) were observed; immunohistochemistry results confirmed overexpression of GSTP1, NQO1, and DUOX2, along with higher 8-OHdG (marker of oxidant stress) in SSL and TSA (Figure 4, B–C). Using data from a Gene Expression Omnibus (GEO) external dataset comprising 10 tubular adenomas (TAs, also known as conventional adenomas), we found these genes related to oxidant stress were not overexpressed in conventional adenomas, indicating this feature specifically existed in premalignant lesions developed from the serrated pathway (Figure 4, D–E). Moreover, we found pathway enrichment of cellular response to oxidative stress was a distinct feature of BRAF-mutant CRCs, as confirmed by both The Cancer Genome Atlas (TCGA) and GEO external datasets (Figure 4, F–G).
Figure 4Identification of an Epi-SL cluster with dysregulated expression of genes related to redox balance in serrated polyps.A, GO enrichment analysis showed increased oxidative stress-related pathway activity in Epi-SL cluster (scored per cell by GSVA). B, Representative images of GSTP1, NQO1, DUOX2, 8-OHdG immunohistological staining of FFPE tissue from NC, SSL, and TSA. C, Heatmap showing expression of oxidative and anti-oxidative genes for each cluster. D, Expression of oxidative stress-related genes in NC (n = 10), conventional adenomas (n = 10), and SSLs (n = 21). Data was curated from GSE76987, which included RNA-seq data of various right-sided colon tissues. E, GSVA scores of response to oxidative stress pathway in NC, conventional adenomas, and SSL. Data were obtained from GSE76987. F, GSVA scores of response to oxidative stress pathway in BRAF wild-type (n = 634) and BRAF mutant CRCs (n = 64). Data were obtained from TCGA. G, GSVA scores of response to oxidative stress pathway in BRAF wild-type (n = 461) and BRAF mutant CRCs (n = 51). Data were obtained from GSE39582. H, SerpinB6 knockdown resulted in increased ROS levels in RKO cell line. ∗P < .05.
We also identified upregulation of a novel tumor-promoting gene that regulated ROS level, namely SERPINB6 (encoded SerpinB6), in epithelia of SLs (Figure 4, H and Figure 5, A), whose role has not been explored in colon tumorigenesis yet. SerpinB6 was found overexpressed in Epi-SL cluster (Figure 4, C) and was specifically overexpressed in SSL but not in conventional adenomas (Figure 4, D and Figure 5, B). Its overexpression in SLs was confirmed by quantitative real-time polymerase chain reaction (RT-PCR) and immunohistochemistry (Figure 5, C–D). Of note, SerpinB6 showed significantly higher expression in BRAF-mutated CRCs in TCGA dataset (Figure 5, E), suggesting it might play a role in serrated tumorigenesis. Results of EdU and transwell assay in RKO cell line (a colon cancer cell line with BRAF mutation) revealed that knock down of SerpinB6 significantly attenuated the cell growth and migration ability (Figure 5, F–H). Furthermore, we found knock down of SerpinB6 resulted in increased reactive oxygen species (ROS) levels as indicated by 2'-7'dichlorofluorescin diacetate (DCFH-DA) assay, whereas usage of N-acetyl L-cysteine (a ROS scavenger) could partially rescue the proliferation ability in the SerpinB6-knockdown cells, suggesting that SerpinB6 can promote serrated tumorigenesis via modulating ROS levels (Figure 4, H and Figure 5, I).
Figure 5SerpinB6 was linked to serrated route of colon tumorigenesis and promoted colon cancer cell line proliferation.A, Violin plot showed increased expression level of SerpinB6 in epithelia of HP, SSL, SSLD, and TSA, compared with normal tissue. B, Upregulation of SerpinB6 in SSL (n = 21) compared with tubular adenoma (conventional adenoma, n = 10) in a bulk RNA-seq dataset (GSE76987). C, Immunohistochemistry staining of SerpinB6 in normal, conventional adenoma, and SLs. D, Quantitative RT-PCR analysis of SerpinB6 expression in 8 fresh SSL tissues compared with paired normal tissues. E, SerpinB6 expression in BRAF-mutant (n = 64) and wild-type (n = 634) CRCs in TCGA database. F, Knock down (SerpinB6-KD) and re-introduction (ResSerpinB6) of SerpinB6 in RKO cell line, as validated by mRNA (n = 3 replicates/group for quantitative RT-PCR assay) and protein level. G, EdU assay showed knock down of SerpinB6 suppressed cell proliferation of RKO, while re-introduction of SerpinB6 rescued cell growth (n = 3 replicates/group). H, Transwell experiment showed that knock down of SerpinB6 suppressed migration ability of RKO, while re-introduction of SerpinB6 rescued cell migration (n = 3 replicates/group). I, Treatment of 5 mM N-acetyl L-cysteine (NAC) overnight partially rescued the proliferation ability of RKO after knocking down SerpinB6 (n = 3 replicates/group). ∗P < .05; ∗∗P < .01; ∗∗∗P < .001.
Subclustering of Epi-SL Identifies One Cluster of Ki67hi OLFM4hi Epithelial Cells With Activation of Notch Signaling Pathway in SSL, and Another TSA-specific Cluster Producing Lysosome
To further characterize proliferation ability of epithelia in SLs, we performed cell cycle analysis (Figure 6, A). We found the left part cells of Epi-SL on the UMAP plot were in G2/M phase, whereas the majority of cells in the right part were in G1 phase. Considering the high heterogeneity within the Epi-SL cluster, we generated 2 subclusters (Epi-1 and Epi-2) derived by reclustering these cells (Figure 6, B). Among these, thenEpi-1 subcluster was mainly found in SSLs and accounted for the majority of epithelial cells in SSLD, whereas Epi-2 was exclusively found in TSAs (Figure 6, C). Trajectory analysis indicated that transit amplifying cells and intermediate cells possess the ability to differentiate toward Epi-1 and Epi-2 cell population (Figure 6, D). Moreover, we found that the SSL-associated epithelial cluster (Epi-1) showed high expression of a set of genes related to serrated polyposis syndrome germline predisposition as reported by Soares de Lima et al.
(Figure 6, E); this Epi-1 cluster also displayed significantly higher copy number variation level (Figure 6, F–G), indicating its malignant potential and critical role in the formation of serrated polyp.
Figure 6Identification of two distinct epithelial cell populations (Epi-1 and Epi-2) in SSL and TSA.A, Cell cycle analysis of epithelia in SLs. B, Reclustering of Epi-SL resulted in 2 subclusters. C, Proportion of Epi-1 and Epi-2 subclusters among HP, SSL, and TSA. D, Unsupervised trajectory plot inferred by Monocle 2 showed that Epi-1 and Epi-2 cells derived from transit amplifying and intermediate cells. E, Higher expression of genes related to SPS germline predisposition in Epi-1 subcluster. F, CNV score calculated by inferCNV R package among four different cell types (BEST4+, mature colonocytes, Epi-1, and Epi-2) in serrated polyps. G, Chromosomal landscape of inferred large-scale CNVs distinguishing serrated polyp epithelial cells from normal colonic cells. The SL epithelial cells were shown with different cell types (BEST4+, mature colonocytes, Epi-1, and Epi-2) and chromosomal regions (horizontal). Amplifications (red) or deletions (blue) were plotted, with normal sample as the reference (top). A 6-state CNV model that predicts the following CNV levels was used: 0: complete loss; 0.5: loss of 1 copy; 1: neutral; 1.5: addition of 1 copy; 2: addition of 2 copies; 3: addition of more than 2 copies. H, Volcano plot of DEGs between Epi-1 and Epi-2 (adjusted P < .05; |log2FC| > 0.25). Red points represented genes significantly overexpressed in Epi-1 cluster, and blue points showed genes significantly overexpressed in Epi-2 cluster. I, Violin plot showed expression of OLFM4, SOX4, HES1, and JUN in epithelia of NC, HP, SSL, SSLD, and TSA.
Compared with the Epi-2 subcluster, Epi-1 showed higher expression of proliferation marker MKI67, SOX4, stem cell marker OLFM4, HES1, and JUN, along with downregulation of LYZ, CEACAM5, DUOX2, and S100A11 (see Figure 6, H–I, and Supplementary Table 5).
Upregulation of the genes related to proliferative ability in the SSL-specific Epi-1 cell cluster partially explained why SSL can rapidly develop into serrated CRC. As expected, these proliferative genes showed significantly higher expression in the epithelia of SSLD than those of SSL (Supplementary Table 6). Considering its clinical importance for SSL malignant transformation, we focused on pathway alterations in the Epi-1 subcluster (MKI67hi OLFM4hi SSL-specific). Gene set variation analysis (GSVA) indicated several activated oncogenic and metabolic signatures in this subcluster, such as Myc, G2M checkpoint, Notch signaling, E2F targets, MAPK targets, and FOXO mediated transcription of oxidative stress metabolic genes (Figure 7, A). Besides, transit amplifying cells and Epi-1 cells shared commonly activated signatures of glycolysis, citrate cycle, fatty acid metabolism, arginine metabolism, and folate metabolism, suggesting metabolic reprogramming in these cells.
Figure 7Distinct subtypes of epithelial cells in SSL and TSA tissues.A, Comparison of pathway activities (calculated based on GSVA) among different epithelial cell subtypes. The scores of pathways were normalized. Each column represents a single cell and each cell cluster analyzed was color-coded. B, Higher enrichment of pathways upregulated in Epi-1 was also observed in SSL with cytological dysplasia samples (n = 22, data from EMTAB7960), compared with SSL without dysplasia (n = 20; top). Higher Epi-1 signature score was observed in SSL with cytological dysplasia samples (bottom). C, Heatmap showed expression of genes related to Notch signaling pathway in each epithelial cell type. D, Immunohistochemistry staining of JAG1 and MYC in NC and SSL tissue. E, GO pathway enrichment analysis of upregulated genes in Epi-2. F, Heatmap showing expression of annexin and S100 family genes in epithelial subclusters. G, Violin plot showed the highest expression of LYZ in Epi-2 cluster. H, Confocal staining of MUC2 and LYZ in NC, SSL, and TSA tissue.
Based on these findings, we hypothesized that Notch signaling pathway and its downstream target MYC, which represented hallmark pathways in Epi-1, might be involved in the malignant transformation of SSL. To confirm the single-cell findings, we recalculated GSVA scores of hallmark and curated gene sets that were elevated in Epi-1, in an external bulk RNA-seq dataset comprising of SSL and SSLD samples. Strikingly, these Epi-1-elevated pathways were further enhanced in the dysplastic lesions, suggesting these features of Epi-1 were relevant to cancerous change of SSLs (Figure 7, B, top). We also identified the top 50 marker genes of Epi-1 (compared with normal colonocytes) to construct ‘Epi-1 signature’ (Supplementary Table 7); the GSVA score of ‘Epi-1 signature’ was significantly elevated in SSLs with dysplasia than those without (Figure 7, B, bottom). Furthermore, we performed immunohistochemistry in SSL and NC tissues. Consistent with findings of overexpression of Notch ligand, receptor, and target genes found at single-cell level in SSL epithelium (Figure 7, C), we observed increased expressions of JAG1 (a Notch ligand) and MYC (a Notch target) in SSL tissues at the protein level (Figure 7, D).
On the other hand, the TSA-associated epithelial cluster (Epi-2) was enriched in pathways of lysosome secretion, angiogenesis, and signaling by RNF43 mutants (Figure 7, A), which is consistent with previous reports that identified RNF43 mutation in approximately 25% of TSAs.
We also investigated GO and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment of upregulated genes in Epi-2 compared with Epi-1, and found activation of pathways related to lysosome, protein transport, neutrophil activation, calcium-dependent protein binding and S100 protein binding in Epi-2 (Figure 7, E). In short, Epi-2 was characterized by ectopic expression of LYZ and also higher expression of annexin and S100 family genes (Figure 7, F and Supplementary Table 8). Increased expression of LYZ was observed in the Epi-2 subcluster, suggesting a trend of Paneth cell differentiation in TSA. These LYZhi TSA-specific epithelial cells may serve as a diagnostic biomarker for TSA. Confocal staining established the existence of abundant LYZ secretion in the crypt of TSA, but not in SSL or NC (Figure 7, G–H).
Subsequently, we characterized transcriptional factor (TF) regulatory network in Epi-1 and Epi-2 clusters by identifying super-regulons. As illustrated in Figure 8, upregulated TF activities for Epi-1 included IRF1, KLF10, FOSL1, RAD21, MYC, ZNF143, and BHLHE40, whereas elevated regulons for Epi-2 incorporated KLF4, KLF6, HIVEP1, NR2C2, ATF4, KLF13, and NFIL3, suggesting distinct development paths between SSL and TSA.
Figure 8TF target network created from normal and SL cells, organized into super-regulons for Epi-1 (A) and Epi-2 (B) cells. Color in each vertex represented enrichment scores in corresponding regulon, whereas edge opacities were determined by the inferred TF-target weightings between regulons.
Both Immune Active and Immunosuppressive T Cell Components are Expanded Within Serrated Adenoma Microenvironment
Next, T cell components were subjected to clustering and UMAP visualization (Figure 9, A). Based on known markers, 6 main T cell types were identified (Figure 9, B and Supplementary Table 9). We observed a decrease in the proportion of naïve T cells in serrated neoplasms, but an increase of both CD8+ and CD4+ T cells, as well as regulatory T cells (Treg), compared with normal colon samples (Figure 9, C). Of note, we observed further accumulation of Tregs in SSLD, which suggested that it may facilitate carcinogenic transformation by creating an immunosuppressive microenvironment. The clinical significance of FOXP3+ Tregs was validated by immunohistochemistry analysis of expanded clinical samples of SLs (Figure 9, D–E), which showed a significant increase in the number of FOXP3+ cells in SSL with dysplasia samples (n = 8; P = .01), compared with SSLs without dysplasia (n = 8).
Figure 9Transcriptomic heterogeneity of T cells in the tumor microenvironment of SLs.A, UMAP plot for clustering of T cells, color-coded for 6 clusters. B, Heatmap showed expression of marker genes in each T cell cluster: CD8+ T cells (cluster 1, 6, 8, 11, 12, 16), CD4+ T cells (cluster 0, 7), naïve T cells (cluster 2, 4, 10, 14, with high expression of naïve markers TCF7, SELL, LEF1, CCR7), Treg (cluster 9, high expression of FOXP3), γδ T cells (cluster 3, 5, high expression of TRDC1 and TRGC1), and NK cells (cluster 15, high expression of KLRF1). C, Cell number and proportion of 6 T cell clusters among HP, SSL, and TSA. D, Comparison of numbers of FOXP3+ cells in NC, HP, SSL, TSA, and SSLD tissues. E, Immunohistochemistry staining of regulatory T cells (FOXP3+ cells) in NC and SLs. F, Volcano plot displayed DEGs of CD8+ T cells between SLs and normal tissue. G, GO pathway enrichment analysis of upregulated genes in CD8+ T cells of SLs. H, UMAP plot showing Seurat clusters of CD8+ T cells. I, UMAP plot for subclustering of CD8+ T cells, color-coded for 6 subclusters (top). Cell number and proportion of CD8+ T cell subclusters among HP, SSL, and TSA (bottom). J, Heatmap showed expression of marker genes related to naïve, inhibitory, effector, resident, and cycling function in each CD8+ T cell subcluster. ∗∗P < .01; ∗∗∗P < .001.
Because CD8+ T cells are known to play a central role in anti-tumor immunity, we focused on alterations of CD8+ T cells (including CD8+ naïve T cells, cluster 1, 6, 8, 10, 11, 12, 14, 16 in Figure 9, A). To investigate functional differences of CD8+ T cells in SLs and NC, we listed the DEGs between the 2 groups (Supplementary Table 10 and Figure 9, F). GO pathway analysis of the upregulated genes showed enrichment of pathways related to T cell activation and cytotoxic function (Figure 9, G). Among the DEGs, effector markers (GZMB, GZMK, GZMA, PRF1, TNF, KLRB1, BHLHE40, and NKG7) and resident markers (ITGAE, ITGB1, ITGA1, PRDM1, and CXCR6) of CD8+ T cells were significantly elevated.
Then, CD8+ T cells were reclustered to generate 5 subpopulations by clustering and dimensionality reduction from 11 Seurat clusters (Figure 9, H–I). Two clusters (clusters 4 and 10) expressing high levels of naïve markers were assigned as CD8+-naïve cells, 1 cluster (cluster 11) with high MKI67 was regarded as CD8+ cycling cells. The UMAP plot and heatmap of functional marker genes (Figure 9, J) also revealed 2 distinct subgroups of GZMK+ cells with cytotoxic function. However, one GZMK+ subgroup showed higher granzyme (GZMB, GZMK, GZMH), KLRG1, and interferon-γ (IFNG) levels, along with elevated ITGB2 expression, which was named ‘GZMK+ Effector-1’; another GZMK+ subgroup was expanded in SSL with dysplasia and TSA microenvironment with much lower granzyme activity, which was titled ‘GZMK+ Effector-2.’ Notably, transcriptional traits of GZMK+ Effector-1 cells were similar to a cell population that was previously identified as ‘ITGB2+ tissue resident memory T cells’ in recent studies on human intestine.
Although these cells showed potential of strong cytotoxic function, their proportion was nearly unchanged in HP and SSL.
Another 2 subpopulations of CD8+ T cells (cluster 0, 2, 5, 7, and 8; on the right part of UMAP plot) displayed upregulation of integrin E (ITGAE) expression. Cells in cluster 5 manifested overexpressed immune checkpoint genes (PDCD1, CTLA4, TIGIT, and LAG3) with lower IFNG and granzyme B (GZMB) expression, and thus was named CD8+ Exhausted cluster, which was mainly detected in SSL. Interestingly, the other cell component (cluster 0, 2, 7, and 8) with high expression of ITGAE (encoded CD103), which serves as a marker gene of tissue-resident memory T cells (Trm) with limited recirculation,
was significantly enriched in the microenvironment of HP, SSL, SSL with dysplasia, and TSA (Figure 10, A); we designated these cells as CD103+ Trm cells. These CD103+ Trm cells showed high expressions of IFNG, GZMB, and IL2. Consistent with these findings, using a predefined Trm CD8+ gene signature by Savas et al,
we demonstrated that GSVA scores of Trm CD8+ gene signature were significantly elevated in CD8+ T cells in all 4 kinds of SLs (Figure 10, B) (all P < .001). Cell cycle analysis revealed that G2/M proportion of CD103+ Trm was not increased in SLs, whereas G2/M proportion of naïve CD8+ T was augmented, suggesting that the increase of these cells was probably attributable to differentiation from naïve cells (Figure 10, C).
Figure 10Characteristics of CD103+Trm cells in SLs.A, Violin plot showing the expression of ITGAE/CD103 in CD8+ T cells in NC and SLs. B, GSVA scores of Trm CD8+ gene signature in NC and SLs. C, Cell cycle analysis of different CD8+ T cell subclusters in NC and SLs. The full list of DEGs is summarized in Supplementary Table 10. D, Confocal staining of CD103 (in green) and CD8 (in red) in SSL. E, Representative immunohistochemical staining of CD103 in NC, SLs, conventional adenomas (TAs), and BRAF-mutant CRCs. F, Quantitative of CD103+ cells in NC, SLs, conventional adenomas, and BRAF-mutant CRCs. ∗P < .05; ∗∗P < .01.
In our study, we observed marked enrichment of CD103+ Trm component in SLs; therefore, we further explored the clinical importance and potential regulating factors of these cells. We confirmed co-localization of CD8 and CD103 by confocal immunofluorescence microscopy in SSL slides (Figure 10, D). Next, we evaluated CD103 expression in NC, serrated tumors, and TA (conventional adenoma) by immunohistochemistry (Figure 10, E). In serrated tumors, CD103 was exclusively expressed in lymphocytes; we observed significantly increased number of CD103+ cells in HP (P = .01), SSL (P = .002), SSLD (P = .003), TSA (P = .01), and BRAF-mutant CRC (P = .01), compared with NC (Figure 10, F). In contrast, CD103+ cells were not augmented in conventional adenoma (P = .40).
Considering that previous reports showed that the number of CD8+CD103+ T cells decreases from ascending colon to rectum
and NC location could cause confounding bias, we additionally collected biopsies of colorectal healthy mucosa from ascending colon, transverse colon, descending colon, and rectum (n = 6 for each group), and performed the immunohistochemical staining of CD103 (ITGAE). The results showed that the number of CD103+ cells decreased from ascending colon to rectum (P < .001) (Figure 11, A–B). However, numbers of CD103+ cells in HPs (45.1 ± 15.7), SSLs (58.0 ± 19.7) and TSAs (42.9 ± 17.1) were still higher, even compared with normal mucosa from ascending colon (21.6 ± 2.2). Similar results were observed in the Human Protein Atlas database, which showed that CD103+ cells decreased significantly from small intestine to large intestine (Figure 11, C). However, the difference between colon and rectum was less striking. That is to say, the fluctuation of T cells in different locations of the large intestine do not contribute to our main findings of abundant CD103+ Trm cells in serrated polyps.
Figure 11CD103+Trm in SLs showed cytotoxic activity and might be regulated by retinol metabolism.A, Immunohistochemistry staining of CD103+ cells in colorectal healthy mucosa from ascending colon, transverse colon, descending colon, and rectum. B, Comparison of numbers of CD103+ cells from proximal colon to distal colon tissues. C, Expression of ITGAE/CD103 in histological sections of duodenum, small intestine, colon, and rectum, from Human Protein Atlas (HPA) database (https://proteinatlas.org). D, FCM analysis in activated CD8+ T cells showed increased GZMB in CD103+ cells, and ITGB2+ cells expressed higher levels of GZMK and GZMB. E, FCM analysis showed increased GZMK, GZMB, and KLRG1 in ITGB2+ CD69+ CD8+ T cells (GZMK+ Effector-1 cells) within SSL tissues. F, GSVA scores of Trm CD8+ gene signature in BRAF wild-type (n = 461) and BRAF mutant CRCs (n = 51). Data were obtained from GSE39582. G, GSVA scores of Trm CD8+ gene signature in BRAF wild-type (n = 634) and BRAF mutant CRCs (n = 64). Data were obtained from TCGA. H, Distribution of GSVA scores of Trm CD8+ gene signature in MSI-H BRAF-mutant (n = 71), MSS BRAF-mutant (n = 33), MSI-H BRAF-wildtype (n = 85), MSI-H BRAF-wildtype CRCs (n = 882). Data were integrated from both GSE39582 and TCGA, with batch effect corrected. I, The correlation between estimated CD8+ T subpopulations and recurrence-free survival in 263 patients with BRAF/KRAS-mutant CRC. The hazard ratio and the P value was evaluated by the Cox proportional hazard model with 95% confidence interval (CI). J, Metabolic activities at the single-cell resolution estimated by scMetabolism in each CD8+ T cell subpopulation. ∗P < .05; ∗∗P < .01; ∗∗∗P < .001.
Then, we determined the immunological features of colon CD8+ Trm cells by flow cytometry (FCM) in another 8 fresh SSL tissues. In activated CD8+CD69+ T cells of SSL, we found that CD103+ Trm expressed higher levels of GZMB compared with CD103– cells, whereas ITGB2+ Trm (GZMK+ Effector-1 cluster in our study) expressed higher levels of GZMK, GZMB, and KLRG1 (Figure 11, D–E), consistent with the scRNA-seq data. These results suggested that CD103+ and ITGB2+ Trm cells have higher cytotoxic potential compared with non-Trm cells in the CD8+ T cell responses of SLs.
Further, using 2 external cohorts, we found the Trm CD8+ gene signature
was significantly enriched in BRAF-mutant CRCs in TCGA (both P < .001) (Figure 11, F–G), suggesting its close association with serrated pathway. Subgroup analysis stratified by BRAF and MSI status showed higher Trm signature score in MSI-H BRAF-mutant (n = 71), compared with MSS BRAF-mutant CRCs (n = 33; P = .004), as well as a trend of higher Trm signature score in MSS BRAF-mutant, compared with MSS BRAF-wildtype CRCs (P = .08) (Figure 11, H). We also investigated the prognostic significance of these CD8+ T cell populations in 263 patients with CRC with BRAF or KRAS mutation, using the signature genes for each CD8+ T cell type identified from scRNA-seq (Supplementary Table 11). The CD103+ Trm and CD8+ cycling clusters were associated with favorable outcomes (P = .005 and P = .001, respectively), whereas GZMK+ Effector-1 showed potential protective effect (P = .058) with marginal significance (Figure 11, I).
Immunometabolism was shown to be a promising target for improving anti-tumor immunity.
for quantifying single-cell metabolism of CD8+ T cells. As shown in Figure 11, J, we identified remarkable activation of retinol metabolism in CD103+ Trm cells, but this phenomenon was not observed in the other activated CD8+ T cells. To test the effect of retinol metabolism on differentiation from naïve T cells to CD103+ Trm, we added retinoic acid (RA), the product of retinol metabolism, in a well-established in vitro Trm differentiation system (Figure 12, A; see Materials and Methods section). We found that addition of both 10 nM and 100 nM RA significantly induced more CD69+CD103+ T cells (Figure 12, B–D). We also observed significant upregulation of SDR16C5 (a rate-limiting retinol dehydrogenase in the biosynthesis of RA) in the epithelia of SLs but not in conventional adenomas (Figure 12, E–G), as well as in BRAF-mutant CRCs (Figure 12, H), suggesting a potential role of epithelium-derived RA in the modulation of the TME in SLs.
Figure 12Effects of RA on in-vitro generation of human CD8+CD69+CD103+T cells.A, Diagram of in-vitro Trm generation assay, with addition of 10 nM or 100 nM RA. B, The gating strategy used for flow cytometry analysis of CD103+ Trm cells. C, Effects of RA on in-vitro generation of human CD8+CD69+CD103+ T cells. FCM analysis of sequential exposure to IL-15 and TGF-β induced CD8+ Trm cells was performed for 3 times. RA (10 nM or 100 nM) was added to the culture medium; expansion of CD8+CD69+CD103+ T cells was determined on day 6. A representative scatter plot of an experiment was shown (gating on CD8+ T cells). D, Statistical analysis of the expansion of CD8+CD69+CD103+ T cells in each group. E, Increased SDR16C5 expression in SSL compared with conventional adenomas (left), as well as in Epi-SL compared with other epithelial cell types (right). F, Expression of SDR16C5 in different global cell types showed it mainly expressed in the epithelia. G, Quantitative RT-PCR analysis of SDR16C5 expression in 7 fresh SSL tissues compared with paired normal tissues. H, SDR16C5 expression in BRAF-mutant (n = 64) and wild-type (n = 634) CRCs in TCGA database. ∗P < .05; ∗∗P < .01; ∗∗∗P < .001.
Altered Myeloid and B Cell Composition Across the SL Spectrum
We identified 2 macrophage clusters and 2 dendritic cell (DC) clusters in SLs (Figure 13, A–C). Among 2 CD14+ macrophage subgroups, FCN1+ inflammatory macrophages expressed various proinflammatory genes, including FCN1, VCAN, S100A8, S100A9, and IL1B, whereas another subgroup of anti-inflammatory macrophages expressed typical M2 markers such as MRC1, CD163, C1QA, APOE, SELENOP, and MAF. The 2 DC clusters, namely conventional and plasmacytoid DC, were designated by high expression levels of CD1C, CLEC9A, FCGR1A, and MZB1, respectively. GSVA analysis demonstrated that inflammatory response-related pathways were enriched in FCN1+ inflammatory macrophages, whereas these were downregulated in the anti-inflammatory macrophage group. The 2 DC clusters showed enrichment of pathways related to antigen presentation and hypoxia (Figure 13, D–E and Supplementary Table 12).
Figure 13Myeloid cell subclusters in SLs.A, UMAP plot showing Seurat clusters of myeloid cells. B, UMAP plot for myeloid cells, color-coded for 4 subclusters (left) and proportion of myeloid cell types in NC and SLs (right). C, Feature plot of markers for each myeloid cell subcluster. D, Difference in pathway activity calculated by GSVA per cell between different myeloid cell types. E, Heatmap showing the canonical markers of each myeloid cell cluster. F, Heatmap representing the TF activity in each myeloid cell subtype. The TF activity was scored using AUCell. G, UMAP plot of myeloid cells color coded according to the expression of ETV5 and the area under the recovery curve (AUC) of the estimated regulon activity of ETV5, which corresponded to the degree of expression regulation of its target genes. ∗∗∗P < .001.
Interestingly, expansion of M2-like anti-inflammatory macrophages was observed in SLs, particularly in SSL with dysplasia and TSA. Moreover, the proportion of cDCs was increased in SSL, but decreased in TSA. To further investigate the TFs regulating this myeloid differentiation process, we performed single-cell regulatory network inference and clustering (SCENIC) analysis (Figure 13, F). Several known and novel TFs regulating anti-inflammatory macrophage (M2) differentiation were identified. In accordance with previous reports,
MAF, SPIC, NR1H3 (LXRα) were identified as candidate TFs for anti-inflammatory macrophages. We also identified some TFs whose roles in macrophage differentiation are not known, such as ETV5 (Figure 13, G). These results provided potential targets for inhibiting or reversing the aberrant polarization of immunosuppressive macrophages in serrated TME.
As shown in Figure 14, A–B, we clustered B cell population (CD79A+) into plasma cells (SDC1+) and non-plasma B cells (CD19+MS4A1+). In large intestine, we found that the vast majority of plasma cells showed high expression of IGHA1, and thus we named this cluster IgA+ plasma cells (Figure 14, C and Supplementary Table 13). In SSL and TSA, we found decreased proportion of non-plasma B cells, but higher proportion of IgA+ plasma cells (Figure 14, B). Besides the production of antibodies, tumor-infiltrating B cells also produce cytokines and facilitate antigen presentation and formation of tertiary lymphoid structures via CXCL13-CXCR5 axis.
Compared with non-plasma B cells in normal tissue, naïve mature B cells and memory B cells in SSL and TSA displayed diminished expressions of TNF-α and CXCR5, implying weakened antitumor secretion and recruitment ability of tumor-infiltrating B cells to form tertiary lymphoid structures. We also found higher expression of Absent in melanoma 2 (AIM2), an inflammasome component, in B cells within serrated adenoma. Because a recent study has uncovered the role of AIM2 regulating B-cell differentiation toward antibody-producing plasma cell via PRDM1 and XBP1,
upregulation of AIM2 may contribute to increased proportion of plasma cells in SLs.
Figure 14Dissection and clustering of B cells in SLs.A, UMAP plot showing Seurat clusters of B cells. B, UMAP plot for B cells (left) and proportion of B cell subclusters in NC and SLs (right). According to classic B cell markers, we further identified naïve mature B cells (IGHD+), memory B cells (CD27+), and plasmablasts (CD38+). C, Feature plot of markers for each B cell subcluster (upper) and violin plot showing different expression level of chemokine and cytokine in naïve mature and memory B cells between NC and SL samples (lower). D, Network of potential interactions between IgA+ plasma cells and myeloid or CD8+ T cells. E, Confocal staining of IgA and MDK in B cells in NC and SSL tissue. F, Heatmap representing the TF activity in each B cell subtype. The TF activity is scored using AUCell. G, UMAP plots of B cells color coded according to the expression of ATF5 or the AUC of the estimated regulon activity of ATF5, which corresponded to the degree of expression regulation of its target genes. H, Gating strategy for FCM validation of IgA+ PC and non-plasma B cells in SL. I, Validation of higher ATF5 expression in IgA+ plasma cells in SSL tissues by FCM.
To investigate altered cell-cell communications in SLs, we inferred the cell-cell interaction between plasma cells and CD8+T/myeloid cells in normal tissue and serrated adenomas, respectively, using Cellchat. As shown in Figure 14, D, putative gain of midkine (MDK)-related interactions was observed between IgA+ plasma cells and other immune cells (particularly myeloid cells) in SLs. MDK was reported to promote cell growth and proliferation and promote immunosuppressive macrophage differentiation;
we further validated the increased MDK expression on IgA+ plasma cells in SSL by confocal staining (Figure 14, E). A recent study also reported a cluster of IgG4+ B cell clones that produce MDK and VEGFA to promote angiogenesis.
The cytokine production function of B cell clusters in serrated adenoma and adenocarcinoma warrant further investigation.
To identify underlying regulators for increased IgA+ plasma cells in serrated adenoma, SCENIC analysis among 4 types of B cells was performed (Figure 14, F). For IgA+ plasma cell, several well-established TFs were identified (eg, CREB3, XBP1, CREB3L2, and PRDM1). Of note, ATF5 was also found to be an underlying TF regulating plasma cell differentiation, which has not been reported before. Both ATF5 and 65 ATF5 targets were predominantly expressed in IgA+ plasma cell, as depicted in the UMAP plot (Figure 14, G). To validate the SCENIC results, we isolated B cells and plasma cells from the SSL tissue samples by FCM (see gating strategy in Figure 14, H) and found significantly higher median fluorescence intensity of ATF5 in plasma cells (P = .0009) (Figure 14, I).
PDGFRA+ Fibroblasts Accumulate in Precancerous SLs and BRAF-mutant CRCs
We detected 15 sub-cell types from 5139 stromal cells (Figure 15, A–C), including 5 fibroblast types, 4 endothelial cell (EC) types, and 4 minor stromal cell types (enteric glial cells, pericytes, vascular smooth muscle cells, and endothelial progenitor cells) with distinct marker genes (Figure 15, D and Supplementary Table 14). Fibroblasts (COL3A1+) in SLs were clustered into 5 subpopulations: PDGFRA+ fibroblast, myofibroblast (TAGLN, ACTA2, and MYL9), Stromal-1 (high level of SFRP2, a canonical Wnt inhibitor), Stromal-2, and Stromal-3 (expressed markers of ADAMDEC1 and APOE, suggesting activated EGF signaling pathway; Stromal-3 also manifested high CCL4 and CCL5 expression). Four types of EC (ENG+ and PECAM1+) were assigned: tip-like EC, stalk-like EC, EC-1, EC-2. EC-1 existed predominantly in NC with upregulation of ICAM1, whereas EC-2 was mainly observed in SLs. Stalk-like EC showed lower expression of CD36, a key modulator of fatty acid metabolism.
Figure 15Accumulation of PDGFRA+fibroblasts in SLs.A, UMAP plot of 5,139 stromal cells in serrated lesions, color-coded by sub-cell type (left). Comparison of the average cell number and proportion of sub-cell types in normal mucosa and SLs (right) was shown. B, UMAP plot showing Seurat clusters of stromal cells in SLs. C, UMAP plot showing group distribution of stromal cells in SLs. D, Heatmap of canonical marker genes for fibroblasts and ECs. The color of each square indicates the z-scaled average gene expression (cyan to red). E, Immunohistochemistry staining and quantitative analysis of PDGFRA+ cells in the stromal of normal, HP, SSL, TSA, and SSLD tissues. F, Transcriptome differences between PDGFRA+ fibroblasts and other fibroblasts. ∗∗P < .01; ∗∗∗P < .001.
Among the fibroblast populations, we observed progressive accumulation of PDGFRA+ fibroblasts with the disease spectrum of the serrated pathway, from 2.67% in normal colon, 13.46% in HP, 23.32% of SSL, 28.84% in TSA, to 32.05% in SSL with dysplasia. PDGFRA+ fibroblasts featured high expression of PDGFRA, POSTN, BMP5, WNT5A, and WNT5B, indicating activation of WNT and BMP pathways, which play roles in regulating the proliferation and function of colon stem cells.
Immunohistochemistry analysis confirmed an increased number of PDGFRA+ stromal cells along the ‘NC-HP-SSL-SSLD’ sequence (Figure 15, E).
Interestingly, consistent with these findings, a previous study in HBUS mouse model of serrated polyps demonstrated that PDGFRA+ fibroblasts contributed to serrated polyp development by increasing MMP3 expression and facilitating membrane-bound HBEGF cleavage.
Interleukin 1 beta and matrix metallopeptidase 3 contribute to development of epidermal growth factor receptor-dependent serrated polyps in mouse cecum.
However, in our human samples, we found these fibroblasts did not secrete MMP3 (Figure 15, F and Figure 16, A). We compared the gene expression profile between PDGFRA+ fibroblasts and the other fibroblasts (the full list of DEGs is summarized in Supplementary Table 15) and found that MMP11 was highly expressed in PDGFRA+ fibroblasts (Figure 15, F and Figure 16, A). In other words, different from observations in mouse model, PDGFRA+ fibroblasts secrete MMP11 to promote HBEGF cleavage and SL development in human serrated polyps. PDGFRA+ fibroblasts also displayed high expression of periostin (POSTN), a secreted extracellular matrix protein that was shown to correlate with immunosuppressive premetastatic niche formation and tumor progression.
of three right-sided BRAF-mutant CRCs. As expected, similar PDGFRA+ subpopulation was also found in BRAF-mutant CRCs, which co-expressed high level of MMP11 and POSTN (see Figure 16, B).
Figure 16Stromal cell dynamics in normal colonic tissue and colorectal tumors derived from serrated pathway.A, Feature plot showing expression of MMP3 and MMP11 in stromal cells of SLs. B, Identification of MMP11-secreting PDGFRA+ fibroblasts in the stromal cells of BRAF-mutant CRCs (upper). Feature plot showing expression of PDGFRA, MMP11, and POSTN in stromal cells of BRAF mutated CRCs (bottom). Data were retrieved from Samsung Medical Center cohort reported in a previous study.
C, Dot plot showing cell-cell interaction of CXCL12-CXCR4 axis between stromal cells and immune cells. D, Pathway enrichment analysis of ECs in SLs vs NCs. E, Diagram illustrating early molecular and immune alterations underlying the serrated pathway inferred from scRNA-seq.
Moreover, Stromal-2 and Stromal-3 cell proportions were also increased in the 3 types of serrated adenomas, whereas Stromal-1 cell proportion decreased across the SL spectrum. We estimated the interaction scores of prominent CXCL12-CXCR4 axis between various stromal cells and immune cells and found the strongest interactions with Stromal-1 cluster, which predominantly existed in normal tissues. The interaction was weakened in PDGFRA+ fibroblasts, suggesting attenuated chemotaxis function of fibroblasts in SLs (Figure 16, C).
Because metabolism of EC has recently been recognized as a driver of tumor progression and angiogenesis,
we investigated the metabolic pathway alterations of ECs in SLs. Compared with ECs in normal samples, ECs in SL samples showed significant enrichment of oxidative phosphorylation/TCA cycle and fatty acid metabolism pathways (Figure 16, D). Active citric acid cycle and fatty acid metabolism may provide more energy for biomass synthesis and angiogenesis in proliferating tumor ECs.
Together, early molecular and immune alterations underlying the serrated pathway inferred from our scRNA-seq are illustrated in Figure 16, E.
Discussion
In the scRNA-seq landscape of serrated pathway to CRC, we concentrated on several driving questions, including: (1) what is the major characteristic of SL epithelia compared with normal colonocytes; (2) what is the heterogeneity of epithelia in SSL and TSA; (3) what is about the immune contexture in SL, and are there any metabolic factors that regulate the disease-specific tissue T cells; and (4) do any cell subpopulations expand in the microenvironment as the disease progresses to SSLD?
We found obvious antioxidant responses in the transformed epithelium of premalignant colon serrated neoplasms. Notably, the cellular responses to oxidative stress were markedly higher in SLs and BRAF-mutant CRCs, compared with conventional adenomas and BRAF wild-type CRCs. These findings suggested that redox imbalance might play a distinct role in serrated tumorigenesis. In fact, ROS is a double-edged sword in cancer initiation and progression.
Mutant BRAF induces DNA strand breaks, activates DNA damage response pathway, and up-regulates glucose transporter-1 in nontransformed epithelial cells.
Therefore, CRC originating from the serrated pathway could be treated differently than conventional CRC via modulation of oxidative stress. Additionally, oxidative stress-related genes could serve as biomarkers for tumors originating from the serrated pathway. We also noticed that SerpinB6 was specifically upregulated in epithelia of SLs, which has been shown to serve as protease inhibitor of cathepsin G in neutrophils.
In the in vitro experiments, SerpinB6 was shown to promote proliferation and migration and suppress ROS production; its exact function in colon serrated tumorigenesis remains to be elucidated.
SSL and TSA originate from 2 types of epithelium showing both similar and different molecular alterations, suggesting therapeutic interventions could be developed based on molecular subtypes of SLs. In epithelium of SSL, we observed aberrant overexpression of a direct intracellular marker of proliferation (Ki67) along with Notch signaling activation. Notch signaling, which maintains stemness by blocking differentiation of intestinal stem cells to secretory lineages,
was found associated with SSL formation and SSLD progression in our study. Consistent with previous experimental studies in both BrafV637E/+ and KrasG12D mouse models, Notch signaling activation has been shown to play a vital role in the progression and invasion in serrated colon tumorigenesis.
For serrated CRC developed from SSL, γ-secretase inhibitor may have potential therapeutic value. Different from SSL epithelium which showed gastric metaplasia,
we discovered the TSA epithelium manifested a trend of Paneth cell differentiation. The abundant lysosome secretion in TSA may have complex interaction with mucosal microbiota,
and participate in TSA-associated carcinogenesis. Interestingly, a small fraction of Epi-SL cells was observed in normal tissue; These cells featuring high expression of MIK67 and OLFM4, suggesting their stem-cell like properties and proliferation potential. They might represent precursors of transformed cells as well as stem cells. In fact, a small amount of transformed cells can also exist in healthy tissue, but this population can be restrained and removed by immune system surveillance as well as epithelial defense under healthy conditions.
Our data showed increased CD8+ T cell activation in the entire spectrum of SLs (HP, SSL, SSLD, and TSA), which was confirmed by a recent scRNA-seq research highlighting cytotoxic immunity in SSL (MSS status) prior to MLH1 silencing and neoantigen exposure that usually occurred after the SSLD stage.
Furthermore, our study is the first to point out that the cytotoxic immunity in SLs resulted from increased fraction of CD103+ Trm cells expressing GZMB and IFNG. In CRC, CD103+ Trm was demonstrated to be associated with better prognosis and represented a large portion of CD8+ lymphocytes in cancerous tissue but not in normal mucosa.
Consistent with our findings that CD103+ Trm was more enriched in colonic neoplasms originated from serrated pathway, increased number of CD103+ Trm cells was found in BRAF-mutant tumors.
We also observed that RA metabolism might be involved in CD103+ Trm augmentation in serrated tumors, which could be a target of chemoprevention. Moreover, CD8+ CD103+ T cells were reported to predict response to PD-L1 blockade therapy in several studies,
Patients with serrated CRC, including at least a portion of BRAF-mutant CRC but with MSS, might be more likely to benefit from immunotherapy due to ample infiltration of Trm. This hypothesis was justified by a recent study in which the response rate of combined BRAF inhibitor and immunotherapy reached 42% in 12 patients with BRAF-mutant CRC with MSS.
The mechanistic details of immunometabolism in SLs and the clinical application of immunotherapy in Trm-abundant serrated CRC warrant further well-designed studies.
In previous studies, anti-tumor T cell responses against CRC were shown to increase at the early stage but decrease with tumor progression.
were already increased, even in the premalignant stage. The imbalanced features of the TME may have implications for development of therapies targeting serrated colon tumors by preserving the cytotoxic function of Trm while suppressing the immunosuppressive cells.
Stromal components can also play vital roles in serrated tumor formation.
Interleukin 1 beta and matrix metallopeptidase 3 contribute to development of epidermal growth factor receptor-dependent serrated polyps in mouse cecum.
In this study, we demonstrated PDGFRA+ fibroblasts enrichment in the microenvironment of all kinds of human serrated tumors, most evident in SSLD. Thus, PDGFRA expression in the stromal may serve as a biomarker for progression from SSL to SSLD. These findings were in line with previous experimental evidence of pathological PDGFRA+ fibroblasts in a mouse model of HP. These pathological PDGFRA+ stromal cells specifically expressed MMP11 and POSTN, which might activate EGFR signaling by promoting cleavage of HBEGF and promote tumor metastasis by facilitating immunosuppressive microenvironment, which represents a candidate intervention target for serrated CRC.
Some limitations of our study should be considered while interpreting the results. First, SSLD samples are uncommon in clinical practice, which is partly attributable to the rapid transformation of SSL to invasive cancer; therefore, we could only obtain 2 samples in this group, although we validated the findings in more slides of SSLD by immunohistochemistry. Second, considering limited tissue material of SL available (<1 cm in size), mutation status was not tested in our study, although a number of studies have shown 60% to 80% of SSL harbor BRAF mutation and approximately 50% of TSA harbor KRAS mutation.
Finally, in-depth mechanistic studies need to be carried out in preclinical models of SL to validate and transform our discoveries.
In summary, our study unveiled the early molecular traits and immune responses across the 4 classes of SLs at single-cell resolution. Our findings help enhance the understanding of the pathogenesis of serrated neoplasia pathway and provide insights for further investigation of novel therapeutic targets.
Materials and Methods
Patient Samples
We obtained 19 human colorectal samples from Renji Hospital, School of Medicine, Shanghai Jiao Tong University for scRNA-seq. Patient clinicopathological data and endoscopy reports of lesion size were recorded. Each participant has provided written informed consent. Ethical approval was obtained from Shanghai Renji Hospital Ethics Committee (KY2021-004).
Single-cell Suspensions for SL Samples
Fresh serrated polyp tissues resected from colonoscopy were promptly washed with phosphate buffered saline (PBS) and divided into 2 parts. One part of specimen was fixed in formalin and sent for pathological evaluation. Formalin-fixed paraffin-embedded (FFPE) biopsy tissue with hematoxylin-eosin staining were interpreted by an experienced gastrointestinal pathologist who was blinded to patients' clinical and endoscopic data. The pathological diagnosis of HP, SSL, SSLD, and TSA were according to the 2019 World Health Organization classification for serrated lesions (fifth edition).
Another part of tissue was stored in MACS Tissue Storage Solution (cat#130-100-008, Miltenyi Biotec) at 4 °C, which was then used to generation of single-cell suspension within 24 hours.
To disassociate the serrated polyp tissue, we collected the tissue in a 6-cm plate added with RPMI 1640 medium supplemented with 10% fetal bovine serum (FBS) and kept on ice. Then the tissue was cut into small pieces of about 1 mm2 and transferred into a 1.5-mL EP tube. We added complete media (supplemented with 10% FBS) and resuspended to 1 mL, followed by adding 20 μg/mL DNaseI (D8071-100mg, Solarbo), and 200 μg/mL TL enzyme (Liberase TL Research Grade, cat#5401020001, Sigma). The tube was placed on a preheated rocking bed, shaking at 250 rpm at 37 °C, and digested for 1 hour. We turned the tube upside down several times every 15 minutes. The resulting suspension was added with an additional 200 uL FBS and passed through a 70-μm sieve. Then, it was subjected to centrifugation at 500g, 4 °C for 5 minutes. The supernatant was discarded, and then the precipitate was lysed the with 1 mL Red Blood Cell Lysis Buffer (cat#00-4300-54, Invitrogen) for 10 minutes at room temperature. After washing cells with PBS twice (centrifugation at 600g, 4 °C for 5 minutes), we collected cells form the bottom and resuspended to 500 uL with PBS. Then, single-cell suspension was filtered with a 40-μm sieve and counted using a cell count plate.
Library Preparation and Sequencing
Cells were counted, and cell density was adjusted to that recommended for 10× Genomics Chromium single-cell 3' v3 processing and library preparation. Sequencing was performed on an Illumina platform (NovaSeq 6000), by GENERGY BIO (Shanghai, China), at a sequencing depth of approximately 50,000 reads per single cell. The sequencing data were deposited at the Genome Sequence Archive database under the accession number HRA002611.
Single-cell Gene Expression Quantification
Single-cell 3' libraries were constructed using a commercial 10× Genomics platform (10× Genomics, Pleasanton, CA). Single-cell transcriptomic amplification and complementary deoxyribonucleic acid (cDNA) library preparation were performed by GENERGY BIO (Shanghai, China) using 10× Genomics Chromium single-cell 3' v3 according to the manufacturer’s instructions. The cDNA libraries were sequenced on NovaSeq 6000 System (Illumina, San Diego, CA), at a sequencing depth of about 50,000 reads per single cell.
Pre-processing of scRNA-seq sequencing data was conducted with Cell Ranger (version 3.0.2, 10× Genomics) using default parameters. The raw sequencing data were demultiplexed from raw base call (BCL) files and converted into fastq files using the ‘mkfastq’ command. The fastq files were aligned to the Genome Reference Consortium Human Build 38 (GRCh38). Counting all confidently mapped reads, filtered gene-barcode matrices were generated using the ‘count’ command for each sample, with low unique molecular identifiers (UMIs) droplets removed for minimizing the number of empty droplets. Low-quality cells were removed if cells had either fewer than 1001 UMIs, over 4000 or less than 1001 expressed genes, less than 5% UMIs from ribosomal genes or over 25% UMIs derived from the mitochondrial genome. We filtered genes that were detected in less than 3 cells. Canonical correlation analysis was used to remove the batch effect from 18 patients in which every patient was assigned to one batch. We annotated cell components using a combination of reference based (R package ‘SingleR’, version 0.99.10) and manual annotation.
Clustering of Single-cell Data Matrix
We performed the clustering of various cell types using R package ‘Seurat’ (version 4.0.3). Data normalization was completed via ‘NormalizeData’ function, with the default scaling parameter of 10,000 and log normalization method. ‘FindVariableGenes’ function was employed to find 4000 genes of the highest variance. We standardized the data with the ‘ScaleData’ function. After performing principal component analysis using highly variable genes, the top 20 principal components and a resolution of 0.5 were selected for the following cluster analysis and visual dimensionality reduction by UMAP for dimension reduction. We used the ‘FindAllMarkers’ or ‘FindMarkers’ function to determine the marker genes of each cluster compared to all other clusters or to specific cluster(s). ‘FeaturePlot’ or ‘VlnPlot’ function was used for presentation of gene expression levels. We labeled the obtained clusters as epithelial cells, T cells, B cells, myeloid cells, or stromal cells through known classic markers (epithelial cells: EPCAM; T cells: CD3D; B cells: MS4A1, CD79A; myeloid cells: CD14, CD163; stromal cells: COL1A1, PECAM1) and finally analyzed each of these clusters separately to identify the finer clusters by repeating the above operations. Cell-cycle phases were predicted using a function included in Seurat that scores each cell based on expression of canonical marker genes for S and G2/M phases.
Copy Number Variation Estimation
The ‘InferCNV’ R package (version 1.8.1) was used to estimate the copy number variations (CNVs) in epithelial cells based on transcriptome profiles (inferCNV of the Trinity CTAT Project, https://github.com/broadinstitute/inferCNV).
The CNV scores of epithelial cells from the normal samples was as the reference cells in the final estimation of CNVs. For the inferCNV analysis, the default parameters were used: ‘True’ for denoise choice, ‘True’ for hidden Markov model setting, and ‘0.1’ for signal cutoff. To reduce the possibility of false-positive CNV calls, the default Bayesian latent mixture model was implemented to identify the posterior probabilities of alterations in each cell. Low-probability CNVs were filtered using the default value of ‘0.5’ for the threshold. To determine the clonal CNV changes in each SL sample, the ‘subcluster’ method was utilized on the CNVs generated by the hidden Markov model. GRCh38 cytoband information was used to convert each CNV to a p- or q- arm level change for simplification based on its location. Each CNV was annotated to be either a gain or a loss. The CNV score of each cell was calculated as quadratic sum of CNV region.
Differentially Expressed Gene and Pathway Analysis
To gain a biological understanding of the identified DEGs, we conducted GO and KEGG analyses, which were conducted through the ‘clusterProfiler’ R package (version 3.18.1) with the gene set background from Molecular Signatures Database (version 7.4).
Volcano plot for DEGs were generated by ‘ggplot2’ R package (version 3.3.3). In single-cell sequencing analysis, DEGs between cell clusters were defined as those with adjusted P < .05 and |log2FC| > 0.25.
Gene Set Variation Analysis Among Different Cell Types
A subset-enriched term of biological relevance was represented as enrichment score using the corresponding genes or markers and was calculated by the gene set variation analysis approach through using ‘GSVA’ R package (version 1.36.2).
Heatmap for biological pathways of interest were generated by ‘ComplexHeatmap’ R package (version 2.5.5). Evaluation of metabolic activity per single cell was performed using the recently-proposed ‘scMetabolism’ R package, based on VISION, AUCell, and ssGSEA methods.
We explored the epithelial cell transition from transit amplifying cells and intermediate cells to SL-specific cells (Epi-1 and Epi-2) using R package ‘Monocle’ (version 2.20.0). The standard workflow was employed to identify high variance gene between these epithelial cell subtypes using ‘dispersionTable’ function. Next, the function ‘reduceDimension’ with the ‘DDRTree’ method was employed to reduce the dimensions, and the cells were then ordered based on the predicted pseudotime via ‘orderCells’ function.
SCENIC Analysis
pySCENIC (version 0.11.0) was used to identify gene regulatory networks and key TFs in scRNA-seq datasets using default parameters. In the first step, given a predefined list of TFs (hs_hgnc_tfs.txt, https://github.com/aertslab/pySCENIC/tree/master/resources), co-expression modules between these TFs and putative target genes were inferred using Arboreto algorithm (GRNBoost) from normalized expression, followed by pruning the insignificant co-expression pair. Secondly, the function ‘ctx’ will measure the enrichment of the TFs motif in putative regulatory regions of the target gene based on the cisTarget databases and the motif annotation tables (https://resources.aertslab.org/cistarget/). The TF and its predicted targets will be retained if the motif of the TF is significantly enriched in one of its modules. Finally, the activity score of the predicted regulons in the individual cells was calculated using AUCell algorithm, which measures the enrichment of target genes in whole-genome ranking.
To measure the specificity of TFs for clusters, the regulon specificity score was calculated as the Jensen-Shannon divergence between the distribution of AUCell score and clusters. TFs with top regulon specificity score will be kept as cluster-specific TFs.
We predicted and visualized regulon networks in epithelial cell subclusters, as reported in a recent study.
The cisTarget step in SCENIC produced a list of the TF and all its target genes in each regulon, with the prediction importance from GRNBoost as the corresponding weight of each TF-target pair. The list was transformed to a target-by-TF weight matrix via ‘recast’ function, which characterized each regulon based on its transcriptional weights on all its targets. Based on the Euclidean distance of weight matrix, k nearest neighbors of each regulon were found, with k selected as the square root of the number of regulons. TFs with similar transcription regulation profile (target genes and corresponding weight) was grouped together using laiden algorithm at the resolution_parameter of 2. In each laiden cluster, the adjacency between each pair of regulon was set as the weight of the TF-target pair if there is a TF-target relationship between them, generating a transcription regulatory network. The network is visualized using R package ‘igraph’ (version 1.2.10), with the color of vertices as mean scaled AUCell score and the width of edges as adjacency between regulons in this cluster.
Receptor-ligand Pairing Analysis
To qualify the cell-cell communication differences, we applied R package ‘CellChat’ (version 1.1.0) to infer and analyze the intercellular communication network on normal and SL samples separately. Briefly, a CellChat object was created from raw count matrix and cell-type annotation. Inferred interactions were filtered using ‘subsetCommunication’ function with the default parameter. The strength of cell interactions was quantitated as the mean log-normalized count of the ligand gene in the source cell type and receptor gene(s) in the target cell type.
Public Database Analysis
We obtained external bulk RNA-seq data from GSE76987 (right-sided colon tissues of 10 normal, 10 conventional polyps, 21 SSLs), EMTAB7960 (22 SSL tissues and 20 SSL with dysplasia), and GSE39582 (51 patients with BRAF-mutant CRC, 217 patients with KRAS-mutant CRC, and 257 patients with wild-type CRC). TCGA (http://cancergenome.nih.gov) CRC data were obtained by R package ‘TCGAbiolinks’ (version 2.20.1).
For integrated analysis of TCGA and GEO, the transcriptome profiles of patients with COAD and READ from the TCGA cohort were downloaded using R package TCGAbiolink and transformed as logarithm ‘Transcripts Per Kilobase Million (TPM)’. The expression profiling by microarray from GSE39582 was log-transformed. The samples were removed if they did not have matched information about BRAF mutation or microsatellite instability. And then an integrated analysis of the 2 cohorts was conducted with batch effect removed by the function ‘combat’ in R package ‘sva’ (version 3.40.0). We accessed the scRNA-seq data of 3 right-sided BRAF-mutant CRCs (SMC-03, 10, 17) from Samsung Medical Center dataset in a previous study.
FFPE tissue blocks of 4 normal colon tissue, 4 HPs, 8 SSLs, 5 TSAs, 4 TAs, 9 SSLs with cytological dysplasia, and 7 BRAF-mutant right-sided CRCs from Shanghai Renji Hospital were used for immunohistochemistry staining. All FFPE colonic samples were stained with hematoxylin and eosin for histopathological evaluation at the Shanghai Renji Hospital. For immunohistochemistry staining, HP, SSL, and SSLD samples were collected from proximal colon, NCs were collected from ascending (n = 1), transverse (n = 1), and descending colon (n = 2), whereas TSA samples were collected from distant colon or rectum. In brief, after antigen retrieval, sections of colonic samples were blocked with goat serum for 30 minutes. Then, the sections were incubated overnight with a primary antibody against GSTP1 (1:250, ab138491, Abcam), NQO1 (1:1000, ab28947, Abcam), DUOX2 (1:1000, NB110-6157655, NOVUS), 8-OHdG (1:100, ab48508, Abcam), SERPINB6 (1:200, 14962-1-AP, Proteintech), JAG1 (1:200, ab109536, Abcam), MYC (1:100, ab32072, Abcam), FOXP3 (1:200, ab4728, Abcam), ITGAE (1:500, ab224202, Abcam), and PDGFRA (1:500, ab134123, Abcam) at 4 °C followed by incubation with a horseradish-peroxidase-conjugated secondary antibody at room temperature for 30 minutes. Then, 3,3-diaminodbenzidine substrate was added. Five random fields were analyzed under a light microscope.
Immunofluorescence staining was performed on human SL tissues. Briefly, NC, SSL, and TSA tissues were fixed in 4% paraformaldehyde, paraffin-embedded, and sectioned. Slides were washed 3 times in PBS, followed by blocking with 10% goat serum for 20 minutes. Sections were then boiled in sodium citrate solution to retrieve antigens. Then, sections were incubated overnight with primary antibodies against MUC2 (1:200, ab231427, Abcam) and LYZ (1:250, ab108508, Abcam); CD8 (1:50, ab33786, Abcam) and ITGAE (1:500, ab224202, Abcam); IgA (1:100, ab124716, Abcam) and MDK (1:50, sc-46701, Santa Cruz) at 4 °C. After washing 3 times in PBS, slides were incubated with fluorescent dye-labeled secondary antibodies (1:500; Invitrogen, Carlsbad, CA) for 30 minutes at room temperature. Slides were then stained with DAPI (Southern Biotech, Birmingham, AL) after washing 3 times with PBS. Images were captured under a confocal fluorescence microscope (LSM-710; Carl Zeiss, Jena, Germany).
Cell Line Experiments
RKO cells were obtained from American Type Culture Collection and were cultured in RPMI1640 medium (Gibco) supplemented with 10% FBS. The cell line was verified by short tandem repeat profile analysis. For siRNA transfection, RKO cells were transfected with 50 pmol siRNA (siNC or siSERPINB6; TSINGKE Biotechnology, Shanghai, China) for 6 hours, then used for different assays after 48 hours. The sense and antisense sequences of the SERPINB6 siRNA are as follows: 5'-CAAAUCUUGGUGCUUCCAUAU -3' and 5'-AUAUGGAAGCACCAAGAUUUG -3', respectively. For SerpinB6 rescue experiment, RKO cells were transfected with FuGENE HD transfection reagent (Promega, Madison, WI) after 4 hours of siSERPINB6 transfection. We transfected 0.6 ug SerpinB6 plasmid DNA with 0.6 uL FuGENE Reagent per 6-well plate for another 2 hours. Transfection medium was replaced by cell growth medium 6 hours after transfection. Proteins from cells were extracted as previously described.
The follow primary antibodies were used for Western blotting: anti-SERPINB6 (14962-1-AP, Proteintech, IL) and anti-GAPDH (KC-5G5, Kangchen Company, Shanghai, China).
Cell proliferation was measured by EdU assay kit (C10338-1, RiboBio, Guangdong, China) according to the manufacturer’s instructions. For transwell assay, a total of 2 × 105 transfected RKO cells were seeded into the upper chamber of a 24-well polycarbonate transwell filter (8 μm pore size, Corning Incorporated, Corning, NY). RPMI1640 media containing 30% FBS was placed in the lower chamber. After 48 hours of incubation, the migration cells adhering to the lower surface were fixed with 4% paraformaldehyde and stained with 0.1% crystal violet. The area of positive staining was calculated in 5 random fields using Image-Pro Plus software (version 6.0, Media Cybernetics). ROS levels were detected using Reactive Oxygen Species Assay Kit (S0033, Beyotime, China) with DCFH-DA probe. Briefly, DCFH-DA was diluted at 1:1000 with serum-free RPMI1640 medium and incubated at 37 °C for 20 minutes. Then, cells were washed in serum-free RPMI1640 medium and collected for FCM analysis. N-acetyl L-cysteine (ST1546, Beyotime, China) was used as a ROS scavenger.
Quantitative Real-time Polymerase Chain Reaction
Total RNA was extracted from 8 SSL samples with paired normal tissues collected from Shanghai Renji Hospital, using TRIzol reagent (Invitrogen, Carlsbad, CA). The cDNA was synthesized using PrimeScript RT reagent Kit (Takara, Japan). Quantitative RT-PCR was performed using SYBR-Green Master mix (Takara, Japan) according to the manufacturer’s instructions. The expression levels of genes were analyzed with GAPDH serving as an internal control. The primer sequences used were as follows: SERPINB6 (forward primer: AGGGAAACACCGCTGCACAGAT; reserve primer: GTGCCAGTCTTGTTCACTTCGG); SDR16C5 (forward primer: TGCACGCCTATACCTGCGATTG; reserve primer: GGCATTGTTGATTAGGATGGAAAC); GAPDH (forward primer: GTCTCCTCTGACTTCAACAGCG; reserve primer: ACCACCCTGTTGCTGTAGCCAA).
Immune Cell Isolation of Human Colon Adenomas
To validate the findings of T and B cells in the SL microenvironment, we obtained single cell suspension of another 8 fresh SSL tissues resected from colonoscopy resection, as described in the first part of the Supplementary Methods. The histological diagnosis was confirmed by 2 experienced pathologists blinded to the endoscopic imaging independently. To isolate immune cells, after filtered through a 70-μm sieve and centrifugation, the cells at the bottom were resuspended and then overlaid onto 33% Percoll (GE Healthcare, Chicago, IL). The suspension was then centrifuged at 900 g for 30 minutes at room temperature with the brake off. After the spin, red blood cells were lysed as previously described. The isolated leukocytes were then washed twice before staining.
In Vitro Generation of Human CD8+CD69+CD103+ T Cells
We used an in-vitro model to generate CD103+ Trm cells from human human peripheral blood mononuclear cells (PBMCs, as previously reported.
In brief, PBMCs were cultured in 24-well plates at 5 × 105 cells /mL in RPMI1640 medium supplemented with 10% heat-inactivated FBS, 1% penicillin/streptomycin (Gibco), and 50 mM 2-mercaptoethanol (Invitrogen). To enrich CD103+ Trm cells, we added 2 ng/mL rhIL-2 (PeproTech, Cranbury, NJ) for the whole course, and added rhIL-15 (50 ng/mL; R&D Systems, Minneapolis, MN) for 3 days, followed by a further 3 days of rhTGF-β (50 ng/mL; PeproTech). The medium was refreshed on day 3. RA (HY-14649, MedChemExpress) at the concentration of 10 nM and 100 nM was used for investigating the effects on the expansion of CD8+CD69+CD103+ T cells. The cells were collected for FCM analysis on day 6 (see Figure 12A).
Flow Cytometric Analysis
Cells were harvested and washed with FACS buffer. For surface staining, the 1 × 106 cells in 100-μL single-cell suspension were stained with antibodies for 30 minutes in the dark at 4 °C. The following antibodies were used: live and dead cell stain FVS (BD Biosciences), anti-CD3, anti-CD4, anti-CD8, anti-CD103/ITGAE, anti-CD19, anti-CD27, anti-IgD, anti-CD38, and anti-CD138 (BD Biosciences, San Jose, CA); and anti-CD69, anti-CD18/ITGB2, anti-KLRG1 (BioLegend, San Diego, CA). For intracellular cytokine staining, immune cells were incubated in complete RPIM1640 containing 10% FBS and leukocytes activation cocktail with GolgiPlug (BD Biosciences) at 37 °C for 5 hours, and then surface-stained cells were fixed and permeabilized with Cytofix/Cytoperm solution (BD Biosciences) for 20 minutes at 4 °C. Subsequently, cells were stained with anti-GZMB, anti-GZMK, and anti-IgA (BD Biosciences) for 30 minutes at 4 °C. For staining TFs, cells were fixed and permeabilized with Transcription Factor Buffer Set (BD Biosciences). Briefly, after surface staining, cells were fixed and permeabilized for 50 minutes at 4 °C, then washed and stained with ATF5 (Abcam, Cambridge, UK) for 50 minutes at 4 °C. FCM was carried out using Celesta (BD Biosciences) and analyzed using FlowJo software (version 10.6.2, Tree Star).
Survival Analysis
We explored the prognostic relevance of CD8+ T cell types, using signature genes derived from our scRNA-seq data. The patient data used in survival analysis were obtained from GSE39582, incorporating 263 patients with CRC with BRAF or KRAS mutation and available survival information. The selected criteria of signature genes for one cell type were DEGs detected more than 25% of the target cell type, with adjusted P value < .05 and the threshold of log2(fold change) at 0.8, according to previous studies.
The survival analysis was performed by Cox regression model, using recurrence-free survival as the outcome. The analysis was done by R package ‘survival’ (version 3.2-3) and ‘survminer’ (version 0.4.8).
Statistical Analysis
All statistical analyses were performed using R version 4.1.0 (https://www.r-project.org/), Python version 3.7.10, and Graphpad version 7.0. For all statistical tests, a 2-sided P value less than .05 was regarded as statistically significant. Comparisons between 2 groups were performed by the unpaired Student t test. Data were expressed as mean ± standard deviations. Scripts and code for the analysis are available at Github (https://github.com/xlucpu/scRNAseq-SL).
Acknowledgment
The authors thank Dr Qiwei Qian and Dr Nana Cui from the Shanghai Institute of Digestive Disease for their kindly technical assistance.
Traditional serrated pathway-associated colorectal carcinoma: morphologic reappraisal of serrated morphology, tumor budding, and identification of frequent PTEN alterations.
Interleukin 1 beta and matrix metallopeptidase 3 contribute to development of epidermal growth factor receptor-dependent serrated polyps in mouse cecum.
Mutant BRAF induces DNA strand breaks, activates DNA damage response pathway, and up-regulates glucose transporter-1 in nontransformed epithelial cells.
Conflicts of interest The authors disclose no conflicts.
Funding This work was supported by the National Natural Science Foundation of China (81973145, 82273735), the Active Components of Natural Medicines and Independent Research Projects of the State Key Laboratory in 2020 (SKLNMZZ202016), the Key R&D Program of Jiangsu Province [Social Development] (BE2020694), Health Technology Project of Pudong New District Health Commission (PW2020D-12), Science and Technology Commission of Shanghai Municipality (Grant No.19411951606), and the Program for Promoting Advanced Appropriate Technology of Shanghai Health Commission (2019SY003).
00212-0&id=gr1.jpg){kind=link}
00212-0&id=gr2.jpg){kind=link}
00212-0&id=gr3.jpg){kind=link}
00212-0&id=gr4.jpg){kind=link}
00212-0&id=gr5.jpg){kind=link}
00212-0&id=gr6.jpg){kind=link}
00212-0&id=gr7.jpg){kind=link}
00212-0&id=gr8.jpg){kind=link}
00212-0&id=gr9.jpg){kind=link}
00212-0&id=gr10.jpg){kind=link}
00212-0&id=gr11.jpg){kind=link}
00212-0&id=gr12.jpg){kind=link}
00212-0&id=gr13.jpg){kind=link}
00212-0&id=gr14.jpg){kind=link}
00212-0&id=gr15.jpg){kind=link}
00212-0&id=gr16.jpg){kind=link}