A trajectory prediction method based on contrastive learning feature dimensionality reduction

Abstract

The present invention discloses a trajectory prediction method based on contrast learning feature dimensionality reduction. The method first obtains the positive sample cell set and the negative sample cell set of each cell according to the original matrix X of the cell to be reduced in dimension. Secondly, the original matrix X of the cell to be reduced in dimension is input into a fully connected neural network to obtain the reduced dimension feature Y of the cell, which is then input into the cell-level pairing module and the cluster-level contrast module respectively, and the loss of the feature is calculated. Then, according to the loss of the feature, the output of the cell-level pairing module and the cluster-level contrast module is optimized. Finally, after the optimization is completed, the feature Y is used to perform trajectory inference to obtain the predicted trajectory. The present invention improves the accuracy of the learned cell features and makes the learned cell features more conducive to trajectory inference.

Description

Technical Field

[0001] The present invention relates to the field of cell trajectory prediction in single-cell sequencing analysis technology, and in particular to a trajectory prediction method based on comparative learning feature dimensionality reduction. Background Art

[0002] Single-cell sequencing (scRNA-seq) technology can analyze gene expression within the transcriptome at the single-cell level, providing great help for people to study cell heterogeneity and biological relationships. Usually, scRNA-seq data analysis includes two main steps: the first step is upstream analysis, which includes denoising the acquired data matrix and clustering cells to annotate cell types or cell states; the second step is downstream analysis based on the results of the first step analysis to explain biological mechanisms, including tasks such as cell trajectory prediction and gene analysis. As one of the important downstream tasks in scRNA-seq data analysis, cell trajectory analysis can model cell data and estimate pseudo-time series to reveal some dynamic processes in biological development, such as cell differentiation processes and key nodes of cell development.

[0003] At present, scholars at home and abroad have made many valuable research results in the field of cell trajectory prediction. These trajectory prediction methods are divided into two categories according to their inference methods: the first category is based on the inference method at the cell level, which uses principal component analysis (PCA) or unified manifold approximation and projection (UMAP) and other dimensionality reduction methods to reduce the gene dimension of each cell, and then construct the development path between cells; the second category is based on the partition level, first using PCA or UMAP to learn the characteristics of the cell, and then divide the cell into multiple independent regions and regard them as different "cell states", and construct the trajectory based on these cell states. However, most of the existing two types of trajectory prediction methods use PCA or UMAP and other methods to reduce the dimension before inference. These dimensionality reduction methods are not specially designed for scRNA-seq cell data. In addition, these methods cannot use the results of upstream analysis (such as cell clustering) to learn cell characteristics, resulting in poor trajectory prediction results. Therefore, it is necessary to design a dimensionality reduction method specifically for trajectory prediction methods. Summary of the invention

[0004] The present invention aims at the deficiency that the dimensionality reduction methods such as PCA or UMAP used in the existing trajectory prediction methods cannot make full use of the prior information, and proposes a scRNA-seq cell data dimensionality reduction method based on contrastive learning for trajectory prediction. The present invention adopts two feature learning components, a cell-level pairing module and a cluster-level contrast module, to learn cell features that help infer cell trajectories. First, a two-layer fully connected neural network is used to map the gene expression data (original features) of scRNA-seq to the dimensionality reduction space to obtain cell representations (dimensionality reduction features), and then the cell-level and cluster-level modules simultaneously optimize the cell representation. Specifically, in each round of training, the method selects a positive sample and a negative sample for each cell sample based on the Euclidean distance of the original cell features, and the cell-level pairing module minimizes the distance between the cell and the positive sample in the dimensionality reduction space while maximizing the distance of the corresponding negative sample. In addition, two enhanced clusters are constructed for each cell type based on the cell type label obtained by the clustering task in the upstream analysis, and the cluster-level module performs contrastive learning in the enhanced clustering to maximize the similarity of the positive cluster pairs and minimize the similarity of all other negative cluster pairs. In short, the cell-level pairing module enables cells to maintain the positional relationship between cells in the original feature space in the dimensionality reduction space, while the cluster-level contrast module prevents cells from being over-dispersed. These two modules ensure that the model learns accurate cell features. Finally, the cell lineage is predicted based on the learned features.

[0005] The specific steps of the method of the present invention include:

[0006] Step (1). Calculate the distance between each cell in the original matrix X of the cells to be reduced, and select the closest k distances for each cell based on the calculated distances. p The cells are taken as the positive sample cell set, and the farthest k n cells as negative sample cells.

[0007] Step (2). Input the original matrix X of the cell to be reduced in dimension into a fully connected neural network to obtain the reduced dimension feature Y of the cell.

[0008] Step (3). Input the dimension reduction feature Y into the cell-level pairwise module and the cluster-level comparison module respectively.

[0009] Step (3.1). In the cell-level module, each round of training is for each cell y i Randomly select a cell from the positive sample set To form a positive sample pair And select a cell from the negative sample set To form a negative sample pair

[0010] Step (3.2). In the clustering module, each round of training is performed on each class c. iAll cells (the category to which each cell belongs obtained by the clustering results in the upstream analysis) are randomly downsampled twice to form two enhanced cell categories and in Represents enhanced class The tth cell in |c i | indicates class c i The number of cells included, γ represents the downsampling rate. As a positive sample pair, all cells generated from different types are classified as or are considered as negative sample pairs.

[0011] Step (4). Calculate the feature loss based on the outputs of the cell-level pairwise module and the cluster-level contrast module, and perform optimization.

[0012] Step (4.1). By minimizing the cell distance of the positive sample pair And maximize the distance of negative sample pairs To make cells maintain the relationship of high-dimensional space in low-dimensional space. Distance between cell sample pairs and The calculation formula is as follows:

[0013]

[0014]

[0015] Where d represents the dimension of the low-dimensional features of the cell. Then the cross entropy loss function Loss is used i To minimize the number of positive samples of cell i and maximize the distance of negative sample pairs The distance is calculated as:

[0016]

[0017] The loss function L of the cell-level pairwise module is obtained by weighted average of the losses of all cells cell :

[0018]

[0019] Step (4.2). By maximizing the cell class c i The similarity between enhanced cell classes is minimized and the similarity of all other negative sample enhanced classes is minimized to prevent the same cells from being over-dispersed in the low-dimensional space. First, each enhanced class is calculated by weighted average The center is considered as its characteristic

[0020]

[0021] Then use the cosine similarity function to calculate the center r of any two enhanced classes i and j i and r j The similarity between s(r i ,r j ), the cosine similarity function is calculated as follows, where ||·|| represents the modulus length of the variable:

[0022]

[0023] By maximizing the similarity of the features of the two enhanced classes, the cells of the two can be brought closer. i Generated enhanced positive sample pairs The similarity of can make cells of the same class closer, and minimizing the similarity of different classes can make cells of different classes more separated in space. Using the cross entropy loss function To calculate c i The enhanced sample a center With its positive sample All its negative sample pairs and The losses are as follows:

[0024]

[0025]

[0026] The loss function L of the cluster-level pairwise module is obtained by weighted averaging the losses of all cell classes. cluster :

[0027]

[0028] Finally, the cell-level loss function L cell And the cluster-level loss function L cluster The overall loss function Loss of the model is obtained by adding them together, and its formula is:

[0029] Loss = L cell +L cluster

[0030] Step (5). After the optimization in step 4 is completed, the dimensionality reduction feature Y learned in step 2 is used for trajectory prediction. First, the Euclidean distance between cells is calculated based on the dimensionality reduction feature Y, and then k is selected for each cell. cThe nearest cells are connected. Then the connections between cells of different classes are summed up, and based on the number of connections between different classes, the minimum spanning tree is constructed between cell classes with the largest number of connections. Finally, a cell class is defined as the initial class, and the cell class from the initial class to any leaf node in the spanning tree is regarded as a cell development trajectory to complete the trajectory prediction.

[0031] The present invention uses a fully connected neural network to obtain low-dimensional features of cells, and together with the clustering or annotation results of the upstream task in the single-cell sequencing analysis process as the input of the model, combines the cell-level paired feature learning module and the cluster-level contrast feature learning module to construct the final cell dimensionality reduction model. The present invention has the following beneficial effects: 1. The analysis results of the upstream task are introduced as prior information to improve the accuracy of the learned cell features; 2. The constructed cell-level and cluster-level modules enable the cell dimensionality reduction features to maintain the cell-cell position relationship in the original feature space while preventing the cells in the same class from being overly dispersed, so that the learned cell features are more conducive to trajectory prediction. BRIEF DESCRIPTION OF THE DRAWINGS

[0032] Figure 1 Network architecture diagram of the dimensionality reduction method for feature learning based on contrastive learning;

[0033] Figure 2 Flowchart of the cell-level pairwise feature learning module;

[0034] Figure 3 Flowchart of cluster-level contrast feature learning module;

[0035] Figure 4 Flowchart of trajectory prediction based on cell dimensionality reduction features. DETAILED DESCRIPTION

[0036] The present invention proposes a dimensionality reduction method for representation learning based on contrastive learning. Figure 1 As shown in Figure 1, the overall network architecture consists of a fully connected neural network, a cell-level pairwise feature learning module, and a cluster-level contrast feature learning module. First, in the fully connected neural network, a two-layer feedforward fully connected neural network is used to map the gene expression data of scRNA-seq to a low-dimensional space to obtain cell features.

[0037] Then, the cell-level paired feature learning module is used to construct a positive sample and a negative sample for each cell for paired learning, so that the cell can maintain the cell-cell position relationship in the original feature space in the reduced dimensional space. At the same time, the cluster-level contrast feature learning module constructs two enhanced cell classes for each cell class and performs contrast learning, so that cells from the same class are more clustered. Finally, the learned cell features are obtained by optimizing the network architecture, and trajectory prediction is performed based on the cell features.

[0038] In order to make the purpose, technical solutions and advantages of the present invention more clearly understood, the implementation steps are described in detail below with reference to the accompanying drawings.

[0039] Step (1). Cell-based gene expression matrix X n×m Calculate the distance between each cell, where n is the number of cells and m is the number of genes contained in each cell. The calculation formula is:

[0040]

[0041] Among them, d i,j Represents the distance between the i-th cell and the j-th cell. According to the calculated distance, the 20 cells closest to each cell i are selected as the positive sample set P i , the farthest n / 5 cells are taken as negative sample set N i .

[0042] Step (2). Convert the original matrix X n×m Input into a two-layer fully connected neural network for dimensionality reduction to obtain the cell's dimensionality reduction feature Y n×d The parameters of the input layer, hidden layer, and output layer of the neural network are [m, 128, d], where m is the number of genes and d is the feature size, which is 16.

[0043] Step (3). Input the dimension reduction feature Y and its enhanced feature into the cell-level pairwise module and cluster-level contrast module respectively to obtain the network loss value:

[0044] Step (3.1). In the cell-level module, Figure 2 As shown, for each cell i, from its positive sample set P i A cell is randomly selected from the From its negative sample set N i A cell is randomly selected from the negative sample pair Then use the following formula to calculate the Euclidean distance between two cells in the positive sample pair and the negative sample pair respectively. and

[0045]

[0046]

[0047] Step (3.2). In the clustering module, Figure 3 As shown, from each cell type A certain number of cells are randomly selected to form two enhanced cell classes and where y i represents cells belonging to the i-th cell class, y i,a Indicates that the i-th cell is selected as the enhanced cell type The cells in The same cell class c i Two enhanced cell classes generated Considered as positive sample pairs, all enhanced classes between different classes is considered as a negative sample pair, where i≠j. Then the weighted average of the cell features in the cell class is used to represent the center of the class The calculation is as follows:

[0048]

[0049] Finally, the cosine similarity function is used to measure the center r of two cell classes i and j i and r j The similarity between s(r i ,r j ), where ||·|| represents the modulus of the variable:

[0050]

[0051] Step (4). Calculate the feature loss based on the outputs of the cell-level pairwise module and the cluster-level contrast module, and perform optimization.

[0052] Step (4.1). In the cell-level module, in order to keep the cell-cell position relationship in the original feature space in the dimensionality reduction space, the cross entropy loss function is used to minimize the distance between the positive sample pairs of each cell i. While maximizing the distance between negative sample pairs To optimize the model. Calculate the loss loss i The formula is:

[0053]

[0054] The loss function L of the cell-level pairwise module is obtained by weighted average of the losses of all cells cell :

[0055]

[0056] Step (4.2). In the clustering module, in order to prevent cells from being overly dispersed, the cell clustering results in the prior information are used to make the distribution of cells in the same class closer. The cross entropy loss function is used To maximize the cell class c i The enhanced sample a center With its positive sample Similarity of center And minimize all its negative sample pairs and Similarity and The loss is calculated as follows:

[0057]

[0058]

[0059] The loss function L of the cluster-level pairwise module is obtained by weighted averaging the losses of all cell classes cluster :

[0060]

[0061] Step (4.2). Finally, the cell-level and cluster-level loss functions are added together to obtain the overall loss function of the model, Loss, whose formula is:

[0062] Loss = L cell +L cluster

[0063] Step (5). After the model optimization is completed, the cell feature Y output by the fully connected neural network is used for trajectory prediction. The following introduces a commonly used trajectory prediction method based on inter-class cell distance.

[0064] Step (5.1). First, use the Euclidean distance to calculate the distance between each cell after dimensionality reduction. The formula is as follows:

[0065]

[0066] where dis i,j represents the distance between cell j and cell i.

[0067] Step (5.2). Secondly, construct an adjacency matrix A n×n , then select the k cells closest to each cell i and set the adjacency matrix to A i,j =1, where j is one of the k selected cells. And all other unselected cell values ​​are set to A i,others =0, k is 20 in the present invention. Then the connection probability between classes is calculated based on the constructed adjacency matrix:

[0068] B i,j =∑A a,b ,

[0069] where a∈c i and b∈c j , B i,j represents the connection probability between cell classes i and j.

[0070] Step (5.3). Figure 4 As shown, according to the calculated inter-class connection rate B i,j A minimum spanning tree is constructed between all classes, and then a starting class is defined. The path from the starting class to any root node on the tree is considered a developmental path.

[0071] In order to facilitate the comparison of the effects of the present invention, the present invention method (scCRT for short) was compared with several traditional prediction methods in this embodiment. The prediction methods used as controls are: (1) scTite: a method for calculating the transition state between cell classes based on entropy; (2) Slingshot: constructing developmental trajectories based on the center distance of each cell class after PCA dimensionality reduction; (3) scShaper: constructing trajectories between cell classes using multiple k-means algorithms and merging them to generate developmental trajectories; (4) Monocle3: constructing connection paths between cells after UMAP dimensionality reduction; (5) PAGA: calculating the connectivity probability between cell classes after PCA dimensionality reduction to indicate the possibility of the existence of a path. (6) TSCAN: connecting each cell class to generate distances based on the cell features after dimensionality reduction. This example uses the Hamming-Ipsen-Mikhailov distance (HIM), the branch F1 score (F1branches) and the trajectory node F1 score (F1milestones) as prediction accuracy indicators, where HIM measures the similarity between the predicted trajectory and the real trajectory topology; F1branches measures the accuracy of the cell mapping to the predicted path and the real path after dimensionality reduction; F1milestones measures the similarity between the cell distance to the nearest predicted trajectory node and the real trajectory node after dimensionality reduction. Table 1 is a comparison of the prediction results:

[0072] Table 1

[0073]

[0074] The experimental results are shown in the table. The results show that the method proposed in the present invention is significantly better than the currently widely used trajectory prediction methods in terms of trajectory prediction accuracy (HIM) and the results of cell distribution after dimensionality reduction (F1branches and F1milestones), which proves the effectiveness of the method of the present invention.

Claims

1. A trajectory prediction method based on contrastive learning feature dimensionality reduction, characterized in that: The following steps are involved: Step 1. According to the original matrix X of the cells to be reduced in dimension, obtain the positive sample cell set and the negative sample cell set of each cell; Step 2. Input the original matrix X of the cell to be reduced into a fully connected neural network to obtain the reduced dimension feature Y of the cell; Step 3. Input the dimension reduction feature Y into the cell-level pairwise module and the cluster-level contrast module respectively. The specific process is as follows: Step 3.

1. In the cell-level module, each round of training is for each cell y i Randomly select a cell from the positive sample set To form a positive sample pair And select a cell from the negative sample set To form a negative sample pair Step 3.

2. In the cluster-level module, each round of training is performed on each class c. i All cells of are randomly downsampled twice to form two enhanced cell classes and in Represents enhanced class The tth cell in |c i | indicates class c i The number of cells included, γ represents the downsampling rate; two enhanced cell classes generated from the same class As a positive sample pair, all cells generated from different types are classified as or Considered as negative sample pairs; Step 4. Calculate the feature loss based on the output of the cell-level pairwise module and the cluster-level contrast module, and perform optimization; The specific process of calculating feature loss and optimizing is as follows: Step 4.

1. Minimize the cell distance of positive sample pairs And maximize the distance of negative sample pairs Distance between cell sample pairs and The calculation formula is as follows: Where d represents the dimension of the low-dimensional features of the cell, and then the cross entropy loss function loss is used i Minimize the number of positive pairs of cell i and maximize the distance of negative sample pairs The distance is calculated as: The loss function L of the cell-level pairwise module is obtained by weighted average of the losses of all cells cell : Step 4.

2. By maximizing the cell class c i The similarity between the enhanced cell classes and minimize the similarity of all other negative sample enhanced classes; First, each enhanced class is calculated by weighted average The center is considered as its characteristic Then use the cosine similarity function to calculate the center r of any two enhanced classes i and j i and r j The similarity between s(r i ,r j ), the cosine similarity function is calculated as follows, where ‖·‖ represents the modulus length of the variable: Using the cross entropy loss function Calculate c i The enhanced sample a center With its positive sample All its negative sample pairs and The losses are as follows: The loss function L of the cluster-level pairwise module is obtained by weighted averaging the losses of all cell classes. cluster : Finally, the cell-level loss function L cell and the cluster-level loss function L cluster The overall loss function Loss of the model is obtained by adding them together, and its formula is: Loss=L cell +L cluster Step 5. After the optimization in step 4 is completed, the dimensionality reduction feature Y learned in step 2 is used to predict the trajectory and obtain the predicted trajectory.

2. The trajectory prediction method based on contrastive learning feature dimensionality reduction according to claim 1, characterized in that: The specific process of step 1 is: calculate the distance between each cell in the original matrix X of the cells to be reduced, and select the closest k distance for each cell based on the calculated distance p The cells are taken as the positive sample cell set, and the farthest k n cells as negative sample cells.

3. The trajectory prediction method based on contrastive learning feature dimensionality reduction according to claim 2, characterized in that: The specific process of trajectory prediction described in step 5 is as follows: First, the Euclidean distance between cells is calculated based on the dimension reduction feature Y; Next, select k for each cell. c The nearest cell is used to build a connection; Then the connections existing between cells of different classes are summed, and based on the number of connections between different classes, the minimum spanning tree is constructed between the cell classes with the condition of the largest number of connections; Finally, a cell class is defined as the initial class, and the cell class from the initial class to any leaf node in the spanning tree is regarded as a cell development trajectory to complete the trajectory prediction.

Images