A Method for Constructing a Graph Neural Network Dataset Based on Gene Expression

Abstract

The present invention belongs to the field of gene expression data analysis, and particularly relates to a method for constructing a graph neural network data set based on gene expression; it includes obtaining gene expression profile data of a disease and extracting RNA data, processing the RNA data to generate an RNA expression matrix; performing standardization processing and differential analysis on the generated RNA expression matrix to obtain a differential expression matrix; using WGCNA to analyze the differential expression matrix and construct an overlapping topological matrix; exporting all edge information and node information in the overlapping topological matrix, and encoding the exported information; screening the encoded information, using the screened information to construct a graph neural network data set, and identifying biomarkers according to the graph neural network data set; the present invention performs WGCNA analysis on differentially expressed genes, pays attention to the relationship between genes, and constructs a graph neural network data set by screening node and edge information, which is beneficial to discovering biomarkers of tumors.

Description

Technical Field

[0001] The present invention belongs to the field of gene expression data analysis, and in particular relates to a method for constructing a graph neural network dataset based on gene expression. Background Art

[0002] Cancer is a major global public health issue. Biomarkers are the most direct, rapid, and effective diagnostic tool in cancer treatment. The screening and acquisition of biomarkers plays a vital role in disease diagnosis, treatment, and efficacy monitoring, and they are also important targets for drug development. Developing innovative biomarker discovery and identification methods, and thereby identifying and identifying valuable biomarkers, has become a key research hotspot in the biomedical field.

[0003] With the rapid development of high-throughput sequencing technology, multi-molecular data, particularly gene expression profile data, are becoming increasingly abundant, from the genomic to the phenotypic level. This opens the door to the use of computational methods to discover biomarkers for complex diseases such as tumors. However, the use of appropriate algorithms to process these massive datasets is a significant challenge facing current bioinformatics and computational biology research. The field of machine learning has continued to develop in recent years, and deep learning is an emerging technology within machine learning algorithms. Many scholars have conducted extensive research on biomarker discovery using deep learning methods, but some challenges remain in practical analysis:

[0004] (1) Most studies use supervised learning methods to study gene expression profile data, but obtaining labeled data through bioinformatics analysis is expensive, and the existing labeled data cannot meet the needs;

[0005] (2) Public databases only record the correspondence between samples and gene expression levels, but lack the gene correlation between different samples. Compared with focusing only on differentially expressed genes, discussing the relationship between genes is more important for studying gene expression regulatory networks;

[0006] (3) Graph neural networks are widely used in social networks, knowledge graphs, recommendation systems, and life sciences due to their excellent ability to process unstructured data, especially in image data. However, it is obviously unreasonable to convert large amounts of gene expression profile data in public databases into images. Summary of the Invention

[0007] To solve the above problems, the present invention provides a method for constructing a graph neural network dataset based on gene expression, comprising the following steps:

[0008] S1. Obtain gene expression profile data of the disease;

[0009] S2. Extract RNA data from the gene expression profile data, and process the RNA data to generate an RNA expression matrix;

[0010] S3. Calculate a differential expression matrix based on the generated RNA expression matrix;

[0011] S4. Use WGCNA to construct an overlapping topological matrix through the differential expression matrix;

[0012] S5. Export all the edge information and node information in the overlapping topological matrix, and encode the exported information;

[0013] S6. Screen the encoded information, construct a graph neural network dataset using the screened information, and identify biomarkers based on the graph neural network dataset.

[0014] Furthermore, the process of obtaining gene expression profile data in step S1 is as follows:

[0015] Download disease sample data from a public database;

[0016] Obtain a human reference genome annotation file;

[0017] Sort out the disease sample data according to the human reference genome annotation file to obtain gene expression profile data.

[0018] Furthermore, the process of processing RNA data includes:

[0019] Perform ID conversion on the RNA data, converting the ensemble ID of the RNA data to symbol ID;

[0020] Set a low expression threshold, and filter out the RNA data with a symbol ID that is less than the low expression threshold.

[0021] Furthermore, the process of constructing an overlapping topological matrix using WGCNA includes:

[0022] S11. Calculate the degree of linear relationship between genes in the differential expression matrix according to the Pearson correlation coefficient, and generate a relationship matrix based on the linear relationship;

[0023] S12. Through a weighted network, convert the relationship matrix into an adjacency matrix;

[0024] S13. Construct an overlapping topological matrix according to the generated adjacency matrix.

[0025] Furthermore, the relationship matrix is expressed as: S = [s ij = [|cor(i, j)|];

[0026] The adjacency matrix is expressed as: A = [aij = [|cor(x i , x j )| β ;

[0027] The overlapping topological matrix is expressed as:

[0028] where i and j are genes, s ij represents the correlation coefficient between gene i and gene j, S is the relationship matrix, cor() is the correlation coefficient function of the R package, β is the soft threshold for the network to conform to the scale-free criterion, x i is the gene expression level of gene i, A is the adjacency matrix, u is the traversal subscript, a uj is the gene expression level of the j-th column in the u-th row, a iu is the gene expression level of the u-th column in the i-th row, a ij is the gene expression level of the j-th column in the i-th row;

[0029] Furthermore, after encoding the edge information and node information derived from the overlapping topological matrix, the scale function is used for standardization, which is expressed as:

[0030]

[0031] where n is the number of disease samples downloaded from the public database, X is the gene expression level of the disease samples, Z is the standard score, and the standard score is the solution process, i.e., the standardization process of the data.

[0032] Furthermore, the information after standardization using the scale function is screened, including:

[0033] First, a correlation threshold is set, and the edge information after standardization is screened according to the correlation threshold;

[0034] After the edge information screening is completed, the node information is screened to ensure that the final number of nodes is more than 90% of the original number of nodes.

[0035] Advantages of the present invention:

[0036] The present invention performs WGCNA analysis on differentially expressed genes, focuses on the relationships between genes, thereby facilitating the discovery of tumor biomarkers; all the edge information and nodes in the overlapping topological matrix are exported, all the exported information is encoded, and then the encoded information is screened using an edge filter and a node filter, and the screened information is used to construct the input data set of the graph neural network algorithm; enabling the use of the graph neural network to identify biomarkers, and the screened data set results in a better identification effect, improving the identification efficiency and accuracy. Brief Description of the Drawings

[0037] Figure 1 This is the flowchart of the construction method of the graph neural network dataset based on gene expression data of the present invention;

[0038] Figure 2 This is the screening result of differential genes of the present invention;

[0039] Figure 3 This is the screening result of soft threshold of the present invention;

[0040] Figure 4 This is the visualization heat map of module characteristic genes of the present invention. Detailed implementation manners

[0041] Next, the technical solutions in the embodiments of the present invention will be clearly and completely described in conjunction with the accompanying drawings in the embodiments of the present invention. Obviously, the described embodiments are only a part of the embodiments of the present invention, rather than all the embodiments. All other embodiments obtained by those of ordinary skill in the art based on the embodiments of the present invention without creative efforts shall fall within the protection scope of the present invention.

[0042] A method for constructing a graph neural network dataset based on gene expression, as Figure 1 shown, includes the following steps:

[0043] S1. Obtain gene expression profile data of a disease;

[0044] S2. Extract RNA data from the gene expression profile data, process the RNA data, and generate an RNA expression matrix;

[0045] S3. Calculate a differential expression matrix based on the generated RNA expression matrix;

[0046] S4. Use WGCNA to construct an overlapping topological matrix through the differential expression matrix;

[0047] S5. Export all edge information and node information in the overlapping topological matrix, and encode the exported information;

[0048] S6. Screen the encoded information, construct a graph neural network dataset using the screened information, and identify biomarkers based on the graph neural network dataset.

[0049] In one embodiment, the gene expression profile data of glioma in the TCGA database was used for analysis. The downloaded gene expression profile data of glioma was annotated with the human CRCH38 reference genome annotation file, and the RNA data in the gene expression profile data of glioma was extracted. The RNA data was subjected to ID conversion, that is, the ensemble ID of the RNA data was converted to the symbol ID. The ensemble ID is the ID identifier used by the database, and the present invention analyzes the relationship between gene expression levels, so it needs to be converted to the symbol ID, and the symbol ID is the corresponding gene name;

[0050] Specifically, the RNA data includes mRNA and lncRNA data;

[0051] A low expression threshold was set. After ID conversion, according to the low expression threshold, the RNA data with low expression levels was screened out to obtain an RNA expression matrix. In this embodiment, the low expression threshold is 40, that is, the gene is not expressed in 95% of the disease samples. A gene with a total expression less than 40 in all disease samples was filtered out. The purpose of doing this is to remove the influence of genes with low expression levels on the overall gene expression level.

[0052] Normalization is to scale the data proportionally so that it falls into a relatively small specific range,

[0053] Preferably, the obtained RNA expression matrix was subjected to normalization processing. Here, the normalization processing was mainly carried out using the quantile normalization function normalizeBetweenArrays() in R language. After the normalization processing was completed, a log transformation was performed to obtain a normalized expression matrix. The normalized expression matrix was grouped, and malignant glioma was used as the reference group and low-grade glioma was used as the control group for differential analysis. The results of differential gene screening are as Figure 2 shown, Figure 2 It is a volcano plot for differential analysis. The abscissa is log2(Fold Change), which represents the logarithm of the difference multiple of the expression level of a gene in two samples. The larger the absolute value, the greater the difference multiple of the expression level between the two samples. The ordinate is -log10(P value), which represents the negative logarithm of FDR. The larger the value, the more significant the differential expression, and the more reliable the differentially expressed genes screened. Each point in the volcano plot for differential analysis represents a gene;

[0054] Preferably, the differential analysis in this embodiment was completed using the limma package and edgeR package in R language, including:

[0055] 101. Use the calcNormFactors() function to calculate the normalization factor of the expression matrix to align the columns of the count matrix;

[0056] 102. Prepare for linear modeling by converting the expression matrix data using the voom() function;

[0057] 103. Fit a linear model for each gene using the lmFit() function;

[0058] 104. Calculate the estimated coefficients and standard errors for the given contrast groups according to the constructed contrast matrix using the contrasts.fit() function;

[0059] 105. Perform Bayesian statistics for differential expression using the eBayes() function.

[0060] In one embodiment, the WGCNA is adopted to analyze the differential expression matrix and construct an overlapping topological matrix. The main processes include:

[0061] S11. Calculate the degree of linear relationship between genes in the differential expression matrix according to the Pearson correlation coefficient, and generate a relationship matrix based on the linear relationship;

[0062] Specifically, the relationship matrix is expressed as: S = [s ij = [|cor(i, j)|];

[0063] S12. Convert the relationship matrix into an adjacency matrix through a weighted network; Since the points in an unweighted network only have a relationship of 0 and 1, the information that can be translated is less, and the strength of the relationship between points cannot be reflected. Therefore, in this embodiment, a weighted network is adopted. By generating a sequence of numbers containing different thresholds and traversing the given thresholds, a soft threshold that meets the scale-free network criterion is found. As Figure 3 shown, the line in the left figure represents the value of the scale-free fitting index selected subjectively, such as 0.9 in the figure, and the relationship matrix is converted into an adjacency matrix;

[0064] Specifically, the adjacency matrix is expressed as: A = [a ij = [|cor(x i , x j )| β ;

[0065] S13. Construct an overlapping topological matrix according to the generated adjacency matrix;

[0066] Specifically, the overlapping topological matrix is expressed as:

[0067]

[0068] Among them, i and j are genes, s ij represents the correlation coefficient between gene i and gene j, S is the relationship matrix, cor() is the correlation coefficient function of the R package, β is the soft threshold for the network to meet the scale-free criterion, x iis the gene expression level of gene i, A is the adjacency matrix, u is the traversal sub, and a uj is the gene expression level of the j-th column in the u-th row, a iu is the gene expression level of the u-th column in the i-th row, a ij is the gene expression level of the j-th column in the i-th row, TOM is the overlapping topological matrix; the more similar the connection situation between gene i and gene j in the network, the TOM ij the larger the value, the more likely gene i and gene j are in the same expression module, that is, in the same classification group, Figure 4 is the visualization graph of module characteristic genes.

[0069] In one embodiment, WGCNA is used to analyze the differential expression matrix and construct an overlapping topological matrix. All the edge information and node information in the overlapping topological matrix are exported. After encoding the exported information, the scale function is used for standardization, which is expressed as:

[0070]

[0071] where n is the number of disease samples downloaded from the public database, X is the gene expression level of the disease samples, and Z is the standard score.

[0072] Preferably, the information after standardization using the scale function is screened, including:

[0073] First, the edge information after standardization is screened: a correlation threshold is set. If the correlation of an edge is greater than the correlation threshold, the edge is retained; if the correlation of an edge is less than the correlation threshold, the edge is removed;

[0074] In this embodiment, according to the weight information standardized by the scale function, the correlation threshold is set to 0.5;

[0075] After the edge information screening is completed, the node information is screened to ensure that the final number of nodes is more than 90% of the original number of nodes.

[0076] where the correlation is the weight value between nodes included in the TOM matrix obtained by WGCNA analysis.

[0077] Since the information of the edges and nodes exported from the TOM matrix obtained by WGCNA is a complete graph, and the complete graph is useless for subsequent classification and prediction of the dataset using the graph neural network algorithm, and the complete graph does not conform to the reality in gene expression, it is necessary to prune the complete graph by screening the edges and nodes to obtain the correlation between genes. The two filters for screening are the screening of edges and the screening of nodes.

[0078] The edge information screening filter is: at least one edge corresponds to each node, and the average degree of the edges is optimal;

[0079] The node information screening filter is: the node ratio is greater than 90%;

[0080] Specifically, the screening of edges follows the rule that each node in this graph has at least one edge connecting to another node. Since the connection correlation between edges represents the correlation between genes, the fewer the remaining edges, the more relevant the nodes connected by these edges. The calculation of the average degree of edges is the total number of edges in the graph divided by the number of nodes. The smaller the degree, the more relevant the nodes in the graph.

[0081] Use the screened information to construct a graph neural network dataset. As shown in Table 1, according to the above process, 7 graph neural network datasets are constructed, as shown in Table 2. Use Python language coding to convert the 7 graph neural network datasets into numbers; in this embodiment, a total of 7 graph embedding algorithms (HOPE, GF, GraRep, SDNE, DeepWalk, LINE, Node2Vec) are used to evaluate the generated graph neural network datasets. The evaluation metrics are micro (micro-average) and macro (macro-average), and they are compared with the graph embedding public dataset. Identify biomarkers based on the evaluation results of the graph neural network datasets.

[0082] Table 1 Graph Neural Network Dataset

[0083]

[0084] Table 2 Basic Encoding for Conversion of Graph Neural Network Datasets

[0085]

[0086] Although the embodiments of the present invention have been shown and described, for those of ordinary skill in the art, it can be understood that various changes, modifications, substitutions, and variations can be made to these embodiments without departing from the principles and spirit of the present invention. The scope of the present invention is defined by the appended claims and their equivalents.

Claims

1. A method for constructing a graph neural network dataset based on gene expression, characterized in that, It includes the following steps: S1. Obtain the gene expression profile data of the disease; S2. Extract RNA data from the gene expression profile data, and process the RNA data to generate an RNA expression matrix; S3. Calculate the differential expression matrix based on the generated RNA expression matrix; S4. Through the differential expression matrix, use WGCNA to construct an overlapping topological matrix; The process of using WGCNA to construct an overlapping topological matrix includes: S11. Calculate the degree of linear relationship between genes in the differential expression matrix according to the Pearson correlation coefficient, and generate a relationship matrix according to the linear relationship; S12. Through a weighted network, convert the relationship matrix into an adjacency matrix; S13. Construct an overlapping topological matrix according to the converted adjacency matrix; The relationship matrix is expressed as: S = [s ij = [|cor(i, j)|]; The adjacency matrix is represented as: A = [a ij = [|cor(x i , x j )| β ; The overlapping topology matrix is represented as: Among them, i and j are genes, and s ij represents the correlation coefficient between gene i and gene j, S is the relationship matrix, cor() is the correlation coefficient function of the R package, β is the soft threshold for the network to conform to the scale-free criterion, and x i is the gene expression level of gene i, A is the adjacency matrix, u is the traversal sub, and a uj is the gene expression level of the j-th column in the u-th row, and a iu is the gene expression level of the u-th column in the i-th row, and a ij is the gene expression level of the j-th column in the i-th row; S5. Export all the edge information and node information in the overlapping topological matrix, and encode the exported information; S6. Screen the encoded information, use the screened information to construct a graph neural network dataset, and identify biomarkers according to the graph neural network dataset.

2. The method for constructing a graph neural network dataset based on gene expression according to claim 1, wherein The process of obtaining gene expression profile data in step S1 is: Download disease sample data from a public database; Obtain the human reference genome annotation file; Sort out the disease sample data according to the human reference genome annotation file to obtain the gene expression profile data.

3. A method for constructing a graph neural network dataset based on gene expression according to claim 1, characterized in that, The process of processing the RNA data includes: Perform ID conversion on the RNA data, and convert the ensemble ID of the RNA data into a symbol ID; Set a low expression threshold, and filter out the RNA data after conversion to symbol ID that is less than the low expression threshold.

4. A method for constructing a graph neural network dataset based on gene expression according to claim 1, characterized in that, After encoding the edge information and node information exported from the overlapping topological matrix, perform standardization processing using the scale function, which is expressed as: Among them, n is the number of disease samples downloaded from the public database, X is the gene expression level of the disease samples, and Z is the standard score.

5. A method for constructing a graph neural network dataset based on gene expression according to claim 4, characterized in that The screening of the information after standardization processing using the scale function includes: First, set a correlation threshold, and screen the edge information after standardization processing according to the correlation threshold; After the edge information screening is completed, screen the node information to ensure that the final number of nodes is more than 90% of the original number of nodes.

Images