Method and system for improving prognosis and determining treatment sensitivity of NBL patient

Abstract

The invention belongs to the technical field of neuroblastoma, and particularly discloses a prognosis improvement and treatment sensitivity determination method and system for NBL patients, and the method comprises the steps: collecting and preprocessing data, constructing a gene co-expression network, screening immune-related modules in the gene co-expression network, and obtaining immune-related prognosis features. According to the method, neuroblastoma patients are subjected to immunotherapy, a novel immune-related risk stratification model NIRSM is constructed, prognosis, immunotherapy and chemotherapy reactions of the neuroblastoma patients are predicted, the novel immune-related risk stratification model NIRSM shows that the patients in a low-risk group show better immunotherapy reactions, and differences of reactions of the patients in different groups to certain chemotherapy drugs are predicted. Cyclophosphamide and topotecan are more sensitive in a low-risk group, and cis-platinum and docetaxel are more sensitive in a high-risk group. By adopting the technical scheme, the optimal clinical prognosis and treatment strategy are promoted, and the progress of neuroblastoma patients with different risks can be known.

Description

Technical Field

[0001] The invention belongs to the technical field of neuroblastoma, and relates to a method and a system for improving the prognosis and determining the treatment sensitivity of NBL patients. Background Art

[0002] Neuroblastoma (NBL) is the main extracranial solid tumor in children, originating from the sympathetic nervous system. NBL accounts for 7%-8% of childhood malignancies and about 15% of childhood cancer-related mortality. According to clinical and biological factors, such as age at diagnosis, disease stage, and MYCN status, NBL patients are traditionally divided into low-risk, intermediate-risk, and high-risk groups. INSS stage 1 / 2 patients without MYCN amplification are considered low-risk NBL, INSS stage 2 / 3 / 4 / 4 patients younger than 18 months and without MYCN amplification are considered intermediate-risk NBL, and INSS stage 4 / 4s patients older than 18 months or patients older than INSS stage 1 with MYCN amplification are considered high-risk NBL. The 5-year survival rate of patients with low-risk or intermediate-risk NBL is more than 90%, while the 5-year survival rate of patients with high-risk NBL is less than 50%.

[0003] The main treatment strategies for NBL patients are surgery, chemotherapy, and radiotherapy. Immunotherapy is a promising treatment method and has also been used to combat NBL. In particular, anti-GD2 immunotherapy has been widely used in NBL and has improved the overall survival (OS) of HR-NBL patients. In addition, Heczey et al. demonstrated that chimeric antigen receptor (CAR)-based immunotherapy, such as L1-cell adhesion molecule (L1-CAM) and bisphosphonoside-specific CAR (GD2-CAR), is safe and has anti-tumor efficacy in heavily pretreated relapsed / refractory NBL patients.

[0004] However, due to the heterogeneity of NBL, not all patients can benefit from immunotherapy. Considering that the risk stratification system is the main treatment guideline for NBL patients, there is an urgent need to construct a risk stratification model to improve the prognosis of NBL patients and optimize the treatment method. Summary of the invention

[0005] The purpose of the present invention is to provide a method and system for improving the prognosis and determining the treatment sensitivity of NBL patients, promote the best clinical prognosis and treatment strategy, and help understand the progress of neuroblastoma patients with different risks.

[0006] In order to achieve the above object, the basic scheme of the present invention is: a method for improving the prognosis and determining the treatment sensitivity of NBL patients, comprising the following steps:

[0007] S1, extract mRNA based on Gencode GTF annotation file and perform preprocessing;

[0008] S2, using the lmQCM algorithm, construct a gene co-expression network and calculate the characteristic gene values ​​of the gene co-expression network module;

[0009] S3, using the ESTIMATE algorithm to screen immune-related modules in the gene co-expression network;

[0010] S4, through Spearman correlation analysis, using gene significance GS, the correlation between module characteristic genes and single module genes MM, and protein interaction network PPI, to screen immune key genes in immune-related modules;

[0011] S5, Gene Ontology enrichment analysis of immune-related modules was performed using the “clusterProfiler” R package, and the p-value of gene set enrichment was calculated;

[0012] S6, univariate Cox proportional hazard regression analysis was performed on the immune key genes in the immune-related module to evaluate their relationship with patient prognosis, and genes with p-values ​​less than 0.01 were screened as immune-related prognostic genes;

[0013] Lasso regression analysis was used to screen the genes with the most immune-related prognostic characteristics. Genes with non-zero coefficients were selected as the screened immune-related prognostic characteristics, and a new immune-related risk stratification model NIRSM was constructed.

[0014] S7, based on the novel immune-related risk stratification model NIRSM, the risk score of each patient was calculated, and the median was used as the cutoff point to divide all patients into high-risk and low-risk groups;

[0015] S8, using the "limma" R package, compare the expression patterns of immune signature genes related to the new immune-related risk stratification model NIRSM and immune-related co-expression network module genes in the low-risk group and the high-risk group, calculate the false discovery rate BH FDR using the Limma package, and then calculate the q value of multiple verification. If the q value is less than 0.05 and |Fold Change| (the standard for screening differentially expressed genes) is greater than 1.2, it is considered to be significantly different;

[0016] S9, Chi-squared test was used to compare the age, INSS stage, MYCN status and traditional stratification clinical characteristics between the low-risk group and the high-risk group to test the relationship between the new immune-related risk stratification model NIRSM and clinical characteristics;

[0017] S10, analyzing the immune microenvironment of the novel immune-related risk stratification model NIRSM, using the Wilcoxon test to compare the immune scores of immune cells between the low-risk group and the high-risk group;

[0018] The CIBERSORT algorithm was used to evaluate the proportion of 22 tumor-infiltrating immune cells in each patient, analyze the relationship between the risk stratification model and the immune microenvironment, and evaluate the immune checkpoint inhibitors ICIs and human leukocyte antigen HLA between the low-risk and high-risk groups;

[0019] S11, predicting immunotherapy response, using the tumor immune dysfunction and rejection TIDE model to predict each patient's immunotherapy response and obtain each patient's TIDE score. A high TIDE score is associated with a poor immunotherapy response.

[0020] S12, predict chemotherapy sensitivity, obtain drug sensitivity data from the Cancer Drug Sensitivity Genomics GDSC database, use the "oncoPredict" R package to calculate the half-inhibitory concentration IC50 of 5 neuroblastoma chemotherapy drugs, and determine the sensitivity of high-risk and low-risk groups of patients to different neuroblastoma chemotherapy drugs.

[0021] The working principle and beneficial effects of this basic scheme are: the new immune-related risk stratification model NIRSM used in this technical scheme is associated with poor prognosis and can be used as an independent prognostic prediction indicator. NIRSM can be used to predict the immunotherapy response of NBL patients, and NIRSM can be used to predict the chemotherapy response of NBL patients.

[0022] There are differences in the responses of patients in different groups to certain chemotherapy drugs. Cyclophosphamide and topotecan are more sensitive in the low-risk group, and cisplatin and docetaxel are more sensitive in the high-risk group. Patients in the low-risk group show better response to immunotherapy. NIRSM is associated with poor prognosis. NIRSM is related to immunotherapy and chemotherapy. Patients in the low-risk group are more sensitive to immunotherapy and cyclophosphamide and topotecan, and patients in the high-risk group are more sensitive to cisplatin and docetaxel. NIRSM not only promotes the best clinical prognosis and treatment strategy, but also helps us understand the progression of NBL patients with different risks.

[0023] Further, the method for preprocessing the mRNA gene data is:

[0024] The average expression value of all samples was calculated and used as the expression value of the replicate genes;

[0025] Genes with zero expression values ​​in 80% of the samples were excluded, neuroblastoma patients without follow-up information were removed, and outlier samples were excluded.

[0026] Preprocess the data for subsequent use.

[0027] Furthermore, in step S2, a gene co-expression network is constructed, and the method of using the lmQCM algorithm to calculate the characteristic gene values ​​of the gene co-expression network module is as follows:

[0028] Calculate the median absolute deviation of each gene and sort them from large to small to select the top 8000 genes;

[0029] The expression matrix was filtered using the “fastFilter” function, with the thresholds set to: lowest percentile mean = 0.2, lowest percentile variance = 0.2, and variance function = var;

[0030] The gene co-expression network was constructed using the “lmQCM” function with the following parameters: γ = 0.7, λ = 1.0, t = 1.0, β = 0.4, minimum cluster size = 10, correlation coefficient calculation = “spearman”. The “lmQCM” function automatically calculated the characteristic gene values ​​of each gene co-expression network module.

[0031] The characteristic gene value is the first principal component of the gene co-expression module expression matrix and can represent the overall expression level of the gene co-expression module.

[0032] Furthermore, in step S3, the ESTIMATE algorithm is used to screen the immune-related modules in the gene co-expression network:

[0033] Based on the gene expression profiles of 141 stromal cell characteristic genes and 141 immune cell characteristic genes, the single-sample gene set enrichment analysis algorithm was used to calculate the base score and immune score, respectively, and finally the two scores were combined to generate the ESTIMATEscore;

[0034] Based on the gene expression profile of each neuroblastoma sample, the ESTIMATE algorithm was used to predict the immune score of immune cells in tumor tissue;

[0035] The Spearman correlation coefficient ρ between the characteristic genes of each gene co-expression network module and the immune score was calculated, and the gene co-expression network module with the highest correlation was selected as the immune-related module.

[0036] Screen immune-related modules in gene co-expression networks and mine immune-related information.

[0037] Further, the steps of screening immune key genes in the immune-related module in step S4 are as follows:

[0038] Through Spearman correlation analysis, strict thresholds were set: genes with GS>0.8 and MM>0.8 were considered as potential key genes;

[0039] The protein-protein interaction network PPI was constructed for the genes in the immune-related module using the STRING database, and the genes ranked in the top 10% of the node degree were considered as key genes;

[0040] The key genes screened by the two methods were intersected to finally obtain the key immune genes.

[0041] Identify key immune-related genes for subsequent use.

[0042] Further, in step S7, the risk score of each patient is calculated based on the novel immune-related risk stratification model NIRSM, and the specific steps are as follows:

[0043] The immune-related prognostic gene dataset was randomly divided into a training set and a validation set in a ratio of 7:3. In the training set, the immune-related prognostic genes were subjected to Lasso regression analysis using the "glmnet" R package to screen the genes with the most immune-related prognostic characteristics;

[0044] In the training set, genes with non-zero coefficients were selected as the immune-related prognostic features for screening: Bruton's tyrosine kinase BTK, CD247 molecule, interleukin 7 receptor IL7R, and interferon regulatory factor 8 IRF8, to construct a new immune-related risk stratification model NIRSM;

[0045] The risk score is the weighted sum of the Lasso-Cox regression coefficient and the corresponding characteristic gene expression value. The risk score model is:

[0046] (0.20)*BTK gene expression value+(0.14)*CD247 gene expression value+

[0047] (-0.18)*IL7R gene expression value+(-0.43)*IRF8 gene expression value

[0048] A risk score was measured for each patient.

[0049] Measuring the risk score for each patient is helpful to predict the immunotherapy response of NBL patients.

[0050] Furthermore, the specific steps of analyzing the immune microenvironment of the novel immune-related risk stratification model NIRSM in step S10 are as follows:

[0051] Based on the principle of linear support vector regression (LSVR), the expression matrix of immune cell subtypes is deconvoluted to estimate the abundance of immune cells. CIBERSORT provides LM22 expression data of 22 common immune infiltrating cells, including immune cells of different cell types and functional states.

[0052] The Wilcoxon test was used to compare the distribution of 22 infiltrating immune cells between the low-risk group and the high-risk group, and immune checkpoint inhibitors (ICIs) and human leukocyte antigen (HLA) were evaluated between the low-risk group and the high-risk group. A p value less than 0.05 was considered a significant difference.

[0053] Analyze the immune microenvironment of the risk stratification model and reveal the immune basis of the risk stratification model.

[0054] The present invention also provides a system for improving the prognosis and determining the treatment sensitivity of NBL patients, comprising a processing unit, wherein the processing unit executes the method of the present invention to improve the prognosis and optimize the treatment of NBL patients.

[0055] This system uses processing units to construct gene co-expression network analysis, analyze the association network between gene expressions in different physiological states or different disease states, and thus reveal the interaction and regulatory mechanism between genes. A risk stratification model is constructed as a prognostic predictor for patients, improving the prognosis of NBL patients and optimizing treatment strategies. BRIEF DESCRIPTION OF THE DRAWINGS

[0056] Figure 1 It is a schematic flow chart of the method for improving the prognosis and determining the treatment sensitivity of NBL patients of the present invention. DETAILED DESCRIPTION

[0057] Embodiments of the present invention are described in detail below, examples of which are shown in the accompanying drawings, wherein the same or similar reference numerals throughout represent the same or similar elements or elements having the same or similar functions. The embodiments described below with reference to the accompanying drawings are exemplary and are only used to explain the present invention, and cannot be understood as limiting the present invention.

[0058] In the description of the present invention, it is necessary to understand that the terms "longitudinal", "lateral", "up", "down", "front", "back", "left", "right", "vertical", "horizontal", "top", "bottom", "inside", "outside", etc., indicating the orientation or position relationship, are based on the orientation or position relationship shown in the accompanying drawings, and are only for the convenience of describing the present invention and simplifying the description, rather than indicating or implying that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation, and therefore cannot be understood as a limitation on the present invention.

[0059] In the description of the present invention, unless otherwise specified and limited, it should be noted that the terms "installed", "connected" and "connected" should be understood in a broad sense. For example, it can be a mechanical connection or an electrical connection, or it can be the internal connection between two components. It can be a direct connection or an indirect connection through an intermediate medium. For ordinary technicians in this field, the specific meanings of the above terms can be understood according to the specific circumstances.

[0060] The present invention discloses a method for improving the prognosis and determining the treatment sensitivity of NBL patients, such as Figure 1 As shown, the following steps are included:

[0061] S1, mRNA was extracted based on the Gencode GTF annotation file and preprocessed, and four outlier samples (GSM1205655, GSM1205656, GSM1205657, and GSM1205658) were excluded according to the expression data sample density map;

[0062] S2, using the lmQCM (local maximum quasi-cluster merging) algorithm, construct the gene co-expression network (GCN) and calculate the characteristic gene values ​​of the gene co-expression network module;

[0063] S3, using the ESTIMATE algorithm to screen immune-related modules in the gene co-expression network;

[0064] S4, through Spearman correlation analysis, using gene significance GS (GS indicates the correlation between module genes and immune characteristics), the correlation MM between module characteristic genes and single module genes, and protein interaction network PPI, to screen immune key genes in immune-related modules;

[0065] S5. In order to explore the biological functions of immune-related modules, the “clusterProfiler” R package was used to perform gene ontology (GO) enrichment analysis on immune-related modules and calculate the p-value of gene set enrichment. The hypergeometric test was used to calculate the p-value of gene set enrichment, and a p-value less than 0.05 was considered statistically significant. The enrichment analysis covered three aspects: biological process (BP), cellular component (CC), and molecular function (MF). The top 10 enriched pathways in each were selected for display.

[0066] S6, univariate Cox proportional hazard regression analysis was performed on the immune key genes in the immune-related module to evaluate their relationship with patient prognosis, and genes with p-values ​​less than 0.01 were screened as immune-related prognostic genes;

[0067] Lasso regression analysis was used to screen the genes with the most immune-related prognostic characteristics. Genes with non-zero coefficients were selected as the screened immune-related prognostic characteristics, and a new immune-related risk stratification model NIRSM was constructed; the new immune-related risk stratification model NIRSM can be used as an independent prognostic factor to predict the prognosis of neuroblastoma patients and optimize treatment methods. The model showed that patients in the low-risk group showed better immunotherapy responses, and also predicted that different groups of patients had different responses to certain chemotherapy drugs. Cyclophosphamide and topotecan were more sensitive in the low-risk group, and cisplatin and docetaxel were more sensitive in the high-risk group. The risk stratification model promotes the best clinical prognosis and treatment strategies, and helps us understand the progress of neuroblastoma patients with different risks.

[0068] S7, based on the novel immune-related risk stratification model NIRSM, the risk score of each patient was calculated, and the median was used as the cutoff point to divide all patients into high-risk and low-risk groups;

[0069] S8, test the predictive performance of the risk stratification model, and compare whether the predictive performance of the risk stratification model is better than other clinical predictive indicators of neuroblastoma, such as age, INSS stage, and MYCN status, by time-dependent receiver operating characteristic curve (ROC) analysis;

[0070] S9, using the "limma" R package, the expression patterns of immune signature genes related to the new immune-related risk stratification model NIRSM and immune-related co-expression network module genes in the low-risk group and the high-risk group were compared. The false discovery rate BH FDR was calculated using the Limma package, and then the q value of multiple verification was calculated. If the q value was less than 0.05 and |Fold Change| (the standard for screening differentially expressed genes) was greater than 1.2, it was considered to be significantly different; the molecular basis of the risk stratification model was explored, and differential expression analysis was performed to compare the expression patterns of genes related to the risk stratification model and module 1 genes in the low-risk group and the high-risk group. Compared with the high-risk group, the expression levels of all genes and module 1 genes in the low-risk group were higher (q value < 0.05, |Fold change|> 1.2). There were different molecular expression patterns between different risk groups. The downregulation of key immune genes affected the expression patterns and biological functions of module genes, thereby changing the downstream related biological functions or pathways, and promoting the occurrence and development of neuroblastoma.

[0071] S10, test the relationship between the risk stratification model and clinical characteristics, compare the age, INSS stage, MYCN status and traditional stratification clinical characteristics between the low-risk group and the high-risk group by Chi-squared test, and test the relationship between the new immune-related risk stratification model NIRSM and clinical characteristics;

[0072] S11, analyzing the immune microenvironment of the novel immune-related risk stratification model NIRSM, using the Wilcoxon test to compare the immune scores of immune cells between the low-risk group and the high-risk group;

[0073] The CIBERSORT algorithm was used to evaluate the proportion of 22 tumor-infiltrating immune cells in each patient (to obtain the immune cell score, NIRSM is closely related to immune-related functions and may be associated with immunotherapy), analyze the relationship between the risk stratification model and the immune microenvironment, and evaluate the immune checkpoint inhibitors ICIs and human leukocyte antigen HLA between the low-risk group and the high-risk group;

[0074] S12, predicting immunotherapy response, using the TIDE (Tumor Immune Dysfunction and Exclusion) model to predict each patient's immunotherapy response and obtain each patient's TIDE score. A high TIDE score is associated with a poorer immunotherapy response.

[0075] The TIDE model is a computational tool for evaluating tumor response to immunotherapy. It predicts a patient's potential response to immunotherapy by analyzing specific immune and tumor features in the tumor microenvironment. A high TIDE score is generally associated with a poorer response to immunotherapy because a high score means there are more immunosuppressive factors in the tumor microenvironment. The Wilcoxon test was then used to compare the TIDE scores between the low-risk and high-risk groups. A P value less than 0.05 was considered statistically significant.

[0076] Patients in the low-risk group were significantly associated with lower TIDE scores (p value = 0.028). Since lower TIDE scores were associated with better immunotherapy response and outcomes, it was suggested that NIRSM was associated with immunotherapy and low-risk patients might benefit from immunotherapy.

[0077] S13, predict chemotherapy sensitivity, obtain drug sensitivity data from the Genomics of Drug Sensitivity in Cancer (GDSC) database, use the “oncoPredict” R package to calculate the half-inhibitory concentration IC50 of five neuroblastoma chemotherapy drugs (cisplatin, cyclophosphamide, docetaxel, topotecan, vinorelbine), and determine the sensitivity of high-risk and low-risk groups to different neuroblastoma chemotherapy drugs;

[0078] IC50 is an important parameter in drug sensitivity testing, which indicates the drug concentration required to inhibit the target biological activity (such as cell proliferation, enzyme activity, etc.) by 50% under specific conditions. Lower IC50 values ​​are generally associated with higher drug sensitivity. Subsequently, the IC50 between low-risk and high-risk subgroups was compared using the Wilcoxon test, and a p value less than 0.05 was considered statistically significant.

[0079] By comparing the IC50 of the two groups of chemotherapy drugs, it was found that the anti-tumor activity of 4 out of 5 chemotherapy drugs (cisplatin, cyclophosphamide, docetaxel, and topotecan) was significantly different between the two groups (p value < 0.05). Among them, patients in the high-risk group were more sensitive to cisplatin and docetaxel, and patients in the low-risk group were more sensitive to cyclophosphamide and topotecan. The results suggested a link between NIRSM and chemotherapy, indicating that patients in different risk groups may respond differently to specific chemotherapy drugs.

[0080] In a preferred embodiment of the present invention, the method for preprocessing mRNA gene data is:

[0081] The average expression value of all samples was calculated and used as the expression value of the replicate genes;

[0082] Genes with zero expression values ​​in 80% of the samples were excluded, neuroblastoma patients without follow-up information were removed, and outlier samples were excluded.

[0083] In a preferred embodiment of the present invention, the method for constructing a gene co-expression network in step S2 and using the lmQCM algorithm to calculate the characteristic gene values ​​of the gene co-expression network module is as follows:

[0084] Calculate the median absolute deviation of each gene and sort them from large to small to select the top 8000 genes;

[0085] The expression matrix was filtered using the “fastFilter” function, with the thresholds set to: lowest percentile mean = 0.2, lowest percentile variance = 0.2, and variance function = var;

[0086] The “lmQCM” function was used to construct a gene co-expression network with the following parameters: γ = 0.7, λ = 1.0, t = 1.0, β = 0.4, minimum cluster size = 10, correlation coefficient calculation = “spearman”. The “lmQCM” function (lmQCM algorithm explored gene co-expression modules and identified 15 gene co-expression modules in total) automatically calculated the characteristic gene value of each gene co-expression network module. The characteristic gene value is the first principal component of the expression matrix of the gene co-expression module, which can represent the overall expression level of the gene co-expression module.

[0087] In a preferred embodiment of the present invention, the ESTIMATE algorithm is used in step S3 (ESTIMATE is a tool that uses gene expression features to infer the proportion of stromal cells and immune cells in tumor tissue and the purity of tumor), and the steps of screening the immune-related modules in the gene co-expression network are:

[0088] Based on the gene expression profiles of 141 stromal cell characteristic genes and 141 immune cell characteristic genes, the single-sample gene set enrichment analysis algorithm was used to calculate the base score and immune score, respectively, and finally the two scores were combined to generate the ESTIMATEscore;

[0089] Based on the gene expression profile of each neuroblastoma sample, the ESTIMATE algorithm was used to predict the immune score (infiltration level) of immune cells in tumor tissue;

[0090] The Spearman correlation coefficient ρ between the characteristic genes of each gene co-expression network module and the immune score was calculated, and the gene co-expression network module with the highest correlation was selected as the immune-related module.

[0091] In a preferred embodiment of the present invention, the step of screening the key immune genes in the immune-related module in step S4 is as follows:

[0092] Gene significance (GS) and correlation between module genes (MM) were used to evaluate the relationship between module genes and immune scores. GS indicated the correlation between module genes and immune features, and MM indicated the correlation between module feature genes and individual module genes. In this example, strict thresholds were set by Spearman correlation analysis: genes with GS>0.8 and MM>0.8 were considered as potential key genes;

[0093] The protein-protein interaction network PPI was constructed for the genes in the immune-related module using the STRING database (https: / / string-db.org / ), and the top 10% genes in terms of node degree (the number of edges directly connected to the node, directly calculated by the STRING database) were considered key genes;

[0094] The key genes screened by the above two methods were intersected to finally obtain the key immune genes for further analysis.

[0095] In a preferred embodiment of the present invention, in step S7, the risk score of each patient is calculated based on the novel immune-related risk stratification model NIRSM, and the specific steps are:

[0096] The immune-related prognostic gene dataset was randomly divided into a training set and a validation set in a ratio of 7:3. In the training set, the immune-related prognostic genes were subjected to Lasso regression analysis (least absolute shrinkage and selection operator) using the "glmnet" R package to screen the genes with the most immune-related prognostic characteristics. Lasso regression is a linear model that achieves shrinkage and selection of coefficients by introducing an L1 regularization term, thereby reducing the complexity of the model while maintaining the model's predictive performance and improving the model's interpretability.

[0097] In the training set, genes with non-zero coefficients were selected as the screened immune-related prognostic features: Bruton's tyrosine kinase BTK, CD247 molecule, interleukin 7 receptor IL7R and interferon regulatory factor 8 IRF8, to construct a new immune-related risk stratification model NIRSM.

[0098] The risk score is the weighted sum of the Lasso-Cox regression coefficient and the corresponding characteristic gene expression value. The risk score model is:

[0099] (0.20)*BTK gene expression value+(0.14)*CD247 gene expression value+

[0100] (-0.18)*IL7R gene expression value + (-0.43)*IRF8 gene expression value,

[0101] The risk score of each patient was measured. The Lasso-Cox regression coefficient is in the brackets. The neuroblastoma dataset was randomly divided into a training dataset and a validation dataset in a ratio of 7:3. In the training dataset, Lasso regression analysis was used to select the best characteristic genes for the above genes. Ten-fold cross validation was used and the optimal lambda (λ) value was selected as 0.00875. Finally, four characteristic genes, BTK, CD247, IL7R and IRF8, were obtained, and a novel immune-related risk stratification model was further constructed. The model is a linear combination model, which is the weighted sum of the expression value of each characteristic gene and its corresponding Lasso-Cox regression coefficient.

[0102] In a preferred embodiment of the present invention, the specific steps of analyzing the immune microenvironment of the novel immune-related risk stratification model NIRSM in step S10 are as follows:

[0103] Based on the principle of linear support vector regression (LSVR), the expression matrix of immune cell subtypes is deconvoluted to estimate the abundance of immune cells. CIBERSORT provides LM22 expression data of 22 common immune infiltrating cells, including immune cells of different cell types and functional states.

[0104] The Wilcoxon test was used to compare the distribution of 22 infiltrating immune cells between the low-risk group and the high-risk group. The expression patterns of immune checkpoint inhibitors ICIs and human leukocyte antigens HLA, such as CD274, CTLA4, HAVCR2, LAG3, PDCD1, PDCD1LG2, TNFRSF18, and human leukocyte antigen (HLA) genes (HLA-A, HLA-B, HLA-C, HLA-DMA, HLA-DMB, HLA-DOA, HLA-DOB, HLA-DPA1, HLA-DPB1, HLA-DQA1, HLA-DQA2, HLA-DQB1, HLA-DQB2, HLA-DRA, HLA-DRB1, HLA-DRB5, HLA-E, HLA-F, and HLA-G) were evaluated between the low-risk group and the high-risk group. A p value less than 0.05 was considered to be a significant difference.

[0105] The present embodiment specifically downloads the gene expression profile of neuroblastoma of GSE49710 data set from GEO database.In order to mine the gene co-expression network module, the top 8000 genes with the highest median absolute deviation are selected.The GCN module is constructed using the lmQCM algorithm to obtain 15 GCN modules.The correlation between the GCN module and the level of immune cell infiltration (immune score) is estimated, and module 1 is very prominent.At the same time, biological function analysis confirms that module 1 is related to immune-related functions.Secondly, the key genes of module 1 are analyzed, and a novel immune-related risk stratification model (NIRSM) based on these key genes is constructed using Lasso regression analysis.

[0106] GS analysis was used to calculate the correlation between module 1 genes and immune scores, and 93 key genes were mined (ρ>0.8). MM analysis was used to calculate the correlation between module 1 characteristic genes and module 1 genes, and 288 key genes were mined (ρ>0.8). Next, PPI analysis was performed on the genes of module 1, and the top 10% of genes with the highest node degree were selected. In order to explore common IRGs, we took the intersection of the key genes mined by the above methods and obtained 15 IRGs (BTK, CD2, CD247, CD48, CD74, GZMA, IKZF1, IL7R, IRF8, LCP2, PLEK, PTPRC, SYK, VAV1, ZAP70).

[0107] The relationship between NIRSM and patient survival was studied and it was observed that NIRSM was associated with poor prognosis and could serve as an independent predictor. At the same time, the survival prediction ability of NIRSM was also evaluated and confirmed in the test and full datasets.

[0108] Further comparison of the predictive performance of NIRSM with other clinical factors (such as age, INSS stage, and MYCN status) showed that NIRSM had better predictive performance. The molecular basis, clinical characteristics, and tumor immune microenvironment (TIM) of NIRSM were studied, and the relationship between NIRSM and immunotherapy was evaluated. It was found that patients in the low-risk group responded better to immunotherapy.

[0109] By comparing the clinical characteristics of the high- and low-risk groups, such as age, INSS stage, MYCN status, and traditional stratification, it was found that the two groups of patients had significant differences in these clinical characteristics. MYCN amplification, age > 18 months, INSS stage (stage 4 or 4s), clinical high risk, and worsening are mostly concentrated in the high-risk group, which is consistent with the characteristics of the model. There are different molecular expression patterns between different risk groups. The downregulation of key immune genes affects the expression pattern and biological function of module genes, thereby changing the downstream related biological functions or pathways, which in turn promotes the occurrence and development of neuroblastoma.

[0110] Comparing the differences in immune cell infiltration levels (immune scores) between the low-risk group and the high-risk group, patients in the low-risk group generally had higher immune scores. Comparing the relative abundance of 22 infiltrating immune cells between the low-risk group and the high-risk group, 59.09% (13 / 22) of the 22 infiltrating immune cells showed significant differences between the high-risk group and the low-risk group. Compared with the high-risk group, the number of immature B cells, macrophages M0, macrophages M1, memory-activated CD4 T cells, follicular helper T cells, and regulatory T cells (Tregs) in the low-risk group increased, while the number of activated mast cells, activated natural killer cells, plasma cells, and CD8 T cells decreased.

[0111] Examining immune checkpoint inhibitors (ICIs) and human leukocyte antigen (HLA) genes between the low-risk and high-risk groups, the high expression patterns of ICIs and HLA genes were more obvious in the low-risk group compared with the high-risk group. NIRSM is closely related to immune-related functions and may be associated with immunotherapy.

[0112] A new immune-related risk stratification model (NIRSM) was successfully constructed, which can predict the prognosis, immunotherapy and chemotherapy response of NBL patients. NIRSM plays an important role in improving the prognosis prediction and optimizing the treatment strategy of NBL patients.

[0113] The present invention also provides a system for improving the prognosis and determining the treatment sensitivity of NBL patients, comprising a processing unit, and the processing unit executes the method of the present invention to improve the prognosis and optimize the treatment of NBL patients.

[0114] The present invention relates NIRSM to poor prognosis and can be used as an independent prognostic predictor. NIRSM can be used to predict the immunotherapy response of NBL patients. NIRSM can be used to predict the chemotherapy response of NBL patients. However, there are differences in the responses of different groups of patients to certain chemotherapy drugs. Cyclophosphamide and topotecan are more sensitive in the low-risk group, and cisplatin and docetaxel are more sensitive in the high-risk group. Patients in the low-risk group show better immunotherapy response. NIRSM is associated with poor prognosis. NIRSM is related to immunotherapy and chemotherapy. Patients in the low-risk group are more sensitive to immunotherapy and cyclophosphamide and topotecan, and patients in the high-risk group are more sensitive to cisplatin and docetaxel. NIRSM not only promotes the best clinical prognosis and treatment strategy, but also helps us understand the progress of NBL patients with different risks.

[0115] In the description of this specification, the description with reference to the terms "one embodiment", "some embodiments", "examples", "specific examples", or "some examples" means that the specific features, structures, materials or characteristics described in conjunction with the embodiment or example are included in at least one embodiment or example of the present invention. In this specification, the schematic representation of the above terms does not necessarily refer to the same embodiment or example. Moreover, the specific features, structures, materials or characteristics described may be combined in any one or more embodiments or examples in a suitable manner.

[0116] Although the embodiments of the present invention have been shown and described, it will be appreciated by those skilled in the art that various changes, modifications, substitutions and variations may be made to the embodiments without departing from the principles and spirit of the present invention, and that the scope of the present invention is defined by the claims and their equivalents.

Claims

1. A method for improving the prognosis and determining the treatment sensitivity of NBL patients, characterized in that: The steps include: S1, extract mRNA based on Gencode GTF annotation file and perform preprocessing; S2, using the lmQCM algorithm, construct a gene co-expression network and calculate the characteristic gene values ​​of the gene co-expression network module; S3, using the ESTIMATE algorithm to screen immune-related modules in the gene co-expression network; S4, through Spearman correlation analysis, using gene significance GS, the correlation between module characteristic genes and single module genes MM, and protein interaction network PPI, to screen immune key genes in immune-related modules; S5, Gene Ontology enrichment analysis of immune-related modules was performed using the "clusterProfiler" R package, and the p-value of gene set enrichment was calculated; S6, univariate Cox proportional hazard regression analysis was performed on the immune key genes in the immune-related module to evaluate their relationship with patient prognosis, and genes with p-values ​​less than 0.01 were screened as immune-related prognostic genes; Lasso regression analysis was used to screen the genes with the most immune-related prognostic characteristics. Genes with non-zero coefficients were selected as the screened immune-related prognostic characteristics, and a new immune-related risk stratification model NIRSM was constructed. S7, based on the novel immune-related risk stratification model NIRSM, the risk score of each patient was calculated, and the median was used as the cutoff point to divide all patients into high-risk and low-risk groups; S8, using the "limma" R package, compare the expression patterns of immune signature genes related to the new immune-related risk stratification model NIRSM and immune-related co-expression network module genes in the low-risk group and the high-risk group, calculate the false discovery rate BH FDR using the Limma package, and then calculate the q value of multiple verification. If the q value is less than 0.05 and |Fold Change| is greater than 1.2, it is considered to be significantly different; S9, chi-square test was used to compare age, INSS stage, MYCN status and traditional stratification clinical characteristics between low-risk group and high-risk group to test the relationship between the new immune-related risk stratification model NIRSM and clinical characteristics; S10, analyzing the immune microenvironment of the novel immune-related risk stratification model NIRSM, using the Wilcoxon test to compare the immune scores of immune cells between the low-risk group and the high-risk group; The CIBERSORT algorithm was used to evaluate the proportion of 22 tumor-infiltrating immune cells in each patient, analyze the relationship between the risk stratification model and the immune microenvironment, and evaluate the immune checkpoint inhibitors ICIs and human leukocyte antigen HLA between the low-risk and high-risk groups; S11, predicting immunotherapy response, using the tumor immune dysfunction and rejection TIDE model to predict each patient's immunotherapy response and obtain each patient's TIDE score. A high TIDE score is associated with a poor immunotherapy response. S12, predict chemotherapy sensitivity, obtain drug sensitivity data from the Cancer Drug Sensitivity Genomics GDSC database, use the "oncoPredict" R package to calculate the half-inhibitory concentration IC50 of 5 neuroblastoma chemotherapy drugs, and determine the sensitivity of high-risk and low-risk groups of patients to different neuroblastoma chemotherapy drugs.

2. The method for improving the prognosis and determining the treatment sensitivity of NBL patients according to claim 1, characterized in that: The method for preprocessing mRNA gene data is: The average expression value of all samples was calculated and used as the expression value of the replicate genes; Genes with zero expression values ​​in 80% of the samples were excluded, neuroblastoma patients without follow-up information were removed, and outlier samples were excluded.

3. The method for improving the prognosis and determining the treatment sensitivity of NBL patients according to claim 1, characterized in that: In step S2, a gene co-expression network is constructed and the method of using the lmQCM algorithm to calculate the characteristic gene values ​​of the gene co-expression network module is as follows: Calculate the median absolute deviation of each gene and sort them from large to small to select the top 8000 genes; The expression matrix was filtered using the "fastFilter" function, with the thresholds set to: lowest percentile mean = 0.2, lowest percentile variance = 0.2, and variance function = var; The "lmQCM" function was used to construct the gene co-expression network, with the parameters set as: γ = 0.7, λ = 1.0, t = 1.0, β = 0.4, minimum cluster size = 10, correlation coefficient calculation = "spearman", and the "lmQCM" function automatically calculated the characteristic gene value of each gene co-expression network module.

4. The method for improving the prognosis and determining the treatment sensitivity of NBL patients according to claim 1, characterized in that: In step S3, the ESTIMATE algorithm is used to screen the immune-related modules in the gene co-expression network: Based on the gene expression profiles of 141 stromal cell characteristic genes and 141 immune cell characteristic genes, the single-sample gene set enrichment analysis algorithm was used to calculate the base score and immune score, respectively, and finally the two scores were combined to generate the ESTIMATEscore; Based on the gene expression profile of each neuroblastoma sample, the ESTIMATE algorithm was used to predict the immune score of immune cells in tumor tissue; The Spearman correlation coefficient ρ between the characteristic genes of each gene co-expression network module and the immune score was calculated, and the gene co-expression network module with the highest correlation was selected as the immune-related module.

5. The method for improving the prognosis and determining the treatment sensitivity of NBL patients according to claim 1, characterized in that: The steps for screening immune key genes in the immune-related module in step S4 are as follows: Through Spearman correlation analysis, strict thresholds were set: genes with GS>0.8 and MM>0.8 were considered as potential key genes; The protein-protein interaction network PPI was constructed for the genes in the immune-related module using the STRING database, and the genes ranked in the top 10% of the node degree were considered as key genes; The key genes screened by the two methods were intersected to finally obtain the key immune genes.

6. The method for improving the prognosis and determining the treatment sensitivity of NBL patients according to claim 1, characterized in that: In step S7, the risk score of each patient is calculated based on the novel immune-related risk stratification model NIRSM. The specific steps are as follows: The immune-related prognostic gene dataset was randomly divided into a training set and a validation set in a ratio of 7:

3. In the training set, the "glmnet" R package was used to perform Lasso regression analysis on the immune-related prognostic genes to screen the genes with the most immune-related prognostic characteristics; In the training set, genes with non-zero coefficients were selected as the immune-related prognostic features for screening: Bruton's tyrosine kinase BTK, CD247 molecule, interleukin 7 receptor IL7R, and interferon regulatory factor 8 IRF8, to construct a new immune-related risk stratification model NIRSM; The risk score is the weighted sum of the Lasso-Cox regression coefficient and the corresponding characteristic gene expression value. The risk score model is: (0.20)*BTK gene expression value+(0.14)*CD247 gene expression value+ (-0.18)*IL7R gene expression value + (-0.43)*IRF8 gene expression value, A risk score was measured for each patient.

7. The method for improving the prognosis and determining the treatment sensitivity of NBL patients according to claim 1, characterized in that: The specific steps of analyzing the immune microenvironment of the novel immune-related risk stratification model NIRSM in step S10 are as follows: Based on the principle of linear support vector regression (LSVR), the expression matrix of immune cell subtypes is deconvoluted to estimate the abundance of immune cells. CIBERSORT provides LM22 expression data of 22 common immune infiltrating cells, including immune cells of different cell types and functional states. The Wilcoxon test was used to compare the distribution of 22 infiltrating immune cells between the low-risk group and the high-risk group, and immune checkpoint inhibitors (ICIs) and human leukocyte antigen (HLA) were evaluated between the low-risk group and the high-risk group. A p value less than 0.05 was considered a significant difference.

8. A system for improving the prognosis and determining the treatment sensitivity of NBL patients, characterized in that: The method comprises a processing unit, wherein the processing unit executes the method according to any one of claims 1 to 7 to improve the prognosis and optimize the treatment of NBL patients.

Images