Automatic sketching and functional classification model of three-level lymphatic structure related ecological niche based on space transcriptome and application

Abstract

The invention belongs to the technical field of biomedicine, and discloses an automatic sketching and functional classification model of three-level lymphatic structure related ecological niche based on a space transcriptome and application. Aiming at the problem that traditional TLS typing depends on maturity and neglects spatial localization and function association, by integrating spatial transcriptome and pathological region analysis, TLS is divided into three functional subtypes: microbial response type TLS-G1, localization in a mucous membrane region, high expression of immunoglobulin A genes and mucous membrane defense genes; the inflammation-related bystander TLS-G2 is distributed in a para-tumor tissue and is used for enriching memory T cells and mast cells; and the anti-tumor reaction type TLS-G3 is positioned in a tumor core region and is used for highly expressing cytotoxic genes and immune checkpoint molecules. The method comprises the following steps: defining a TLS-flare feature gene set; the score of the gene set can predict immune checkpoint inhibitor treatment response, and is used for developing an immunotherapy prediction kit or optimizing a targeted therapy strategy.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of biomedical and bioinformatics, and particularly relates to an automatic delineation and functional classification model of a tertiary lymphoid structure (TLS) related niche based on a spatial transcriptome and an application thereof. BACKGROUND

[0002] The tertiary lymphoid structure (TLS) is an immune aggregate in the tumor microenvironment (TME) that is structurally and functionally similar to a secondary lymphoid organ. It participates in the anti-tumor immune response by recruiting and activating immune cells. The presence of TLS is closely related to the response and prognosis of immunotherapy in tumor patients. Esophageal squamous cell carcinoma (ESCC) is an ideal model for studying human multi-stage malignant transformation, and has clear pathological staging characteristics from normal tissue to metastatic cancer. Distinguishing different inflammatory environments and mucosal immune states in the stages of precancerous lesions, carcinoma in situ and invasive carcinoma helps to analyze the unique functional niche they shape. Correspondingly, as a key immune hub, the tertiary lymphoid structure (TLS) can respond to both chronic inflammation and tumor occurrence, but its specific function may depend on the local microenvironment it is in. Therefore, developing a method that can quickly and accurately identify TLS regions in tissue sections and use spatial transcriptomic characteristics to identify and functionally classify these TLSs has important application value and will provide a new perspective for designing stage-adaptive immune regulation strategies.

[0003] In recent years, the development of single-cell and spatial transcriptome technologies has provided tools for analyzing the heterogeneity of TLS, but existing technologies still have the following limitations: (1) TLS identification method is not accurate: Traditional TLS identification mainly relies on pathological morphology observation or limited gene set scoring, lacking integrated analysis of spatial coordinate information and immune cell composition, resulting in blurred definition of the core and peripheral regions of TLS, affecting the accuracy of subsequent functional analysis.

[0004] (2) TLS functional classification is not clear: Existing researches mainly focus on the mature state of TLS (such as mature / immature), without fully considering the spatial niche differences in different stages of tumor evolution, thus unable to effectively distinguish different TLS functional subtypes (such as microbe response type, chronic inflammation related type or anti-tumor oriented type).

[0005] (3) The response prediction efficiency of immunotherapy is insufficient: The currently reported TLS-related gene sets (such as TLS_9, TLS_12, TLS_50) have limited predictive ability for the response of immune checkpoint inhibitors (such as anti-PD-1 / CTLA-4 antibodies) treatment, and lack of extensive verification across cancer types, making it difficult to meet the needs of clinical precise screening of patients who benefit from treatment.

[0006] Therefore, it is urgent to establish a TLS functional zoning model integrating spatial transcriptomic features, cell composition and functional phenotypes, and to develop a feature gene set with high predictive performance, to provide new tools for the precise implementation of tumor immunotherapy. SUMMARY

[0007] The application provides an automatic delineation and functional classification model of a spatial transcriptome-based tertiary lymphoid structure (TLS) related niche and application.

[0008] The application is implemented by the following technical solutions: an automatic delineation and functional classification model of a spatial transcriptome-based tertiary lymphoid structure (TLS) related niche, which divides the TLS into three functional subtypes corresponding to different niche zoning: Microbial response type (TLS-G1): It is located in the mucosa and responds to microbial exposure as a mucosal immune sentinel, revealing the niche zoning of TLS in mucosal surveillance; it highly expresses immunoglobulin A related genes (IGHA1, JCHAIN, PIGR) and mucosal defense genes (MUC5B, DEFB1), and lowly expresses interferon stimulated genes, cytotoxic effect molecules and immune checkpoint molecules; its core and peripheral zones are enriched with plasma cells, proliferative MKI67+B cells, CCR7+initial B cells, GPR183+resident B cells and memory T cells, CCR7+initial T cells are highly expressed in the peripheral zone and have the highest number of BCR unique sequences; Inflammation-related bystander type (TLS-G2): It is located in the peritumoral area, reflecting the "bystander" state related to chronic inflammation, revealing the niche zoning of TLS in chronic inflammation; it has the highest inflammation response score but the expression of classic inflammatory genes is not significantly increased; it is enriched with memory CD8 T cells (TCF7+ CD8 cells), memory CD4+ T cells (ITGB1+ CD4 T, GPR183+ CD4 T cells) and mast cells, showing the characteristics of memory cell retention driven by chronic inflammation; the number of BCR unique sequences remains high; Anti-tumor response type (TLS-G3): It is located in the intratumoral area, represents the effector hub involved in tumor antigen-specific immunity, and reveals the niche partitioning of TLS in tumor-specific immunity; cytotoxic activity, T cell exhaustion, phagocytosis and graft rejection related genes are significantly up-regulated, with high expression of cytotoxic effector genes (GNLY, GZMB, PRF1), immune checkpoint molecules (PDCD1, CTLA4), interferon stimulated genes (STAT1, ISG15), and chemokines (CXCL10, CXCL11, CCL21, CXCL13); its core and peripheral regions are significantly enriched with ISG15+B / T cells, IFNG+CD4 T cells, neutrophils and CXCL10+inflammatory macrophages; the number of unique sequences of its BCR is significantly reduced, accompanied by clonal expansion of tumor-specific B cells.

[0009] The TLS-G3 subtype presents a "flare" like radial distribution in space, which is named TLS-flare and constitutes an anti-tumor immune niche with spatial organization; the TLS-flare characteristic gene set contains 75 genes, and the gene set consists of: CXCL10, CHI3L1, HMOX1, GBP5, ADAMDEC1, DUSP4, SH2D1A, MMP19, MX2, FAM20A, IL2RG, PLEK, CST7, SRGN, GIMAP4, LY75, CCL19, ADAM19, AIF1, IL32, SLAMF1, IFI44L, NLRC3, SAMD9L, GBP2, CD2, IL27RA, RASSF4, TGFB3, IKZF1, SELPLG, GPRIN3, LCK, ISG15, SAMSN1, RAB8B, IKZF3, PAG1, SLA, SPN, CD68, HCK, CD48, SLC2A6, CD96, GMFG, PSTPIP1, ST8SIA4, SOCS1, ADAM8, HCLS1, P2RY8, MPP1, SCIMP, TTYH2, DOK2, CPNE5, LST1, PREX1, CCL4, CD40, LMCD1, PSMB9, JAK1, SASH3, ADAMTS4, ZAP70, STK10, SP100, LTF, DOK3, TBC1D2B, PSMB8, PEAK1, TRIM22.

[0010] The application also provides an automatic delineation and functional typing method of a three-level lymphoid structure (TLS) related niche based on spatial transcriptome, comprising the following steps: (a) Obtain esophageal cancer tissue spatial transcriptome data, score the data based on the published TLS_50 gene set, which contains the following genes: FDCSP, CR2, CXCL13, LTF, CD52, MS4A1, CCL19, LINC00926, LTB, CORO1A, CD79B, TXNIP, CD19, LIMD2, CD37, ARHGAP45, BLK, TMC8, CCL21, PTPN6, ATP2A3, IGHM, SPIB, TMSB4X, CXCR4, NCF1, CD79A, ARHGAP9, DEF6, EVL, TBC1D10C, RASAL3, INPP5D, RNASET2, RASGRP2, TNFRSF13C, RAC2, CD22, ARHGEF1, AC103591.3, TRAF3IP3, HLADQB1, CD53, ARHGAP4, TRBC2, POU2AF1, TRAF5, OGA, FCRL3, HLA-DQA1; determine the sites with scores exceeding 90% quantile of all spatial coordinates sites as candidate TLS positions; (b) Based on the results of step (a), define the area containing more than or equal to 5 spatially continuous candidate detection points as the TLS core region; and extend the distance of 5 points outward from the edge point of the core region, define the spatial sites within this range as the TLS peripheral region; (c) Perform unsupervised clustering on the gene expression profile of the TLS core region, and divide the three functional subtypes according to the molecular characteristics: TLS-G1, TLS-G2, and TLS-G3 according to claim 1; (d) Verify the functional characteristics of each subtype obtained in step (c) by analyzing cell composition through RCTD deconvolution, and combining tissue morphological characteristics and differential gene expression analysis; (e) Extract the differentially expressed genes of TLS-G3 and other subtypes, and screen the genes that meet the following conditions to define them as the TLS-flare characteristic gene set: (e.1) The expression amount of the gene in the TLS-G3 subtype is significantly higher than that in other subtypes, and the specific standard is that the log2 ratio of the differential expression fold of the gene between the TLS-G3 subtype and other subtypes is greater than 0.5, and the gene is expressed in not less than 25% of the TLS-G3 spatial sites; (e.2) There is a potential difference in expression between the immunotherapy response group and the non-response group, and the Wilcoxon rank sum test p value is less than 0.2; The gene set that meets the above conditions is the TLS-flare characteristic gene set.

[0011] The application also provides a method for constructing an immunotherapy response prediction model using the model, and the TLS-flare gene set is used to construct the model by the following steps: (a) Based on the tumor patient's transcriptome sequencing data, calculate the GSVA enrichment score of the TLS-flare feature gene set; (b) According to the median value of the GSVA score of all patients, the patients are divided into high-score group and low-score group; wherein the high-score group corresponds to the patients with good response to clinical immunotherapy (CR+PR). The application also provides the use of the immunotherapy response prediction model constructed by the method in the preparation of a product for predicting the treatment response and survival period of esophageal squamous cell carcinoma (ESCC) patients to PD-1 / CTLA-4 inhibitors.

[0012] The application also provides an immunotherapy response prediction system, comprising: The spatial positioning module: performs steps (a)-(b) of the automatic delineation and functional partitioning method of the spatial transcriptome-based tertiary lymphoid structure TLS-related niche, for determining the TLS core area and the peripheral area; The gene set generation module: performs steps (c)-(e) of the automatic delineation and functional partitioning method of the spatial transcriptome-based tertiary lymphoid structure TLS-related niche, for dividing the TLS functional subtype and extracting the TLS-flare gene set; The scoring prediction module: performs steps (a)-(b) of the method for constructing an immunotherapy response prediction model, for calculating the GSVA score and predicting the immunotherapy response.

[0013] Clinical value of TLS-flare score: it can not only predict the efficacy of immunotherapy, but also provide new insights for TLS targeting therapy strategy of mucosal cancer. The "bystander" type TLS with chronic inflammation characteristics also has predictive potential, and the mucosal flora-induced immune aggregates may have higher organ conservation.

[0014] TLS functional framework based on spatial ecological context (rather than static mature stage): By comparing mucosal (TLS-G1), paraneoplastic (TLS-G2) and intratumoral (TLS-G3) TLS, we revealed their niche partitioning in mucosal surveillance, chronic inflammation and tumor-specific immunity. This spatial diversity suggests that TLS is not a fixed immune hub, but a structure dynamically shaped by antigen availability, local microenvironment and chronic tissue signals. Another point of interest is TLS-G2, which is distributed in the peritumoral and distal regions, showing chronic inflammation characteristics but lacking clear tumor-oriented activation. This "bystander" state may reflect immune readiness but insufficient antigen stimulation, suggesting that it can be reprogrammed into a neoantigen-responsive anti-tumor TLS through treatment, thus becoming an early biomarker of high-risk lesions and a potential target for tumor antigen vaccines, providing new ideas for chemoprevention and postoperative adjuvant therapy. BRIEF DESCRIPTION OF DRAWINGS

[0015] Figure 1 For sampling and experimental flowchart; Figure 2 For esophageal squamous cell carcinoma (ESCC) multi-modal atlas constructed by single-cell transcriptome sequencing (scRNA-seq), single-cell TCR sequencing (scBCR-seq), and spatial transcriptome sequencing; (a) Distribution of each global cluster cell subpopulation in the scRNA-seq dataset; (b) Unsupervised clustering of spatial transcriptome data based on ISCHIA algorithm; (c) Bar chart showing the composition of niche types at different disease stages; (d) Representative images of hematoxylin-eosin staining (H&E) and spatial deconvolution; (e) Representative images of cell subpopulations and niche types in spatial sections displayed by PhenoCycler-Fusion protein detection technology; Figure 3 For spatial microenvironment type analysis based on transcriptome data; (a) Spatial pathological staging of esophageal squamous cell carcinoma development process; NOR (normal epithelium); BCH (basal cell hyperplasia); mD (mild dysplasia); MD (moderate dysplasia); SD&CIS (severe dysplasia / carcinoma in situ); ICA (invasive carcinoma); MCA (metastatic carcinoma); Scale bar: 100 μm in the top panel; 2 mm in the rest; (b) Microenvironment composition in epithelial / tumor-dominated regions at different disease stages; (c) Disease stage-specific composition corresponding to the microenvironment region; (d) Heatmap showing the cell composition difference of 15 CCs; (e) Microenvironment type definition of cell composition category (CC); Figure 4The image shows grade 3 lymphoid structures (TLSs) in esophageal squamous cell carcinoma (ESCC). In the figure: (a) H&E staining images of grade 3 lymphoid structures (TLSs) in different disease stages; red dashed circles indicate TLSs; Top image: tissue section from a patient with high-grade dysplasia; Middle image: tissue section from a patient with TLSs in adjacent tissue; Bottom image: tissue section from a patient without TLSs in adjacent tissue; Scale bar: the first image in each figure is 3 mm, and the remaining images are 300 μm; (b) Representative images of matched H&E staining and spatial deconvolution results, indicating TLSs; dashed lines and asterisks ( (c) Box plots show the CC8 niche types in different disease stages, with statistical significance determined by the Wilcoxon rank-sum test; (d) Representative images of PhenoCycler-Fusion multiplexed panoramic tissue markers, showing TLSs in matched slices of homologous samples, corresponding to extended data plot b; the same numerical markers represent different but matched homologous sample slices. Solid lines indicate TLSs, and boxes of the same color represent the same specimen sample. The first image on the left is a panoramic view, and the right image is a magnified view of it (or them); Figure 5 The authentication process for TLS in a three-level lymphatic structure; Figure 6 (a) A bar chart showing the number of TLS for each patient and the proportion of TLS in the three groups. (b, c) Representative H&E images containing TLS. b, panoramic view; c, magnified view. The same numbers indicate TLS of the same type. Figure 7 This heatmap illustrates the functional atlas of the three groups of tertiary lymphoid structures and their peripheral regions in esophageal squamous cell carcinoma. It displays the artificial annotation types (aTLS, iTLS, nTLS, pTLS) for the three TLS_G1 / G2 / G3 lesions; patient origin; tissue origin (N: normal esophagus, Adj: tumor margin, T: tumor); B-cell receptor diversity and mutational characteristics (BCR); and immune-related cell pathway scores. Figure 8Fig. 6. Spatial distribution of TLS-G1 / G2 / G3 and matched gene sets scores in esophageal squamous cell carcinoma (ESCC) sections. (a) Spatial localization distribution of TLS-G1 / G2 / G3 and matched FCER2, CXCL13, CXCL9, CR1, TAP1, STAT1, ISG15, SAA1 gene and TLS-flare gene set scores in esophageal squamous cell carcinoma (ESCC) sections. (b) Gene signature scores of TLS-G3 (named as TLS-flare) in each specified subgroup of spatial sites. (c) Functional differences of TLS-G3 and TLS-G2. (d) Spatial localization map of specified niche types (CCs). (e) Venn diagram showing the number of BCR unique sequences in tumor tissue, peritumoral tissue and distal normal tissue. Spatial: BCR sequences obtained by TRUST4 algorithm in spatial transcriptomic data. scBCR: BCR data from single-cell BCR sequencing. Figure 9 Fig. 7. Heatmap of expression levels of specific genes in three groups of tertiary lymphoid structures (TLS) and their surrounding areas. Figure 10 Fig. 8. Boxplot of abundance distribution of specific cell subpopulations in three groups of TLS areas (a) and surrounding areas (b). Figure 11 Fig. 9. Boxplot of sequence species number and mutation frequency of BCR three chains in three groups of TLS areas (a) and surrounding areas (b). Figure 12 Fig. 10. Differential gene expression and pathway enrichment analysis of three groups of TLS. Mfuzz algorithm was used to analyze the dynamic expression pattern of different TLS and their corresponding peripheral areas, and GSVA algorithm was used for pathway enrichment analysis. Figure 13Clinical significance of TLS-flare signature in esophageal squamous cell carcinoma PD-1 immunotherapy cohort: (Based on bulk RNA-seq data, analyze gene set difference between responders and non-responders, predictive performance and survival association) (a) Boxplot showing the difference in GSVA score of TLS-flare signature between response group and non-response group (Wilcoxon rank-sum test). (b) ROC curve of TLS-flare gene score predicting immunotherapy response (AUC=0.793). (c) Survival curve of patients in high / low TLS-flare group (difference in progression-free survival PFS, log-rank test). (d) Survival curve of patients in high / low TLS-flare group (difference in overall survival OS, log-rank test). (e) Boxplot of GSVA score of Non-TLS, TLS_G2, TLS_G1 signature groups (Wilcoxon rank-sum test). (f) ROC curve of TLS_G2 and TLS_G1 gene score (AUC=0.737 / 0.717). (g) Survival analysis of TLS_G2 and TLS_G1 grouping (difference in progression-free survival PFS, log-rank test). (h) Survival analysis of TLS_G2 and TLS_G1 grouping (difference in overall survival OS, log-rank test). Figure 14 Clinical significance of TLS_9, TLS_12 and TLS_50 signature in esophageal squamous cell carcinoma PD-1 immunotherapy cohort: (Based on bulk RNA-seq data, analyze gene set difference between responders and non-responders, predictive performance) (a) Boxplot showing the difference in GSVA score of TLS_9, TLS_12, TLS_50 signature between response group and non-response group (Wilcoxon rank-sum test). (b) ROC curve of TLS_9, TLS_12 and TLS_50 gene score predicting immunotherapy response (AUC=0.599 / 0.589 / 0.717). Figure 15Clinical significance of TLS-flare signature in melanoma cohort (PRJEB23709): (based on bulk RNA-seq data, analyzing gene set difference between responders and non-responders, predictive performance and survival association) (a) Boxplot showing the difference in GSVA score of TLS-flare signature between response group vs. non-response group (Wilcoxon rank-sum test). (b) ROC curve of TLS-flare gene score predicting immunotherapy response (AUC = 0.807). (c) Survival curve of patients in high / low TLS-flare group (overall survival OS difference, log-rank test). (d) Boxplot of GSVA score of Non-TLS, TLS_G2, TLS_G1 signature groups (Wilcoxon rank-sum test). (e) ROC curve of TLS_G2 and TLS_G1 gene score (AUC = 0.772 / 0.491). (f) Survival analysis of TLS_G2 and TLS_G1 grouping (overall survival OS difference, log-rank test). DETAILED DESCRIPTION

[0016] In order to make the objects, technical solutions and advantages of the embodiments of the present application clearer, the technical solutions in the embodiments of the present application will be described clearly and completely below. Obviously, the described embodiments are part of, rather than all of the embodiments of the present application; based on the embodiments in the present application, all other embodiments obtained by those of ordinary skill in the art without creative work fall within the protection scope of the present application.

[0017] Unless otherwise defined, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this application belongs, and the materials cited and referenced herein are incorporated by reference in their entirety. Equivalents to the described procedures and techniques described herein can be made and used by those of ordinary skill in the art without departing from the scope of the present application. The experimental methods in the following examples are routine methods unless otherwise specified. The instruments and equipment used in the following examples are routine laboratory instruments and equipment unless otherwise specified. The experimental materials used in the following examples are purchased from routine biochemical reagent stores unless otherwise specified.

[0018] I. Experimental Methods 1. Ethical Statement: This study was approved by the Institutional Review Board of Shanxi Tumor Hospital (Shanxi Tumor Institute) (Approval No. 2019LL245). All patients have signed a written informed consent form, agreeing to sample collection and data analysis.

[0019] 2. Tissue sample collection: Fresh tissue samples for single-cell sequencing and paraffin-embedded sections for multiplexed immunofluorescence staining and TLS / TIL (tumor-infiltrating lymphocyte) evaluation in this study were obtained from Shanxi Tumor Hospital (Shanxi Institute of Oncology) with approval from the institutional review board of the hospital. Informed consent was obtained from all patients. Four milliliters of peripheral blood was collected from each patient before surgery for PBMC isolation. The following samples were obtained after surgical resection: metastasis-positive lymph nodes (pLNs), metastasis-negative lymph nodes (nLNs), peritumoral tissue (adjacent tissue ≤1 cm from the tumor margin), tumor margin tissue (T1), tumor center tissue (T3), transition tissue between tumor margin and tumor center (T2), and distal normal tissue >5 cm from the tumor. The above tissues were divided into two parts: one part was used for tissue digestion to prepare cell suspension within 30 minutes after surgery; the other part was snap-frozen in liquid nitrogen and OCT-embedded for frozen section preparation (histopathological diagnosis) and 10x Genomics Visium spatial transcriptome sequencing. The cell suspension was divided into two parts: paired single-cell transcriptome sequencing and single-cell VDJ sequencing were performed using a 5' library kit; 8,000 live cells were taken from each sample for single-cell transcriptome sequencing (scRNA-seq), and human T cell receptor (TCR) and B cell receptor (BCR) amplification were performed according to the manufacturer's protocol using the 10X Genomics Chromium Next GEM Single-cell 5' Kit (version 2, item number: PN-1000263, PN-1000190, PN-1000256, PN-1000252, PN-1000253). Sequencing was performed on the Illumina NovaSeq-6000 platform.

[0020] 3. Single-cell RNA-seq data analysis S1, Quality control and alignment: The Cell Ranger toolkit (v2.1.0) provided by 10x Genomics was used to align the sequencing reads and generate gene-cell UMI matrices based on the GRCh38 reference genome. High-quality cell criteria were retained: the number of detected genes was ≥200, and the UMI count was between 500 and 50,000. Cells with mitochondrial gene expression accounted for >10% were excluded. The R package Harmony (v1.0) was used for batch correction of different patient samples.

[0021] S2, Unsupervised clustering, dimension reduction and identification of signature genes: Unsupervised graph clustering algorithm of Seurat v4 (version 4.1.1) was applied. 3,000 highly variable genes (excluding mitochondrial genes and immunoglobulin genes) were selected by FindVariableFeatures, and principal component analysis was performed by RunPCA. Based on the first 20 Harmony corrected components, the nearest neighbor graph was constructed using FindNeighbours, and the initial clustering was performed by FindClusters using Louvain algorithm (resolution = 1), and the clusters were annotated according to known markers. Finally, 10 main cell types were identified: 5 types of immune cells (NK / T cells, B cells, plasma cells, mast cells, myeloid cells) and 5 types of non-immune cells (epithelial / cancer cells, endothelial cells, fibroblasts, platelets, neurons). UMAP dimension reduction was performed by RunUMAP function to visualize the clustering results. The FindMarkers function of Wilcoxon rank-sum test was used to identify differentially expressed genes (DEGs), with the following parameters: logfc.threshold = 1, min.pct = 0.3, only.pos = T, max.cells.per.ident = 1000.

[0022] S3, Sub-cluster clustering of cell subpopulations: Secondary clustering was performed for different main cell types (parameter adjustments as follows): CD4+ T cells / B cells: first 20 Harmony components, resolution = 0.3; CD8+ T cells: first 20 Harmony components, resolution = 0.5; Endothelial cells / fibroblasts: first 20 Harmony components, resolution = 0.1; Epithelial / cancer cells: first 25 Harmony components, resolution = 0.8; Myeloid cells: first 25 Harmony components, resolution = 0.6 Secondary clustering DEG identification parameters: logfc.threshold = 0.6, min.pct = 0.2, only.pos = T, max.cells.per.ident = 500.

[0023] 4. Spatial transcriptome sequencing: OCT-embedded tissue blocks were sectioned on a Leica CM3050 microtome, and brightfield images were acquired at 20x resolution by a Leica Aperio Versa8 whole slide scanner. Agilent 2100 Bioanalyzer was used to assess the RNA quality of OCT-embedded blocks (RNA integrity number RIN > 7 were used for subsequent analysis). Visium Spatial Tissue Optimization Slide & Reagent Kit (Cat. No.: PN-1000193, PN-1000194) was used to perform tissue optimization (TO) according to the standard protocol. The tissue permeation time was optimized as follows: 16 minutes for distal normal tissue, 18 minutes for paracancerous tissue, 30 minutes for tumor tissue and lymph node. Visium Spatial Gene Expression Slide & Reagent Kit (Cat. No.: PN-1000187) was used to prepare spatial transcriptome library. Sequencing was performed on the NovaSeq PE150 platform (Illumina).

[0024] Spatial sample processing and visualization: All raw FASTQ files were processed using the spatial transcriptome analysis software SpaceRanger (version 1.1.0, 10x Genomics) with default parameters, and the transcripts were aligned to the GRCh38 human reference genome. The filtering criteria for spots were as follows: (1) total UMI number > 1000; (2) number of detected genes > 500; (3) mitochondrial gene proportion < 15%; (4) hemoglobin gene proportion < 5%. Data processing and visualization were completed by Seurat software (version 4.1.1)56. The raw counts of each section were normalized using the SCTransform method (Poisson distribution model) to show gene expression levels.

[0025] 5. Calculate feature scores in spatial data: The Seurat function AddModuleScore (default parameters) was used to calculate the scores. This function calculates the average expression level of the provided feature genes at the single spot level, subtracting the aggregate expression of the control feature set. All analyzed features were binned according to the average expression, and control features were randomly selected from each bin. For feature scores of cell modules, the top five genes with the highest differential expression were collected from each cell subpopulation as a gene set representing the features of each cell module. The Seurat function AddModuleScore was applied to calculate the average expression level of each program at the single spot level.

[0026] The signature genes for the main cell types are as follows: T cells (PTPRC, CD3D, CD3E, CD4, CD8A), NK cells (KLRB1, NCR1, FGFBP2, CX3CR1, NCAM1, NKG7, GNLY), B cells (CD79A, SLAMF7, BLNK, FCRL5, IGHG1, MZB1, SDC1), fibroblasts (DCN, SFRP4, LUM, COL1A1, PCLAF), endothelial cells (PECAM1, VWF), myeloid cells (MMP19, CD163, CD14, FCGR3A, FCGR3B, CD68, C1QA, C1QB, TPSAB1, TPSB2, MNDA, CSF3R, S100A9, XCR1, CLEC9A, CD1C, LAMP3, LILRA4).

[0027] 6. TLS gene set scores: The three published TLS gene set scores were calculated for each checkpoint using AddModuleScore, containing genes as follows: TLS_9: PTGDS, RBP5, EIF1AY, CETP, SKAP1, LAT, CCR6, CD1D, CD79B; TLS_12: CCL2, CCL3, CCL4, CCL5, CCL8, CCL18, CCL19, CCL21, CXCL9, CXCL10, CXCL11, CXCL13; TLS_50: FDCSP, CR2, CXCL13, LTF, CD52, MS4A1, CCL19, LINC00926, LTB, CORO1A, CD79B, TXNIP, CD19, LIMD2, CD37, ARHGAP45, BLK, TMC8, CCL21, PTPN6, ATP2A3, IGHM, SPIB, TMSB4X, CXCR4, NCF1, CD79A, ARHGAP9, DEF6, EVL, TBC1D10C, RASAL3, INPP5D, RNASET2, RASGRP2, TNFRSF13C, RAC2, CD22, ARHGEF1, AC103591.3, TRAF3IP3, HLA-DQB1, CD53, ARHGAP4, TRBC2, POU2AF1, TRAF5, OGA, FCRL3, HLA-DQA1.

[0028] TLS_9 and TLS_12 are from Ma Y et al. published in Integrating tertiary lymphoid structure-associated genes into computational models to evaluate prognostication and immune infiltration in pancreatic cancer. J Leukoc Biol. 2024 Sep 2;116(3):589-600., and TLS_50 gene set is from Wu R et al. published in Comprehensive analysis of spatial architecture in primary liver cancer. Sci Adv. 2021 Dec 17;7(51):eabg3750.

[0029] 7. Spatial cell type deconvolution analysis: To determine the cell composition of each spot, the RCTD algorithm in Cable DM et al. published in Robust decomposition of cell type mixtures in spatial transcriptomics. Nat Biotechnol. 2022 Apr;40(4):517-526. was applied to the ST dataset (spatial transcriptome dataset), which is implemented in the R package spacexr (v2.2.0). The RCTD full model was applied to each spatial transcriptome sample. For spatial transcriptome samples from different tissues, the clustering data of single-cell subpopulations from the corresponding tissue were used as reference. Then the cell deconvolution scores of each spot were normalized using the normalize_weights function.

[0030] 8. Identification of niche composition in spatial transcriptomic data: To resolve the cell composition in spatial transcriptomic data, the R language implementation of the ISCHIA algorithm in the published literature [Identifying Spatial Co-occurrence in Healthy and InflAmed tissues (ISCHIA). Mol Syst Biol. 2024 Feb;20(2):98-119.] by Lafzi A et al. was applied. By combining the pathological structure characteristics with the average clustering variance plot (k value range 2-20, determined according to the elbow rule), the optimal cluster number (k = 15) was determined using the Composition.cluster.k function. Then the Composition.cluster function was applied to 134,181 detected points, which were divided into 15 spatial ecotypes representing different cell compositions. The division of these ecotypes is based on the principle that "functionally related cell types exhibit specific spatial co-occurrence patterns in the tumor microenvironment" and is defined by analyzing the cell composition of the detected points under each composition state.

[0031] 9. Automatic delineation of spatial TLS core and peripheral regions: To analyze the TLS features at different stages of disease progression, a spatial transcriptomic-based TLS identification process (Figure Figure 5 ) was developed. It includes the following 5 steps: Tissue sample preparation: Obtain tumor and associated tissue samples (including distal tissue, tumor adjacent tissue, tumor tissue, metastatic lymph node, etc.) from different patients (e.g., Patient 1204, Patient 1119), and determine the spatial location and type of the samples to lay the foundation for subsequent analysis.

[0032] Cell clustering and typing: Spatial co-occurrence analysis using ISCHIA reveals the interaction network of cell types and ligand-receptor (LR) pairs in the tissue microenvironment, and then unsupervised clustering of spatial sites is achieved.

[0033] Calculation of spatial site TLS feature scores: According to the three published gene sets TLS_9, TLS_12 and TLS_50 described in Experimental Method 6, each site in the spatial transcriptome is scored using Addmodulecore, and the top 10% quantile of TLS_50 is used for quantitative display to evaluate the strength and distribution of TLS-related features in the tissue and determine the TLS core and peripheral regions.

[0034] Key cell and component validation: Analysis of key cell (including B cell, plasma cell, follicular dendritic cell, follicular helper T cell, vascular endothelial cell and fibroblast) signature score (B cell score, Plasma score, Follicular DC score, Tfh cell score, Endothelial cell score, CAF score) and validation of key TLS component (key TLS cell component) to confirm the TLS related cell composition and functional status from the perspective of specific cell population, and to assist in determining the TLS characteristics.

[0035] Based on the automatic delineation of TLS_50 score, the tertiary lymphatic structure and its surrounding area that meet the conditions on each spatial transcriptome slice and the scattered TLS_50 score high points are determined. For the areas that meet the TLS judgment standard, the tertiary lymphatic structures are numbered from large to small according to the area of the tertiary lymphatic structure. Combined with multi-dimensional data such as cell signature score and CAF (cancer associated fibroblast) score, the cell distribution, clustering type and signature score result of the tissue sample are comprehensively organized to identify and confirm the existence, type, function and role in the tumor microenvironment of the TLS.

[0036] The specific process integrates the histopathological features of H&E staining, unsupervised spatial clustering based on ISCHIA and published TLS gene features (TLS_9, TLS_12, TLS_50). The sliding quantile window method of TLS_50 signature score is used to select the top 10% (quantile_90) as the candidate TLS area. The spatial positioning of key TLS related immune and stromal components (including B cells, plasma cells, follicular dendritic cells, follicular helper T cells, T cells, endothelial cells and fibroblasts) by RCTD algorithm deconvolution is used for secondary screening of the candidate area. The aggregation with high composite score and containing ≥5 spatial transcriptome detection points is defined as TLS, and finally 80 TLSs are identified.

[0037] The automatic workflow of TLS core recognition and surrounding positioning is as follows: S1, preprocessing of spatial transcriptomic gene score and coordinate data: TLS feature thresholding: calculate TLS_50 score data for each site in each tissue section (img_name), mark the detection points with TLS_50 score ≥ the 90th quantile (quantile_90) as TLS_mark=TRUE (candidate core points); coordinate integration: extract spatial coordinates (imagerow, imagecol) by the GetTissueCoordinates function of the Seurat package in R language and merge with the metadata.

[0038] S2, neighborhood construction and spatial connectivity: define the identification conditions of adjacent points in the spatial hexagonal lattice, as follows: a. Hexagonal distance calculation: calculate the median of the minimum Euclidean distance between all spatial detection points (hex_distance), reflecting the hexagonal grid geometry of the Visium chip.

[0039] b. Adaptive radius neighborhood search: define the neighborhood radius as radius_multiplier × hex_distance (default radius_multiplier = 1.2). This parameter setting ensures: 1. Cover all nearest neighbors (distance = hex_distance), considering floating point rounding errors; 2. Exclude the second nearest neighbor (distance = √3 × hex_distance ≈ 1.732 × hex_distance) to prevent over-connection; 3. Identify neighborhood points through radius search (RANN::nn2).

[0040] S3, TLS core region identification: the purpose of this step is to find multiple spatial regions that meet the core point marking conditions in S1 in the spatial hexagonal lattice, and exclude scattered high-frequency regional sites, so as to ensure the accuracy of TLS region identification. The process is as follows: a. Graph-based clustering: construct an adjacency matrix by connecting TLS-marked positive points within ≤1.2 × hex_distance. Use igraph::graph.adjacency and igraph::components functions to identify connected components.

[0041] b. TLS core naming: define connected components with ≥5 detection points as TLS cores. Name them in descending order of the number of detection points as [img_name]-TLS[size_rank] (e.g. TF02ST021022-TLS1). Smaller clusters (<5 points) are classified as scattered TLS points and are not considered as target TLS.

[0042] S4, Expanding the perimeter region of each TLS outward from the spatial location of the detected TLS core region, steps as follows: a. Edge point identification: If a core point is adjacent to any non-TLS point (TLS_mark = FALSE) within its neighborhood radius (defined in step 3), mark it as an edge point (is_edge = TRUE).

[0043] b. Perimeter region localization: Expanding 5 times hex_distance distance outward from each edge point to capture the perimeter non-TLS points.

[0044] c. Shared perimeter region assignment: Assign the captured perimeter points to the nearest TLS core according to hex_distance. When there is a multiple core equidistance conflict, preferentially assign to the core with the most number of detection points.

[0045] S5, TLS spatial auto-delineation visualization: Visualize the core and perimeter points using ggplot2. Rotate the Cartesian coordinates by 180° to match the Visium orientation, and draw the cluster boundary by concaveman::concaveman(). TLS cores (labeled as TLS1, TLS2, etc. in descending order of size) are assigned unique colors, and the corresponding perimeter regions inherit the core color. The results are cross-validated with H&E morphology and RCTD cell type deconvolution results (Appendix Figure 5 ).

[0046] 10. Molecular typing with TLS gene expression: To resolve the functional heterogeneity of TLS clusters, perform unsupervised molecular typing based on the transcriptomic profile. First, process the spatial transcriptomic data within the TLS core points, and identify 2000 highly variable genes (HVGs) using the variance stabilization method (FindVariableFeatures). Select the top 200 HVGs to construct the TLS aggregate expression matrix (average gene expression within each spatial TLS cluster), and perform Z-score standardization before hierarchical clustering. Calculate the Euclidean distance between clusters using the Ward minimum variance method (hclust, method = "ward.D2"), and divide the dendrogram into three biological subtypes (TLS-G1 / G2 / G3, corresponding to the heatmap from left to right position) by dynamic tree cutting (cutree, k = 3). These molecular subtypes are mapped to the spatial coordinates by different colors, and their anatomical distribution in the TLS core region and related perimeter regions are visualized (Appendix Figure 7 ).

[0047] 11. Multidimensional feature profiling of TLS typing a. Gene set expression analysis: Calculate the average expression level of each group of detection points in TLS-G1-3 and its surrounding area in the gene set related to cytotoxicity, immune checkpoint, TLR receptor, B cell maturation, chemokine, chemokine receptor, TLS maturation, and anti-tumor, anti-viral mucosal immune activity, etc.

[0048] b. Mfuzz and enrichment analysis: Apply Mfuzz to detect the consensus expression trend genes of TLS-G1-3 and its surrounding area detection points. By observing the within-group sum of squares of different cluster numbers, 9 clusters are selected. The known functional genes in each cluster are annotated, and the fgsea package is used to perform gene set enrichment analysis based on the msigdbr database. The ClusterGVis (v0.1.1) can be used to visualize the clustering and enrichment results (Appendix Figure 12 ).

[0049] c. TLS subtype comparison characterization: Based on spatial resolution characteristics, the functional heterogeneity is quantified by the following dimensions: spatial characteristics (core / surrounding markers); BCR dynamics (clonal diversity / mutation load); pathway activity (region-specific averaging) (Appendix Figures 7-9 ).

[0050] d. Integrate these features into the annotation of hierarchical clustering heatmaps, associate molecular subtypes with functional landscapes. By comparing the cell type composition of RCTD deconvolution, analyze the proportion difference of immune cells of each TLS subtype and its matching surrounding area, and use box plots to show the spatially layered cell abundance patterns (Appendix Figure 10 ).

[0051] e. Spatial and single-cell data BCR profiling: Use the TRUST4 algorithm to perform BCR sequence analysis on the Visium spatial transcriptome, and define the same nucleotide sequence as a unique BCR sequence. Quantify the BCR indicators (IGH / IGK / IGL chain unique sequence number, mutation frequency) of each spatial detection point and integrate them into the TLS molecular typing heat map (Appendix Figure 7 , 11).

[0052] 12. Identification of TLS-flare gene set: Spatial regions were classified into seven categories (Non_TLS, TLS-G1 / G2 / G3, periTLS-G1 / G2 / G3) excluding scattered TLS positive points. Region-specific molecular signatures were identified by differential expression analysis between Non_TLS and TLS core regions (G1-G3). Ten rounds of independent analysis (1000 randomly selected points per group) were performed to screen for differentially expressed genes with log2fc > 0.5 and detected in > 25% of cells in a group (pct. 0.25 criteria). Genes detected in > 6 rounds of repeated analysis were defined as robust signatures, and the top 500 differentially expressed genes per group were selected by fold change. To enhance clinical relevance, candidate genes were evaluated for association with treatment response in independent ESCC immunotherapy cohorts, and genes with nominal P-value < 0.2 were retained to balance stringency and biological findings. From the top 500 spatially derived signature genes, only genes associated with treatment were selected, and genes with clinical relevance were retained, while genes with no clear association with treatment efficacy were excluded.

[0053] A representative gene set for TLS-G1 / G2 / G3 was finally selected, and TLS-G3 was named as TLS-flare gene set, with the corresponding genes as follows: TLS-G3 (TLS-flare): CXCL10, CHI3L1, HMOX1, GBP5, ADAMDEC1, DUSP4, SH2D1A, MMP19, MX2, FAM20A, IL2RG, PLEK, CST7, SRGN, GIMAP4, LY75, CCL19, ADAM19, AIF1, IL32, SLAMF1, IFI44L, NLRC3, SAMD9L, GBP2, CD2, IL27RA, RASSF4, TGFB3, IKZF1, SELPLG, GPRIN3, LCK, ISG15, SAMSN1, RAB8B, IKZF3, PAG1, SLA, SPN, CD68, HCK, CD48, SLC2A6, CD96, GMFG, PSTPIP1, ST8SIA4, SOCS1, ADAM8, HCLS1, P2RY8, MPP1, SCIMP, TTYH2, DOK2, CPNE5, LST1, PREX1, CCL4, CD40, LMCD1, PSMB9, JAK1, SASH3, ADAMTS4, ZAP70, STK10, SP100, LTF, DOK3, TBC1D2B, PSMB8, PEAK1, TRIM22; TLS-G1: SAA4, VTCN1, PIGR, MUC5B, AZGP1, SAA2, LTF, DEFB1, BPIFB1, PROM1, WFDC2, CHST9, SLC34A2, TFF3, TFF1, PPP1R9A, SPDEF, DMBT1, TF, MUC16, LOXL4, KRT23, SLC6A14, GOLM1, DEPTOR, CX3CL1, CLIC6, FOXQ1, LZTS3, ARHGEF37, HID1, FOXC1, TLE2, OVOL2, IRX1, CXCL17, HLA-DQA2, SYT8, SFRP1, PALMD, SELENBP1, PPM1L, MUC1, NFIB, EFNA5, CGNL1, IRX3, EVA1C, IGKC, SEMA3B, CREB3L2, TRPM4, IL20RA, CD1C, ERN1, PPP1R26, WNK2, ENPP4, POU2AF1, GPRC5A, ZSWIM4, CFH, MEG3, SH3BP4, ELF3, ACBD4, PRSS8, SCNN1D, STARD10, TCF7L1, COL9A3, SARM1, SCD5, FOXO3, DAPK2, SULT1A1, UNC13B, DHRSX, PTPRN2, PHYH, XRRA1, FGFR2, SLC17A9, COL9A2, NEIL1, ARSD, SMCO4, GLT8D1, FZD7, KLF9, CCNG1, RNF14, FUT3, MANSC1, PTGES; TLS-G2: TCF7, KIF21B, WDFY4, GIMAP7, NCF4, ZAP70, PASK, SASH3, P2RY8, FLT3LG, CD48, C1orf56, CCL19, KLHL6, POU2AF1, DOK3, BCAS4, LCK, GMFG, IKZF3, IGHG1, IFF01, PREX1, IKZF1, CPNE5, SLA, TMC6, ABTB1, PBXIP1, CD2, DENND6B, DOK2, STK10, SELPLG, HCLS1, SAMSN1, AN09, MZFl, PSTPIP1, 1-Mar, SLAMF1, DEF6, CST7, ADAM8, FCGR2B, C5orf56, USP20, LST1, RABEP2, STMN3, SBF1, TNFRSF18, SH2B3, 9-Mar, RALGPS2, GRK6, CCL4, CD40, CNPPD1, MICAL1, KCTD17, IL2RG, RBM38, ZNF335, IGKC, GALK1, GIMAP4, FLYWCH1, ST8SIA4, RASSF4, TRIM22, CNTRL.

[0054] 13. Validation of TLS-flare and other two TLS signature gene sets in clinical cohorts: To evaluate the clinical utility of spatially resolved signatures, a comprehensive validation was performed in two independent immunotherapy cohorts (ESCC and melanoma PRJEB23709). The expression matrix was converted into gene x sample format suitable for gene set variation analysis (GSVA), and signature scores for each spatial region (Non_TLS, TLS-G1, TLS-G2, TLS-G3) were calculated.

[0055] This study validated the clinical value of TLS signature gene sets (including TLS-flare, TLS-G2, TLS-G1 gene sets and published TLS_9 / 12 / 50 three gene sets) in two independent immunotherapy cohorts of esophageal squamous cell carcinoma (n=43) and metastatic melanoma (n=91). Specifically, it includes: 1. Data acquisition: The ESCC cohort was derived from an anti-PD-1 monotherapy study (DOI: 10.1016 / j.ccell.2022.12.004), and the melanoma cohort was derived from an anti-PD-1 ± anti-CTLA-4 combination therapy study (DOI: 10.1016 / j.ccell.2019.01.003), both of which contain RNA-seq data and records of treatment response (responders / non-responders), survival (PFS / OS).

[0056] 2. Feature scoring: The gene expression matrix obtained from RNA-seq sequencing of each cohort was converted into gene × sample format. Gene set variation analysis (GSVA) ​​was used to calculate different TLS feature scores for each sample, which correspond to the spatially defined TLS regions (Non_TLS, TLS-G1, TLS-G2, TLS-G3).

[0057] 3. Efficacy association: Nonparametric tests (Wilcoxon rank-sum test) were used to compare the differences in scores between responders and non-responders; predictive performance was assessed by receiver operating characteristic (ROC) curves and the area under the curve (AUC) was reported.

[0058] 4. Survival validation: The sample was divided into high / low subgroups based on the median score. The Kaplan-Meier method was used to analyze the survival differences, and the Log-rank test was used to determine statistical significance (the endpoints for the ESCC cohort were PFS and OS).

[0059] This validation process confirms that the TLS-flare gene set has consistent predictive value for efficacy across cancer-type immunotherapy cohorts, providing supporting examples for the claims.

[0060] II. Experimental Results 1. Multidimensional Analysis of Esophageal Squamous Cell Carcinoma Based on Single-Cell RNA, BCR, and Spatial Multiomics: To systematically characterize the cellular heterogeneity of esophageal squamous cell carcinoma (ESCC) and its dynamic changes during tumor progression, this study performed single-cell RNA / BCR sequencing on 407,437 cells from 92 samples from 14 patients (see attached). Figure 1 Based on artificial annotation of classical lineage markers, 11 major cell types were identified: myeloid cells, B cells, NK / T cells, plasma cells, mast cells, fibroblasts, endothelial cells, epithelial / cancer cells, platelets, and neurons. To achieve high-resolution analysis of the cell population microenvironment, we further subdivided each major cell type into subpopulations, ultimately identifying 61 significantly heterogeneous cell subpopulations (see appendix). Figure 2 a).

[0061] This study performed spatial transcriptome sequencing on 68 samples from 15 patients. To comprehensively elucidate the mechanisms of esophageal carcinogenesis, the WHO three-tier classification of esophageal squamous epithelial dysplasia (ESD) (mild, moderate, and severe / carcinoma in situ) was adopted, distinct from the clinical two-tier system that only distinguishes between low / high grade dysplasia (LGD / HGD). This classification covers seven consecutive pathological stages (see appendix). Figure 3 a): Normal epithelium (NOR), basal cell hyperplasia (BCH), mild dysplasia (mD), moderate dysplasia (MD), severe dysplasia / carcinoma in situ (SD&CIS), invasive carcinoma (ICA), and lymph node metastasis (MCA). Three main regions are defined (see appendix). Figure 3a): epithelial basal layer / tumor dominant zone, lamina propria microenvironment (ME). Due to the limited integrity of the mucosal muscularis / submucosal layer structure in frozen samples, the latter was excluded from microenvironment analysis.

[0062] By reconstructing spatial co-localization and ligand-receptor networks through ISCHIA framework, 15 cell composition categories (CCs) were identified based on spatial 10x Visium data (Fig. 1 Figure 2 d; Fig. 1 Figure 3 b-e), including: epithelial dominant (CC9 / 11 / 12), tumor dominant (CC3 / 5 / 10), lamina propria-like (CC15), tertiary lymphoid structure (TLS; CC8), lymphocyte aggregation zone (CC1 / CC6), outer muscularis-like (CC2), fibroblast-enriched (CC4 / 7 / 14), and LYVE1+ macrophage-enriched (CC13). These CCs constitute a highly conserved spatial architecture from pre-cancerous lesions (NOR to MD) to invasive cancer (Fig. 1 Figure 2 e, Fig. 1 Figure 3 b,c). Spatial distribution analysis showed that the spatial architecture was consistent with the anatomical structure at the pre-cancerous stage (basal layer: CC9; suprabasal layer: CC12 / 11), but the epithelial structure was disrupted and the microenvironment was mixed at the SD / CIS stage (Fig. 1 Figure 2 d). Key niches such as TLS, lamina propria, and glands were verified by PhenoCycler-Fusion (PCF, formerly CODEX) technology, and it was confirmed that the destruction of epithelial structure in cancer tissue was accompanied by stromal / immune infiltration (Fig. 1 Figure 2 e).

[0063] 2. Spatial zoning and functional plasticity of tertiary lymphoid structures in esophageal squamous cell carcinoma: Tertiary lymphoid structures (TLS) have been shown to play a crucial role in tumor immunity as immune hubs. Although TLS has been systematically characterized from the perspective of maturation status, new insights suggest that its microanatomical localization and context-dependent antigen / inflammatory signals may fundamentally determine its functional programming. Therefore, it is reasonable to speculate that esophageal tissue has a unique immune geographic niche of stage-specific inflammatory environment and mucosal immune landscape during carcinogenesis, which may drive the formation of functionally heterogeneous TLS subgroups.

[0064] In the analysis of esophageal tissue of carcinoma in situ (CIS) and invasive squamous cell carcinoma, we not only observed TLS in the peritumoral and intratumoral regions, but also found its presence in the tumor-adjacent tissue and the distal tissue that was histologically normal. Spatial transcriptome analysis verified this distribution feature (Fig. 1 Figure 4 b, top), and further analysis of CC8, a TLS niche based on ISCHIA algorithm deconvolution analysis (Fig. 1 Figure 4b). Notably, TLSs were more abundant in early pre-invasive stages such as basal cell hyperplasia (BCH) and severe dysplasia / CIS (SD&CIS) (Fig. 1c), suggesting the existence of inflammation-associated and tumor-associated TLSs. By PhenoCycler-Fusion detection on adjacent tissue sections matched with the samples, we verified the ubiquitous presence of TLSs in different disease stages (Fig. 1e). Figure 4 Figure 2

[0065] To systematically profile TLSs in disease progression, we developed a spatial transcriptome-based TLS identification pipeline (Fig. 2a). This method integrated H&E staining histopathological features, ISCHIA unsupervised spatial clustering, and established TLS gene signatures (TLS_9, TLS_12, TLS_50), adopted a sliding quantile window method to filter the top 10% (quantile_90) of TLS_50 signature scores as candidate regions, and further refined them by spatial localization of key components such as B cells, plasma cells, follicular dendritic cells (FDCs), follicular helper T cells (Tfh), T cells, endothelial cells, and fibroblasts. Aggregates with high composite scores and containing ≥5 spatial transcriptome spots were defined as TLSs (Fig. 2b), and finally, 80 TLSs were identified. Figure 4 Figure 4 We speculated that the TLS surrounding microenvironment might reflect the dynamic cell and pathway influx and efflux states of its function, thus the 5 nearest neighbor spots around each TLS were designated as the peripheral zone (Fig. 2c). In 15 patients with spatial transcriptome data, TLSs were detected in 14 esophageal squamous cell carcinoma (ESCC) samples, distributed in distal normal, paracancer, and cancer tissues, with 1-12 TLSs per sample (Fig. 2d). The only case without TLS detection was a precancerous lesion with sparse microenvironment components. Notably, TLSs were also observed around the tumor bed of metastatic lymph nodes (Fig. 2e).

[0066] Unsupervised clustering based on transcriptomic features divided TLSs into three functional groups (TLS-G1, G2, G3; Fig. 3a). TLS-G1 (n=7) was mainly located in the paracancer region, TLS-G2 (n=42) was widely distributed in distal, paracancer, and cancer tissues, and TLS-G3 (n=31) was enriched in cancer tissues (corresponding to pTLS and iTLS). H&E and matched PCF images showed that TLS-G1 was located near esophageal glands to form lymphatic aggregates (Fig. 3b). TLS-G2 was mainly located in the paracancer region and cancer tissues (Fig. 3c). TLS-G3 was mainly located in cancer tissues (Fig. 3d). Figure 5 Figure 7 Figure 5

[0067] Figure 6 Figure 7 Figure 7 ​​​​​​​​​TLS-G2 / G3 showed size and morphology heterogeneity, which might contain both mature and immature TLS, suggesting that TLS grouping was more relevant to anatomical location than maturity.

[0068] TLS maturation-related genes (FCER2, CXCL13, CXCL9, CR119) were differentially expressed in all three groups, and their levels were correlated with TLS size but not with grouping Figure 8 a). Functional marker analysis showed that antigen presentation gene (TAP1) and type I interferon response genes (STAT1, ISG15) were significantly upregulated in TLS-G3 and showed a radial expression pattern (Fig. Figure 8 a, Fig. Figure 9 ); acute phase reactant SAA1 was highly expressed in G1 / G3. These findings suggest that the functional status of TLS and its impact on the immune microenvironment are more dependent on spatial positioning than structural maturity.

[0069] By analyzing immune activation status, cell composition, BCR somatic hypermutation (SHM) and clonal diversity, we found that TLS-G1 highly expressed immunoglobulin A-related genes (IGHA1, JCHAIN, PIGR) and mucosal defense genes (MUC5B, DEFB1), but lower expression of interferon-stimulated genes, cytotoxic effectors and immune checkpoint molecules (Fig. Figure 9 ). Its core and peripheral zones were enriched with plasma cells, proliferative MKI67+ B cells, CCR7+ naive B cells, GPR183+ resident B cells and memory T cells, while CCR7+ naive T cells were overexpressed in the peripheral zone of G1 (Fig. Figure 10 ). Although the SHM frequency was comparable to other groups, G1 had the highest number of BCR unique sequences (Fig. Figure 11 ), reflecting its polyclonal response characteristics to mucosal microbiota.

[0070] TLS-G2 showed the highest inflammatory response score (Fig. Figure 6 ), but the expression of classical inflammatory genes was not significantly elevated (Fig. Figure 9 ). Its memory CD8 T cells (CD8_C3_TCF7), memory CD4 T cells (CD4_C8_ITGB1, CD4_C4_GPR183) and mast cells were the most enriched among the three groups (Fig. Figure 6 , Fig. Figure 9 ), suggesting a chronic inflammation-driven memory cell retention feature. Although the SHM rate was comparable to G3, the number of BCR unique sequences remained high (Fig. Figure 11 ), supporting the hypothesis of a polyclonal response triggered by tissue stress or microbiota dysbiosis.

[0071] TLS-G3 exhibited the strongest features of effector immune activation: significant upregulation of cytotoxic activity, T cell exhaustion, phagocytosis, and transplant rejection-related genes (Supplemental Figure 9 ), high expression of cytotoxic effector genes (GNLY, GZMB, PRF1), immune checkpoint molecules (PDCD1, CTLA4), ISGs, and chemokines (Supplemental Figure 7 ). ISG15+B / T cells, IFNG+CD4 T cells, neutrophils, and CXCL10+inflammatory macrophages were significantly enriched in G3 and its peripheral zone (Supplemental Figure 10 ). While SHM frequency was comparable to G2, the number of BCR unique sequences was significantly reduced (Supplemental Figure 11 ), suggesting clonal expansion of tumor-specific B cells. These features defined TLS-G3 as an anti-tumor immune niche with spatial organization.

[0072] Interestingly, tumor antigens might shape the immune niche in a spatially restricted manner. As shown in Figure 8 a, adjacent G2 / G3 TLSs had distinct microenvironmental features, with more tumor niche CC5 enriched around G3 (Supplemental Figure 8 d). Spatial transcriptome and single-cell BCR analysis showed little overlap of BCR sequences across different organizational zones (Supplemental Figure 8 e), suggesting spatially compartmentalized antigen exposure driving local B cell clonal selection.

[0073] Through Mfuzz clustering analysis, we found that TLS-G3 and its peripheral zone specifically enriched adaptive / natural immune activation pathways (clusters C4 / C9); while G1 uniquely upregulated MYD88-TLR4 pathway (cluster C5), possibly reflecting microorganism / DAMP-driven microenvironmental features (Supplemental Figure 12 ). These findings support a functional zonation model of TLS: G1 as a mucosal immune sentry responding to microbial exposure; G2 reflecting a chronic inflammation-associated “bystander” state; and G3 representing an effector hub involved in tumor antigen-specific immunity.

[0074] 3. Clinical relevance of TLS zonation: To explore the association of TLS zonation with immunotherapy efficacy, we compared the differential gene features across TLS core zones. We found that TLS-G3 features (named TLS-flare) were progressively elevated from G1 to G3, and spatially presented a “solar flare”-like radial distribution (Supplemental Figure 8 a, b). Bulk mRNA analysis of an anti-PD-1 treated ESCC cohort showed that responders had significantly higher TLS-flare scores (p=0.0016), and TLS-flare score effectively predicted treatment response (AUC=0.793) and survival benefit (PFS p=0.013; OS p=0.075) (Supplemental Figure 13a-d). Notably, TLS-G2 signature performed comparably to TLS-flare in predicting efficacy (Supp Figure 13 e-h), while the three TLS gene set originally used for spatial localization, TLS_9 / 12 / 50 signature, had no predictive value (Supp Figure 14 e,f). In melanoma cohort (PRJEB23709) validation, TLS-flare also showed excellent predictive performance (p=1.6e-07, AUC=0.807) (Supp Figure 15 ).

[0075] In summary, by dissecting TLS functional partitioning across spatial stages, we revealed the clinical value of TLS-flare score (Supp Figure 8 c): both in predicting immunotherapy efficacy and providing new insights into TLS targeting strategies for mucosal cancers. The "bystander" TLS with chronic inflammation signature also has predictive potential, while mucosal microbiota-induced immune aggregates might have higher organ conservation.

[0076] Although the critical role of tertiary lymphoid structures (TLS) in tumorigenesis has been recognized, the functional dynamics of TLS under different inflammatory and antigenic stress conditions remain unclear. As a semi-open organ system, the esophagus is highly influenced by mucosal immune responses and chronic inflammation, both of which are crucial for maintaining tissue homeostasis. This unique background provides an ideal model for studying how endogenous immune responses, chronic inflammation signals, and tumor neoantigens collectively shape TLS functions.

[0077] The present invention proposes a TLS functional framework based on spatial ecological context (rather than static maturation stages): by comparing mucosal (TLS-G1), paratumoral (TLS-G2), and intratumoral (TLS-G3) TLS, we revealed their niche partitioning in mucosal surveillance, chronic inflammation, and tumor-specific immunity. This spatial diversity suggests that TLS is not a fixed immune hub, but a structure dynamically shaped by antigen availability, local microenvironment, and chronic tissue signals. Of particular interest is TLS-G2, which is distributed in the peritumoral and distal regions, exhibiting chronic inflammation characteristics but lacking clear tumor-directed activation. This "bystander" state might reflect immune readiness but insufficient antigen stimulation, suggesting that it can be reprogrammed into a neoantigen-responsive anti-tumor TLS through treatment, thereby serving as an early biomarker for high-risk lesions and a potential target for tumor antigen vaccines, providing new ideas for chemoprevention and postoperative adjuvant therapy.

[0078] Notably, we found that the peripheral microenvironment of different TLS subtypes has unique features, and the representative features show irregular flare-like radiation patterns. These observations suggest that the spatial composition of the TLS-peripheral zone, which defines the process of cell and cytokine entry and exit, may dynamically record the overall functional status of the TLS. Importantly, the TLS-flare feature derived from the TLS-G3 adjacent microenvironment not only captures the anti-tumor immune activity, but is even superior to TLS-G3 itself in predicting the efficacy of immunotherapy. These findings highlight the extreme importance of characterizing the TLS-related niche for understanding its formation mechanism, functional status, and therapeutic significance.

[0079] Finally, it should be noted that: the above embodiments are only used to illustrate the technical solutions of the present application, and not to limit it; although the present application has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand: it can still modify the technical solutions recorded in the foregoing embodiments, or make equivalent replacement for part or all of the technical features; and these modifications or replacements do not make the essence of the corresponding technical solutions deviate from the scope of the technical solutions of the embodiments of the present application.

Claims

1. An automatic delineation and functional classification model for TLS-related niches based on the three-level lymphoid structure of spatial transcriptomics, characterized in that: The model classifies TLS into three functional subtypes, each corresponding to a different niche: The microbial responsive type, TLS-G1, is located in the mucosa and acts as a mucosal immune sentinel to respond to microbial exposure, revealing the niche partitioning of TLS in mucosal surveillance. It highly expresses immunoglobulin A-related genes IGHA1, JCHAIN, PIGR and mucosal defense genes MUC5B, DEFB1, and low expression of interferon-stimulated genes, cytotoxic effector molecules and immune checkpoint molecules. Its core and peripheral regions are enriched with plasma cells, proliferating MKI67+ B cells, CCR7+ naive B cells, GPR183+ settled B cells and memory T cells. CCR7+ naive T cells are highly expressed in its peripheral region and have the highest number of unique BCR sequences. Inflammation-associated bystander TLS-G2: It is located in the peritumoral region, reflecting the "bystander" status associated with chronic inflammation, revealing the niche zoning of TLS in chronic inflammation; it has the highest inflammatory response score but no significant increase in the expression of classical inflammatory genes; it is enriched with memory CD8 T cells (TCF7+ CD8 cells), memory CD4+ T cells (ITGB1+ CD4 T cells), GPR183+ CD4 T cells, and mast cells, exhibiting the characteristics of memory cell retention driven by chronic inflammation; the number of unique BCR sequences remains high; Antitumor-responsive TLS-G3: Located within the tumor region, it represents the effector hub involved in tumor antigen-specific immunity, revealing the niche division of TLS in tumor-specific immunity; cytotoxic activity, T cell exhaustion, phagocytosis, and transplant rejection-related genes are significantly upregulated, with high expression of cytotoxic effector genes GNLY, GZMB, and PRF1, immune checkpoint molecules PDCD1 and CTLA4, interferon-stimulating genes STAT1 and ISG15, and chemokines CXCL10, CXCL11, CCL21, and CXCL13; its core and peripheral regions are significantly enriched with ISG15+ B / T cells, IFNG+ CD4 T cells, neutrophils, and CXCL10+ inflammatory macrophages; Its BCR unique sequence number was significantly reduced, accompanied by clonal expansion of tumor-specific B cells.

2. The automatic delineation and functional classification model for TLS-related niches based on the three-level lymphoid structure of spatial transcriptomics according to claim 1, characterized in that: The TLS-G3 subtype exhibits a spatially distributed "immune flare" pattern, named TLS-flare, and constitutes a spatially organized anti-tumor immune niche. The TLS-flare characteristic gene set comprises 75 genes, consisting of the following genes: CXCL10, CHI3L1, HMOX1, GBP5, ADAMDEC1, DUSP4, SH2D1A, MMP19, MX2, FAM20A, IL2RG, PLEK, CST7, SRGN, GIMAP4, LY75, CCL19, ADAM19, AIF1, IL32, SLAMF1, IFI44L, NLRC3, SAMD9L, GBP2, CD2, IL27RA, and RASS. F4, TGFB3, IKZF1, SELPLG, GPRIN3, LCK, ISG15, SAMSN1, RAB8B, IKZF3, PAG1, SLA, SPN, CD68, HCK, CD48, SLC2A6, CD96, GMFG, PSTPIP1, ST8SIA4, SOCS1, ADAM8, HCLS1, P2 RY8, MPP1, SCIMP, TTYH2, DOK2, CPNE5, LST1, PREX1, CCL4, CD40, LMCD1, PSMB9, JAK1, SASH3, ADAMTS4, ZAP70, STK10, SP100, LTF, DOK3, TBC1D2B, PSMB8, PEAK1, TRIM22.

3. An automatic delineation and functional typing method for TLS-related niches based on the tertiary lymphoid structure of spatial transcriptomics, characterized in that: Includes the following steps: (a) Spatial transcriptome data of esophageal cancer tissues were obtained and scored based on the publicly available TLS_50 gene set, which includes the following genes: FDCSP, CR2, CXCL13, LTF, CD52, MS4A1, CCL19, LINC00926, LTB, CORO1A, CD79B, TXNIP, CD19, LIMD2, CD37, ARHGAP45, BLK, TMC8, CCL21, PTPN6, ATP2A3, IGHM, SPIB, TMSB4X, CXCR4, NCF1, CD79A, ARHGAP9, DEF6, EVL, TB. C1D10C, RASAL3, INPP5D, RNASET2, RASGRP2, TNFRSF13C, RAC2, CD22, ARHGEF1, AC103591.3, TRAF3IP3, HLADQB1, CD53, ARHGAP4, TRBC2, POU2AF1, TRAF5, OGA, FCRL3, HLA-DQA1; sites with scores exceeding the 90th percentile of all spatial coordinate sites were identified as candidate TLS locations; (b) Based on the results of step (a), the region containing 5 or more spatially consecutive candidate detection points is defined as the TLS core region; and the spatial points within this range are defined as the TLS peripheral region, starting from the edge point of the core region and extending outward by 5 points. (c) Unsupervised clustering of gene expression profiles in the TLS core region was performed, and the three functional subtypes were identified based on their molecular characteristics: TLS-G1, TLS-G2, and TLS-G3. (d) Cell composition was analyzed by RCTD deconvolution, and the functional characteristics of each subtype obtained in step (c) were verified by combining tissue morphology characteristics and differential gene expression analysis. (e) Extract differentially expressed genes between TLS-G3 and other subtypes, and screen genes that meet the following criteria to define them as the TLS-flare characteristic gene set: (e.1) The spatial expression level of the gene in TLS-G3 is significantly higher than that in other subtypes. The specific criteria are: the log2 ratio of the differential expression fold between the gene within and outside the TLS-G3 subtype is greater than 0.5, and the gene is expressed in no less than 25% of the TLS-G3 spatial sites; (e.2) There was a potential difference in expression between the immunotherapy response group and the non-response group, with a p-value of less than 0.2 according to the Wilcoxon rank-sum test; The set of genes that meets the above conditions is the TLS-flare characteristic gene set.

4. The method for constructing the model according to claim 2, characterized in that, Includes the following steps: (a) Calculate the GSVA enrichment score of the TLS-flare feature gene set based on transcriptome sequencing data from cancer patients; b) Based on the median GSVA score of all patients, patients were divided into a high-score group and a low-score group; the high-score group corresponds to patients with CR+PR who have a good response to clinical immunotherapy.

5. The application of the immunotherapy response prediction model constructed by the method of claim 4 in the preparation of the product, wherein the product is used to predict the treatment response and survival of esophageal squamous cell carcinoma (ESCC) patients to PD-1 / CTLA-4 inhibitors.

6. An immunotherapy response prediction system, characterized in that: include: A spatial positioning module is configured to perform steps (a) to (b) of claims 3 for determining the core area and peripheral area of ​​TLS; The subtype segmentation and gene set generation module is configured to perform steps (c) to (e) of claims 3, for segmenting TLS functional subtypes and extracting TLS-flare feature gene sets; The scoring and prediction module is configured to perform steps (a) to (b) of claims 4 for calculating the GSVA enrichment score and predicting the immunotherapy response.

7. The application of the automatic delineation and functional classification model of the TLS-related niches based on spatial transcriptomics as described in claim 1 in the preparation of products for the diagnosis and treatment of esophageal squamous cell carcinoma (ESCC), wherein the application includes: (a) Locating tumor immune-active regions based on the spatial distribution of TLS-G3 subtypes to guide the selection of immunotherapy target areas; (b) Based on the characteristics of the TLS-G2 subtype, develop formulations for immune reprogramming and convert them into antitumor TLS.

Images