Single cell level dynamic gene regulatory network inference method and system based on graph neural network and comparative learning

Abstract

The invention provides a single cell level dynamic gene regulation network inference method based on a graph neural network and comparative learning, and aims to solve the problem that dynamic change of a gene regulation relation under single cell resolution is difficult to reveal due to the fact that cell heterogeneity cannot be described by an existing monomics method. The method comprises the following steps: firstly, pre-training to generate a cell token and a gene token through comparative learning and graph representation learning; then, carrying out multi-modal fusion on the two genes and the interpolated gene expression data, and constructing context specific representation of each gene in each cell; further, a graph neural network with a cross attention mechanism is adopted to model cell specific regulation interaction between a transcription factor and a target gene, and dynamic gene regulation network inference of a single cell level is realized; and finally, carrying out joint optimization on model parameters through comparative learning and expression reconstruction. According to the method, the dynamic regulation and control relation of cell state dependence can be accurately deduced, and an important calculation tool is provided for development analysis and research.

Description

Technical Field

[0001] This invention belongs to the fields of bioinformatics and computational biology, and relates to a method and system for inferring dynamic gene regulation networks in single cells based on graph neural networks and contrastive learning. Background Technology

[0002] Gene regulatory networks are central to understanding cellular life activities and functional mechanisms. Traditional single-omics gene regulatory network inference methods are mostly based on population cell data, which can only obtain averaged regulatory relationships and cannot reveal heterogeneity at the single-cell level. Single-cell RNA sequencing (scRNA-seq) technology is an important tool in current genomics research. It can analyze the transcriptome of a single cell through high-throughput sequencing and reveal differences in gene expression in different developmental, differentiation, and pathological states. However, accurately inferring dynamic gene regulatory networks at single-cell resolution still faces significant challenges: cellular heterogeneity makes it difficult to capture dynamic changes in regulatory relationships; single-cell data is characterized by high dimensionality, sparsity, and high noise, which seriously affects the accuracy of GRN inference; and there is a lack of computational models that can simultaneously capture cellular heterogeneity and dynamic regulatory relationships. Summary of the Invention

[0003] The purpose of this invention is to overcome the shortcomings of the prior art and provide a method and system that can accurately infer the dynamic gene regulatory network of a single cell.

[0004] To achieve the above objectives, the present invention adopts the following technical solution:

[0005] A method for inferring dynamic gene regulation networks at the single-cell level based on graph neural networks and contrastive learning includes:

[0006] Data preprocessing: The raw scRNA-seq data were preprocessed to obtain the gene expression matrix;

[0007] Pre-trained cell tokens: Based on the gene expression matrix, the neighborhood relationships of cells are determined to construct positive samples; cell embedding representation learning is performed on cells to obtain cell tokens;

[0008] Pre-trained gene tokens: Gene tokens are extracted by pre-training an inductive graph neural network model guided by a priori gene regulatory network.

[0009] Gene expression matrix imputation: Based on the preprocessed gene expression matrix, the gene expression of cells is weighted and fused with the average gene expression of its positive sample neighborhood to obtain the imputed gene expression matrix;

[0010] Feature fusion: Multimodal fusion of cell tokens, gene tokens, and imputed gene expression data to generate a specific representation vector for each gene in each cell;

[0011] Network inference and cell embedding generation: The specific representation vector of cell genes is input into a graph neural network model. Based on the multi-head cross attention mechanism, the regulatory relationship between transcription factors and target genes in the cell is calculated to output the gene regulatory network of the cell. At the same time, the transcription factor representation is weighted and integrated based on the regulatory relationship, and then processed with the initial transformation features of the target gene through network layers to generate the cell embedding vector.

[0012] Expression reconstruction: The cell embedding vector is passed through a feedforward network decoder to reconstruct the cell's gene expression;

[0013] Model training: The parameters of the graph neural network model are trained by jointly optimizing the comparative learning loss function and the weighted reconstruction loss function;

[0014] Regulatory network inference: Obtain the specific representation vectors of all genes in the target cell, input them into the trained graph neural network, and infer the gene regulatory network of the target cell.

[0015] Furthermore, the preprocessing of the raw scRNA-seq data includes: removing cells with fewer than a preset number of genes, removing genes expressed in fewer than a preset number of cells, standardizing the gene expression value of each cell to the same level, and screening for highly variable genes.

[0016] And / or, the K-nearest neighbor algorithm is used to determine the neighborhood relationships of cells;

[0017] And / or, use a momentum contrastive learning framework to learn the embedding representation of cells.

[0018] Furthermore, the inductive graph neural network model includes two SAGEConv convolutional layers, and a batch normalization layer, a PReLU activation function, and a Dropout layer are provided between the SAGEConv convolutional layers; the output of the last SAGEConv convolutional layer is used as the gene token.

[0019] The inductive graph neural network model guided by a priori gene regulation network includes: defining positive sample edges based on the edge index of the priori gene regulation network; generating negative sample edges through a structured negative sampling method; calculating the dot product of the gene tokens at both ends of the positive and negative sample edges and using the sigmoid function as the score of the corresponding edge; training the inductive graph neural network model by maximizing the log-likelihood of the positive sample edges and minimizing the log-likelihood of the negative sample edges; the structured negative sampling refers to: for each positive sample edge... By fixing the head node gene And randomly sample a tail node gene To construct negative sample edges ,and , Represents the relationship between genes in the priori gene regulatory network. Genes that are connected by edges.

[0020] Furthermore, the feature fusion generates the specificity representation vector using the following formula:

[0021]

[0022] in, Representing the The first cell Cell-gene-specific expression vectors of individual genes. Represents hidden dimensions; Representing the A cell token for each cell; Representing the Gene token for each gene; Representing the The first cell Gene expression data after imputation of gene expression values.

[0023] Furthermore, the regulatory relationships between transcription factors and target genes in cells are calculated based on the multi-head cross-attention mechanism, specifically including:

[0024] The specific representation vectors of transcription factors and target genes are subjected to nonlinear transformations using different multilayer perceptrons to obtain the query matrix. Key matrix Sum matrix ;

[0025]

[0026] in, and This represents two different multilayer perceptrons. and Generate query matrix respectively The trainable weight matrix and bias vector used; and Generate the key matrix respectively The trainable weight matrix and bias vector used; and These are the generated value matrices. The trainable weight matrix and bias vector used; represent The fusion feature vector of a transcription factor in n cells. represent Fusion feature vectors of target genes in n cells The dimension of the fused feature vector; , It is the number of multiple heads; Represents transcription factor characterization, Representative target gene characterization;

[0027] Constructing a mask matrix using prior gene regulatory networks This is used to filter out non-existent regulatory relationships in attention calculations; based on the query matrix... Key matrix and mask matrix The attention score between transcription factors and target genes is calculated through a multi-head cross-attention mechanism;

[0028]

[0029] in, , Representing the first The cell number Transcription factors in the size and target genes The representation; Representing the The cell number Transcription factors in the size With target genes The attention score is the relevance score;

[0030] The calculated attention score is first processed After normalization, the adjacency matrix representing the strength of the regulatory relationship is obtained by averaging across multiple dimensions. That is, the first Gene regulatory networks in individual cells:

[0031]

[0032]

[0033] in, Representing the Transcription factors of individual cells With target genes Regulatory relationships within gene regulatory networks.

[0034] Furthermore, the specific steps for generating the cell embedding vector include:

[0035] The attention scores of the cell's gene regulatory network are used to perform a weighted summation of transcription factor representations to obtain context-aware representations of target genes. The specific calculation formula is as follows:

[0036]

[0037] in It is an attention score; It is a transcription factor characteristic; It is the context-aware representation of the target gene obtained;

[0038] The context-aware representation is nonlinearly transformed using a first multilayer perceptron to obtain the transformed target gene representation.

[0039] The transformed target gene characterization is residually linked with the initial transformed features of the target gene to obtain an enhanced target gene characterization.

[0040] The enhanced target gene characterization is subjected to layer normalization, and then feature extraction is performed using a second multilayer perceptron.

[0041] The outputs of the second multilayer perceptron for all target genes are summed in the gene dimension to generate a cell embedding vector representing the overall cell state.

[0042] Furthermore, the model training steps specifically include:

[0043] Construct a contrastive learning loss function that brings positive samples closer and pushes negative samples further away based on cell embedding vectors;

[0044] A weighted reconstruction loss function is constructed, which assigns higher weights to the reconstruction errors of non-zero values ​​in the gene expression matrix;

[0045] The contrastive learning loss and weighted reconstruction loss are jointly optimized to train the parameters of the graph neural network and decoder, while keeping the pre-trained cell token and gene token parameters fixed.

[0046] Furthermore, the weighted reconstruction loss function is:

[0047]

[0048] in, Represents the weighted reconstruction loss; This represents the calculated mean square error; This represents the portion of the gene expression matrix with values ​​of 0 after preprocessing. This represents the portion of the gene expression matrix after preprocessing that has non-zero values; Representative at The graph neural network model predicts the value at the corresponding position. Representative at The graph neural network model prediction value at the corresponding position.

[0049] Furthermore, the contrastive learning loss function is:

[0050]

[0051] in, Represents contrastive learning loss; For a single target cell, the cell embedding vector. To and The embedding vector of a single positive sample cell that constitutes a positive sample pair; To compare the set of negative sample cell embedding vectors in the dynamic learning queue, To compare the number of negative samples in the dynamic learning queue; The temperature coefficient is a hyperparameter that is greater than zero.

[0052] A single-cell-level dynamic gene regulation network inference system based on graph neural networks and contrastive learning includes a memory and a processor. The memory stores a computer program, which, when executed by the processor, causes the processor to implement the method described above.

[0053] Beneficial effects:

[0054] Through systematic experimental verification, the technical solution provided by this invention has the following advantages:

[0055] 1. High-precision static gene regulatory network inference: This invention constructs a highly integrated inference system strongly constrained by prior biological knowledge. By combining a multimodal fusion strategy to integrate inherent gene attributes, cell state, and expression dynamics, it provides a high-information-density feature foundation for inference. Furthermore, it employs a cross-attention mechanism based on prior network masks to directly simulate TF-TG regulatory interactions, making the inference process both data-driven and biologically plausible. Experiments on multiple benchmark datasets demonstrate that this method significantly outperforms existing methods in the AUPRC and AUROC metrics.

[0056] 2. Transcription Factor-Specific Identification: The fusion feature vector of this invention carries cell state information, making the attention weights dynamically dependent on the specific cellular environment. This allows the same Transcription Factor (TF) to exhibit differentiated regulatory strengths in different cells, naturally leading to cell type-specific regulatory networks. Therefore, this invention can accurately identify cell type-specific regulatory factors and their regulatory networks, providing a new perspective for understanding cellular functional heterogeneity.

[0057] 3. Dynamic analysis: The cell embedding space constructed based on contrastive learning continuously represents the developmental trajectory, which can track the dynamic evolution of gene regulatory networks during development and capture the regulatory reprogramming events of key transcription factors; Attached Figure Description

[0058] Figure 1Schematic diagram of the overall architecture of this invention;

[0059] Figure 2 The static gene regulatory network inference performance comparison chart shows the comparison results of the method of the present invention with the AUPRC and AUROC indices on multiple benchmark simulations and real datasets. Figure 2 (a1)-(d2) respectively show the comparison of AUPRC and AUROC indices on the hematopoietic stem cell dataset, gonadal sex-determining cell dataset, myeloid dendritic cell dataset, and human hepatocyte dataset.

[0060] Figure 3 This demonstrates the ability of the method of the present invention to reconstruct the dynamics of branching regulation during the development of multiple lineages of the mouse pancreas.

[0061] Figure 3 (a) Shows cell distribution based on raw scRNA-seq data. Figure 3 (b) shows the cell distribution based on scDGRN learning, both of which show discrete cell clusters corresponding to the major pancreatic lineages.

[0062] Figure 3 (c) is a hierarchical structure tree diagram of pancreatic cell types based on the similarity of transcriptional regulatory factors, revealing the regulatory relationships between cell types and their consistency with known developmental trajectories.

[0063] Figure 3 (d) is a heatmap of the regulatory activities of key transcription factors in different cell types. The left side shows the results inferred by scDGRN, and the right side shows the results of the comparative method, which verifies the advantages of scDGRN in accurately identifying cell type-specific regulatory factors.

[0064] Figure 3 (e) is a heatmap showing the dynamic changes of regulatory edges during cell differentiation, illustrating the dynamic changes of the regulatory intensity of key transcription factors such as Neurog3, Fev, and Nkx6-1 and their target genes over time.

[0065] Figure 3 (f) is a bar chart showing the enrichment of target genes in the Reactome pathway corresponding to the regulatory edge that significantly increases with time, indicating significant enrichment in pathways related to β-cell development and function. Detailed Implementation

[0066] The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0067] based on Figure 1 The overall architecture shown in this embodiment is based on a single-cell level dynamic gene regulation network inference method using graph neural networks and contrastive learning, and is implemented according to the following steps:

[0068] 1. Data preprocessing: Scanpy was used to remove low-quality cells with fewer than 100 expressed genes and to remove low-abundance genes expressed in fewer than 3 cells. Then, the gene expression values ​​of each cell were standardized to the same level and the variance was stabilized by logarithmic transformation. Finally, highly variable genes were screened to obtain the preprocessed gene expression matrix.

[0069] 2. Pre-trained cell tokens: using the preprocessed gene expression matrix as input. , Represents cell number, This represents the number of genes. The specific steps are as follows: First, the KNN algorithm is used to mine cell neighborhood relationships and determine positive samples for contrastive learning; then, Momentum Contrast (MoCo) is used to train cell tokens. , Represents cell number, Representing the hidden dimension, this embodiment is set as follows: The obtained cell tokens form the basis for subsequent extraction of cell state characterization from dynamic gene regulatory networks.

[0070] The core of this invention, which uses momentum contrastive learning to train cell tokens, lies in the following: An online encoder encodes a batch of cellular gene expression data and compares it with a queue of positive and negative samples generated by a momentum encoder maintained through a momentum update mechanism. Then, by optimizing the contrastive learning loss function, the online encoder is trained to map similar cells to nearby positions in the embedding space and dissimilar cells to distant positions. After model training is complete, the online encoder is fixed, and its output embedding vector serves as the cell token used downstream.

[0071] 3. Gene Token Pre-training: An inductive graph neural network (SNN) serves as the core model for gene embedding learning. It comprises two SAGEConv convolutional layers, with batch normalization, PReLU activation, and Dropout layers between them. The output of the last SAGEConv convolutional layer is used as the gene token. Preprocessed gene expression data and a comprehensive prior gene regulation network proposed by NicheNet are used as input. The two SAGEConv convolutional layers extract local topological features of the genes. A BatchNorm layer is introduced to alleviate gradient vanishing, the PReLU activation function enhances non-linear expression, and the Dropout layer suppresses overfitting. The training phase focuses on "link prediction": based on the edge index of the prior network, negative sample edges are generated through structured_negative_sampling, and the model is trained using positive sample log-likelihood loss and negative sample log-likelihood loss. Finally, the result after the two SAGEConv convolutional layers is used as the gene token. , Represents the number of genes. Represents hidden dimensions.

[0072] 4. Imputation of gene expression matrix: Based on the preprocessed gene expression matrix, the K-nearest neighbor algorithm is used to identify cell neighborhoods as positive samples. The original expression matrix and the average expression matrix of the positive sample neighborhood are imputed by weighted linear fusion, with each weighting 0.5 to achieve balanced integration, and finally the imputed gene expression data is obtained.

[0073] 5. Feature Fusion: The cell tokens, gene tokens, and imputed gene expression data are fused using a multimodal approach to generate a specific representation vector for each gene in each cell. Specifically, to effectively integrate inherent gene attributes, cell state, and expression dynamics, spatial consistency is first ensured through dimensional alignment: gene tokens are broadcast along the cell dimension to... Cell tokens are broadcast along the gene dimension to The gene expression matrix is ​​expanded to [a certain size] through dimensional repetition. Subsequently, multimodal fusion was achieved through element-wise addition, ultimately generating a dimension of... Cell-gene specific representation vectors.

[0074]

[0075] in, Representing the A token for each cell. Representing the A token for each gene. Representing the The gene in the first Expression values ​​in individual cells, Representing the The first cell in the cell A gene-specific representation vector (fusion feature vector).

[0076] 6. Network Inference (Inference of Single-Cell Level Dynamic Gene Regulation Network Based on Graph Neural Network and Momentum Contrast Learning): Accurate inference of single-cell level dynamic gene regulation network is achieved through multi-head cross-attention graph neural network.

[0077] First, the encoder receives the single-cell-specific representation vectors of the fused transcription factor (TF) and target gene (TG). , .

[0078]

[0079] in, and This represents two different neural networks. represent TFs fuse feature vectors in n cells represent Each TG is fused into a feature vector in n cells, and each vector undergoes a layer of MLP nonlinear transformation.

[0080] Then, after another linear transformation, the calculation is performed. , and .

[0081]

[0082] in, , It is the number of multiple heads; Represents TF characterization, This represents the TG characterization.

[0083] We then modeled the specific regulatory interactions between TF and TG within a single cell using a cross-attention mechanism.

[0084]

[0085] in Representing the The cell number Transcription factors in the size With target genes The attention score, or relevance, is the score that corresponds to the relevance score. This represents the prior network (where 1 represents the existence of a regulatory relationship, and -1e9 represents the absence of a regulatory relationship).

[0086]

[0087]

[0088] Here, softmax normalization along the transcription factor dimension can effectively highlight its regulatory specificity, and then averaging across multiple factors yields... , representing the A gene regulatory network in each cell.

[0089] 7. Generate cell embedding vectors: using transcription factors for characterization After reconstructing the representation of each target gene in each cell from the gene regulatory network, a layer-norm residual connection is used to accelerate convergence. Then, a multilayer perceptron is used to generate the cell's feature embedding, denoted as the cell embedding vector. Specifically, this includes:

[0090] The attention weight matrix of the gene regulation network is used to perform a weighted summation of the value matrix V to obtain the context-aware representation of the target gene.

[0091] The context-aware representation is nonlinearly transformed using a first multilayer perceptron to obtain the transformed target gene representation.

[0092] The transformed target gene characterization is residually linked with the initial transformed features of the target gene to obtain an enhanced target gene characterization.

[0093] The enhanced target gene characterization is subjected to layer normalization, and then feature extraction is performed using a second multilayer perceptron.

[0094] The outputs of the second multilayer perceptron for all target genes are summed in the gene dimension to generate a cell embedding vector representing the overall cell state.

[0095]

[0096] Where norm represents layer_norm. Representative cells The embedding representation vector.

[0097] 8. Expression remodeling: Cellular embedding to The scRNA-seq data were reconstructed using a feedforward network (a simple non-linear mapping).

[0098]

[0099] in This represents reconstructed scRNA-seq data.

[0100] 9. Model Training: The parameters of the graph neural network model are trained by jointly optimizing the comparative learning loss function and the weighted reconstruction loss function.

[0101] (1) Construct a contrastive learning loss function to bring positive samples closer and push negative samples further away based on cell embedding vectors. Specifically, a negative sample queue is maintained by a momentum encoder to maximize the similarity of cell embeddings of positive samples (K nearest neighbors are used as positive samples on the preprocessed gene expression matrix) while minimizing the similarity of cell embeddings of negative samples, so as to enhance the discriminativeness of cell state representation.

[0102]

[0103] in, Represents contrastive learning loss; For a single target cell, the cell embedding vector. To and The embedding vector of a single positive sample cell that constitutes a positive sample pair; To compare the set of negative sample cell embedding vectors in the dynamic learning queue, To compare the number of negative samples in the dynamic learning queue; The temperature coefficient is a hyperparameter greater than zero used to adjust the sharpness of the similarity distribution, thereby controlling the model's ability to distinguish difficult negative samples.

[0104] The momentum encoder updates parameters through gradient backpropagation; the momentum encoder has the same neural network structure as the encoder, but instead of updating parameters directly through gradient backpropagation, it updates them through the exponential moving average (EMA) of the encoder parameters. (2) Construct a weighted reconstruction loss function to give higher weight to the reconstruction error of non-zero values ​​in the gene expression matrix. Specifically, the decoder reconstructs the cell embedding vector into a gene expression matrix in reverse and uses a weighted reconstruction loss function to give higher weight to the reconstruction error of non-zero values ​​in the gene expression matrix, ensuring that GRN inference not only captures the regulatory topology but also remains consistent with the statistical characteristics of the original transcriptome data.

[0105]

[0106] in This represents the portion of the gene expression matrix with values ​​of 0 after preprocessing. This represents the portion of the gene expression matrix after preprocessing that has non-zero values; This represents the portion of the pre-reconstructed gene expression matrix where the value is 0. This represents the portion of the reconstructed gene expression matrix where the values ​​are not zero. The zero portion represents the loss. As a sparse regularization term, it is used to penalize the model for making blind predictions on genes that shouldn't be expressed, preventing the model from overfitting to noise and guiding it to learn the inherent, large number of zero-value patterns in gene expression data. This is crucial for processing data with significant "missing" phenomena, such as single-cell RNA sequencing. (Non-zero partial loss) The core task of the model is to directly measure the accuracy of the model's reconstruction of existing gene expression levels and minimize this loss to force the model to learn how to accurately recover known biological signals.

[0107] (3) Jointly optimize the contrastive learning loss and the weighted reconstruction loss to train the parameters of the graph neural network and the decoder, while keeping the pre-trained cell token and gene token parameters fixed.

[0108] 10. Inferring the Gene Regulatory Network: Obtain the specific representation vectors of all genes in the target cell, input them into the trained graph neural network, and infer the gene regulatory network (GRN) of the target cell.

[0109] To verify the effectiveness of the proposed scDGRN method, a systematic comparison was conducted between the scDGRN method and existing methods on simulated datasets (HSC and GSD) and real transcriptome datasets (mDC and hHep), using the AUROC and AUPRC metrics. Specifically, let TP represent the number of true positives, TN represent true negatives, FP represent false positives, and FN represent false negatives. The formulas for calculating the True Positive Rate (TPR) and False Positive Rate (FPR) are as follows:

[0110]

[0111]

[0112] Precision refers to the proportion of samples that are actually positive out of all samples predicted as positive by the model; it measures the model's accuracy in predicting the positive class. Recall refers to the proportion of samples that are actually positive out of all samples predicted as positive by the model; it measures the model's ability to identify the positive class. Therefore, the formulas for calculating precision and recall are:

[0113]

[0114]

[0115] AUROC is the area under the ROC curve, which plots the model's true positive rate and false positive rate on the same graph to measure the model's ability to distinguish between positive and negative classes. AUPR is the area under the PR curve, which plots the model's precision and recall on the same graph to measure the model's ability to identify the positive class.

[0116] During implementation, the GRNs inferred from all cells are first averaged, and then the edge weights of each target gene are normalized by min-max scaling. Figure 2 The results show that the gene regulatory network GRN inferred by the scDGRN method of this invention achieves the highest average AUPRC and AUROC on simulated data, and exhibits the smallest performance fluctuation with increasing dropout rates. In real-world data evaluation, scDGRN also outperforms all comparative algorithms, confirming its ability to reconstruct accurate GRNs from complex transcriptome data. Figure 2 (a) Presents the performance evaluation results of various gene regulatory network prediction algorithms on the hematopoietic stem cell (HSC) dataset, with performance metrics including AUPRC and AUROC. Here, 0%, 50%, and 70% represent the random masking ratios applied to the dataset, used to simulate different levels of technical noise in single-cell RNA-seq data. Ten independent datasets were generated for each masking ratio to evaluate robustness. Box plots show that the scDGRN method of this invention has significantly higher box positions than the comparative algorithms at all noise levels, indicating that it has the best and most stable gene regulatory network prediction performance. Figure 2 (b1)(b2) are the evaluation of the AUPRC and AUROC metrics on the gonadal sex determination (GSD) dataset. Figure 2 (c1)(c2) are the evaluation of AUPRC and AUROC metrics on the Myeloid Dendritic Cells (mDC) dataset. Figure 2 (d1)(d2) represent the AUPRC and AUROC metrics evaluated on the human hepatocyte (hHep) dataset.

[0117] based on Figure 3 The accompanying figures illustrate the application of this invention in analyzing complex developmental processes:

[0118] 1. Cell embedding visualization analysis (corresponding to) Figure 3(a)(b): UMAP dimensionality reduction was performed based on the original gene expression matrix and the cell embeddings generated by scDGRN, respectively. Both revealed a clear pancreatic lineage structure, including major populations such as α cells, β cells, δ cells, ε cells, and endocrine progenitor cells. The continuous transitions in the cell embedding space reflect developmental continuity, indicating that the model captures not only static cell types but also the dynamic processes of cell state transitions.

[0119] 2. Construction of cell type-specific regulatory networks (corresponding) Figure 3 (c): A hierarchical clustering dendrogram was constructed based on the similarity of transcriptional regulatory factors. The specific implementation included: aggregating single-cell GRNs for each cell type to obtain a consensus network; retaining the top 10% of the regulatory transcription factors for each target gene; selecting the top 20% of TFs based on total regulatory strength; and performing hierarchical clustering based on Jaccard similarity. The results showed that Ngn3-low-expressing endocrine progenitor cells and Ngn3-high-expressing cells had the highest similarity, forming the most closely related branch, consistent with their shared progenitor identity.

[0120] 3. Validation of key transcription factor activity (corresponding to) Figure 3 (d): Systematic comparison of transcription factor activity predicted by the model with known cell type-specific regulators reported in the literature. Key functional transference factors (TFs) including Sox9 (a hypo-endocrine precursor of Ngn3), Neurog3 (a hyper-endocrine precursor of Ngn3), Fev (pre-endocrine cells), Pax6 (α cells), Nkx6-1 and Pdx1 (β cells), Arx (ε cells), and Hhex (δ cells) were selected. The importance score of each TF for each cell type was calculated from the top 50% of regulatory edges. The results showed that the importance scores assigned by scDGRN were highly consistent with the expected cell type specificity of these TFs.

[0121] 4. Developmental dynamic regulation analysis (corresponding) Figure 3 (e), 3(f)): Changes in the intensity of TF-target gene regulatory edges predicted along the pseudo-time axis of endocrine differentiation. Transcription factors enriched in cell types such as Neurog3, Fev, and Nkx6-1 were selected, and the Spearman correlation between the intensity of each regulatory edge and pseudo-time was calculated. Reactome pathway enrichment analysis was performed on target genes corresponding to regulatory edges that significantly increased with pseudo-time. The results showed significant enrichment in pathways such as "β-cell development regulation" and "β-cell gene expression regulation," reflecting key processes accompanying endocrine maturation.

Claims

1. A method for inferring dynamic gene regulation networks at the single-cell level based on graph neural networks and contrastive learning, characterized in that, include: Data preprocessing: The raw scRNA-seq data were preprocessed to obtain the gene expression matrix; Pre-trained cell tokens: Based on the gene expression matrix, the neighborhood relationships of cells are determined to construct positive samples; Cell embedding representation learning is performed to obtain cell tokens; Pre-trained gene tokens: Gene tokens are extracted by pre-training an inductive graph neural network guided by a priori gene regulatory network. Gene expression matrix imputation: Based on the preprocessed gene expression matrix, the gene expression of cells is weighted and fused with the average gene expression of its positive sample neighborhood to obtain the imputed gene expression matrix; Feature fusion: Multimodal fusion of cell tokens, gene tokens, and imputed gene expression data to generate a specific representation vector for each gene in each cell; Network inference and cell embedding generation: The specific representation vector of cell genes is input into a graph neural network model. Based on the multi-head cross attention mechanism, the regulatory relationship between transcription factors and target genes in the cell is calculated to output the gene regulatory network of the cell. At the same time, the transcription factor representation is weighted and integrated based on the regulatory relationship, and then processed with the initial transformation features of the target gene through network layers to generate the cell embedding vector. Expression reconstruction: The cell embedding vector is passed through a feedforward network decoder to reconstruct the cell's gene expression; Model training: The parameters of the graph neural network model are trained by jointly optimizing the comparative learning loss function and the weighted reconstruction loss function; Regulatory network inference: Obtain the specific representation vectors of all genes in the target cell, input them into the trained graph neural network, and infer the gene regulatory network of the target cell.

2. The method according to claim 1, characterized in that, Preprocessing of raw scRNA-seq data includes: removing cells with fewer than a preset number of genes, removing genes expressed in fewer than a preset number of cells, standardizing the gene expression value of each cell to the same level, and screening for highly variable genes. And / or, the K-nearest neighbor algorithm is used to determine the neighborhood relationships of cells; And / or, use a momentum contrastive learning framework to learn the embedding representation of cells.

3. The method according to claim 1, characterized in that, The inductive graph neural network model includes two SAGEConv convolutional layers, and between the SAGEConv convolutional layers are batch normalization layers, PReLU activation functions, and Dropout layers; the output of the last SAGEConv convolutional layer is used as the gene token. The inductive graph neural network model guided by a priori gene regulation network includes: defining positive sample edges based on the edge index of the priori gene regulation network; generating negative sample edges through a structured negative sampling method; calculating the dot product of the gene tokens at both ends of the positive and negative sample edges and using the sigmoid function as the score of the corresponding edge; training the inductive graph neural network model by maximizing the log-likelihood of the positive sample edges and minimizing the log-likelihood of the negative sample edges; the structured negative sampling refers to: for each positive sample edge... By fixing the head node gene And randomly sample a tail node gene To construct negative sample edges ,and , Represents the relationship between genes in the priori gene regulatory network. Genes that are connected by edges.

4. The method according to claim 1, characterized in that, The feature fusion generates the specificity representation vector using the following formula: ； in, Representing the The first cell Cell-gene-specific expression vectors of individual genes. Represents hidden dimensions; Representing the A cell token for each cell; Representing the Gene token for each gene; Representing the The first cell Gene expression data after imputation of gene expression values.

5. The method according to claim 1, characterized in that, The calculation of the regulatory relationship between transcription factors and target genes in cells based on the multi-head cross-attention mechanism specifically includes: The specific representation vectors of transcription factors and target genes are subjected to nonlinear transformations using different multilayer perceptrons to obtain the query matrix. Key matrix Sum matrix ; ； in, and This represents two different multilayer perceptrons. and Generate query matrix respectively The trainable weight matrix and bias vector used; and Generate the key matrix respectively The trainable weight matrix and bias vector used; and These are the generated value matrices. The trainable weight matrix and bias vector used; represent The fusion feature vector of a transcription factor in n cells. represent Fusion feature vectors of target genes in n cells The dimension of the fused feature vector; , It is the number of multiple heads; Represents transcription factor characterization, Representative target gene characterization; Constructing a mask matrix using prior gene regulatory networks This is used to filter out non-existent regulatory relationships in attention calculations; based on the query matrix... Key matrix and mask matrix The attention score between transcription factors and target genes is calculated through a multi-head cross-attention mechanism; ； in, , Representing the first The cell number Transcription factors in the size and target genes The representation; Representing the The cell number Transcription factors in the size With target genes The attention score is the relevance score; The calculated attention score is first processed After normalization, the adjacency matrix representing the strength of the regulatory relationship is obtained by averaging across multiple dimensions. That is, the first Gene regulatory networks in individual cells: ； ； in, Representing the Transcription factors of individual cells With target genes Regulatory relationships within gene regulatory networks.

6. The method according to claim 1, characterized in that, The specific steps for generating cell embedding vectors include: The attention scores of the cell's gene regulatory network are used to perform a weighted summation of transcription factor representations to obtain context-aware representations of target genes. The specific calculation formula is as follows: ； in It is an attention score; It is a transcription factor characteristic; It is the context-aware representation of the target gene obtained; The context-aware representation is nonlinearly transformed using a first multilayer perceptron to obtain the transformed target gene representation. The transformed target gene characterization is residually linked with the initial transformed features of the target gene to obtain an enhanced target gene characterization. The enhanced target gene characterization is subjected to layer normalization, and then feature extraction is performed using a second multilayer perceptron. The outputs of the second multilayer perceptron for all target genes are summed in the gene dimension to generate a cell embedding vector representing the overall cell state.

7. The method according to claim 1, characterized in that, The specific steps of model training include: Construct a contrastive learning loss function that brings positive samples closer and pushes negative samples further away based on cell embedding vectors; A weighted reconstruction loss function is constructed, which assigns higher weights to the reconstruction errors of non-zero values ​​in the gene expression matrix; The contrastive learning loss and weighted reconstruction loss are jointly optimized to train the parameters of the graph neural network and decoder, while keeping the pre-trained cell token and gene token parameters fixed.

8. The method according to claim 7, characterized in that, The weighted reconstruction loss function is: ； in, Represents the weighted reconstruction loss; This represents the calculated mean square error; This represents the portion of the gene expression matrix with values ​​of 0 after preprocessing. This represents the portion of the gene expression matrix after preprocessing that has non-zero values; Representative at The graph neural network model predicts the value at the corresponding position. Representative at The graph neural network model prediction value at the corresponding position.

9. The method according to claim 7, characterized in that, The contrastive learning loss function is: ； in, Represents contrastive learning loss; For a single target cell, the cell embedding vector. To and The embedding vector of a single positive sample cell that constitutes a positive sample pair; To compare the set of negative sample cell embedding vectors in the dynamic learning queue, To compare the number of negative samples in the dynamic learning queue; The temperature coefficient is a hyperparameter that is greater than zero.

10. A single-cell-level dynamic gene regulation network inference system based on graph neural networks and contrastive learning, comprising a memory and a processor, wherein the memory stores a computer program, characterized in that, When the computer program is executed by the processor, the processor causes the processor to implement the method as described in any one of claims 1 to 9.

Images