A spatial transcriptome data clustering method

Abstract

The application discloses a spatial transcriptome data clustering method, comprising the following steps: firstly, preprocessing the original data, and constructing a spot space neighbor network; secondly, further learning spatial information and low-dimensional potential representation of gene expression through a graph attention automatic encoder, learning self-expression coefficient matrixes of different layers in the encoder through a multi-scale self-expression module and fusing the same together, adopting spectral clustering in a deep subspace clustering module, and feeding back the clustering label to a self-supervision module; the self-supervision module guides the learning of the potential representation by constructing a self-supervision path returning to the encoder; finally, biological analysis is carried out; the application learns low-dimensional potential embedding by integrating spatial information and gene expression profiles, and improves the clustering performance of spatial resolution transcriptome by using STMSGAL.

Description

Technical Field

[0001] This invention relates to the field of systems bioinformatics, and in particular to a method for clustering spatial transcriptome data. Background Technology

[0002] Recent advances in spatial transcriptomics (ST) have enabled gene expression analysis and location information within tissues. ST refers to transcriptomics techniques that preserve spatial information and gene expression profiles of tissue samples. Spatial transcriptomics (ST) has made significant progress in recent years. Depending on the data generation method, ST techniques can be categorized into NGS-based (Next-Generation Sequencing) and image-based methods. NGS-based ST techniques obtain spatially resolved data by attaching spatially located barcodes to tissue slices. Therefore, each spot captured in an NGS-based ST dataset typically contains multiple cells. Many NGS-based ST methods have been developed, including Visium from 10XGenomics, GeoMx from NanoString, Slide-Seq, Slide-SeqV2, and Stereo-Seq. These methods obtain RNA transcripts through in situ sequencing or in situ hybridization and preserve cellular spatial information through images of stained tissue samples. Image-based ST techniques, such as STARMap, merFISH, and seqFISH+, typically achieve single-cell or subcellular resolution. Typically, an ST dataset consists of a gene expression matrix, where each row represents a gene, each column represents a spot or cell, and a spatial location matrix, which records the spatial coordinates of the spots or cells.

[0003] Clustering analysis is an indispensable step in transcriptome data analysis. Clustering is crucial for annotating cell types, understanding tissue structure, identifying co-expressed gene modules, and many downstream analyses, such as trajectory inference and intercellular communication. Most existing clustering methods do not effectively utilize available spatial information. These non-spatial methods can be broadly categorized into two types. The first type uses traditional clustering methods, such as K-means and the Louvain algorithm. Depending on the resolution of the spatial transcriptome technique, these methods are limited to a small number of spots, and the clustering results may be discontinuous in tissue slices. The second type uses cell type features defined by single-cell RNA sequences to deconvolve the spots, but these methods are not suitable for resolution at the cellular or subcellular level.

[0004] In fact, despite the remarkable performance achieved in clustering, existing methods still have several issues to address. Many of them use only shallow models to extract features, which may fail to capture the complex underlying structure of large-scale datasets, while others use deep models but still neglect to fully utilize feature information embedded at different layers. Furthermore, some methods do not employ cluster-oriented loss functions, and some fail to effectively leverage cluster labels to guide node representation learning; all of these problems lead to suboptimal clustering results.

[0005] Patent CN114091603A discloses a spatial transcriptome cell clustering method, including the following steps: preprocessing the gene expression of each cell in the spatial transit group; generating an adjacency matrix A based on the cell coordinates of the spatial transit group to obtain a graph structure representation of the spatial transcriptome cells, with cell feature matrix X representing cell gene expression; inputting the adjacency matrix A and cell feature matrix X into a trained graph convolutional neural network model DGI; the graph convolutional neural network model DGI outputs node feature representations with spatial information; and after processing the node feature representations using dimensionality reduction and clustering algorithms, identifying and obtaining the spatial transcriptome cell types. This patent does not employ a clustering-oriented loss function, resulting in suboptimal clustering results. Summary of the Invention

[0006] To address the problems mentioned in the background section, this invention provides a spatial transcriptome data clustering method (STMSGAL) to solve the problems of the prior art.

[0007] The technical solution adopted in this invention is:

[0008] A spatial transcriptome data clustering method includes the following steps:

[0009] S1 preprocesses the raw data: selects highly variable genes, normalizes, and log-normalizes to correct the sequencing depth.

[0010] The S2 spatial neighbor network is constructed by building a spatial neighbor network based on a predefined radius, and then pruning it according to the pre-clustering of gene expression to construct another spatial neighbor network.

[0011] The integration of S3 gene expression and spatial location information is achieved through a graph attention autoencoder, which further learns low-dimensional latent representations of spatial information and gene expression. A multi-scale self-expression module learns and fuses the self-expression coefficient matrices of different layers in the encoder. In the deep subspace clustering module, spectral clustering is employed, and cluster labels are fed back to the self-supervised module. The self-supervised module guides the learning of latent representations by constructing a self-supervised path back to the encoder, calculating the total loss function. The effectiveness of the loss function directly affects whether the desired results can be achieved; a good loss function allows the model to learn in the expected direction. With a reasonable loss function, the model can better extract the main features of the data during training.

[0012] S4 was used for biological analysis.

[0013] Further, the specific steps of step S2 are as follows: Convert the spatial information into an undirected neighbor network according to a predefined radius r, let A be the adjacency matrix of the spatial neighbor network, and let A be the adjacency matrix of the spatial neighbor network if and only if the Euclidean distance between spot i and spot j is less than r. ij =1, for 10x Visium data, the network is set to contain the six nearest neighbors of each point, and then the spatial neighbor network is pruned based on the pre-clustering of gene expression using a cell type-aware module.

[0014] Furthermore, the pre-clustering of gene expression is performed using the Louvain algorithm, which has a small resolution value on PCA embeddings and prunes the edges belonging to different clusters.

[0015] Furthermore, the graph attention autoencoder includes: taking the gene expression from step S1 as input, generating a spot embedding by incorporating neighboring information, where x... i Here, L is the normalized expression for spot i, and L is the number of encoder layers. This is achieved by treating the expression spectrum as the initial spot embedding. The embedding of spot i generated by the k-th encoder layer (k∈{1, 2, ..., L-1}) in layer k is given by the formula:

[0016]

[0017] Among them, W k It is a trainable weight matrix, σ is a non-linear activation function, and S i It is the set of neighbors of spot i in the spatial neighbor network, including spot i itself. It is the edge weight between spot i and spot j in the output of the k-th graph attention layer. The encoder vector at the (k-1)th layer;

[0018] The final encoder layer L does not use an attention layer; the formula is as follows:

[0019]

[0020] Among them: W L It is the trainable weight matrix of the Lth layer. It is the vector of spot i in the encoder at layer L-1;

[0021] The encoder's output is considered the final spot embedding, and and Equivalence means wanting to convert the original rice... Save it, and make a copy. In contrast, the decoder inverts the potential embeddings back into the reconstructed normalized expression spectrum. It is input into the decoder network to reconstruct the embedding, where C is the self-expression matrix. It is the vector representation of spot i in the Lth layer, which is the input of the decoder;

[0022] The formula for reconstructing spot i in layer k-1 of the k-th decoder layer (k∈{2,…,L-1,L}) is:

[0023]

[0024] in: It is the edge weight between spot i and spot j in the (k-1)th graph attention layer decoder. It is the trainable weight matrix of the k-layer decoder. It is the vector of spotj in the k-th layer decoder;

[0025] Similar to the encoder, the final layer formula for the decoder is:

[0026]

[0027] The decoder output is considered a normalized expression for reconstruction. To avoid overfitting, settings are configured separately. and in, It's transpose, att (k) It represents the weights of all edges in the k-th layer.

[0028] Using the sum of squared residuals to minimize the reconstruction loss of node attributes, it can be calculated as follows:

[0029]

[0030] Where: F is the F-norm of the vector;

[0031] Weight decay is equivalent to penalizing the l2 norm of the weights, thereby minimizing the regularization loss. Therefore, the total loss of the graph attention autoencoder can be calculated as follows:

[0032]

[0033] Furthermore, the graph attention autoencoder also employs a self-attention mechanism widely used in graph neural networks. Simply put, the attention mechanism involves dynamic weights, which are parameterized. In the k-th encoder layer, the edge weight from spot i to its neighbor spot j is calculated as follows:

[0034]

[0035] in, and It is a trainable weight vector, and Sigmoid is the activation function;

[0036] Normalize it using the softmax function to obtain Normalization makes spatial similarity weights comparable.

[0037]

[0038] When using a cell type-aware module, both types of spatial neighbor networks employ self-attention mechanisms, allowing... and These represent the spatial similarities learned from the spatial neighbor network and the cell type-aware spatial neighbor network, respectively. The final spatial similarity is the linear sum of these two values, att. ij It is the total spatial similarity (i.e., the normalized edge weight) from spot i to its neighbor spot j in one layer of encoder. It is the linear sum of the similarity based on the spatial neighbor network and the similarity of the spatial neighbor network after using the cell type awareness module.

[0039]

[0040] Here, α is a hyperparameter representing the cell type-aware spatial neighbor network, which is set to 0.5 by default.

[0041] Furthermore, the multi-scale self-expression module integrates all spot embeddings. The final self-expression coefficient matrix is ​​obtained by using individual self-expression coefficients. The self-expression layer is a fully connected layer without any activation or bias, and its weights correspond to the coefficient matrix C. This coefficient matrix is ​​further used to construct the affinity matrix Λ for spectral clustering. Assuming the normalized gene expression encoding after the i-th layer is Z... (k)(k = 1, 2, ..., L), where L represents the total number of layers in the encoder. Each layer of the encoder can be computed to produce a self-expression coefficient matrix C. (k) By integrating different self-expression matrices to fully utilize the information embedded in each layer, the loss of the multi-scale self-expression module can be expressed as:

[0042]

[0043] Where: C (i) Z is the self-expression coefficient matrix of the i-th layer encoder, obtained during model training; (i) Z is the gene expression vector of the i-th encoder layer. (i) =H (i) The matrix diag is a diagonal matrix;

[0044] To prevent the self-expression coefficient matrix C of each layer from being... (i) This becomes redundant, so a regularization loss is introduced, the formula of which is:

[0045]

[0046] Where P is the P-norm of the matrix;

[0047] The total loss of the multi-scale self-expression module can be expressed as:

[0048]

[0049] Furthermore, the deep subspace clustering module includes: using an adaptive method to cluster the self-expression matrices C of different layers. (k) The final self-expression matrix C is obtained by fusion. F This method can be described as

[0050]

[0051] The deep subspace clustering module can then use it to construct an affinity matrix Λ for spectral clustering.

[0052]

[0053] Then, spectral clustering is used on the affinity matrix Λ to obtain the clustering result Y. clu .

[0054] Furthermore, the self-supervised module uses fully connected layers for classification, calculating the cross-entropy loss between the classification result P and the clustering result Q. The motivation for this module is to use the current clustering labels to self-supervised learn the node representation H. (k) Furthermore, the node representation incorporates the current clustering information. Since the result of spectral clustering is an n-dimensional vector representing the relationship between n subspaces (i.e., clusters), the dimension of the fully connected layer is designed to be {d}.L ×D1×D2×D3×k}, where d L The extracted node representation Z (L) The dimensions of the fully connected layer are D1, D2, and D3, respectively, and P is the pseudo-label, which is the classification label of the fully connected layer. n is the number of nodes, m is the number of clusters, and the clustering labels of the deep clustering module are... As the learning objective of P, this module can be constrained using cross-entropy loss, setting iteration T0 and temporary labels. Use Q instead of Y clu As a self-supervised learning objective, Y clu Considering Q as the initial value, Y is used every T0 periods. clu Update Q to make the learning objective Q stable over T0 periods. The loss function of the self-supervised module is:

[0055]

[0056] Furthermore, the total loss function can be expressed as:

[0057]

[0058] Here, λ1 and λ2 are two trade-off parameters of the loss function mentioned above.

[0059] Furthermore, the biological analysis described in step S4 includes the following steps:

[0060] S4.1 Identify the spatial domain using Leiden clustering or Louvain clustering based on the obtained spot embedding feature matrix;

[0061] S4.2 uses the t-test from the SCANPY package to perform differential expression analysis. Differential expression analysis is a major downstream analysis method for transcriptome data. It helps identify biomarkers for new cell types or detect gene signatures of cellular heterogeneity, and further provides data for other secondary analyses (such as pathway analysis and network analysis).

[0062] S4.3 Finally, trajectory inference is performed, which allows us to better study the underlying dynamics of biological processes such as cell development, differentiation, and immune responses.

[0063] Compared with the prior art, the beneficial effects of the present invention are:

[0064] The advantage of this invention lies in the fact that STMSGAL uses a graph attention mechanism to consider spatial neighbor information. It adaptively learns the edge weights of the SNN and further uses them to update the point representation by collectively aggregating information from its neighbors.

[0065] The advantage of this invention lies in using a multi-scale self-representation module to learn the self-representation coefficient matrices of different layers in the encoder and fuse them together. Some methods, such as SEDR and CCST, only utilize the representation of the deepest hidden layer of the encoder for clustering tasks, which wastes a large amount of useful feature information embedded in other layers of the encoder. The multi-scale self-representation module, however, fully utilizes the node representation relationships of each layer of the encoder, helping to learn a more discriminative self-representation coefficient matrix using embedded multi-scale information, and mapping to more accurate subspace clustering.

[0066] The advantage of this invention lies in employing a deep subspace clustering module to obtain cluster labels with a clustering loss function, and introducing a self-supervised module to effectively guide the learning of latent representations. Therefore, the learned node representations will contain clustering information, which will greatly help achieve better clustering performance. Attached Figure Description

[0067] Figure 1 STMSGAL flowchart;

[0068] Figure 2 Clustering results of STMSGAL on the dorsolateral prefrontal cortex dataset;

[0069] Figure 3 Clustering accuracy of STMSGAL under different hyperparameters in the multi-scale self-expression module;

[0070] Figure 4 Clustering results of STMSGAL on the human breast cancer (Block ASection 1) dataset;

[0071] Figure 5 Clustering results of STMSGAL on the human breast cancer (Ductal Carcinoma InSitu) dataset;

[0072] Figure 6 Clustering results of STMSGAL on the Adult Mouse Brain (FFPE) human breast cancer dataset; Detailed Implementation

[0073] To clearly illustrate the technical features of this solution, the present invention will be described in detail below through specific embodiments and in conjunction with the accompanying drawings. Many specific details are set forth in the following description to provide a thorough understanding of this application; however, this application may also be implemented in other ways different from those described herein. Therefore, the scope of protection of this application is not limited to the specific embodiments disclosed below. Unless otherwise defined, all technical terms used herein have the same meaning as commonly understood by those skilled in the art. The technical terms used herein are for the purpose of describing specific embodiments only and are not intended to limit the scope of protection of the present invention.

[0074] Example 1

[0075] Please refer to Figure 1 The present invention provides an embodiment of a spatial transcriptome data clustering method (STMSGAL), comprising the following steps:

[0076] S1 preprocesses the raw data.

[0077] Four 10x Visium datasets were collected: human dorsolateral prefrontal cortex (DLPFC) tissue, mouse brain, and two human breast cancer datasets. Maynard et al. manually labeled the cortical layers (L1-L6) and white matter (WM) of 12 DLPFC slices using marker genes and cellular structures. This is a publicly available 10xVisium ST benchmark dataset, with the number of spots per slice ranging from 3498 to 4789. The DLPFC datasets can be downloaded from the spatialLIBD database (http: / / spatial.libd.org / spatialLIBD).

[0078] The Adult Mouse Brain (FFPE) dataset contains 2,264 spots and 19,465 genes. The Human Breast Cancer (Block ASection 1) dataset contains 3,798 spots and 36,601 genes, and the Human Breast Cancer (Ductal Carcinoma In Situ) dataset contains 2,518 spots and 19,743 genes. Both can be downloaded from the 10x Genomics website (https: / / support.10xgenomics.com / spatial-gene-expression / datasets).

[0079] Table 1

[0080]

[0081] As shown in Table 1, spots outside the main tissue regions were first removed from all datasets. Then, the original gene expression data were logarithmically transformed using the SCANPY package and normalized according to the library size. Finally, an appropriate number of hypervariable genes were selected as input.

[0082] Construction of S2 Spatial Neighbor Network

[0083] To incorporate the similarity of a given spot's neighboring spots, spatial information is transformed into an undirected neighbor network based on a predefined radius r. Let A be the adjacency matrix of the spatial neighbor network (SNN). A is defined as follows: A = ... ij =1. For 10x Visium data, the network is set up to contain the six nearest neighbors of each point. STMSGAL then employs a cell type-aware module to prune the SNN based on pre-clustering of gene expression. Specifically, the pre-clustering of gene expression is performed using the Louvain algorithm with a small resolution value (default setting 0.2) on PCA embeddings, pruning edges belonging to different clusters.

[0084] Integration of S3 gene expression with spatial location information

[0085] S3.1 uses a graph attention autoencoder to further learn low-dimensional latent representations of spatial information and gene expression.

[0086] The encoder in the graph attention autoencoder takes normalized gene expression as input and generates spot embeddings by incorporating neighboring information. Let x i Let L be the normalized expression for spot i, and L be the number of encoder layers. This is achieved by treating the expression spectrum as the initial spot embedding. The embedding of spot i generated by the k-th encoder layer (k∈{1, 2, ..., L-1}) in layer k is given by the formula:

[0087]

[0088] Among them W k It is a trainable weight matrix, σ is a non-linear activation function, and S i It is the set of neighbors of spot i in the SNN (including spot i itself). It is the edge weight between spot i and spot j in the output of the k-th graph attention layer. It is the vector of spotj in the (k-1)th layer, which is the encoder layer before the kth layer.

[0089] Layer L is the last encoder layer, and its output is used for clustering. Writing it as layer L makes the formula easier to visualize. The Lth encoder layer does not use an attention layer, and the formula is:

[0090]

[0091] Among them: W L It is the trainable weight matrix of the Lth layer. It is the vector of spot i in the (L-1)th layer encoder, which is the penultimate layer encoder.

[0092] The encoder's output is considered the final spot embedding, and In contrast, the decoder reverses the latent embeddings back to the reconstructed normalized representation spectrum. The decoder and encoder are symmetric and have basically the same formula structure, with only the vector superscripts being different. In general, a vector in a certain layer is calculated based on the vectors of the previous layer encoder or decoder through weights and activation functions. It is input into the decoder network to reconstruct the embedding, where C is the self-expression matrix. It is the vector representation of spot i in the Lth layer, which is the input of the decoder. Adding a hat symbol to the vector represents the encoder vector.

[0093] The formula for reconstructing spot i in the k-1 layer of the k-th decoder layer (k∈{2,…,L-1,L}) is:

[0094]

[0095] in: It is the edge weight between spot i and spot j in the (k-1)th graph attention layer decoder. It is the trainable weight matrix of the k-th layer decoder. It is the vector of spotj in the k-th layer decoder.

[0096] Similar to the encoder, the final layer of the decoder is formulated as follows:

[0097]

[0098] The decoder's output is considered a normalized expression of the reconstruction. To avoid overfitting, STMSGAL sets... and in, Yes, W k Transpose The att value is the edge weight between spot i and spot j in the output of the k-th graph attention layer. (k) This refers to the weights of all edges in the k-th layer.

[0099] Using the sum of squared residuals to minimize the reconstruction loss of node attributes, it can be calculated as follows:

[0100]

[0101] Where: F is the F-norm of the vector;

[0102] Weight decay is equivalent to penalizing the l2 norm of the weights, thereby minimizing the regularization loss. Therefore, the total loss of the graph attention autoencoder can be calculated as follows:

[0103]

[0104] The graph attention autoencoder also employs a self-attention mechanism widely used in graph neural networks. In the encoder, the number of hidden layer neurons is much lower than that in the input layer, so fewer features (neurons) can be used to represent the input data, thus achieving dimensionality reduction. The decoder can restore the dimensionality-reduced data to data of the same dimension as the initial input data. Generally, formula (5), i.e., the sum of squared residuals, is used to calculate the loss of the encoder input data and the decoder output data. The purpose of model training is to make its loss smaller and smaller, i.e., the input and output are closer and closer, so that the low-dimensional data obtained by the encoder can better represent the original data. This low-dimensional data saves a lot of computational cost when clustering. For the attention mechanism, simply put, it is a weight, and it is a dynamic weight. The weight is parameterized and learned in the network, while the weight of the graph encoder without attention is static and determined by the graph structure of the input data. Therefore, the attention model has great flexibility, and great flexibility means that the model has great potential for performance improvement. The attention mechanism is simply a dynamic weight, and the weight is parameterized. In the k-th encoder layer, the edge weight from spot i to its neighbor spot j is calculated as follows:

[0105]

[0106] in, and It is a trainable weight vector, and Sigmoid is the activation function;

[0107] Normalize it using the softmax function:

[0108]

[0109] When using a cell type-aware module, both types of spatial neighbor networks employ self-attention mechanisms, allowing... and These represent the spatial similarities learned from the spatial neighbor network and the cell type-aware spatial neighbor network, respectively. The final spatial similarity is the linear sum of these two values.

[0110]

[0111] Here, α is a hyperparameter representing the cell type-aware spatial neighbor network, which is set to 0.5 by default.

[0112] S3.2 learns and merges the self-expression coefficient matrices of different layers in the encoder through a multi-scale self-expression module.

[0113] Excellent subspace clustering methods rely on the self-representative properties of the data, meaning that each data point in a subspace union can be represented as a linear combination of the other data. The multi-scale self-representative module integrates all spot embeddings. The final self-expression coefficient matrix is ​​obtained by using individual self-expression coefficients. The self-expression layer is a fully connected layer without any activation or bias, and its weights correspond to the coefficient matrix C, which is further used to construct the affinity matrix Λ for spectral clustering. The original objective function in the deep subspace clustering network consists of regularization and self-expression. Since the deep subspace clustering network ignores the multi-scale information from other layers of the encoder, the original self-expression method is developed by considering the multi-scale information from different layers. Assume that the normalized gene expression encoding after the i-th layer is Z. (k) (k = 1, 2, ..., L), where L represents the total number of layers in the encoder. Each layer of the encoder can be computed to produce a self-expression coefficient matrix C. (k) For example, C (1) Let represent the self-representation matrix of layer 1 in the encoder. The information embedded in each layer can be fully utilized by integrating different self-representation matrices; the loss of the multi-scale self-representation module can be expressed as:

[0114]

[0115] Where: C (i) Z is the self-expression coefficient matrix of the i-th layer encoder, obtained during model training; (i) Z is the gene expression vector of the i-th encoder layer. (i) =H (i) The matrix diag is a diagonal matrix;

[0116] To prevent the self-expression coefficient matrix C of each layer from being... (i) This becomes redundant, so a regularization loss is introduced, the formula of which is:

[0117]

[0118] Where P is the P-norm of the matrix;

[0119] The total loss of the multi-scale self-expression module can be expressed as:

[0120]

[0121] In S3.3, the deep subspace clustering module employs spectral clustering and feeds back the clustering labels to the self-supervised module. The deep subspace clustering module uses an adaptive method to combine the self-expression matrices C of different layers. (k) The final self-expression matrix C is obtained by fusion. F This method can be described as

[0122]

[0123] The deep subspace clustering module can then use it to construct an affinity matrix Λ for spectral clustering.

[0124]

[0125] Then, spectral clustering is used on the affinity matrix Λ to obtain the clustering result Y. clu Spectral clustering is a fundamental clustering method. It primarily involves calculating and decomposing the normalized Laplacian matrix, followed by k-means clustering. Like other basic activation functions, its formulas are not presented in detail; they are simply called directly. Spectral clustering is similar; therefore, detailed explanations and formulas are not provided here. The clustering result Y... clu The self-supervised module mentions that the clustering label Y... clu Used as the learning objective of P, but with Y clu Updating P is risky because it can hinder model learning and convergence. Therefore, a training trick is employed: setting iteration T0 and temporary labels Q. Q is used instead of Y. clu As a self-supervised learning objective, and Y clu Consider it as the initial value of Q.

[0126] The S3.4 self-supervised module guides the learning of latent representations by constructing a self-supervised path that returns to the encoder.

[0127] The self-supervised module uses fully connected layers for classification, calculating the cross-entropy loss between the classification result P and the clustering result Q. The module is motivated by using the current clustering labels to self-supervise the learning of node representation H. (k)The node representation incorporates the current clustering information. Although the current spectral clustering result may not represent the correct category labels for all samples, it still retains much useful information. Once a better node representation is implemented, clustering performance will improve, and more accurate cluster labels can be obtained for self-supervised training in the next iteration. Since the result of spectral clustering is an n-dimensional vector representing the relationships among n subspaces (i.e., clusters), the dimension of the fully connected layer is designed to be {d}. L ×D1×D2×D3×k}, where d L The extracted node representation Z (L) The dimensions of the fully connected layer are D1, D2, and D3, respectively. P is the pseudo-label, which is the classification label of the fully connected layer. n is the number of nodes, and m is the number of clusters. The cluster labels of the deep clustering module are... As the learning objective of P, this module can be constrained using cross-entropy loss. Note the cluster label Y. clu It changes in each epoch, represented by Y. clu Updating P is risky because it can hinder model learning and convergence. To address this, a training trick is employed, setting iteration T0 and temporary labels. Use Q instead of Y clu As a self-supervised learning objective, Y is used. clu Considering Q as the initial value, Y is used every T0 periods. clu Update Q to make the learning objective Q stable over T0 periods. The loss function of the self-supervised module is:

[0128]

[0129] The total loss function of S3.5STMSGAL can be expressed as:

[0130]

[0131] S4 was used for biological analysis.

[0132] S4.1 Identify the spatial domain using Leiden clustering or Louvain clustering based on the obtained spot embedding feature matrix;

[0133] S4.2 uses the t-test from the SCANPY package to perform differential expression analysis. Differential expression analysis is a major downstream analysis method for transcriptome data. It helps identify biomarkers for new cell types or detect gene signatures of cellular heterogeneity, and further provides data for other secondary analyses (such as pathway analysis and network analysis).

[0134] S4.3 Finally, trajectory inference is performed, which allows us to better study the underlying dynamics of biological processes such as cell development, differentiation, and immune responses.

[0135] Example 2

[0136] This invention provides an embodiment of a spatial transcriptome data clustering method (STMSGAL), including all the steps described in the embodiment, and using five metrics to evaluate the performance of the STMSGAL algorithm: Adjusted Land coefficient (ARI), DB coefficient (DB), CH coefficient (CH), silhouette coefficient (SC), and S_dbw coefficient. ARI calculates the similarity between the cluster labels predicted by the algorithm and the true cluster labels; the closer the value is to 1, the more similar the clusters. If there are no true labels, four commonly used clustering metrics are compared: SC, DB, CH, and S_dbw. SC is calculated using the average intra-cluster distance *a* and the average nearest cluster distance *b* for each sample; the closer the value is to 1, the more similar the clusters. DB is defined as the average similarity measure between each cluster and its most similar cluster, where similarity is the ratio of intra-cluster distance to inter-cluster distance. The lowest score is zero, and lower values ​​indicate better clusters. CH evaluates the effectiveness of clusters based on the average of the sum of squares between and within clusters. The silhouette coefficient score for a sample is (ba) / max(a,b), with the best value being 1 and the worst value being -1. The basic idea behind the S_DBW exponent is that for every pair of cluster centers, at least one cluster center should have a density greater than the density of its midpoint. The minimum value of S_DBW signifies the optimal number of clusters.

[0137] STMSGAL is compared with two non-spatial algorithms (SCANPY and Seurat) and four recently published spatial clustering algorithms (stLearn, SEDR, CCST, STAGATE, and DeepST), with parameter settings shown in Table 2.

[0138] Table 2

[0139]

[0140]

[0141] Figure 2The results show the clustering results of STMSGAL on the dorsolateral prefrontal cortex dataset. The results demonstrate that the spatial domain identified by STMSGAL is consistent with the hand-annotated DLPFC data and the definition of cortical stratification in neuroscience, and STMSGAL can effectively identify the expected cortical structures. Overall, the median ARI of the method is 0.9 percentage points higher than the second-ranked method. In DLPFC slice 151508, STMSGAL clearly delineates layer boundaries, achieving the best clustering accuracy (ARI = 0.490). In contrast, the non-spatial method SCANPY's cluster assignments can roughly follow the expected stratification pattern, but its cluster boundaries are discontinuous, with a large number of outliers, affecting clustering accuracy. Notably, the algorithm utilizing spatial information significantly outperforms the non-spatial clustering method SCANPY. These results illustrate the superiority of STMSGAL in spatial domain identification and the necessity of utilizing spatial information.

[0142] like Figure 3 The robustness of STMSGAL was further tested by comparing the clustering accuracy of STMSGAL under different hyperparameters, and it was found that STMSGAL is not sensitive to the trade-off parameters of the multi-scale self-expression module.

[0143] To illustrate the generalization ability of the method to cancer tissue, it was first applied to the 10xVisium human breast cancer dataset (Block A Section 1). Table 3 shows the comparison results of each method, with the best value for each metric highlighted in bold. The method's ARI value is 6 percentage points higher than the second-ranked method. Compared with the regions identified by other spatial algorithms, the structural domains obtained by the method are highly consistent with manually annotated regions. Figure 4 The results show that STMSGAL discovered regions with greater regional continuity and achieved the best clustering accuracy (ARI = 0.605). As expected, STMSGAL identified some subclusters within the tumor region, such as domains 4 and 13. Meanwhile, regions with lower heterogeneity, such as healthy regions, showed high consistency with manual annotations, demonstrating the advantages of STMSGAL.

[0144] Table 3

[0145]

[0146]

[0147] The study primarily used differential expression analysis to compare the transcriptional differences within tumors in domain 1 (ductal / lobular carcinoma, DCIS / LCIS), domains 4, and 13 (invasive ductal carcinoma, IDC). In domain 1, CPB1 significantly distinguished DCIS from other subtypes of breast cancer. COX6X may play a crucial role in differentiating between ER-positive (estrogen receptor) and ER-negative subtypes. Notably, IL6ST expression was associated with a lower risk of invasion, metastasis, and recurrence, and thus with better prognosis. On the other hand, IGFBP5 and CRISP3 were observed in domains 4 and 13, which may be related to breast cancer treatment. The study showed that downregulation of CRISP3 not only significantly inhibited the migration and invasion of breast cancer cells but was also associated with the inhibition of the ERK1 / 2MAPK signaling pathway. This implies that CRISP3 could serve as an indicator for assessing the clinical prognosis of breast cancer patients. Furthermore, IGFBP5 is a potential therapy for controlling tamoxifen resistance in breast cancer. In summary, this demonstrates that STMSGAL can identify more refined regions with different biological functions.

[0148] As shown in Table 4 and Figure 5 As shown, another ST dataset for human breast cancer (Ductal Carcinoma In Situ) was also examined. The STMSGAL domain was smoother and more continuous than other spatial algorithms, reflecting STMSGAL's ability to perform finer segmentation of complex tissues. AZGP1 levels determined the histological grade of breast cancer tumors in Region 1, while CD24 was a key marker for triple-negative breast cancer. In particular, ERBB2 overexpression classified breast cancer into an ERBB2 / HER2-positive subclass, accounting for approximately 20%–30% of all breast cancer types and generally associated with poor prognosis. Accordingly, ERBB2 targeting has proven to be an effective approach against ERBB2-positive breast cancer.

[0149] Table 4

[0150]

[0151]

[0152] As shown in Table 4, STMSGAL outperforms previous methods on mouse brain datasets. Figure 6As shown, the clustering results identified by SCANPY broadly delineate tissue structures containing different cell types, but lack the ability to identify small spatial domains. For example, SCANPY's clustering assignment failed to identify the "cord-like" structure—Angus Amon—and the "arrow-like" structure—Dentate gyrus—within the hippocampus. Furthermore, SEDR only smooths domain boundaries but cannot depict smaller spatial domains. The direct application of STMSGAL brings some improvements in spatial domain identification, enabling the identification of the granular cell layer within the Dentate gyrus.

[0153] Example 3

[0154] One embodiment of the present invention is a spatial transcriptome data clustering method (STMSGAL), which includes all the steps in Example 1, and further analyzes the multi-scale self-expression module and the self-supervised module.

[0155] As shown in Table 5, an ablation study was conducted on 10x Visium data to further analyze the importance of the multi-scale self-expression module and the self-supervised module in STMSGAL. Each loss corresponds to a proposed module. It is the reconstruction loss of the normalized expression. It is the loss of the multi-scale self-expressive module. This is a loss of the self-supervised module.

[0156] Table 5

[0157]

[0158]

[0159] Table 5 shows that the multi-scale self-expression module and the self-supervised module can significantly contribute to clustering performance, and they work well with the graph attention autoencoder. The attempt to use the self-supervised module is meaningful; it utilizes the learning of cluster labels from supervised point embeddings and helps learn more current clustering information. Compared to the self-supervised module, the multi-scale self-expression module exhibits a more significant impact, implying that an appropriate clustering-oriented loss function plays a crucial role in clustering tasks.

[0160] It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above, and that the invention can be implemented in other specific forms without departing from its spirit or essential characteristics. Therefore, the embodiments should be considered illustrative and non-limiting in all respects, and the scope of the invention is defined by the appended claims rather than the foregoing description. Thus, all variations falling within the meaning and scope of equivalents of the claims are intended to be included within the present invention, and no reference numerals in the claims should be construed as limiting the scope of the claims.

Claims

1. A method of spatial transcriptome data clustering, the method comprising: Includes the following steps: S1. Preprocessing of raw data: Selecting highly variable genes, normalization, and log normalization to correct the sequencing depth of the data; S2. Construction of spatial neighbor network: A spatial neighbor network is constructed based on a predefined radius, and then pruned according to the pre-clustering of gene expression to construct another spatial neighbor network; S3. Integration of gene expression and spatial location information; S31. Further learn low-dimensional latent representations of spatial information and gene expression through graph attention autoencoders; S32. Learn and fuse the self-expression coefficient matrices of different layers in the encoder through a multi-scale self-expression module; The multi-scale self-expression module obtains a final self-expression coefficient matrix by integrating all the spot embeddings in a single self-expression coefficient The normalized gene expression code after the i-th layer is denoted as , where N represents the total number of layers of the encoder, and each layer of the encoder calculates a self-expression coefficient matrix The loss of the multi-scale self-expression module is represented as: in: It is the self-expression coefficient matrix of the i-th layer encoder, which is obtained during training in the model; It is the gene expression vector of the i-th layer encoder. = The matrix diag is a diagonal matrix; S33. In the deep subspace clustering module, spectral clustering is used, and the clustering labels are fed back to the self-supervised module; The original objective function in the deep subspace clustering network consists of regularization and self-expression. The deep subspace clustering module uses an adaptive method to combine the self-expression matrices of different layers. The final self-expression matrix is ​​obtained by fusion. The method is described as The deep subspace clustering module uses it to construct an affinity matrix. Perform spectral clustering: ； Then, in the affinity matrix Clustering results obtained using spectral clustering ; S34. The self-supervised module guides the learning of latent representations by constructing a self-supervised path that returns to the encoder; The self-supervised module uses fully connected layers for classification and calculates the classification results. and clustering results The cross-entropy loss between them, where , It is the number of nodes. It is the cluster number, which represents the cluster labels of the deep clustering module. As The learning objective is to constrain this module using cross-entropy loss and set the iteration. and temporary tags ,use Instead As a self-supervised learning objective, Considered The initial value, every Used in this period renew , to make learning objectives exist The loss function of the self-supervised module is expressed as follows: (The loss function is stable over a certain period.) S35. Calculate the total loss function of STMSGAL, expressed as: ； S4 was used for biological analysis.

2. The spatial transcriptome data clustering method according to claim 1, characterized in that, The specific steps of step S2 are as follows: Convert the spatial information into an undirected neighbor network according to the predefined radius r, and let... Let the adjacency matrix of the spatial neighbor network be such that the spot and spot When the Euclidean distance between them is less than r, For 10x Visium data, the network is set up to contain the six nearest neighbors of each point, and then a cell type-aware module is used to prune the spatial neighbor network based on the pre-clustering of gene expression.

3. The spatial transcriptome data clustering method according to claim 2, characterized in that, The pre-clustering of gene expression was performed using the Louvain algorithm, which uses PCA embeddings with resolution values ​​to prune the edges belonging to different clusters.

4. The spatial transcriptome data clustering method according to claim 1, characterized in that, The graph attention autoencoder includes: taking the gene expression from step S1 as input, generating spot embeddings by incorporating neighboring information, and setting... It is spot The normalized expression, It is the number of layers in the encoder, by treating the expression spectrum as the initial spot embedding. , No. indivual The spot generated by the encoder layer exist The embedding formula in the layer is: in, σ is a trainable weight matrix, and σ is a non-linear activation function. It is a spot in a spatial neighbor network The neighborhood set that includes spot itself, It is the first Note the spot in the layer output of the image. and spot Edge weights between them The encoder vector of spotj at layer k-1; Last layer The encoder layer does not use an attention layer, and the formula is: in: It is the trainable weight matrix of the Lth layer. It is the vector of spot i in the encoder at layer L-1; The encoder's output is considered the final spot embedding, and In contrast, the decoder reverses the potential embeddings back into the reconstructed normalized expression spectrum. It is input into the decoder network to reconstruct the embedding, where ， It is a self-expression matrix. It is the vector representation of spot i in the Lth layer, which is the input of the decoder; No. indivual ( ) The spot reconstructed by the decoder layer exist The formula for the layer is: in: It is the first -1 image attention layer decoder spot and spot The boundary rights between them It is the trainable weight matrix of the k-layer decoder. It is the vector of spot j in the decoder at the k-th layer; Similar to the encoder, the final layer formula for the decoder is: The decoder output is considered a normalized expression for reconstruction. To avoid overfitting, settings are configured separately. and ,in, It is transpose. These are the weights of all edges in the k-th layer; The reconstruction loss of node attributes is minimized using the sum of squared residuals, calculated as follows: in: Vector F norm; Weight decay is equivalent to the weight's... The norm imposes a penalty, thereby minimizing the regularization loss. Therefore, the total loss of the graph attention autoencoder is calculated as follows:

5. The spatial transcriptome data clustering method according to claim 1, characterized in that, The graph attention autoencoder also employs a self-attention mechanism for graph neural networks. The attention mechanism uses dynamic weights, which are parameterized. In each encoder layer, from spot to its neighbor spot The edge weights are calculated as follows: in, and It is a trainable weight vector. It is an activation function; Normalize it using the softmax function: When using a cell type-aware module, both types of spatial neighbor networks employ self-attention mechanisms, allowing... and These represent the spatial similarities learned from the spatial neighbor network and the cell type-aware spatial neighbor network, respectively. The final spatial similarity is the linear sum of these two values. It is a spot in a layer encoder. to its neighbor spot The total spatial similarity, i.e., the normalized edge weights, is a linear sum of the similarity based on the spatial neighbor network and the similarity of the spatial neighbor network after using the cell type-aware module. in, It is a hyperparameter representing the cell type-aware spatial neighbor network, with a default setting of 0.

5.

6. The spatial transcriptome data clustering method according to claim 1, characterized in that, The self-expression layer is a fully connected layer with no activations or biases; its weights correspond to the coefficient matrix. This coefficient matrix is ​​further used to construct the affinity matrix. Perform spectral clustering and analyze the self-expression coefficient matrix at each layer. Introducing regularization loss, its formula is: Where P is the P-norm of the matrix; The total loss of the multi-scale self-expression module is expressed as: 。 7. The spatial transcriptome data clustering method according to claim 1, characterized in that, The dimensional design of the fully connected layer of the self-supervised module is as follows: ,in It is the extracted node representation dimensionality and These are the dimensions of the fully connected layer.

8. The spatial transcriptome data clustering method according to claim 1, characterized in that, The biological analysis described in step S4 includes the following steps: S4.1 Identify the spatial domain using Leiden clustering or Louvain clustering based on the obtained spot embedding feature matrix; S4.2 It uses the t-test in the SCANPY package to perform differential expression analysis; S4.3 Finally, trajectory inference is performed.

Images