A method for circulating tumor cell RNA sequencing based on a graph neural network

Abstract

The application discloses a circulating tumor cell RNA sequencing method based on a graph neural network, and specifically comprises the following steps: pretreating circulating tumor cell RNA sequencing data to obtain a non-zero expression gene group; constructing a graph neural network GNN for RNA sequencing analysis; quantifying genes encoded between different cell states by using an adaptive Gaussian mixture model; reconstructing a gene expression matrix by using iterative calculation of a feature encoder, a graph encoder and a clustering encoder; and training the reconstructed gene expression matrix group by using an interpolation encoder, so that cells with strong cell gene expression of lower quality are simulated, and interpolation of inferior genes is realized. The application starts from RNA sequencing data, constructs a graph neural network for circulating tumor cell RNA sequencing analysis, realizes interpolation of inferior genes by designing an adaptive Gaussian mixture model and iterative calculation and training between the encoders, has high sequencing and recognition efficiency, and has a wide application prospect in the field of RNA sequencing.

Description

TECHNICAL FIELD

[0001] The present application belongs to the field of machine learning, and particularly relates to a circulating tumor cell RNA sequencing method based on a graph neural network (GNN). BACKGROUND

[0002] There are no two identical cells in the world. For multicellular organisms, there are also differences between cells, and the differences between different groups of cells are not the same. Such differences not only manifest in macroscopic morphology, such as the huge differences in cell types, states and interactions between various tissues of an organism; but also manifest in genetic information, such as differences in genomic information and gene expression levels. Circulating tumor cell RNA sequencing technology can analyze the transcriptional state of a single cell or a group of cells, and obtain the gene expression status of different cells in the same tissue, thereby more finely and accurately reflecting the state of the biological tissue. Circulating tumor cell RNA sequencing mainly includes three processes: first, the separation process of circulating tumor cells. In order to obtain the transcriptome of circulating tumor cells, the key step is to separate single cells from the cell population. Cells can be isolated from a cell suspension or a tissue sample. There are currently many methods for separating cells, such as flow cytometry sorting, enzyme digestion, microfluidic or other emerging separation technologies. If the sample is already in a suspended state (such as circulating tumor cells) and the content is relatively abundant, flow cytometry sorting works well. If the sample is a solid tissue, collagen and other extracellular proteins can be digested by enzymes; second, the RNA isolation and extraction process, mainly including guanidine isothiocyanate centrifugation, hot phenol and other methods. Finally, the transcription sequencing process. After extracting RNA, reverse transcription of RNA to cDNA is performed by using reverse transcriptase, and the amplified cDNA is sequenced using next-generation sequencing technology to obtain a large amount of cell expression data. However, direct analysis of it will face many difficulties, including the sparsity of sequencing and the complex differences in gene expression, etc. These factors will affect the accuracy of prediction. Therefore, how to process the cell gene expression data to realize the interpolation of genes and the specific clustering of cells is a difficult point in the process of circulating tumor cell RNA sequencing.

[0003] As a machine learning algorithm, GNN has become a tool to develop features and patterns in graph-structured data and meet the needs of graph learning tasks such as classification, prediction, segmentation, and generation. In the biomedical field, GNN is applied to the direction of circulating tumor cell RNA sequencing, non-coding RNA regulation, etc. In September 2021, Xin et al. designed a weighted GNN, taking circulating tumor cells and genes as the nodes of the graph neural network, and realized the automatic annotation of the sequencing cell type of circulating tumor cells; in the same year, Viplove et al. used GNN to predict the RNA-protein interaction of circulating tumor cells, and the prediction accuracy was 90.6% by using the whole network information instead of predicting a single protein target. SUMMARY

[0004] The purpose of the present application is to provide a circulating tumor cell RNA sequencing method based on a graph neural network. The present application starts from cell RNA sequencing data preprocessing, constructs a graph neural network for circulating tumor cell RNA sequencing analysis, quantifies cell expressed genes by designing an adaptive Gaussian mixture model, and reconstructs the gene expression matrix through iterative calculation of the feature encoder, graph encoder, and clustering encoder, thereby achieving the interpolation of inferior genes and the clustering of cells, with high clustering interpolation efficiency and accuracy.

[0005] The technical solution for achieving the purpose of the present application is: a circulating tumor cell RNA sequencing method based on a graph neural network, comprising the following steps:

[0006] S1, preprocessing the collected circulating tumor cell RNA sequencing data, removing low-quality gene expression cells, and obtaining a screened non-zero expression gene group.

[0007] S2, constructing a graph neural network GNN for circulating tumor cell RNA sequencing analysis, wherein the GNN includes an adaptive Gaussian mixture model, a feature encoder, a graph encoder, a clustering encoder, and an interpolation encoder.

[0008] S3, quantifying the non-zero expression genes encoded between different cell states in the screened non-zero expression gene group by the adaptive Gaussian mixture model, and obtaining a quantized non-zero expression gene group.

[0009] S4, inputting the quantized screened non-zero expression gene group into the feature encoder to construct and prune the gene unit graph, then inputting the gene unit graph into the graph encoder to realize the clustering of gene expression, and finally using the clustering encoder to reconstruct the gene expression matrix, and iteratively training until convergence, to obtain a reconstructed gene expression matrix group.

[0010] S5, the reconstructed gene expression matrix group after training iteration is regularized and trained by using an interpolation encoder, an interpolated gene expression matrix group is obtained, so that cells with strong cell gene expression of lower quality are simulated, and interpolation of inferior genes is realized.

[0011] Compared with the prior art, the present application has the following advantages:

[0012] (1) The present application quantifies cell expression genes by designing an adaptive Gaussian mixture model, the retention rate of genes is higher, the input gene expression data is converted into discretized gene regulation signals, and the subsequent regularization processing of feature encoding is facilitated.

[0013] (2) The present application reconstructs the gene expression matrix through iterative calculation of the feature encoder, the graph encoder and the clustering encoder in the graph neural network, and then trains the reconstructed gene expression matrix through the interpolation encoder, so as to realize the interpolation of inferior genes and the clustering of cells, and has high clustering interpolation efficiency and accuracy. BRIEF DESCRIPTION OF DRAWINGS

[0014] Figure 1 is a flowchart of the present application based on the graph neural network circulating tumor cell RNA sequencing method.

[0015] Figure 2 is a schematic diagram of the present application based on the graph neural network circulating tumor cell RNA sequencing method.

[0016] Figure 3 is a schematic diagram of the adaptive Gaussian mixture model of the present application.

[0017] Figure 4 is a schematic diagram of the feature encoder design in the graph neural network of the present application.

[0018] Figure 5 is a schematic diagram of the graph encoder design in the graph neural network of the present application.

[0019] Figure 6 is a schematic diagram of the clustering encoder design in the graph neural network of the present application. DETAILED DESCRIPTION

[0020] The present application will be described in detail below with reference to the accompanying drawings.

[0021] In combination with Figure 1The graph neural network-based circulating tumor cell RNA sequencing method provided by the application is constructed from cell RNA sequencing data preprocessing, constructs a graph neural network for circulating tumor cell RNA sequencing analysis, quantifies cell expression genes by designing an adaptive Gaussian mixture model, reconstructs a gene expression matrix through iterative calculation of a feature encoder, a graph encoder and a clustering encoder, and trains the reconstructed gene expression matrix through an interpolation encoder to realize interpolation of inferior genes and clustering of cells, and comprises the following steps:

[0022] S1, the collected circulating tumor cell RNA sequencing data is preprocessed, the low-quality gene expression cells are removed, and the screened non-zero expression gene group is obtained, and the specific process is as follows:

[0023] Since the loss rate of RNA sequencing gene expression data is high, that is, only 1% of the cell gene expression is in a non-zero state; therefore, in order to retain the non-zero expression genes, the gene expression values are sorted by using a standard deviation, and the 1000 gene expression values with the smallest Euclidean distance are taken for subsequent construction of the quantized non-zero expression gene group, wherein the gene standard deviation is determined by the gene expression sample data, the gene expression sample mean and the sample number.

[0024] S2, a graph neural network GNN for circulating tumor cell RNA sequencing analysis is constructed, the GNN comprises an adaptive Gaussian mixture model, a feature encoder, a graph encoder, a clustering encoder and an interpolation encoder, and the network structure and specific functions are as follows:

[0025] The graph neural network of the application comprises an adaptive Gaussian mixture model, a feature encoder, a graph encoder, a clustering encoder and an interpolation encoder, and a network schematic diagram is as shown in Figure 2 . The adaptive Gaussian mixture model can convert the input gene expression data into discretized gene regulation signals, facilitating subsequent regularization processing of feature encoding; in combination with Figure 4 , the feature encoder comprises a feature encoding layer and a feature decoding layer, and the sizes are 512*128 and 128*512 respectively, and the cell gene unit graph is constructed and pruned by using the feature encoding and decoding mode; in combination with Figure 5 , the graph encoder comprises a graph encoding layer and a graph decoding layer, the graph encoding layer comprises two embedding layers, responsible for embedding the pruned cell gene unit graph, and the graph decoding layer is responsible for inner product operation on the gene embedding graph to reconstruct the adjacency matrix of the cell gene unit graph, and then the cell expression genes are clustered to obtain a plurality of cell clusters; in combination with Figure 6The clustering encoder is distributed in multiple cell cluster units, all cell cluster clustering results are cascaded by constructing a clustering encoder on each cell cluster unit, and the reconstructed gene expression matrix group is obtained through iterative training; the interpolation encoder has the same structure as the feature encoder, and can obtain a large number of cells with strong gene expression on the basis of the high-quality gene expression matrix by combining adaptive mixed Gaussian signals and regularizing training on the reconstructed gene expression matrix, so as to realize the interpolation of poor genes;

[0026] S3, quantize the non-zero expression genes encoded between different cell states in the screened non-zero expression gene group through an adaptive mixed Gaussian model to obtain a quantized non-zero expression gene group, and the specific process is as follows:

[0027] For N cells, the normalized expression value set X of the genes is represented as:

[0028] X={x1...x j ...x N} (1)

[0029] Wherein, the expression value x j ∈X of the jth gene; if X follows a Gaussian mixture distribution, then for k gene regulation signals, the density function P(X; Θ) of X is:

[0030]

[0031] Wherein, N is the number of cells, k is the number of gene regulation signals, α i is a mixed Gaussian weight, p(x j ; Θ i ) represents a weight probability, Θ i represents the gene expression range of the ith Gaussian distribution, μ i represents the mean value of the ith Gaussian distribution, σ i represents the standard deviation of the ith Gaussian distribution, L(Θ; X) is a likelihood function, the number of Gaussian components is selected by the Bayesian criterion, and each cell original gene expression distribution value can be marked down, and for any cell gene expression value x j , the probability that it belongs to the Gaussian distribution i is represented by the following formula:

[0032]

[0033] Through the adaptive Gaussian mixture model, the expression of multiple modes in a single cell gene is inferred, the non-zero expression genes encoded between different cell states in the screened non-zero expression gene group are quantized to obtain a quantized non-zero expression gene group, so as to convert the initial gene expression signal into a discrete gene regulation signal, and the captured signal has a better signal-to-noise ratio. Figure 3For adaptive Gaussian mixture diagram, ① is the Gaussian mixture distribution of the whole data, and ② and ③ are the single Gaussian model distribution constituting ①. For different non-zero gene expression Gaussian model distribution of cells, the adaptive Gaussian mixture model can quantify the non-zero expression genes, facilitate the subsequent regularization training process, and relative to the ordinary high-pass or low-pass Gaussian mixture model, it is more conducive to training to simulate high-quality cells with strong gene expression.

[0034] S4, input the quantized and screened non-zero expression gene group into the feature encoder to construct and prune the gene unit graph, then input the gene unit graph into the graph encoder to realize the clustering of gene expression, and finally use the clustering encoder to reconstruct the gene expression matrix. The specific process is as follows:

[0035] S41, input the quantized non-zero expression gene group into the feature encoder to construct and prune the gene unit graph. As known from S3, the normalized expression value set X of the gene, the gene expression can be reconstructed by feature encoding and decoding. The reconstructed expression is denoted as X1, and X and X1 satisfy:

[0036]

[0037] Wherein, M is the number of input gene expression matrix, M1 is the dimension of training, and Z is the training coefficient. During the training of the feature encoder, a judgment standard is needed to realize the integrity of the constructed cell gene unit graph, that is, to minimize the loss function, which is represented as:

[0038] ∑(X-X1) 2 (5)

[0039] It is the main item used as the mean square error between X and X1, and the mean square error is represented as:

[0040] α∑(X-X1) 2 W (6)

[0041] Wherein, α is a regularization strength parameter, and α∈[0, 1], W represents the loss weight. According to the above analysis, the loss function Loss of the feature encoder is: W

[0042] Loss W =(1-α)∑(X-X1) 2 +α∑(X-X1) 2 W (7)

[0043] ​S42, input the cell gene unit graph into a graph encoder to realize clustering of gene expression, then reconstruct the gene expression matrix by using a clustering encoder, iteratively train until convergence, after pruning the cell gene unit graph, the graph encoder is used for embedding the pruned cell gene unit graph to obtain a gene embedding graph, then inner product operation is performed on the gene embedding graph to reconstruct the adjacency matrix of the cell gene unit graph, the gene embedding graph is denoted as Q, and the reconstructed adjacency matrix of the cell gene unit graph is denoted as P, then a Sigmoid activation function is used to obtain:

[0044]

[0045] After reconstructing the adjacency matrix of the cell gene unit graph, the k-means algorithm is used for clustering the cell genes, then the gene expression matrix is reconstructed by using the clustering encoder, and iterative training is performed until convergence, and during the iterative training, the adjacency matrix of the cell gene unit graph is denoted as P0, and the iterative process is represented as:

[0046]

[0047] wherein V0 is the normalized adjacency matrix of the initial gene unit graph, and λ is a convergence speed parameter; the iteration is stopped when convergence is achieved, and a reconstructed gene expression matrix group is obtained.

[0048] S5, the reconstructed gene expression matrix group after training iteration is subjected to regularized training by using an interpolation encoder, and an interpolated gene expression matrix group is obtained, so that cells with strong gene expression are simulated, and the interpolation of inferior genes is realized. The interpolation encoder has the same structure as the feature encoder, and can combine an adaptive mixed Gaussian signal to perform regularized training on the reconstructed gene expression matrix group; the interpolated gene expression matrix group includes the screened non-zero expression genes and the high-quality cell expression genes simulated by training, and finally the interpolation of inferior genes is realized.

[0049] In summary, the circulating tumor cell RNA sequencing method based on the graph neural network disclosed in the present application starts from cell RNA sequencing data preprocessing, constructs a graph neural network for circulating tumor cell RNA sequencing analysis, quantifies cell expression genes by designing an adaptive mixed Gaussian model, reconstructs a gene expression matrix through iterative calculation of a feature encoder, a graph encoder and a clustering encoder, and then performs regularized training on the reconstructed gene expression matrix by using an interpolation encoder, so that the interpolation of inferior genes and the clustering of cells are realized.

Claims

1. A method for sequencing circulating tumor cell RNA based on graph neural networks, characterized in that, Includes the following steps: S1. The collected circulating tumor cell RNA sequencing data were preprocessed to obtain the screened non-zero expressed genome, as follows: Because RNA sequencing gene expression data has a loss rate, only 1% of cells have non-zero gene expression. In order to retain non-zero expression genes, the gene expression values ​​are sorted by standard deviation, and the 1000 gene expression values ​​with the smallest Euclidean distance are selected for subsequent construction of a quantified non-zero expression genome. S2. Construct a GNN for RNA sequencing analysis of circulating tumor cells. The GNN includes an adaptive Gaussian mixture model, a feature encoder, a graph encoder, a cluster encoder, and an interpolation encoder. The adaptive Gaussian mixture model converts the input non-zero expressed genome into discrete regulatory signals, which facilitates the regularization of subsequent feature encoding. The feature encoder consists of a feature encoding layer and a feature decoding layer. It uses feature encoding and decoding to construct and prune the cell gene unit graph. The graph encoder consists of a graph encoding layer and a graph decoding layer. The graph encoding layer includes two embedding layers responsible for embedding the pruned gene unit graph to obtain a gene embedding graph. The graph decoding layer is responsible for performing inner product operations on the gene embedding graph to reconstruct the adjacency matrix of the cell gene unit graph. Then, the k-means algorithm is used to cluster the cell-expressed genes to obtain multiple cell clusters. The clustering encoder is distributed in multiple cell cluster units. By constructing a clustering encoder on each cell cluster unit, the clustering results of all cell clusters are concatenated and iteratively trained to finally obtain the reconstructed gene expression matrix. The interpolation encoder has the same structure as the feature encoder. By combining an adaptive Gaussian mixture signal, the reconstructed gene expression matrix is ​​regularized and trained to achieve the interpolation of inferior genes. S3. The non-zero expression genes encoded by different cell states in the non-zero expression genome after screening are quantified by an adaptive Gaussian mixture model to obtain the quantified non-zero expression genome; S4. Input the non-zero expression genome after quantitative screening into the feature encoder to construct and prune the gene unit graph, then input the gene unit graph into the graph encoder to realize gene expression clustering, and finally use the clustering encoder to reconstruct the gene expression matrix. Iterate and train until convergence to obtain the reconstructed gene expression matrix set. S5. The reconstructed gene expression matrix after training iteration is regularized using the interpolation encoder to obtain the interpolated gene expression matrix, thereby realizing the interpolation of inferior genes.

2. The method for sequencing circulating tumor cell RNA based on graph neural networks according to claim 1, characterized in that, In S3, non-zero expression genes encoded by different cell states in the non-zero expression genome after screening are quantified using an adaptive Gaussian mixture model, resulting in a quantified non-zero expression genome, as follows: For N cells, the set of normalized gene expression values ​​X is represented as: X={x1…x j …x N } (1) Wherein, the expression value x of the j-th gene j For any x ∈ X, if X follows a Gaussian mixture distribution, then for k gene regulatory signals, the density function P(X; Θ) of X is: Where N is the number of cells, k is the number of gene regulatory signals, and α i For Gaussian mixture weights, p(x) j ;Θ i ) represents the weighted probability, Θ i μ represents the gene expression range of the i-th Gaussian distribution. i Let σ represent the mean of the i-th Gaussian distribution. i Let x represent the standard deviation of the i-th Gaussian distribution, and L(Θ; X) be the likelihood function; for any cell gene expression value x j The probability that it belongs to the Gaussian distribution i is expressed by the following formula: By inferring multiple modes of gene expression in a single cell using an adaptive Gaussian mixture model, the non-zero expressed genes encoded by different cell states in the screened non-zero expressed genome are quantified, and the initial gene expression signal is converted into a discrete gene regulatory signal, thereby obtaining a quantified non-zero expressed genome.

3. The method for sequencing circulating tumor cell RNA based on graph neural networks according to claim 1, characterized in that, In S4, the quantized non-zero expressed genome is input into the feature encoder of the GNN to construct and prune the gene unit graph. Then, the gene unit graph is input into the graph encoder of the GNN to achieve gene expression clustering. Finally, the clustering encoder in the GNN is used to reconstruct the gene expression matrix. Iterative training continues until convergence, resulting in a reconstructed gene expression matrix set, as detailed below: S41. Input the quantified non-zero expression genome into the feature encoder to construct and prune the gene unit map. The set of gene normalized expression values ​​X obtained from S3 is used to reconstruct the gene expression through feature encoding and decoding. The reconstructed expression is denoted as X1. At this time, X and X1 satisfy: Where M is the number of input gene expression matrices, M1 is the training dimension, and Z is the training coefficient. During feature encoder training, an evaluation criterion is needed to ensure the completeness of the constructed cellular gene unit map, i.e., minimizing the loss function, expressed as: ∑(X-X1) 2 (5) It is the main term used as the mean square error between X and X1, and the mean square error is expressed as: αΣ(X-X1) 2 W(6) where α is the regularization strength parameter, α∈[0,1], and W represents the loss weight. From the above analysis, the loss function of the feature encoder is Loss W Represented as: Loss W (1-α)∑(X-X1) 2 +α∑(X-X1) 2 W (7) S42. Input the cell gene unit graph into the graph encoder to achieve gene expression clustering, and then use the cluster encoder to reconstruct the gene expression matrix. Iterate until convergence. After pruning the cell gene unit graph, the graph encoder embeds the pruned cell gene unit graph to obtain the gene embedding graph. Then, perform an inner product operation on the gene embedding graph to reconstruct the adjacency matrix of the cell gene unit graph. Let the gene embedding graph be Q, and the reconstructed adjacency matrix of the cell gene unit graph be P. Then, using the Sigmoid activation function, we get: After reconstructing the adjacency matrix of the cell gene unit graph, the k-means algorithm is used to cluster the cell genes, and then the gene expression matrix is ​​reconstructed using a cluster encoder. This process is iteratively trained until convergence. During iterative training, the adjacency matrix of the cell gene unit graph is denoted as P0. The iterative process is then represented as: Where V0 is the normalized adjacency matrix of the initial gene unit graph, and λ is the convergence rate parameter; the iteration stops when convergence is achieved, and the reconstructed gene expression matrix is ​​obtained.

4. The method for sequencing circulating tumor cell RNA based on graph neural networks according to claim 1, characterized in that: In step S5, the interpolation encoder has the same structure as the feature encoder. It combines an adaptive Gaussian mixture signal to perform regularization training on the reconstructed gene expression matrix group, thereby achieving the interpolation of inferior genes.

Images