Breast cancer cell subpopulation characteristic gene set based on multi-source single cell transcriptome data and screening method and application thereof

Abstract

The present application relates to breast cancer cell subpopulation characteristic gene set based on multi-source single cell transcriptome data and screening method and application thereof. The present application designs a breast cancer cell subpopulation characteristic gene set screening method based on multi-source data, cross-verification of at least three sets of independent source single cell transcriptome data, and takes the intersection of characteristic genes to screen robust characteristic genes, solves the sample heterogeneity of single data set and the corresponding problem of cell population across data, and improves the cross-data set repeatability of characteristic genes. It is proposed to calculate the Euclidean distance of subpopulation center point based on UMAP space coordinates, to realize objective quantitative judgment of cross-data set subpopulation correspondence by adopting hierarchical clustering and adaptive rationality judgment standard, and to avoid the subjectivity and irreproducibility of traditional naked eye observation. Three quantitative indexes are introduced, and the rationality of clustering division is judged by double judgment of separation degree and position relationship, which is not dependent on fixed absolute value and can be adapted to different data sets.

Description

Technical Field

[0001] This invention belongs to the field of biomedical technology and relates to a set of characteristic genes of breast cancer cell subpopulations based on multi-source single-cell transcriptome data, as well as its screening method and application. Background Technology

[0002] Breast cancer is one of the most heterogeneous malignant tumors, with significant differences in gene expression, proliferation capacity, metastatic potential, and treatment response among different patients and even different cancer cells within the same patient. This heterogeneity is driven by different functional cancer cell subpopulations. Therefore, identifying and characterizing these cancer cell subpopulations and searching for characteristic genes that reflect the status of these subpopulations is key to achieving precision diagnosis and treatment of breast cancer.

[0003] The development of single-cell transcriptome sequencing technology has provided a powerful tool for analyzing the cellular composition and heterogeneity of tumors at the single-cell level. In recent years, researchers have used this technology to accumulate a large amount of single-cell transcriptome data on breast cancer, providing valuable resources for a deeper understanding of the cellular heterogeneity of breast cancer.

[0004] Several breast cancer gene sets have been developed for clinical diagnostic support and molecular subtyping. For example, the PAM50 gene set classifies breast cancer into four intrinsic molecular subtypes (Luminal A, Luminal B, HER2-enriched, and Basal-like) by detecting 50 characteristic genes. Furthermore, various studies have reported molecular subtyping systems for breast cancer based on gene expression profiling data, further refining the molecular subtyping of breast cancer and providing more references for precision treatment. In addition, there are various breast cancer testing kits on the market based on gene expression profiling (such as Oncotype DX and MammaPrint) for assessing patient recurrence risk and assisting in chemotherapy decisions.

[0005] However, the aforementioned existing technologies share the following common drawbacks: First, existing gene sets are not derived from single-cell transcriptome data. Gene sets acquired based on gene expression data, PAM50, and various commercial reagent kits are all developed based on traditional bulk transcriptome data (i.e., mixed tissue samples). Bulk data reflects the average gene expression level of the entire tissue and cannot distinguish between different subpopulations of cancer cells and their respective expression characteristics. Therefore, these gene sets are difficult to truly reflect the intrinsic functional subpopulation status of cancer cells, and the signals they detect may be mixed with the influence of other cell types in the tumor microenvironment.

[0006] Second, existing gene sets are insufficient to reflect the intrinsic pathogenesis of breast cancer. Due to a lack of single-cell-level analytical capabilities, existing gene sets cannot characterize the dynamic changes of different functional cancer cell subpopulations within breast cancer. For example, crucial information such as which subpopulations are eliminated during treatment, and which subpopulations develop drug resistance or acquire metastatic ability is difficult to obtain from bulk-derived gene sets. Therefore, existing gene sets have inherent limitations in revealing the biological mechanisms of breast cancer metastasis, recurrence, and drug resistance.

[0007] Third, the clinical applications of existing technologies are relatively limited. PAM50 is mainly used for molecular subtyping and prognostic assessment, and Oncotype DX is mainly used to predict chemotherapy benefits; the functions of each gene set are relatively limited. There is currently no gene set that can be used simultaneously for multiple clinical scenarios such as treatment response assessment, metastasis status determination, relapse risk prediction, and molecular subtyping, lacking a multifunctional, integrated molecular tool.

[0008] In addition, existing single-cell transcriptome data analysis for breast cancer has the following limitations: The research mainly focuses on the subpopulation of different cell types (such as immune cells and fibroblasts) in the tumor microenvironment, while the characterization of the subpopulation features of cancer cells themselves is relatively weak. It mostly stays at the level of basic subpopulation classification and lacks systematic identification of functional characteristic genes of each cancer cell subpopulation.

[0009] The lack of summarization and integration of common characteristics across different data sources has resulted in the underutilization of data value. Breast cancer single-cell transcriptome data comes from diverse sources, involving different laboratories and patients with varying clinical backgrounds. Current technologies lack effective integration of common characteristics across different data sources, leading to unclear subgroup correspondences across datasets and hindering the full development of the potential value of these data resources.

[0010] Due to the aforementioned issues, the application value of breast cancer cell subpopulation characteristic genes obtained based on single-cell data in key clinical issues (such as treatment response assessment, metastasis status assessment, recurrence risk prediction, molecular subtyping, etc.) has not yet been systematically verified. Summary of the Invention

[0011] To address the shortcomings of existing technologies and practical needs, this invention provides a set of characteristic genes for breast cancer cell subpopulations based on multi-source single-cell transcriptome data, along with its screening method and applications. The aim is to solve the following technical problems: First, to provide a robust set of characteristic genes for breast cancer cell subpopulations across datasets and its screening method; second, to provide an objective and quantifiable method for determining subpopulation correspondences across datasets; and third, to provide a breast cancer molecular tool applicable to multiple scenarios.

[0012] To achieve the above objectives, the present invention adopts the following technical solution: In a first aspect, the present invention provides a method for screening a set of characteristic genes of breast cancer cell subpopulations based on multi-source single-cell transcriptome data, the screening method comprising the following steps: S1. Obtain at least three independent breast cancer single-cell transcriptome datasets; S2. Perform quality control on each dataset obtained in step S1, extract the epithelial cell markers, and divide cancer cells that do not express immune cell markers, endothelial cell markers, and mesenchymal cell markers into subpopulations; S3. Screen for highly expressed genes in each subgroup divided in step S2. The screening criteria are: average log2 fold change in expression > 0, and corrected p-value ≤ 0.01. S4. Perform log-normalization on each dataset and select 2000-5000 (e.g., 2000, 2500, 3000, 3500, 4000, 4500, or 5000) highly variable genes; project each dataset onto a shared low-dimensional space; identify the nearest neighbor cell pairs in this space as candidate anchor points; filter and score the anchor points based on neighborhood structural similarity. After obtaining the anchor points, the datasets are integrated into a unified space, and then UMAP dimensionality reduction is performed. S5. After dimensionality reduction using UMAP, extract the coordinates (UMAP1, UMAP2) of each cell in two-dimensional space. For each cell subpopulation k, calculate its center point. : UMAP1 coordinates UMAP1 = the arithmetic mean of the UMAP1 coordinates of all cells within this subpopulation; UMAP2 coordinates UMAP2 = the arithmetic mean of the UMAP2 coordinates of all cells within this subpopulation; For any two cell subpopulations i and j from different datasets, calculate the Euclidean distance between their centroids. : ; S6. Use a hierarchical clustering algorithm to cluster the distance matrix and generate a clustering tree. Select cut points on the clustering tree to divide the subgroups into several clusters. Subgroups within the same cluster are identified as the same large cell class. For a given clustering result, calculate the following three quantitative indicators: Average distance within cluster : The arithmetic mean of the Euclidean distances between all pairs of subgroups within the same cluster; Average inter-cluster distance : The arithmetic mean of the Euclidean distances between all pairs of subgroups in different clusters; average cluster distance : The total arithmetic mean of the Euclidean distances within and between all clusters; A cluster partition is considered a large cell group when it meets the following conditions: (1) Resolution condition: ; in, This is the separation threshold, which is used to quantify the degree to which the inter-cluster distance of any given cluster is significantly greater than the intra-cluster distance. (2) Positional relationship conditions: ; S7. For each major cell class after merging in step S6, take the union of all highly expressed genes of cell subpopulations contained in each dataset within that major cell class to obtain the set of characteristic genes corresponding to that class for each dataset. For each class of characteristic genes, if it is a characteristic gene of that class in at least three datasets, it is considered a robust characteristic gene belonging to that cell class, and finally the characteristic gene set of different cell subclasses of breast cancer is obtained.

[0013] This invention designs a method for screening characteristic gene sets of breast cancer cell subpopulations based on multi-source data. It requires cross-validation using at least three sets of independent single-cell transcriptome data. Robust characteristic genes are screened by taking the intersection of characteristic genes, effectively solving the problems of heterogeneity in single datasets and cross-dataset correspondence of cell populations, significantly improving the reproducibility of characteristic genes across datasets. A method is proposed to calculate the Euclidean distance of subpopulation centroids based on UMAP spatial coordinates, and to achieve objective quantitative determination of cross-dataset subpopulation correspondence using hierarchical clustering and adaptive rationality judgment criteria, avoiding the subjectivity and non-reproducibility of traditional visual observation. Three quantitative indicators are introduced: intra-cluster average distance, inter-cluster average distance, and global average distance, through separation condition (…). ) and positional relationship conditions ( The criteria for determining the rationality of clustering partitions are based on a dual-evaluation method. This method does not rely on fixed absolute values ​​and can adapt to different datasets.

[0014] Optionally, the dataset in step S1 may come from different laboratories or different patient cohorts.

[0015] The feature gene set of this invention is derived from single-cell transcriptome data, which can reflect the intrinsic cellular subpopulation status of breast cancer cells and capture the intrinsic mechanism differences in metastasis, treatment response, recurrence, and different molecular subtypes from the essence of the disease.

[0016] Optionally, the epithelial markers in step S2 include at least one of EPCAM, GATA3, CDH1, KRT8, KRT18, or KRT19.

[0017] Optionally, the immune cell marker includes PTPRC.

[0018] Optionally, the mesenchymal cell markers include VIM.

[0019] Optionally, the endothelial cell markers include at least one of PECAM1, VWF, or CDH5.

[0020] Optionally, the subgroup division in step S2 specifically includes: performing PCA dimensionality reduction and selecting the top 20 principal components; constructing a cell similarity graph based on the shared nearest neighbor algorithm; performing graph clustering using the Louvain or Leiden algorithm, with the resolution parameter set to 0.1-1, and selecting the optimal resolution based on the expression of cell type markers; and obtaining the subgroup division results of cancer cells in each dataset.

[0021] Optionally, the screening in step S3 includes using the FindAllMarkers function in the open-source tool Seurat package in the R working environment to calculate the highly expressed genes in each subpopulation.

[0022] Optionally, the method for projecting each dataset into a shared low-dimensional space in step S4 includes canonical correlation analysis; Optionally, the method of integrating the datasets into a unified space includes integrating them using the IntegrateData function in the Seurat package within the R working environment; Optionally, the parameters for the UMAP dimensionality reduction visualization are set as follows: reduction = "pca", dims = 1:20.

[0023] Optionally, the steps described in step S5 The value range is 1.5 to 3.0.

[0024] Secondly, a set of characteristic genes of breast cancer cell subpopulations based on multi-source single-cell transcriptome data, wherein the set of characteristic genes of breast cancer cell subpopulations is obtained by the screening method described in the first aspect.

[0025] Optionally, the gene set includes three characteristic gene sets, denoted as the first characteristic gene set (GS-A), the second characteristic gene set (GS-B), and the third characteristic gene set (GS-C), respectively. Thirdly, the present invention provides the use of the breast cancer cell subpopulation characteristic gene set and / or its detection reagent based on multi-source single-cell transcriptome data as described in the second aspect in the preparation of any of the following products: (1) Products for dynamic monitoring and evaluation of breast cancer treatment response; (2) Products for assessing and constructing metastatic invasion status of different tissue types of breast cancer; (3) Products for subpopulation analysis and risk stratification of recurrent breast cancer cells; (4) Products for analyzing the cell characteristics and molecular subtyping of different subtypes of breast cancer.

[0026] Fourthly, the present invention provides a multi-scenario assessment system for breast cancer, the assessment system comprising: (1) Input module: used to receive gene expression profile data of the sample to be tested and to standardize the data; (2) Storage module: stores the set of characteristic genes of breast cancer cell subpopulations based on multi-source single-cell transcriptome data as described in the second aspect; (3) Enrichment score calculation module: used to calculate the enrichment score of each sample in each characteristic gene set using gene set enrichment analysis method; (4) Scene determination module: used to evaluate any of the following scenes based on the enrichment score: Dynamic monitoring and assessment of breast cancer treatment response: comparing changes in enrichment scores before and after treatment, and identifying response and resistance subgroups; Assessment and construction of metastatic invasion status of different tissue types in breast cancer: NMF clustering is performed on the enrichment score matrix to determine and construct the tissue type and invasion status of the samples; Breast cancer recurrence cell subpopulation analysis and risk stratification: Based on the multi-gene set enrichment scores, cell subpopulations associated with recurrence were determined and stratified into high / medium / low recurrence risk. Analysis of cellular characteristics and molecular subtyping of different subtypes of breast cancer: Analysis of enrichment fraction distribution to determine the molecular subtype of the sample and the corresponding cellular characteristics of different subtypes; (5) Output module: used to output assessment results and clinical decision recommendations.

[0027] This invention systematically validates the application value of characteristic gene sets in four independent clinical scenarios: treatment response assessment, metastasis and invasion assessment, recurrent cell characteristics and risk stratification, molecular subtyping, and cell characteristics. This invention comprises multiple characteristic gene sets, and by analyzing the enrichment patterns of each gene set in different clinical application scenarios, the status and dynamic changes of different functional cell subpopulations can be simultaneously assessed. This multi-gene set joint analysis strategy overcomes the limitations of traditional single biomarkers or single gene sets in revealing tumor heterogeneity, and can more comprehensively and fundamentally reflect the molecular characteristics and biological behavior of breast cancer at the cell subpopulation level.

[0028] Optionally, the dynamic monitoring and evaluation of breast cancer treatment response specifically includes: 1) Plot the dynamic response curves of the enrichment fraction of each gene set as a function of treatment time; 2) Calculate the effect size Cohen's d, Cohen's d = (M1 – M2) / SD_pooled, where: M1 is the arithmetic mean of the enrichment scores of the characteristic gene set in the treated population; M2 is the arithmetic mean of the enrichment scores of the characteristic gene set in the untreated population; SD_pooled is the combined standard deviation of the two groups; 3) Make judgments based on the curve change pattern and effect size: Targeted response subpopulation identification: If the enrichment score of a certain characteristic gene set decreases or continues to decrease after treatment, it is determined that the treatment has an inhibitory effect on the cell subpopulation marked by the gene set, and the subpopulation is the targeted response subpopulation of the treatment. Acquired drug resistance warning: If the enrichment score of a certain characteristic gene set increases after treatment, or rebounds after an initial decrease, the cell subpopulation marked by that gene set is determined to be a potential acquired drug resistance subpopulation, indicating that the treatment plan needs to be adjusted. Primary drug resistance identification: If the enrichment score of a certain characteristic gene set does not change significantly before and after treatment, the treatment is deemed ineffective for the cell subpopulation marked by that gene set. This cell subpopulation is a potential primary drug-resistant subpopulation, suggesting that the treatment regimen needs to be adjusted.

[0029] Optionally, the assessment and construction of the metastatic invasion status of different tissue types of breast cancer includes: 1) Using the non-negative matrix factorization algorithm, by setting different rank ranges, the non-negative matrix factorization algorithm is repeatedly run for each rank value to calculate the cophenetic correlation coefficient corresponding to each rank value. The rank value corresponding to the first time the cophenetic correlation coefficient reaches the highest value is selected as the optimal number of clusters. If there are multiple local maxima, the rank value with the smallest value is selected. 2) Based on the cluster categories obtained above, perform non-negative matrix factorization clustering using the Brunet algorithm to obtain the sample partitioning; 3) The basismap function is used to calculate and visualize the contribution weight of each feature gene set to the classification of each category, so as to obtain the degree of cancer cell invasion or metastasis in each sample; Optionally, the breast cancer recurrence cell subset analysis and risk stratification includes: 1) For each feature gene enrichment score, the surv_cutpoint function is used to determine the optimal cutoff value for each feature gene set. Patients are divided into high enrichment group and low enrichment group. Then, the Kaplan-Meier method is used to plot the relapse-free survival curves of the high enrichment group and the low enrichment group, and the Log-rank test is performed. With p < 0.05 as the threshold, it is determined whether there is a statistically significant association between the enrichment level of the feature gene set and the risk of relapse, so as to obtain the association between the cell subpopulation marked by the feature gene set and the relapse. 2) For the characteristic gene set that is significantly associated with the risk of relapse, the risk direction is determined according to the survival curve: if the relapse-free survival of the high enrichment group is significantly shorter than that of the low enrichment group, then the high enrichment is determined to be associated with a high risk of relapse; if the relapse-free survival of the low enrichment group is significantly shorter than that of the high enrichment group, then the low enrichment is determined to be associated with a high risk of relapse. 3) For the characteristic gene sets significantly associated with relapse risk, obtain the populations with high and low enrichment for each characteristic gene set: The intersection of each feature set and the population associated with high recurrence risk is used to determine the high recurrence risk population. The intersection of each feature set and the population associated with low recurrence risk is used to determine the population with low recurrence risk. The remaining population was assessed as having a moderate risk of recurrence; 4) Output relapse-related cell subsets and protective cell subsets, and risk stratification results for each patient; Optionally, the analysis of cellular characteristics and molecular subtyping of different breast cancer subtypes includes: The samples were divided according to histological subtypes, and the enrichment fraction distribution of each subtype in each characteristic gene set was calculated to analyze the cell subpopulation characteristics of each molecular subtype. For mixed samples containing multiple histological subtypes or multiple disease stages, nonnegative matrix factorization can be used for automated category identification. For each cluster category, the enrichment score distribution of each feature gene set is calculated, and the cell subpopulation feature enrichment spectrum of that category is output. For longitudinal samples containing different disease stages of the same patient, nonnegative matrix factorization clustering can be used for analysis: Treatment-induced subtype conversion: whether samples before and after treatment belong to different cluster categories; Subtype evolution during metastasis: Do primary lesions and metastatic lesions exhibit different subtype characteristics?

[0030] Fifthly, the present invention provides a computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, performs the functions of the breast cancer multi-scenario assessment system as described in the fourth aspect.

[0031] In a sixth aspect, the present invention provides a computer device comprising the computer-readable storage medium described in the fifth aspect.

[0032] Compared with the prior art, the present invention has at least the following beneficial effects: (1) Cross-dataset robustness The obtained feature gene set has good cross-dataset reproducibility through cross-validation and intersection screening of at least three independent single-cell datasets. (2) Objective quantitative judgment An objective quantification of subgroup correspondence was achieved by adopting a subgroup correspondence determination method based on Euclidean distance clustering and an adaptive rationality determination criterion based on cluster analysis. (3) Revealing the mechanism of subtype heterogeneity The feature gene set screened in this invention comes from different cell subpopulations, which can analyze the activity differences of cell subpopulations in breast cancer under different clinical backgrounds, and provide richer decision-making information for precision treatment; (4) Dynamic classification capability The characteristic gene set of this invention can analyze cell subtype transformation and evolution during treatment and metastasis, overcoming the limitations of traditional static typing; (5) High efficacy in assessing treatment response In paired pre- and post-treatment samples, the characteristic gene set can effectively distinguish between treatment-responsive subgroups and ineffective / potentially resistant subgroups. The characteristic gene set predicts treatment effect AUC of 94.35%, and the effect size reaches the internationally recognized intermediate or higher standard (|d|>1). It can be used to accurately assess the treatment response and anticipate drug resistance in breast cancer patients. (6) Accurate identification of metastatic lesions The characteristic gene set can effectively distinguish metastatic lesions from adjacent tissues, with an AUC of 0.86, which is excellent. It can be used to assist in judging the metastatic status of breast cancer and constructing a metastatic invasion map. (7) Strong recurrence prediction ability Survival analysis showed that the enrichment level of characteristic gene sets was significantly associated with recurrence-free survival in breast cancer patients after treatment (p<0.05). Patients could be stratified into high, intermediate, and low risk levels through the combination of multiple gene sets. Furthermore, it could identify cell subpopulations associated with recurrence, a function that existing recurrence prediction tools cannot achieve. (8) Superior molecular typing efficiency In a dataset containing 310 cases of different subtypes of breast cancer, the AUC of the characteristic gene set distinguishing triple-negative breast cancer from non-triple-negative breast cancer reached 0.82, which was significantly better than the prediction methods based on this data reported in the literature (AUC=0.58). Attached Figure Description

[0033] Figure 1This is a schematic diagram of the roadmap for the screening method of breast cancer cell subpopulation characteristic gene set of the present invention.

[0034] Figure 2 This is a schematic diagram illustrating the application route of the breast cancer cell subpopulation characteristic gene set of this invention.

[0035] Figure 3A This is a UMAP diagram showing the integration of single-cell transcriptome data from three different laboratories in Example 1 of this invention.

[0036] Figure 3B This is a clustering result diagram from Embodiment 1 of the present invention.

[0037] Figure 3C This is a distribution map of cell subpopulations and categories corresponding to the data in Embodiment 1 of the present invention.

[0038] Figure 4 Figure A shows the enrichment distribution of each feature gene set before and after breast cancer treatment and the prediction results of the treatment response in Embodiment 2 of the present invention. Figure B shows the distribution of each feature gene set before and after treatment and the ROC curves of the three feature gene sets used to predict the treatment response.

[0039] Figure 5A This is a graph showing the results of calculating the cophenetic correlation coefficients for each rank value in Example 3.

[0040] Figure 5B This is a diagram showing the sample partitioning results in Example 3.

[0041] Figure 5C This is a graph showing the contribution weight of each feature gene set to the classification of each category in Example 3.

[0042] Figure 5D This is the ROC curve diagram of the characteristic gene set used to distinguish between the two types of tissues in Example 3.

[0043] Figure 6 This is a graph showing the association between the enrichment of each characteristic gene set and relapse in Example 4, as well as the results of the classification of different relapse risk groups.

[0044] Figure 7 Figure A shows the enrichment distribution of each characteristic gene set in different subtypes of breast cancer and the prediction results of subtype types in Example 5. Figure A shows the distribution of each characteristic gene set in the non-triple-negative and triple-negative subtypes, and Figure B shows the ROC curves of the three sets of characteristic gene sets and literature indicators used to distinguish the two subtypes. Detailed Implementation

[0045] The technical solution of the present invention will be further described below with reference to the accompanying drawings and specific embodiments. However, the following embodiments are merely simple examples of the present invention and do not represent or limit the scope of protection of the present invention. The scope of protection of the present invention is determined by the claims.

[0046] Where specific techniques or conditions are not specified in the examples, they shall be performed in accordance with the techniques or conditions described in the literature in this field, or in accordance with the product instructions. Reagents or instruments whose manufacturers are not specified are all conventional products that can be purchased from legitimate channels.

[0047] This invention provides a characteristic gene set for breast cancer cell subpopulations based on multi-source single-cell transcriptome data and a method for screening it. Addressing the problem that existing gene sets, derived from bulk transcriptome data, cannot reflect the intrinsic functional subpopulation status of cancer cells, this invention utilizes at least three independently sourced breast cancer single-cell transcriptome datasets. Through cross-dataset integration and intersection screening, a characteristic gene set that robustly reflects the functional subpopulations of different breast cancer cells is obtained.

[0048] Specifically, the technical roadmap for screening gene sets characteristic of breast cancer cell subsets is as follows: Figure 1 As shown, it includes the following steps: Step 1: Data Acquisition Obtain at least three independent sets of breast cancer single-cell transcriptome data, which may come from different laboratories and different patient cohorts.

[0049] Step 2: Independent subpopulation of each dataset Quality control was performed on each dataset obtained in step 1. Based on the expression of classic epithelial markers (including but not limited to EPCAM, GATA3, CDH1, KRT8, KRT18, and KRT19), and excluding the expression of immune cell markers (including but not limited to PTPRC), stromal cell markers (including but not limited to VIM), and endothelial cell markers (including but not limited to PECAM1, VWF, and CDH5), cancer cell populations were identified and extracted.

[0050] The extracted cancer cells were subgrouped, specifically including: performing PCA dimensionality reduction and selecting the top 20 principal components; constructing a cell similarity graph (k=20) based on the shared nearest neighbor (SNN) algorithm; performing graph clustering using the Louvain or Leiden algorithm with resolution parameters set to 0.1-1, selecting the optimal resolution based on the expression of cell type markers; and obtaining the subgrouping results of cancer cells in each dataset.

[0051] Step 3: Obtaining highly expressed genes from each dataset subpopulation For cancer cells in each dataset, the FindAllMarkers function was used to calculate the highly expressed genes in each subgroup. The screening criteria were: average log2 expression fold change > 0, and corrected p-value ≤ 0.01.

[0052] Step 4: Dataset Integration and Dimensionality Reduction The FindIntegrationAnchors function was used to identify anchor cell pairs across datasets: each dataset was log-normalized, and 2000-5000 hypervariable genes were selected; each dataset was projected to a shared low-dimensional space using canonical correlation analysis (CCA); the nearest neighbor cell pairs in this space were identified as candidate anchors; and the anchors were filtered and scored based on neighborhood structural similarity.

[0053] After obtaining the anchor point, the IntegrateData function is used to integrate multiple sets of data into a unified space, and then UMAP dimensionality reduction visualization is performed with the parameters set as follows: reduction = "pca", dims = 1:20.

[0054] Step 5: Calculation of Euclidean distance between subgroups After dimensionality reduction using UMAP, the coordinates (UMAP1, UMAP2) of each cell in two-dimensional space are extracted. For each cell subpopulation k, its center point is calculated. : The UMAP1 coordinates are the arithmetic mean of the UMAP1 coordinates of all cells within the subpopulation. The UMAP2 coordinates are the arithmetic mean of the UMAP2 coordinates of all cells in this subpopulation.

[0055] For any two cell subpopulations i and j from different datasets, calculate the Euclidean distance between their centroids: Step 6: Determine the correspondence between subgroups across datasets based on Euclidean distance clustering. Hierarchical clustering algorithms (e.g., using the pheatmap package) are employed to cluster the distance matrix, generating a clustering tree. Cut points are selected on the clustering tree to divide the subgroups into several clusters, with subgroups within the same cluster being identified as belonging to the same major cell class.

[0056] For a given clustering result, calculate the following three quantitative indicators: Intra-cluster average distance : The arithmetic mean of the Euclidean distances between all pairs of subgroups within the same cluster; Average inter-cluster distance : The arithmetic mean of the Euclidean distances between all pairs of subgroups in different clusters; average cluster distance : The total arithmetic mean of the Euclidean distances within and between all clusters.

[0057] A clustering partition is considered a reasonable major category partition when it simultaneously meets the following conditions: (1) Resolution condition: in, The separation threshold is used to quantify the degree to which the inter-cluster distance of any cluster is significantly greater than the intra-cluster distance. Those skilled in the art can adjust it according to the actual dataset, and the preferred range is 1.5 to 3.0.

[0058] (2) Positional relationship conditions: This condition ensures the overall structural advantage of the clustering results, namely, the distance between subgroups within the same major class is less than the global average, while the distance between subgroups between different major classes is greater than the global average.

[0059] Step 7: Merging Major Category Characteristic Genes For each major cell category after merging in step 6: The union of the highly expressed genes of all cell subpopulations contained in each dataset within the major category is taken to obtain the set of characteristic genes corresponding to the category for each dataset. For each class of characteristic genes, if a gene is a characteristic gene of that class in at least three sets of data, it is considered a robust characteristic gene of that cell class, thus obtaining the characteristic gene set of different cell subclasses of breast cancer.

[0060] Step 8: Obtain the feature gene set Using the above method, characteristic gene sets of different cell subtypes of breast cancer are finally obtained. In Example 1 of the present invention, three characteristic gene sets are obtained, which are denoted as the first characteristic gene set (GS-A), the second characteristic gene set (GS-B), and the third characteristic gene set (GS-C), respectively.

[0061] The first characteristic gene set (GS-A) includes: B2M, TMSB4X, BTG1, FOS, JUNB, SAT1, S100A6, JUN, NFKBIA, KLF6, HSP90AA1, DUSP1, ZFP36L2, ZFP36L1, IFITM2, MT2A, FTH1, CLU, MT-C O1, MT-CO3, XBP1, IFITM3, CEBPD, S100A11, ANXA2, IER3, HSPB1, IFI27, CCND1, IFI6, SOX4, TNFSF10, CD52, SRGN, S100A4, AREG and HES1.

[0062] The second characteristic gene set (GS-B) includes: TRPS1, NEAT1, AFF3, XIST, SOX5, ZSWIM6, INPP4B, SMYD3, IGF1R, SIK3, CDK14, DPYD, CRIM1, DOCK10, NFKB1, LINC01619, FCHSD2, LPP, ANKRD28, BCL2, NAMPT, TBC1D9, FTX, PDE4D, ZBTB20, WWOX, NALADL2, ETV6, GLCCI1, RAPGEF2, DENND4A, LARGE1, PA RD3B, CCSER1, IMP2L, NFIA, AUTS2, PTPRK, FHIT, CERS6, FNDC3B, PCAT1, TANC2, ALCAM, RUNX1, MBNL1, DLEU2, MACROD2, SAMD12, GNAQ, RAD51B, RORA, ANKRD36C, MAGI3, RAB11FIP1, SMCHD1, RALGPS2, TNIK, RNF19A, MYO1E, ATXN1, LDLRAD4, MALAT1, MT- CYB, MT-CO2, MT-CO1, MT-CO3, UTRN, JMJD1C, PATJ, PTK2, FMN1, CENPP, CD44, TGFBR3, PVT1, LRBA, MACF1, ZFAND3, ITPR2, BTBD9, RABGAP1L, PDE7A, ARID1B, JARID2, FBXL17, GPHN, RBM47, LINC00342, SIPA1L3, LINC01473, SSBP2, BACH2, SIPA1L1, RASG RP3、RAPGEF1、ZHX2、MGAT5、PIP4K2A、TMEM131L、VAV3、FOXP1、SAMD4A、PITPNC1、NR3C1、AGAP1、TRIO、EZR、EXT1、CBLB、ABTB2 、ARHGAP26、FOXO1、CAMK1D、PEAK1、PALLD、CD55、HIVEP2、RAPGEF5、PRKCH、CNOT6L、MAML3、MAP3K8、RFX3、CEMIP2、AC016831.7. SETBP1, PRKCE, SHROOM3, GMDS, ANKS1B, PKP4, DACH1, PBX1, APBB2, ZDHHC14, ARL15, SNED1, MLLT3, THRB, SYT L2, DOCK9, NEDD4L, ANK3, MYO6, MT-ATP8, LRMDA, KAZN, MAGI1, MAST4, TIAM1, STARD13, SNX9, PDE4B, LYST, NCKAP 5. PELI1, PPP3CA, NCAM2, DST, TANC1, NFKBIZ, CADPS2, NRIP1, RALGAPA2, KCNMA1, NEK10, MACC1, THSD4, FAM160 A1, PTPRJ, CMSS1, FAAH2, RASAL2, ESR1, MAML2, TSC22D2, BCAS3, DLGAP1, BRAF, HERC1, RALGAPA1, GAB1 and CAMK2D. .

[0063] The third trait gene set (GS-C) includes: TYMS, TOP2A, MKI67, TUBA1B, STMN1, HMGB2, TUBB, TPX2, UBE2C, CENPF, RRM2, PTTG1, CCNB1, CDC20, S100A4, GAPDH, ACTB, TMSB4X, LGALS1 and IFI6.

[0064] This invention further develops the application of the above-mentioned characteristic gene set in breast cancer diagnosis and prognostic assessment (or the preparation of corresponding products), and the application roadmap is as follows. Figure 2 As shown.

[0065] Application 1: Dynamic monitoring and evaluation of breast cancer treatment response This application can be used for pre-treatment and post-treatment, multiple time points during treatment, or different clinical scenarios for control and treatment groups, including the following steps: (1) Obtain gene expression profile data of biological samples from patients before and after treatment (or at multiple time points during treatment, or between the control group and the treatment group), and list behavioral genes as samples, and standardize the gene expression profile data. (2) Using the characteristic gene sets (GS-A, GS-B, GS-C) described in this invention, the enrichment scores of each gene set at different time points of treatment are calculated by the gene set enrichment analysis (GSEA) method. (3) Plot the dynamic response curves of the enrichment fraction of each gene set as a function of treatment time; (4) Calculate the effect size (Cohen's d): Cohen's d = (M1 – M2) / SD_pooled Where: M1 is the arithmetic mean of the enrichment scores of the characteristic gene set in the treated population; M2 is the arithmetic mean of the enrichment scores of the characteristic gene set in the untreated population; SD_pooled is the pooled standard deviation of the two groups.

[0066] (5) Make a judgment based on the curve change pattern and the effect size.

[0067] The effect size threshold T is set with reference to the internationally accepted Cohen's d effect size standard, Targeted responsive subpopulation identification: If the enrichment score of a characteristic gene set decreases or continuously decreases after treatment (|d|≥T), it is determined that the treatment has an inhibitory effect on the cell subpopulation marked by this gene set, and this subpopulation is the targeted responsive subpopulation of the treatment; Acquired drug resistance warning: If the enrichment score of a characteristic gene set increases after treatment, or rebounds and increases after an initial decrease (|d| ≥T), the cell subpopulation marked by this gene set is determined to be a potential acquired drug resistance subpopulation, suggesting that the treatment regimen needs to be adjusted; Primary drug resistance identification: If the enrichment score of a characteristic gene set does not change significantly before and after treatment (|d| <T), it is determined that the treatment is ineffective on the cell subpopulation marked by this gene set, and this cell subpopulation is a potential primary drug resistance subpopulation, suggesting that the treatment regimen needs to be adjusted.

[0068] (6) In Example 2 of the present invention, the above method is applied to the treatment data of breast cancer patients, and the results show that: The Cohen's d values of the two sets of characteristic gene sets (GS-A and GS-B) for treatment response both exceed 1. According to the internationally accepted standard (d≥0.8 is a large effect size), it indicates that the characteristic gene sets of the present invention can significantly distinguish the treatment-responsive population from the non-responsive population; Among them, the AUC of the GS-A characteristic gene set for predicting treatment response reaches 94.35% (95% CI: 89.88%-98.83%), which further verifies the high accuracy and clinical application value of the characteristic gene set of the present invention in treatment response assessment.

[0069] In addition, the GS-C enrichment score decreases after treatment, while the GS-A and GS-B enrichment scores increase, suggesting that the cell subpopulation marked by GS-C can respond to the treatment (inhibited after treatment), and the cell subpopulations marked by GS-A and GS-B are potential ineffective or drug-resistant subpopulations (enrichment increases after treatment). This result shows that through the combined analysis of multiple gene sets, the present invention can simultaneously identify effective treatment targets and drug resistance-related subpopulations, providing a more comprehensive molecular basis for precise treatment decision-making.

[0070] Application 2: Assessment and construction of metastatic invasion status in different tissue types of breast cancer Clinically, for patients with multiple primary lesions or suspected metastases, accurately determining the presence and extent of cancer cell metastasis in different tissue samples, as well as the molecular relationships between different tissues, is crucial for developing precise treatment plans. Traditional pathological examinations rely on morphological judgment, which is highly subjective, difficult to quantify, and unable to reveal the molecular relationships between cancer cell subpopulations in different tissues.

[0071] This invention provides an objective quantitative assessment tool based on the enrichment score of a feature gene set, which can effectively overcome the aforementioned limitations. This application targets breast cancer that has metastasized or has unknown metastases (initial diagnosis or suspected metastasis), providing an objective tool for assessing cancer cell invasion and constructing invasion status atlases for different tissue types, including different cancer tissues, metastatic lesions, adjacent tissues, and distant tissues. The tool includes the following steps: (1) Obtain gene expression profile data from different tissue samples, including behavioral genes, and list them as samples. Standardize the gene expression profile data. The samples may include, but are not limited to: primary breast cancer lesions, lymph node metastases, distant metastases (such as liver metastases, lung metastases, bone metastases), adjacent tissues of the tumor (<2cm from the tumor edge), and distal normal tissues (>2cm from the tumor edge or from the contralateral breast). (2) Using the feature gene sets (GS-A, GS-B, GS-C) described in this invention, the enrichment score of each sample in each feature gene set is calculated by gene set enrichment analysis method to obtain the enrichment score matrix of sample × feature gene set; (3) Using the NMF algorithm, by setting different rank ranges (preferably 2-15), the NMF algorithm is repeatedly run for each rank value (preferably 50-100 times) to calculate the cophenetic correlation coefficient corresponding to each rank value. The rank value corresponding to the first time the cophenetic correlation coefficient reaches its highest value is selected as the optimal number of clusters. If there are multiple local maxima, the rank value with the smallest value is selected; (4) Based on the clustering categories obtained above, perform NMF clustering using the Brunet clustering algorithm to obtain the sample division.

[0072] (5) The basismap function is used to calculate and visualize the contribution weight of each feature gene set to the classification of each category, so as to obtain the degree of cancer cell invasion or metastasis in each sample.

[0073] (6) In Example 3 of the present invention, the above method was applied to data on breast cancer metastasis and adjacent tissues, and the results showed that: The three characteristic gene sets (GS-A, GS-B, GS-C) can effectively identify and classify metastatic tissues and adjacent tissues, with an AUC of 87.5% (95% CI: 77.44%-97.56%), indicating that the characteristic gene sets of the present invention have excellent diagnostic efficacy in distinguishing metastatic lesions from adjacent tissues.

[0074] The contribution weights of each characteristic gene set to the classification indicate that GS-A and GS-B contribute significantly to the division of adjacent tissues, suggesting that the cell subpopulations they mark may be involved in or drive the invasion of cancer cells into adjacent or distant tissues.

[0075] This embodiment demonstrates that the present invention can not only identify metastatic tissues with high accuracy, but also analyze the functional roles of different cell subpopulations in the metastasis process, providing a molecular basis for understanding the mechanism of breast cancer metastasis and developing intervention strategies.

[0076] Application 3: Subpopulation Analysis and Risk Stratification of Recurrent Breast Cancer Cells Breast cancer patients still face varying degrees of recurrence risk after treatment. Accurately identifying high-risk patients and providing intensive treatment, while avoiding overtreatment of low-risk patients, is a core objective of precision medicine. This application targets breast cancer patients undergoing specific treatment regimens (such as chemotherapy, endocrine therapy, and targeted therapy). By detecting the enrichment of the characteristic gene set described in this invention in pre-treatment biological samples, it identifies high-risk and protective cell subpopulations associated with recurrence, classifying patients into high, intermediate, and low recurrence risk levels, providing a basis for clinical decision support and personalized treatment. The process includes the following steps: (1) For a specific treatment plan, obtain the gene expression profile data of the patient's biological samples before treatment, and list the behavioral genes as samples, and perform standardized processing on the gene expression profile data; (2) Using the feature gene sets (GS-A, GS-B, GS-C) described in this invention, the enrichment score of each sample in each feature gene set is calculated by gene set enrichment analysis method; (3) For each characteristic gene enrichment score, the surv_cutpoint function (from the survminer package) is used to automatically determine the optimal cutoff value for each characteristic gene set, and patients are divided into high enrichment group and low enrichment group. As an alternative, the population can also be divided according to the median or mean of the enrichment score. Then, the Kaplan-Meier method is used to plot the relapse-free survival curves of the high enrichment group and the low enrichment group, and the Log-rank test is performed. With p < 0.05 as the threshold, it is determined whether there is a statistically significant association between the enrichment level of the characteristic gene set and the risk of relapse, so as to obtain the association between the cell subpopulation marked by the characteristic gene set and the relapse. (4) For the characteristic gene set that is significantly associated with the risk of relapse, the risk direction is determined according to the survival curve: if the relapse-free survival of the high enrichment group is significantly shorter than that of the low enrichment group, then the high enrichment is associated with the high risk of relapse; if the relapse-free survival of the low enrichment group is significantly shorter than that of the high enrichment group, then the low enrichment is associated with the high risk of relapse.

[0077] (5) For the characteristic gene sets that are significantly associated with the risk of relapse, obtain the populations with high and low enrichment of each characteristic gene set: The intersection of each feature set and the population associated with high recurrence risk is used to determine the high recurrence risk population. The intersection of each feature set and the population associated with low recurrence risk is used to determine the population with low recurrence risk. The remaining population was assessed as having a moderate risk of recurrence.

[0078] (6) Output relapse-related cell subsets and protective cell subsets, and risk stratification results (high risk / intermediate risk / low risk) for each patient.

[0079] (7) In Example 4 of the present invention, the results of applying the above method to breast cancer treatment data show that: In this treatment regimen, high enrichment of GS-C was significantly associated with a high risk of relapse, classifying them as a high-risk relapse subgroup; while low enrichment of GS-A and GS-B was significantly associated with a high risk of relapse, classifying them as protective subgroups.

[0080] Further, patients were successfully stratified into high, intermediate, and low recurrence risk groups using the above method, with significant differences (P < 0.05). The low-risk group showed a recurrence risk reduction of approximately 90% (HR: 0.1057, 95% CI: (0.03231, 0.3461)), while the intermediate-risk group showed a recurrence risk reduction of approximately 70% (HR: 0.2837, 95% CI: (0.17482, 0.4603)). These embodiments demonstrate that the present invention can effectively identify high-risk and protective cell subpopulations for recurrence, as well as patient groups with different recurrence risk levels, and has clear clinical guiding value.

[0081] Application 4: Cellular Characterization and Molecular Subtyping of Different Subtypes of Breast Cancer Breast cancer is typically classified into several major molecular subtypes, including Luminal A, Luminal B, HER2-enriched, and triple-negative (Basal-like). Clinically, classification often relies on the positivity or negativity of molecules expressed on the surface of cancer cells, such as estrogen receptor (ER), progesterone receptor (PR), and human epidermal growth factor receptor 2 (HER2). However, current classification methods are primarily based on the detection results of the primary tumor. In reality, different types of cancer undergo dynamic changes in their surface receptor molecules during metastasis or after treatment (e.g., ER-positive to ER-negative, or changes in HER2 status), rendering the original classification ineffective. Existing molecular classifications reflect the average expression level of the entire tissue and cannot reveal the composition and dynamic changes of different functional cell subpopulations within the tumor. Furthermore, current classifications lack a unified tool that can be applied simultaneously to primary lesions, metastases, and different treatment scenarios.

[0082] This application can be used for breast cancer samples at different disease stages (including but not limited to primary lesions, metastatic lesions, and different treatment scenarios). By analyzing the enrichment fraction distribution of the characteristic gene set (GS-A, GS-B, GS-C) described in this invention, it can construct cell characteristic maps of different subtypes, thereby achieving objective determination and dynamic analysis of molecular subtypes.

[0083] Includes the following steps: (1) Obtain gene expression profile data of the patient's biological samples before treatment, and list behavioral genes as samples, and standardize the gene expression profile data; (2) Using the feature gene sets (GS-A, GS-B, GS-C) described in this invention, the enrichment score of each sample in each feature gene set is calculated by gene set enrichment analysis method; (3) Samples can be directly divided according to histological subtypes, and the enrichment fraction distribution (mean, median, standard deviation, etc.) of each subtype in each characteristic gene set can be calculated to analyze the cell subpopulation characteristics of each molecular subtype.

[0084] (4) For mixed samples containing multiple histological subtypes or multiple disease stages, the nonnegative matrix factorization (NMF) method in Application 2 can be used for automated category identification. For each cluster category, calculate its enrichment score distribution (mean, median, etc.) in each feature gene set, and output the cell subpopulation feature enrichment spectrum of that category.

[0085] (5) For longitudinal samples containing different disease stages of the same patient (e.g., before treatment, after treatment, and after relapse), NMF clustering can be used to resolve the following: Treatment-induced subtype conversion: whether samples before and after treatment belong to different cluster categories (e.g., from Luminal type to triple-negative sample). Subtype evolution during metastasis: Do primary lesions and metastatic lesions exhibit different subtype characteristics?

[0086] (6) In Example 5 of the present invention, the above method was applied to data on different subtypes of breast cancer. The results showed that the enrichment patterns of GS-B and GS-C in triple-negative and non-triple-negative breast cancer were completely opposite, and the AUC for distinguishing between the two subtypes both exceeded 80%, with GS-B at 82.53% (95% CI: 78.01%-87.06%) and GS-C at 81.9% (95% CI: 77.27%-86.54%). This data has been published, and the original report proposed a set of indicators that could predict the clinical outcome of breast cancer, but its AUC for distinguishing between triple-negative and non-triple-negative subtypes was only 58.71% (95% CI: 53.46%-63.97%). In comparison, the characteristic gene set of the present invention has significantly better distinguishing power, indicating that it better reflects the essential characteristics of breast cancer.

[0087] This invention also provides a multi-scenario assessment system for breast cancer, characterized in that it includes: (1) Input module: used to receive gene expression profile data of the sample to be tested and to standardize the data; (2) Storage module: stores the feature gene set (GS-A, GS-B, GS-C) described in this invention; (3) Enrichment score calculation module: used to calculate the enrichment score of each sample in each characteristic gene set using gene set enrichment analysis method; (4) Scene determination module: used to evaluate any of the following scenes based on the enrichment score: Dynamic monitoring and assessment of breast cancer treatment response: comparing changes in enrichment scores before and after treatment, and identifying response and resistance subgroups; Assessment and construction of metastatic invasion status of different tissue types in breast cancer: NMF clustering is performed on the enrichment score matrix to determine and construct the tissue type and invasion status of the samples; Breast cancer recurrence cell subpopulation analysis and risk stratification: Based on the multi-gene set enrichment scores, cell subpopulations associated with recurrence were determined and stratified into high / medium / low recurrence risk. Analysis of cellular characteristics and molecular subtyping of different breast cancer subtypes: Analysis of enrichment fraction distribution to determine the molecular subtype of the sample and the corresponding cellular characteristics of different subtypes; (5) Output module: used to output assessment results and clinical decision recommendations.

[0088] The system can be integrated into computer equipment, servers, or cloud platforms and implemented through software programs.

[0089] Example 1 This embodiment obtained three sets of single-cell transcriptome data from breast cancer from different laboratories in the GEO and GSA databases, with the search numbers being: GSE180286 (data 1, abbreviated as D1), GSE190870 (data 2, abbreviated as D2), and HRA002051 (data 3, abbreviated as D3). Each set of data was analyzed, including quality control using the open-source tool Seurat in the R environment, extraction of epithelial cell markers, and subpopulation of cells expressing markers for immune cells, endothelial cells, and mesenchymal cells. PCA dimensionality reduction was performed, selecting the top 20 principal components; a cell similarity graph (k=20) was constructed based on the shared nearest neighbor (SNN) algorithm; and graph clustering was performed using either the Louvain or Leiden algorithm. The resolution parameters for the three sets of data were set to 0.5, 0.5, and 0.1, respectively. The FindAllMarkers function in the Seurat package was used, and the results were corrected for a p-value ≤ 0.01 to obtain the highly expressed genes in each subgroup.

[0090] Multiple datasets were integrated using the FindIntegrationAnchors and IntegrateData functions from the Seurat package. The integrated dataset was then visualized using UMAP dimensionality reduction, with parameters set to: reduction = "pca", dims = 1:20. The resulting UMAP plot after integrating the three datasets is shown below. Figure 3A As shown, the results indicate that there is overlap in specific cell subpopulations (starting with 'c', marked with different colors) among the three sets of data (D1, D2, D3).

[0091] Then, based on the scheme, the Euclidean distance between subgroups was calculated, and clustering was performed using the default parameters of pheatmap. The clustering results are as follows: Figure 3B As shown, the results indicate that different data and different cell subpopulations are clustered into the same class (the numbers in the figure represent Euclidean distances).

[0092] For cut point 2, the subgroup is divided into two clusters with the following intra-cluster average distances: Cluster 1, 3.876484; Cluster 2, 1.31148. The inter-cluster average distance is 8.21133. The inter-cluster average distance is 4.549175.

[0093] The separation degree between cluster 1 and cluster 2 is 2.118241, which satisfies the following conditions: average distance within clusters < average distance between clusters < average distance between clusters; The separation degree of cluster 2 to cluster 1 is 6.261116, which satisfies the following conditions: average distance within clusters < average distance between clusters < average distance between clusters.

[0094] For cut point 3, the subgroup is divided into three clusters with the following average intra-cluster distances: Cluster 1, 2.133396; Cluster 2, 1.31148; Cluster 3, 1.212779. The average distance between Cluster 1 and Cluster 2 is 5.9722, between Cluster 1 and Cluster 3 is 6.034874, and between Cluster 2 and Cluster 3 is 10.82365. The average cluster distance is 4.192607.

[0095] The separation degree of cluster 1 to cluster 2 is 2.799387, and the separation degree of cluster 1 to cluster 3 is 2.828765, which satisfies the following conditions: average distance within clusters < average distance between clusters < average distance between clusters; The separation degree of cluster 2 to cluster 1 is 4.553786, and the separation degree of cluster 2 to cluster 3 is 8.253001, which satisfies the following: average distance within clusters < average distance between clusters < average distance between clusters; The separation degree of cluster 3 to cluster 1 is 4.976073, and the separation degree of cluster 3 to cluster 2 is 8.924669, which satisfies the following: average distance within clusters < average distance between clusters < average distance between clusters.

[0096] For cut point 4, the subgroup is divided into four clusters with the following average intra-cluster distances: Cluster 1, 2.133396; Cluster 2, 0; Cluster 3, 1.212779; Cluster 4, 0.5489684. The average distance between clusters 1 and 2 is 4.536378, between clusters 1 and 3 is 6.034874, between clusters 1 and 4 is 6.690111, between clusters 2 and 3 is 9.534946, between clusters 2 and 4 is 2.401862, and between clusters 3 and 4 is 11.468. The average cluster distance is 4.21331.

[0097] The separation degree of cluster 1 to cluster 2 is 2.126365, the separation degree of cluster 1 to cluster 3 is 2.828765, and the separation degree of cluster 4 is 3.135898, which satisfies the following conditions: average distance within clusters < average distance between clusters < average distance between clusters; Cluster 2 has a separation degree of Inf to all other clusters, which does not satisfy the following condition: average distance within clusters < average distance between clusters < average distance between clusters; The separation degree of cluster 3 to cluster 1 is 4.976073, the separation degree of cluster 2 is 7.862067, and the separation degree of cluster 4 is 9.45597, which satisfies the following conditions: average distance within clusters < average distance between clusters < average distance between clusters; The separation degree of cluster 4 to cluster 1 is 12.1867, the separation degree of cluster 2 is 4.375228, and the separation degree of cluster 3 is 20.89009, which does not satisfy the condition: average distance within clusters < average distance between clusters < average distance between clusters.

[0098] The calculated data shows that the minimum separation degree (2.799387) of the cluster partition obtained by cut point 3 is larger than that of cut points 2 and 4, indicating that cut points 2 and 4 either failed to achieve sufficient cluster separation or over-separated the clusters, resulting in poor separation effect. (Separation degree) Setting this value to 2.126365 in this data setting is more consistent with the actual data. Additionally, at cut point 4, the positional relationship between cluster 2 and cluster 4 does not meet the requirements of this invention.

[0099] Therefore, in summary Based on location and other factors, this data can be reasonably divided into three categories, and considering the actual distribution of the data, this three-cluster division is also quite consistent with the data category classification. The distribution of cell subpopulations and categories corresponding to each data point is as follows: Figure 3C As shown, the results indicate that although the classification of cancer cell subpopulations is not exactly the same in each dataset, there are three major categories of these cell subpopulations across the three datasets. These three categories show good separation between classes and good clustering within classes in the UMAP plot, indicating that the classification of the three categories is consistent with the true distribution of the data.

[0100] Finally, according to the technical solution of the present invention, the characteristic gene sets are processed to obtain the characteristic gene sets corresponding to each type of cell: GS-A, GS-B, and GS-C.

[0101] Example 2 Gene expression profiles of 102 breast cancer patients before and after treatment were obtained from the GEO database (retrieval number: GSE147271). This data underwent data standardization, converting it into behavioral genes and listing them in a matrix format. Following the approach in Application 1, the enrichment scores of the three characteristic gene sets obtained in Example 1 were calculated in this data using the GSEA package in the R environment. The distribution of each characteristic gene set before and after treatment is shown below. Figure 4 As shown in Figure A, the results indicate that the enrichment fraction of GS-C decreased after treatment, while that of GS-A and GS-B increased. The Cohen's d effect sizes of GS-A, GS-B, and GS-C were 1.432664, 1.290397, and -0.486555, respectively.

[0102] The `roc` function from the `pROC` package was used to predict treatment response based on treatment information for each patient and the enrichment score of each feature gene set on the gene expression profile of each patient. The results are as follows: Figure 4 As shown in Figure B, the AUCs of the three sets of characteristic genes used to predict treatment response were: GS-A, 94.35% (89.88%-98.83%); GS-B, 89.35% (82.76%-95.95%); and GS-C, 64.68% (53.69%-75.67%).

[0103] If the Cohen's d effect size is set to a threshold T value of 0.8 for a large effect and 0.5 for a moderate effect, then the GS-C labeled cell subpopulation is determined to be a moderately low effective response subpopulation for treatment (suppressed after treatment), while the GS-A and GS-B labeled cell subpopulations are potential ineffective or resistant subpopulations (enriched after treatment), and are strong response subpopulations for treatment.

[0104] Example 3 Gene expression profiles of 42 samples (accession number: GSE184717) from multiple organs of breast cancer metastases and adjacent tissues were obtained from the GEO database. This data underwent data standardization, converting it into behavioral genes and listing them in a matrix format. Following the approach in Application 2, the enrichment scores of the three feature gene sets obtained in Example 1 were calculated in this data using the GSEA package in the R environment. The NMF algorithm was used, with the rank range set to 2-10. The NMF algorithm was run repeatedly for each rank value, 100 times, and the cophenetic correlation coefficients corresponding to each rank value were calculated as follows: Figure 5A As shown, the results indicate that the rank value corresponding to the first peak of the cophenetic correlation coefficient is 2, therefore 2 is the optimal number of clusters. Further NMF clustering is then performed using the Brunet algorithm, resulting in the following sample partitioning: Figure 5B As shown, samples from adjacent tissues and metastatic lesions are clearly distinguished, with class 1 consisting mostly of samples from adjacent tissues and class 2 consisting of samples from metastatic lesions.

[0105] The basismap function was used to calculate and visualize the contribution weight of each feature gene set to the classification of each category. Figure 5C The results showed that GS-C contributed significantly to the segmentation of metastatic tissues, suggesting that its labeled cell subpopulations were more active in metastatic lesion samples; GS-A and GS-B contributed significantly to the segmentation of adjacent tissues, suggesting that their labeled cell subpopulations may participate in or drive the invasion of cancer cells into adjacent or distant tissues. The AUC of the characteristic gene set used to distinguish between the two types of tissues was 87.5% (95% CI: 77.44%-97.56%). Figure 5D ).

[0106] Example 4 Gene expression profiles of 298 breast cancer patients in the pre-treatment cohort were obtained from the GEO database (retrieval number: GSE17705). This data underwent data standardization, converting it into behavioral genes and listing them in a matrix format. Following the approach in Application 3, the enrichment scores of the three characteristic gene sets obtained in Example 1 were calculated in this data using the GSEA package in the R environment.

[0107] For each characteristic gene enrichment score, the optimal cutoff value for each characteristic gene set was automatically determined using the `surv_cutpoint` function (from the `survminer` package), dividing patients into high-enrichment and low-enrichment groups. Kaplan-Meier method was used to plot relapse-free survival curves for the high-enrichment and low-enrichment groups, and Log-rank test was performed. Using p < 0.05 as the threshold, the results showed that high enrichment in GS-C was significantly associated with a high risk of relapse, classifying them as a high-risk relapse subgroup; while low enrichment in GS-A and GS-B was significantly associated with a high risk of relapse, classifying them as protective subgroups. Figure 6 (Figure A in the middle)

[0108] Further, following the scheme in Application 3, the intersection of the high-enrichment GS-C, low-enrichment GS-A, and low-enrichment GS-B groups was used to identify the high-risk group for recurrence; the intersection of the low-enrichment GS-C, high-enrichment GS-A, and high-enrichment GS-B groups was used to identify the low-risk group for recurrence; the remaining group was identified as having an intermediate risk of recurrence. The Kaplan-Meier method was used to plot recurrence-free survival curves for the high, intermediate, and low enrichment groups, and Log-rank tests were performed. The results showed that the recurrence risk was reduced by approximately 90% in the low-risk group (HR: 0.1057, 95% CI: (0.03231, 0.3461)), and by approximately 70% in the intermediate-risk group (HR: 0.2837, 95% CI: (0.17482, 0.4603)). Figure 6 (Figure B in the middle)

[0109] Example 5 Gene expression profiles of 310 breast cancer subtypes were obtained from the GEO database (accession number: GSE25055). This data underwent data standardization, converting it into behavioral genes and listing them in a matrix format. According to the scheme of application four of this invention, the enrichment scores of the three characteristic gene sets obtained in Example 1 were calculated in this data using the GSEA package in the R environment.

[0110] The samples were directly divided into histological subtypes, triple-negative and non-triple-negative breast cancer. The enrichment score distribution showed that GS-C was high in triple-negative breast cancer, while GS-A and GS-B were high in non-triple-negative breast cancer. Figure 7 Figure A shows the activity levels of different cell subpopulations in two breast cancer subtypes. The enrichment scores of three feature gene sets were used to predict breast cancer tissue types. The results showed that GS-B and GS-C could effectively classify different tissue types, with AUCs of 82.53% (95% CI: 78.01%-87.06%) and 81.9% (95% CI: 77.27%-86.54%), respectively. Figure 7Figure B in the middle). This data has been publicly published (Hatzis C, et al. JAMA 2011;305(18):1873-1881). The original report proposed a set of indicators that can predict the clinical outcome of breast cancer and achieved the prediction of clinical outcome for each sample. The AUC of this indicator in distinguishing between triple-negative and non-triple-negative subtypes was only 58.71% (95% CI: 53.46%-63.97%). Figure 7 (See Figure B in the original text). In comparison, the feature gene set of this invention exhibits significantly superior distinguishing power.

[0111] In summary, this invention designs a method for screening characteristic gene sets of breast cancer cell subpopulations based on multi-source data. It cross-validates the data using single-cell transcriptome data from at least three independent sources, selecting robust characteristic genes from the intersection of these genes. This solves the problems of heterogeneity in single datasets and cross-dataset correspondence of cell populations, improving the reproducibility of characteristic genes across datasets. The invention proposes calculating the Euclidean distance between subpopulation centroids based on UMAP spatial coordinates, employing hierarchical clustering and an adaptive rationality criterion to achieve objective quantitative determination of cross-dataset subpopulation correspondence, avoiding the subjectivity and non-reproducibility of traditional visual observation. Three quantitative indicators—intra-cluster average distance, inter-cluster average distance, and global average distance—are introduced. The rationality of clustering is determined by both separation and positional relationship conditions. This criterion does not rely on fixed absolute values ​​and is adaptive to different datasets.

[0112] The applicant declares that the above description is only a specific embodiment of the present invention, but the protection scope of the present invention is not limited thereto. Those skilled in the art should understand that any changes or substitutions that can be easily conceived by those skilled in the art within the technical scope disclosed in the present invention fall within the protection and disclosure scope of the present invention.

Claims

1. A method for screening characteristic gene sets of breast cancer cell subpopulations based on multi-source single-cell transcriptome data, characterized in that, The screening method includes the following steps: S1. Obtain at least three independent breast cancer single-cell transcriptome datasets; S2. Perform quality control on each dataset obtained in step S1, extract the epithelial cell markers, and divide cancer cells that do not express immune cell markers, endothelial cell markers, and mesenchymal cell markers into subpopulations; S3. Screen for highly expressed genes in each subgroup divided in step S2. The screening criteria are: average log2 fold change in expression > 0, and corrected p-value ≤ 0.

01. S4. Perform log-normalization on each dataset and select 2000-5000 hypervariable genes; project each dataset into a shared low-dimensional space; identify the nearest neighbor cell pairs in this space as candidate anchor points; filter and score the anchor points based on neighborhood structural similarity. After obtaining the anchor points, the datasets are integrated into a unified space, and then UMAP dimensionality reduction is performed. S5. After dimensionality reduction using UMAP, extract the coordinates (UMAP1, UMAP2) of each cell in two-dimensional space. For each cell subpopulation k, calculate its center point. : UMAP1 coordinates UMAP1 = the arithmetic mean of the UMAP1 coordinates of all cells within this subpopulation; UMAP2 coordinates UMAP2 = the arithmetic mean of the UMAP2 coordinates of all cells within this subpopulation; For any two cell subpopulations i and j from different datasets, calculate the Euclidean distance between their centroids. : ； S6. Use a hierarchical clustering algorithm to cluster the distance matrix and generate a clustering tree. Select cut points on the clustering tree to divide the subgroups into several clusters. Subgroups within the same cluster are identified as the same large cell class. For a given clustering result, calculate the following three quantitative indicators: Intra-cluster average distance : The arithmetic mean of the Euclidean distances between all pairs of subgroups within the same cluster; Average inter-cluster distance : The arithmetic mean of the Euclidean distances between all pairs of subgroups in different clusters; average cluster distance : The total arithmetic mean of the Euclidean distances within and between all clusters; A cluster partition is considered a large cell group when it meets the following conditions: (1) Resolution condition: ; in, This is the separation threshold, which is used to quantify the degree to which the inter-cluster distance of any given cluster is significantly greater than the intra-cluster distance. (2) Positional relationship conditions: ; S7. For each major cell class after merging in step S6, take the union of all highly expressed genes of cell subpopulations contained in each dataset within that major cell class to obtain the set of characteristic genes corresponding to that class for each dataset. For each class of characteristic genes, if it is a characteristic gene of that class in at least three datasets, it is considered a robust characteristic gene belonging to that cell class, and finally the characteristic gene set of different cell subclasses of breast cancer is obtained.

2. The method for screening characteristic gene sets of breast cancer cell subpopulations based on multi-source single-cell transcriptome data according to claim 1, characterized in that, In step S1, the datasets come from different laboratories or different patient cohorts; Optionally, the epithelial markers in step S2 include at least one of EPCAM, GATA3, CDH1, KRT8, KRT18, or KRT19; Optionally, the immune cell marker includes PTPRC; Optionally, the mesenchymal cell markers include VIM; Optionally, the endothelial cell markers include at least one of PECAM1, VWF, or CDH5; Optionally, the subgroup division in step S2 specifically includes: performing PCA dimensionality reduction and selecting the top 20 principal components; constructing a cell similarity graph based on the shared nearest neighbor algorithm; performing graph clustering using the Louvain or Leiden algorithm, with the resolution parameter set to 0.1-1, and selecting the optimal resolution based on the expression of cell type markers; and obtaining the subgroup division results of cancer cells in each dataset.

3. The method for screening a set of characteristic genes for breast cancer cell subpopulations based on multi-source single-cell transcriptome data according to any one of claims 1-3, characterized in that, The screening described in step S3 includes using the FindAllMarkers function in the open-source tool Seurat package in the R working environment to calculate the highly expressed genes in each subpopulation; Optionally, the method for projecting each dataset into a shared low-dimensional space in step S4 includes canonical correlation analysis; Optionally, the method of integrating the datasets into a unified space includes integrating them using the IntegrateData function in the Seurat package within the R working environment; Optionally, the parameters for the UMAP dimensionality reduction visualization are set as follows: reduction = "pca", dims = 1:

20.

4. The method for screening a set of characteristic genes for breast cancer cell subpopulations based on multi-source single-cell transcriptome data according to any one of claims 1-3, characterized in that, The steps described in step S5 The value range is 1.5 to 3.

0.

5. A set of characteristic genes for breast cancer cell subpopulations based on multi-source single-cell transcriptome data, characterized in that, The set of characteristic genes of the breast cancer cell subpopulation is obtained by the screening method described in any one of claims 1-4.

6. The use of the breast cancer cell subpopulation characteristic gene set and / or its detection reagent based on multi-source single-cell transcriptome data as described in claim 5 in the preparation of any of the following products: (1) Products for dynamic monitoring and evaluation of breast cancer treatment response; (2) Products for assessing and constructing metastatic invasion status of different tissue types of breast cancer; (3) Products for subpopulation analysis and risk stratification of recurrent breast cancer cells; (4) Products for analyzing the cell characteristics and molecular subtyping of different subtypes of breast cancer.

7. A multi-scenario assessment system for breast cancer, characterized in that, The evaluation system includes: (1) Input module: used to receive gene expression profile data of the sample to be tested and to standardize the data; (2) Storage module: stores the set of characteristic genes of breast cancer cell subpopulations based on multi-source single-cell transcriptome data as described in claim 5; (3) Enrichment score calculation module: used to calculate the enrichment score of each sample in each characteristic gene set using gene set enrichment analysis method; (4) Scene determination module: used to evaluate any of the following scenes based on the enrichment score: Dynamic monitoring and assessment of breast cancer treatment response: comparing changes in enrichment scores before and after treatment, and identifying response and resistance subgroups; Assessment and construction of metastatic invasion status of different tissue types in breast cancer: NMF clustering is performed on the enrichment score matrix to determine and construct the tissue type and invasion status of the samples; Breast cancer recurrence cell subpopulation analysis and risk stratification: Based on the multi-gene set enrichment scores, cell subpopulations associated with recurrence were determined and stratified into high / medium / low recurrence risk. Analysis of cellular characteristics and molecular subtyping of different subtypes of breast cancer: Analysis of enrichment fraction distribution to determine the molecular subtype of the sample and the corresponding cellular characteristics of different subtypes; (5) Output module: used to output assessment results and clinical decision recommendations.

8. The breast cancer multi-scenario assessment system according to claim 7, characterized in that, The dynamic monitoring and evaluation of breast cancer treatment response specifically includes: 1) Plot the dynamic response curves of the enrichment fraction of each gene set as a function of treatment time; 2) Calculate the effect size Cohen's d, Cohen's d = (M1 – M2) / SD_pooled, where: M1 is the arithmetic mean of the enrichment scores of the characteristic gene set in the treated population; M2 is the arithmetic mean of the enrichment scores of the characteristic gene set in the untreated population; SD_pooled is the combined standard deviation of the two groups; 3) Make judgments based on the curve change pattern and effect size: Targeted response subpopulation identification: If the enrichment score of a certain characteristic gene set decreases or continues to decrease after treatment, it is determined that the treatment has an inhibitory effect on the cell subpopulation marked by the gene set, and the subpopulation is the targeted response subpopulation of the treatment. Acquired drug resistance warning: If the enrichment score of a certain characteristic gene set increases after treatment, or rebounds after an initial decrease, the cell subpopulation marked by that gene set is determined to be a potential acquired drug resistance subpopulation, indicating that the treatment plan needs to be adjusted. Primary drug resistance identification: If the enrichment score of a certain characteristic gene set does not change significantly before and after treatment, the treatment is deemed ineffective for the cell subpopulation marked by that gene set. This cell subpopulation is a potential primary drug resistance subpopulation, suggesting that the treatment plan needs to be adjusted. Optionally, the assessment and construction of the metastatic invasion status of different tissue types of breast cancer includes: 1) Using the non-negative matrix factorization algorithm, by setting different rank ranges, the non-negative matrix factorization algorithm is repeatedly run for each rank value to calculate the cophenetic correlation coefficient corresponding to each rank value. The rank value corresponding to the first time the cophenetic correlation coefficient reaches the highest value is selected as the optimal number of clusters. If there are multiple local maxima, the rank value with the smallest value is selected. 2) Based on the cluster categories obtained above, perform non-negative matrix factorization clustering using the Brunet algorithm to obtain the sample partitioning; 3) The basismap function is used to calculate and visualize the contribution weight of each feature gene set to the classification of each category, so as to obtain the degree of cancer cell invasion or metastasis in each sample; Optionally, the breast cancer recurrence cell subset analysis and risk stratification includes: 1) For each feature gene enrichment score, the surv_cutpoint function is used to determine the optimal cutoff value for each feature gene set. Patients are divided into high enrichment group and low enrichment group. Then, the Kaplan-Meier method is used to plot the relapse-free survival curves of the high enrichment group and the low enrichment group, and the Log-rank test is performed. With p < 0.05 as the threshold, it is determined whether there is a statistically significant association between the enrichment level of the feature gene set and the risk of relapse, so as to obtain the association between the cell subpopulation marked by the feature gene set and the relapse. 2) For the characteristic gene set that is significantly associated with the risk of relapse, the risk direction is determined according to the survival curve: if the relapse-free survival of the high enrichment group is significantly shorter than that of the low enrichment group, then the high enrichment is determined to be associated with a high risk of relapse; if the relapse-free survival of the low enrichment group is significantly shorter than that of the high enrichment group, then the low enrichment is determined to be associated with a high risk of relapse. 3) For the characteristic gene sets significantly associated with relapse risk, obtain the populations with high and low enrichment for each characteristic gene set: The intersection of each feature set and the population associated with high recurrence risk is used to determine the high recurrence risk population. The intersection of each feature set and the population associated with low recurrence risk is used to determine the population with low recurrence risk. The remaining population was assessed as having a moderate risk of recurrence; 4) Output relapse-related cell subsets and protective cell subsets, and risk stratification results for each patient; Optionally, the analysis of cellular characteristics and molecular subtyping of different breast cancer subtypes includes: The samples were divided according to histological subtypes, and the enrichment fraction distribution of each subtype in each characteristic gene set was calculated to analyze the cell subpopulation characteristics of each molecular subtype. For mixed samples containing multiple histological subtypes or multiple disease stages, nonnegative matrix factorization can be used for automated category identification. For each cluster category, the enrichment score distribution of each feature gene set is calculated, and the cell subpopulation feature enrichment spectrum of that category is output. For longitudinal samples containing different disease stages of the same patient, nonnegative matrix factorization clustering can be used for analysis: Treatment-induced subtype conversion: whether samples before and after treatment belong to different cluster categories; Subtype evolution during metastasis: Do primary lesions and metastatic lesions exhibit different subtype characteristics? 9. A computer-readable storage medium having a computer program stored thereon, characterized in that, When the program is executed by the processor, it implements the functions of the breast cancer multi-scenario assessment system as described in claim 7 or 8.

10. A computer device, characterized in that, The computer device includes the computer-readable storage medium of claim 9.

Images