Figure 1.
K-M survival curves of patients in PD-1 high expression group and low expression group of GSE26939 (A), GSE68465 (B), GSE68571 (C), and GSE72094 (D).
Citation: William A. Toscano, Jr., Hitakshi Sehgal, Emily Yang, Lindsey Spaude, A. Frank Bettmann. Environment-Gene interaction in common complex diseases: New approaches[J]. AIMS Molecular Science, 2014, 1(4): 126-140. doi: 10.3934/molsci.2014.4.126
| [1] | Keqiang Li, Shangjiu Wang . Symmetry of positive solutions of a $ p $-Laplace equation with convex nonlinearites. AIMS Mathematics, 2023, 8(6): 13425-13431. doi: 10.3934/math.2023680 |
| [2] | Keqiang Li, Shangjiu Wang, Shaoyong Li . Symmetry of large solutions for semilinear elliptic equations in a symmetric convex domain. AIMS Mathematics, 2022, 7(6): 10860-10866. doi: 10.3934/math.2022607 |
| [3] | Zihan Yue, Wei Jiang, Boying Wu, Biao Zhang . A meshless method based on the Laplace transform for multi-term time-space fractional diffusion equation. AIMS Mathematics, 2024, 9(3): 7040-7062. doi: 10.3934/math.2024343 |
| [4] | Yuhua Long, Sha Li . Results on homoclinic solutions of a partial difference equation involving the mean curvature operator. AIMS Mathematics, 2025, 10(3): 6429-6447. doi: 10.3934/math.2025293 |
| [5] | Song Wang, Xiao-Bao Shu, Linxin Shu . Existence of solutions to a class of damped random impulsive differential equations under Dirichlet boundary value conditions. AIMS Mathematics, 2022, 7(5): 7685-7705. doi: 10.3934/math.2022431 |
| [6] | Amjad Hussain, Muhammad Khubaib Zia, Kottakkaran Sooppy Nisar, Velusamy Vijayakumar, Ilyas Khan . Lie analysis, conserved vectors, nonlinear self-adjoint classification and exact solutions of generalized $ \left(N+1\right) $-dimensional nonlinear Boussinesq equation. AIMS Mathematics, 2022, 7(7): 13139-13168. doi: 10.3934/math.2022725 |
| [7] | Xiangmei Li, Kamran, Absar Ul Haq, Xiujun Zhang . Numerical solution of the linear time fractional Klein-Gordon equation using transform based localized RBF method and quadrature. AIMS Mathematics, 2020, 5(5): 5287-5308. doi: 10.3934/math.2020339 |
| [8] | Dilara Siraeva . The invariant solution with blow-up for the gas dynamics equations from one-parameter three-dimensional subalgebra consisting of space, Galilean and pressure translations. AIMS Mathematics, 2024, 9(1): 89-98. doi: 10.3934/math.2024006 |
| [9] | Khalid K. Ali, Mohamed S. Mohamed, Weam G. Alharbi, M. Maneea . Solving the time fractional q-deformed tanh-Gordon equation: A theoretical analysis using controlled Picard's transform method. AIMS Mathematics, 2024, 9(9): 24654-24676. doi: 10.3934/math.20241201 |
| [10] | Zhilin Li, Guoping Chen, Weiwei Long, Xinyuan Pan . Variational approach to p-Laplacian fractional differential equations with instantaneous and non-instantaneous impulses. AIMS Mathematics, 2022, 7(9): 16986-17000. doi: 10.3934/math.2022933 |
Lung cancer is one of the most common cancers that endanger human health seriously, and approximately one million people die from it every year [1,2]. Lung adenocarcinoma (LUAD) is widely regarded as the most common histological subtype of lung cancer [3]. It is reported that the 5-year survival rate for LUAD patients in the United States is below 20% [4]. Over the past few decades, despite the progress has been made in the treatment of LUAD, the overall survival rate of this cancer is not evidently upgraded [5]. Recently, it has been proved that tremendous advance in cancer treatment has been obtained through immunotherapy, and immune checkpoint inhibitors are increasingly used to treat LUAD [6].
Evidence has indicated that targeting programmed cell-death protein 1/programmed cell death 1 ligand 1 (PD-1/PD-L1)-related pathways is an emerging concept in cancer immunotherapy [7]. PD-1 (also called PDCD-1) is expressed on both exhausted and activated T cells, which is an immune checkpoint molecular [8]. PD-1 blockade had remarkable efficacy against extensive advanced malignancies [9], and it also serves critical role in lung cancer. For example, Kagamu et al. [10] observed that many advanced lung cancer patients achieved long-term survival without disease progression after treatment with PD-1 blockade, suggesting that PD-1 blockade therapy could evidently improve the outcome of lung cancer patients. In addition, with the development of gene expression profiles, researchers have discovered several gene expression characteristics related to survival outcomes of various cancers. The constructed prognostic model may be helpful for early diagnosis of LUAD and personalized treatment [11]. For instance, Zhang et al. [12] constructed a risk model with four genes that could predict the overall survival of LUAD patients. Moreover, Sun et al. [13] established an immune prognostic model for LUAD, which could independently identify high-risk patients with poor survival. However, the PD-1 associated risk model that was used for predicting the prognosis of LUAD has not been constructed.
In the present study, we aimed to identify PD-1-related gene signature to predict the survival of LUAD patients. The gene expression data and clinical information related to LUAD were downloaded from the Gene Expression Omnibus (GEO) and The Cancer Genome Atlas (TCGA) databases. After a series of bioinformatics analysis, the genes related to both PD-1 expression level and immune cells were identified, and then Lasso Cox regression analysis was performed to obtain key genes to construct prognostic model. Furthermore, the predictive performance of this model was validated by using the data from an independent GEO and TCGA cohorts. Figure S1 shows the analysis process of the present study. This study will provide new insights into the mechanisms of PD-1 immune checkpoints and provide novel therapeutic targets for LUAD patients with poor prognosis.
The National Center for Biotechnology information (NCBI) GEO database [14] was used to search the gene expression profiles of LUAD with "lung adenocarcinoma" as a keyword. The search criteria were as follows: 1) the species was set to "Homo sapiens" 2) the clinical parameters of the samples in dataset must include prognostic information, 3) the sample size of cancer tissues in the dataset was at least 60. Finally, four datasets (GSE26939, GSE68465, GSE72094, and GSE68571) that met the requirements were obtained and included in this analysis.
Concretely, GSE26939 dataset included 116 LUAD samples, among which, one sample without prognostic information was removed. The platform was GPL9053 Agilent-UNC-custom-4X44K. GSE68465 dataset included 443 LUAD samples and 19 normal lung tissue samples. Of which, 442 samples with complete prognostic information were selected for further analysis. The platform was GPL96 [HG-U133A] Affymetrix Human Genome U133A Array. GSE72094 dataset included 442 LUAD samples. After filtering, 398 samples were enrolled for analysis. The platform was GPL15048 Rosetta/Merck Human RSTA Custom Affymetrix 2.0 microarray [HuRSTA_2a520709.CDF]. In addition, GSE68571 dataset included 86 LUAD samples, which were analyzed based on GPL80 [Hu6800] Affymetrix Human Full Length HuGeneFL Array platform.
For these four datasets, standardized probe expression matrix and the corresponding annotation information were downloaded from GEO database. Then, probes that not matched to gene symbol were removed. If different probes mapped to the same gene symbol, the average value of these probes was applied as the final expression value.
Combined with the expression level of PD-1 in each sample and the corresponding survival information, K-M survival analysis was performed by using R package survminer (Version 0.4.3) to determine the optimal cutpoint. According to the optimal cutpoint, samples were divided into PD-1 low expression group (PD-1 expression level in sample < cutpoint) and high expression group (PD-1 expression level in sample > cutpoint). Next, survival analysis and log-rank test were used to compare the prognostic value of low and high expression groups by using the R package survival (Version 2.42-6).
The distribution of clinical factors in each sample of the PD-1 low and high expression groups was counted and plotted using R package ggstatsplot (Version: 0.6.5). In addition, the difference in continuous variables between two groups was compared using Welch's t-test, and the difference in discrete variables between two groups was calculated using Chi-square test. Clinical factors with p-value < 0.05 were regarded as statistically significant.
R package gene set variation analysis (GSVA, Version: 1.36.2) [15] was applied to calculate the enrichment score of each Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway in each sample to obtain a scoring matrix. Then, differential pathway between PD-1 low and high expression groups was identified by using R package limma [16] (Version 3.10.3). Then, Benjamini- Hochberg method was used to perform multiple test correction and obtained adjusted p-value (adj. p-value). The adj. p-value < 0.05 was considered as the cut-off criterion of significantly enriched pathway.
DEGs between PD-1 high expression and low expression groups were screened by using the Bayesian method provided by limma in R package. Further, Benjamini- Hochberg method was used for the adjustment of p-value. |log fold change (FC)| > 0.585 and adj. p-value < 0.05 were considered as statistically significant.
According to the expression matrix of all genes, variance analysis was adopted to find the TOP2000 genes with higher degree of variation among samples. Then, these genes were analyzed by using R package weighted gene co-expression network analysis (WGCNA) [17] (Versions 1.61 and modules closely related to PD-1 were identified. After that, R package clusterprofiler [18] (Version: 3.8.1) was used to perform KEGG pathway enrichment analysis of genes in these modules with the following parameters: pAdjustMethod = "BH" and pvalueCutoff = "0.05". P-value < 0.05 was the threshold value of significantly enriched results.
Subsequently, we narrowed down the range of candidate genes to obtain closely related genes. In brief, the above obtained DEGs and genes in modules were integrated, and then the overlapping genes were selected for further analysis.
The gene sets used to label the type of each immune cell in tumor microenvironment (TME) were obtained from the study by Charoentong et al. [19]. Based on the ssGSEA algorithm in R package GSVA [20], the enrichment score of each immune cell in samples was calculated, which represented the relative abundance of each TME infiltrating cell in each sample. Then, based on the relative abundance, pheatmap package (Version: 1.0.12) was used to perform unsupervised clustering for samples, and the samples were clarified into high and low infiltrating immune groups. After that, DEGs between these two groups were screened using limma in R package. The threshold values of DEGs were |log FC| > 1 and adj. p-value < 0.05, and these genes were considered to be closely related to the immune microenvironment.
Genes related to PD-1 expression as well as immune cells were integrated, and the shared genes were obtained for further analyses. Next, Gene ontology (GO)-biological process (BP) [21] and KEGG [22] pathway analyses of these genes were performed using DAVID [23] (Version 6.8). Count ≥ 2 and p-value < 0.05 were regarded as statistically significant.
Combined with the expression values of genes and the overall survival time of each sample, univariate Cox regression analysis was conducted. Genes with p-value < 0.05 were regarded as prognosis-related genes and included for the next analysis.
Then, based on these prognosis-related genes, a LASSO Cox regression was applied to identify key genes for the construction of prognostic models by using glmnet package [24] (Version 2.0-18). The 50-fold cross-validation was performed to determine the optimal value of Lasso penalty parameter. After that, multivariate Cox regression analysis of these key genes was conducted to construct prognostic risk model using the following formula: Risk score = βgene1*exprgene1 + βgene2*exprgene2 +... +βgeneN*exprgeneN (β represents correlation coefficient of prognostic gene and expr represents the expression value of the corresponding genes). According to the median value of risk score, the samples were divided into low risk group and high risk group. Next, K-M survival analysis of these two groups was performed using log-rank test.
To determine whether the risks core model could be applied as an independent prognostic factor, the prognostic value of the risk score and other clinical parameters, including age, sex, stage, grade, marked lymphocytes, smoking status, TP53 mutation, EGFR mutation, KRAS mutation, and STK11 mutation, was evaluated with univariate and multivariate Cox regression analyses. Variables with p-value < 0.05 were regarded as independent prognostic factors in both analyses.
We used a microarray datasets GSE72094 to validate the stability and repeatability of prognostic model. At first, the expression value of model genes in each sample of GSE72094 was extracted. According to the prognostic coefficient obtained from the previous analysis, the risk score was calculated using the same formula. Then, the samples were divided into low risk group and high risk group based on the median, followed by K-M survival analysis. In addition to the GEO dataset, we also included data from TCGA to validate the prediction performance of the model. In brief, the RNA-seq data of LUAD and survival information of corresponding samples were downloaded from the UCSC Xena platform (https://toil.xenahubs.net). In the previous study, the samples with number "-01A", that is, the primary tumor samples were selected for subsequent analysis. Finally, a total of 497 samples were obtained. Based on the expression value of seven genes in model, the risk score was calculated using the same formula. Similarly, all samples were assigned to low-risk and high-risk groups according to the median value of risk score, followed by a survival analysis.
Moreover, we compared the predictive ability of the seven-gene prognostic model with other molecular biomarkers such as PD-1 and PD-L1 in the GSE26939 dataset. Meanwhile, a 1-, 3-, and 5-year receiver operating characteristic (ROC) curves were plotted using the R package timeROC [25], and then the area under the curve (AUC) value of the survival ROC curves was calculated.
The clinical information of cases in four datasets is presented in Table S1. As shown in Figure 1, results showed that the survival probability of patients in the PD-1 high expression group was higher than that of low expression group. Despite there was no significant difference in survival time between PD-1 high expression group and low expression group in the GSE68465 and GSE72094 datasets, the overall trend was consistent.
According to the analysis results in the previous step, we observed that the number of samples in the PD-1 low expression group of GSE68571 and GSE72094 was significantly lower than that in the high expression group. In order to reduce the impact of the imbalance of samples in the two groups on the results, we only selected samples in the GSE26939 and GSE68465 datasets for further analyses.
The distribution of different clinical factors in the PD-1 high and low expression groups in these two datasets was calculated, and factors with p-value < 0.05 were obtained. In brief, in GSE26939 dataset, the age of the PD-1 high expression group was evidently lower than that of the low expression group (Figure 2A, p value = 0.004); the high expression group had more TP53 wild type (Figure 2B, p-value = 0.037); and the PD-1 high expression group had more Squamoid subtype (Figure 2C, p-value ≤ 0.001). Meanwhile, in the dataset of GSE68465, the age of patients in the PD-1 high expression group was markedly younger than that of the low expression group (Figure 2D, p-value = 0.007); the high expression group had more T1-T2 samples and less T3-T4 samples than the low expression group (Figure 2E, p-value = 0.031); and PD-1 high expression group had fewer white people than low expression group (Figure 2F, p-value ≤ 0.001).
Afterward, GSVA was performed to explore difference in KEGG pathways between the PD-1 high expression and low expression groups from the GSE26939 as well as GSE68465 datasets. In total, there were 72 KEGG pathways enriched by samples in GSE26939 dataset and 137 KEGG pathways in GSE68465 dataset. The top20 pathways are presented in Figure 3. In the GSE26939 dataset, several pathways related to immunity were observed (Figure 3A), such as immunodeficiency, intestinal immune network for IgA production, and primary immunodeficiency. These pathways were increased in high expression group in comparison with that in low expression group. Meanwhile, in the GSE68465 dataset, 12 pathways were significantly higher in the PD-1 high expression group than in the low expression group, such as type II diabetes mellitus pathway and JAK-STAT-signaling pathway (Figure 3B). Notably, JAK-STAT-signaling pathway was evidently increased in the PD-1 high expression group of two GEO datasets.
Combined with the results of the previous analyses, we found that the expression level of PD-1 in GSE26939 was significantly related to prognosis of LUAD, and samples in this dataset were mainly enriched in immune-related pathways. Therefore, GSE26939 dataset was retained for further analysis. According to the threshold (|log FC| > 0.585 and adj. p-value < 0.05), 822 DEGs (709 up-regulated and 113 down-regulated genes) were screened between PD-1 high expression group and low expression group. As shown in Figure S2, DEGs could distinguish the PD-1 low expression group from the high expression group, which proved initially that these DEGs were associated with the expression level of PD-1.
WGCNA was used to identify PD-1 expression related modules. This analysis was performed by a step-by-step method. First, soft threshold power of 8 was selected to calculate the adjacencies because it was the lowest power at which the scale-free topology fit index reached 0.85 (Figure 4A, B). Second, based on clustering and dynamic tree cut methods, the similar modules were clustered. Modules with correlation coefficient > 0.7 were selected and finally seven modules were obtained (Figure 4C). Next, the relationship between module and clinical characteristics (age, TP53 mutation status, PD-1 high- or low-expressed groups, and PD-1 expression level) was analyzed (Figure 4D). For the expression level of PD-1, magenta module, pink module, red module, and purple module represented evidently positive correlations (r > 0 and p-value < 0.05), while yellow module showed significantly negative correlation (r < 0 and p-value < 0.05). In addition, KEGG pathway enrichment analysis for genes in five modules showed that magenta module, pink module, red module, purple module, and yellow module were respectively enriched in 4, 39, 12, 52, and 8 pathways (Figure S3). For example, genes in magenta module, pink module, red module, purple module, and yellow module were mainly enriched in primary immunodeficiency, cytokine-cytokine receptor interaction, protein digestion and absorption, TNF signaling pathway, and cholesterol metabolism, respectively.
Furthermore, genes in the five models and DEGs in PD-1 low expression vs. high expression groups were intersected, and 388 overlapping genes were obtained (Figure S4), which were regarded as candidate genes for further analysis.
As presented in Figure S5A, the samples were divided into low immune infiltration group and high immune infiltration group according to the enrichment fraction of immune cells in each sample. We also observed that the high immune infiltration group was mainly composed of samples with high expression level of PD-1, and the low immune infiltration group was mainly contained samples with low expression level of PD-1. Furthermore, differential analysis was performed between high immune infiltration and low immune infiltration groups. A total of 467 DEGs were identified, including 445 up-regulated and 22 down-regulated genes.
A total of common 237 genes were obtained via the intersection of 467 immune related genes and 388 PD-1 related genes. Thus, these 237 genes were considered as PD-1 related immune genes. To explore the biological function of these genes, GO_BP terms and KEGG pathways analyses were conducted. Results showed that these 237 genes were significantly enriched in 87 GO-BP terms (such as immune response and inflammatory response) and 38 KEGG pathways (such as cytokine-cytokine receptor interaction and rheumatoid arthritis). The top20 results are presented in Figure S5B.
Based on the 237 PD-1 related immune genes, univariate Cox regression analysis was performed. A total of 47 genes that were significantly associated with overall survival of LUAD were obtained. Next, LASSO Cox regression analysis for these 47 genes was conducted to reveal the optimal modeling parameter (λ), and results showed λ = 7 was considered as optimum parameter (Figure 5A). Finally, a seven-gene (DPT, ITGAD, CLECL1, SYT13, DUSP26, AMPD1, and NELL2) signature was used to construct the prognostic model. Risk score was calculated by using the following formula: (-9.322*exprAMPD1) + (-0.169*exprCLECL1) + (-0.176*exprDPT) + (-0.029*exprDUSP26) + (-0.056*exprITGAD) + (-0.063*exprNELL2) + (0.077*exprSYT13). Based on the median of risk score, samples were divided into low risk group and high risk group. Patients in high risk group had significantly worse survival rate than those in low risk group (Figure 5B). Distributions of risk score, follow-up time, and seven-gene expression profiles are shown in Figure 5C. The expression level of SYT13 was up-regulated in the high risk group, and the expression level of the remaining six genes was down-regulated in the high risk group compared with the low risk group. In addition, K-M survival curve showed that the prognosis of patients with high expression of SYT13 was poor, and the prognosis of patients with low expression of the remaining 6 genes was worse (Figure 5D).
Further, univariate and multivariate Cox regression analyses were performed to evaluate whether risk score could be used as an independent prognostic factor for LUAD. Several clinical factors were included in this analysis, including sex, age, stage, grade, marked_lymphocytes, smoking_status, and genetic mutation status (TP53, EGFR, KRAS, STK11). Univariate Cox analysis indicated that sex, age, and risk groups were significantly associated with prognosis of LUAD (Figure 6A). Next, multivariate analysis revealed that the risk score was an independent prognostic predictor (Figure 6B), suggesting that the constructed model could independently predict the prognosis of LUAD.
According to the seven-gene risk model, the risk score of the 398 LUAD patients from GSE72094 was calculated using the same formula. Then, these patients were classified into low and high risk groups according to the median of risk score. K-M survival curve presented that the prognosis of patients in high risk group was significantly worse than those in low risk group (Figure 7A). The risk scores distribution and follow up of the patients in these two groups were analyzed, indicating that the high risk score group had more mortality cases (Figure 7B). Moreover, the expression level of these seven genes in the two groups was consistent with the results of above analysis (Figure 7C). These findings suggested that the constructed prognostic model was reliable and efficient.
To further determine whether the prognostic model was reliable in different populations, we executed the same formula and determined the cut-off value in the TCGA cohort. Next, patients were divided into low- or high-risk groups. Consistent with the results in the GSE72094, patients with high risk score had a lower overall survival than those with low risk score (Figure 7D). The distributions of risk score of LUAD patients as well as the relationship between risk score and follow up were visualized in Figure 7E. In addition, the expression level of SYT13 was up-regulated in the high risk group, and the level of the remaining six genes was down-regulated in the high risk group (Figure 7F), which was consistent with that in the GEO datasets.
Next, we plotted the ROC curves for prognostic model, PD-1 and PD-L1, and then the AUC values were calculated. Results showed that AUC values for prognostic model, PD-1, and PD-L1 were respectively 0.637, 0.541, and 0.504 (Figure 8A–C), indicating that the predictive ability of prognostic model was better than that of PD-1 and PD-L1. In addition, the accuracy of these three biomarkers for survival of LUAD patients 1, 3, and 5 years were evaluated. The AUC values of the ROC curves for seven-gene, PD-1, and PD-L1 predicted the 1-year overall survival were 0.764, 0.363, and 0.344, respectively (Figure 8D); AUC values of ROC predicted the 3-year survival were 0.671, 0.446, and 0.442, respectively (Figure 8E); and AUC values of ROC predicted the 5-year survival were 0.702, 0.423, and 0.418, respectively (Figure 8F). These results further confirmed that the constructed model outperformed other biomarkers with superior predictive performance.
In the present study, a PD-1 expression-related immune prognostic model was constructed, which included seven genes (DPT, ITGAD, CLECL1, SYT13, DUSP26, AMPD1, and NELL2). All samples were divided into high-risk and low-risk groups based on the risk score for this model. Notably, the survival time of the patients in the high risk group was significantly lower than those in the low risk group. Meanwhile, this prognostic model was verified to be effective and reliable by other GEO dataset and TCGA cohort. In addition, survival analysis suggested that the prognosis of patients with higher expression level of SYT13 was poor, and the survival time of patients with lower expression level of the remaining 6 genes was worse. Our constructed PD-1-related prognostic model may be helpful for early diagnosis and treatment of LUAD in clinical practice.
DPT, a noncollagenous extracellular matrix (ECM) protein, is involved in ECM assembly activities and cellular adhesion. Yamatoji et al. [26] indicated that DPT played a significant role in tumor metastasis and invasiveness of oral squamous cell carcinomas (OSCCs), and might be a therapeutic target for the treatment of OSCCs. Guo et al. [27] suggested that DPT could suppress the proliferation of papillary thyroid carcinoma cells via repression of MYC. Some previous studies also suggested that DPT was differentially expressed in some other cancers, such as hepatocellular carcinoma [28]. ITGAD belongs to the beta-2 integrin family of membrane glycoproteins, and is expressed in the tissue and circulating myeloid leukocytes. Evidences demonstrated that ITGAD was associated with metastatic gastric carcinoma [29] and histiocytic sarcomas [30]. However, there was no previous research describing the association between ITGAD and LUAD or lung cancer. Another gene, SYT13, encodes a member of the large synaptotagmin protein family, and it is considered as a target for the treatment of peritoneal metastasis of gastric cancer [31]. Jahn et al. [32] identified SYT13 as a liver tumor suppressor gene. Meanwhile, silencing of SYT13 promoted apoptosis and suppressed proliferation of colorectal cancer cells and thus inhibited tumor growth [33]. The study of Zhang et al. [34] suggested that SYT13 was involved in the growth, invasion, migration, and apoptosis of LUAD. Together, these studies further highlighted that SYT13 was closely associated with the development of cancer. Recently, DUSP26 has been proved to be a potential therapeutic target for human cancer [35]. Bourgonje et al. [36] indicated that lower expression level of DUSP26 was correlated with poor prognosis of gliomas and overexpression of it in E98 glioblastoma cells could reduce tumorigenicity. In addition, the findings of Yu et al. [37] suggested that DUSP26 might be regarded as an oncogene in anaplastic thyroid cancer. For the NELL2, it encodes a glycoprotein containing several von Willebrand factor C domains and epidermal growth factor (EGF)-like domains. Study in mice has shown that the protein plays a role in the growth and differentiation of nerve cells and tumorigenesis [38]. Nakamura et al. [39] indicated that NELL2 played roles in the inhibition of cell migration of renal cell carcinoma. Meanwhile, Kim et al. [40] confirmed that NELL2 was involved in cancer development by the regulation of E2F transcription factor 1. Taken together, these studies suggested that genes such as DPT, ITGAD, SYT13, DUSP26, and NELL2 were involved in pathogenesis of cancers, and were associated with the prognosis of various cancers. Nevertheless, few studies reported the roles of CLECL1 and AMPD1 in cancer especially in the LUAD. Thus, we speculated that these genes also could predict the prognosis of LUAD.
In this analysis, we also explore the relationship between immune cells and genes. It is reported that altered cell adhesion played promoting roles in immune suppression of cancers [41], and DPT was involved in cellular adhesion [26], suggesting that DPT was associated with immunity in cancers. Integrin adhesion is required for aspects of immune function such as antigen presentation and migration in germinal centers, inflammation sites, and lymph nodes [42]. Synaptotagmin (SYT) is a large class of membrane transporters with 13 members, and is considered to be the main calcium sensor in the process of immune cell exocytosis [43]. Among these, SYT16 is reported to be correlated with immune infiltrates and a prognostic biomarker in glioma [44]. The member of dual-specificity phosphatases (DUSPs) family, such as DUSP1, DUSP2 DUSP4 DUSP10, and DUSP16 are important for immune response and metabolic homeostasis [45]. CLECL1 is involved in the regulating immunity [46], and AMPD1 plays a role in the immune function in sepsis patients [47]. These evidences indicated the indirect or direct relationship between genes in the model and immunity. Therefore, these genes may affect the prognosis of LUAD via regulating the immune microenvironment. However, the specific molecular mechanism still needs further experimental verification.
Importantly, a model composed of seven PD-1 related immune genes (DPT, ITGAD, CLECL1, SYT13, DUSP26, AMPD1, and NELL2) was constructed, which could distinguish LUAD patients and predict prognosis. After calculation of the risk score, results of survival analysis showed that the overall survival of patients in the high-risk group was significantly worse than those in the low-risk group. This model could be used as an independent prognostic factor for LUAD patients. According to our research, the models constructed from these seven PD-1 related genes can well predict the prognosis of patients with LUAD. We think that these seven genes are biomarkers for predicting the prognosis of LUAD patients, and may become new research targets. In short, our study provided a new method for evaluating the prognosis of LUAD patients.
Nevertheless, our results have to be interpreted in light of some limitations. First, the molecular functions of these identified genes in LUAD remain largely unknown, which require further studies in vivo and in vitro. Secondly, all genes in the constructed model were derived from the public databases (GEO and TCGA); therefore, more external samples from medical centers are needed to validate the predictive performance of the risk model.
In summary, a PD-1 related seven-gene (DPT, ITGAD, CLECL1, SYT13, DUSP26, AMPD1, and NELL2) prognostic model was constructed, which could predict the overall survival of LUAD patients well. These seven genes help us to further understand the molecular mechanism of PD-1, and they may be potential biomarkers for the treatment of LUAD. However, the lack of experimental verification is a limitation of the present study, thus further researches are needed to verify the present results.
This study was supported by the Shanghai Minhang District Health Commission Research Project (No. 2020MW05).
The authors declare that they have no conflict of interest.
| [1] |
Nakamura J, Mutlu E, Sharma V, et al. (2014) The Endogenous Exposome. DNA Repair (Amst) 19: 3-13. doi: 10.1016/j.dnarep.2014.03.031
|
| [2] |
Wild CP (2012) The exposome: From concept to utility. Int J Epidemiol 41: 24-32. doi: 10.1093/ije/dyr236
|
| [3] | EPA, Persistent Organic Pollutants: A Global Issue, A Global Response. EPA, 2009. Available from: http://www2.epa.gov/international-cooperation/persistent-organic-pollutants-global-issue-global-response#pops (Accessed July 17, 2014). |
| [4] |
Magliano DJ, Loh VH, Harding JL, et al. (2014) Persistent Organic Pollutants and Diabetes: A Review of the Epidemiological Evidence. Diabetes Meab 40: 1-14. doi: 10.1016/j.diabet.2013.09.006
|
| [5] |
Myre M, Imbeault P (2014) Persistent Organic Pllutants Meet Adipose Tissue Hypoxia: Does Cross-talk Contribute to Obesity. Obes Res 15: 19-28. doi: 10.1111/obr.12086
|
| [6] |
Taylor KW, Novak RF, Anderson HA, et al. (2013) Evaluation of the Association Between Persistent Organic pollutants (POPs) and diabetes in Epidemiological Studies: A National Toxicology Program Workshop Review. Environ Health Perspect 121: 774-783. doi: 10.1289/ehp.1205502
|
| [7] |
Gascon M, Morales E, Sunyer J, et al. (2013) Effects of Persistent Organic Pollutants on the Developing Respiratory and Immune Systems: A Systematic Review. Environ Int 52: 51-65. doi: 10.1016/j.envint.2012.11.005
|
| [8] |
Li QQ, Loganath A, Chong YS, et al. (2006) Persistent Organic Pollutants and Adverse Health Effects in Humans. J Toxicol Env Health Part A 69: 1987-2005. doi: 10.1080/15287390600751447
|
| [9] |
Zielinski I, Tsui I-C (1995) Cystic Fibrosis: Genotypic and Phenotypic Variations. Ann, Rev Genet 29: 777-807. doi: 10.1146/annurev.ge.29.120195.004021
|
| [10] | Braun A, Roth R, McGinnis MJ (2003) Technology challenges in screening single gene disorders. Eur J Pediatr 162: S 13-16. |
| [11] | WHO (2014) World Health Statistics 2014. Geneva, Switzerland: World Health Organization, 180. |
| [12] | Russell JY (2010) A Philosophical Framework for an Open and Critical Transdisciplinary Inquuiry. In: Brown VA, Harris JA, Russell JY, editors. Tackling Wicked Problems Through the Transdisciplinary Imagination. NY, NY: Earthscan, 31-60. |
| [13] | Stokols D, Hall KL, Vogel AL (2013) Transdisciplinary Public Health: Definitions, Core Characteristics, and Strategies for Success. In: Haire-Joshu D, McBride TD, editors. Transdisciplinary Public Health: Research, Education, and Practice. San Francisco, CA: Jossey-Bass, 1-30. |
| [14] |
Jacobs JA, Frickel S (2009) Interdisciplinarity: A Critical Assessment. Annu Rev Sociol 35: 43-65. doi: 10.1146/annurev-soc-070308-115954
|
| [15] | Klein JT (2008) Evaluation of Interdisciplinary and Transdisciplinary Research - A Literature Review. Amer J Prev Med Suppl 2 35: S-116-S-123. |
| [16] | Hays TC (2008) The Science of Team Science - Commentary on Measurements of Scientific Readiness. Amer J Prev Med Suppl 2 35: S-193-S195. |
| [17] | DE (1997) Pasteur's Quadrant: Basic Science and Technological Innovation: Brookings Institution, Washington, DC. |
| [18] | NRC (2014) Convergence: Facilitating Transdisciplinary Integration of Life Sciences, Physical Sciences, Engineering and Beyond. Washington, DC: National Academies Press, 137. |
| [19] | Roche M, Glick TF (1944) Thomas Edison and Basic Science, Edison: Myths and Reality. Arbor: 39-50. |
| [20] | Duggan MB (2014) Prevention of childhood malnutrition: immensity of the challenge and variety of strategies. Paediatr Int Child Health 32: 190-203. |
| [21] |
Blackburn GL (2001) Pasteur's Quadrant and Malnutrition. Nature 409: 397-401. doi: 10.1038/35053187
|
| [22] |
Verbeke W (2005) Consumer acceptance of functional foods: Socio-demographic, cognitive and attitudinal determinants. Food Quality and Preference 16: 45-57. doi: 10.1016/j.foodqual.2004.01.001
|
| [23] | Bhutta ZA, Darmstadt GL (2014) A Role for Science Investments in Advancing Newborn Health. Science Translational Medicine 6: 253-258. |
| [24] | Childress DS (1999) Working in Pasteur's Quadrant. J Rehabilitation Res Dev 36: xi-xii. |
| [25] | De Luigi AJ, Cooper RA (2014) Adaptive sports technology and biomechanics: prosthetics. PM R 8 Suppl.: S40-S57. |
| [26] |
Hardbower DM, Peek RM, Jr., Wilson KT (2014) At the Bench: Helicobacter pylori, dysregulated host responses, DNA damage, and gastric cancer. J Leukoc Biol 96: 201-212. doi: 10.1189/jlb.4BT0214-099R
|
| [27] |
Dalal RS, Moss SF (2014) At the Bedside: Helicobacter pylori, dysregulated host responses, DNA damage, and gastric cancer. J Leukoc Biol 96: 213-224. doi: 10.1189/jlb.4BT0214-100R
|
| [28] | Nobel Committee, The Nobel Prize in Physiology or Medicine Prize Announcement. The Nobel Committee for Physiology or Medicine,2005. Available from:http://www.nobelprize.org/nobel_prizes/medicine/laureates/2005/announcement.html. |
| [29] | Konturek SJ, Konturek PC, Brzozowski T, et al. (2005) From nerves and hormones to bacteria in the stomach; Nobel prize for achievements in gastrology during last century. J Physiol Pharmacol 56: 507-536. |
| [30] |
Marshall B (2006) Helicobacter connections. Chem Med Chem 1: 783-802. doi: 10.1002/cmdc.200600153
|
| [31] | Woo KT, Lau YK, Yap HK, et al. (2006) 3rd College of Physicians' lecture--translational research: From bench to bedside and from bedside to bench; incorporating a clinical research journey in IgA nephritis (1976 to 2006). Ann Acad Med Singapore 35: 735-741. |
| [32] |
Hamilton JG, Edwards HM, Khoury MJ, et al. (2014) Cancer screening and genetics: a tale of two paradigms. Cancer Epidemiol Biomarkers Prev 23: 909-916. doi: 10.1158/1055-9965.EPI-13-1016
|
| [33] | Preston R (2000) The Genome Warrier. The New Yorker 76: 66-83. |
| [34] | Tchounwou PB, Toscano WA (2011) Environmental Epidemiology and Human Health-Biomarkers of Disease and Genetic Susceptibility In: Nriagu JO, editor. Encyclopedia of Environmentl Health Sciences. San Diego, CA: Elsevier, 357-366. |
| [35] |
Teo S, Morgan M, Stirling D, et al. (2000) Assessment of the In Vitro and In Vivo Genotoxicity of Thalomid (Thalidomide). Teratog Carcinog Mutagen 20: 301-311. doi: 10.1002/1520-6866(2000)20:5<301::AID-TCM6>3.0.CO;2-2
|
| [36] |
Portin P (2002) Historical Development of the Concept of the Gene. J Med Philosophy 27: 257-286. doi: 10.1076/jmep.27.3.257.2980
|
| [37] | Crick FHC (1958) Central Dogma of Molecular Biology. Symp Soc Exp Biol XII 12: 139-173. |
| [38] |
Crick FHC (1970) Central Dogma of Molecular Biology. Nature (London) 227: 561-563. doi: 10.1038/227561a0
|
| [39] |
Hahn WC, Weinberg RA (2002) Modeling the Molecular Circuitry of Cancer. Nat Rec Cancer 2: 331-341. doi: 10.1038/nrc795
|
| [40] | Hocquette JF (2005) Where are we in genomics? J Physiol Pharmacol Suppl 3: 37-70. |
| [41] |
Sharma S, Kelly TK, Jones PA (2010) Epigenetics in cancer. Carcinogenesis 31: 27-36. doi: 10.1093/carcin/bgp220
|
| [42] |
Long CR, Westhusin ME, Golding MC (2014) Reshaping the transcriptional frontier: epigenetics and somatic cell nuclear transfer. Mol Reprod Dev 81: 183-193. doi: 10.1002/mrd.22271
|
| [43] | Nilsson E, Skinner MK (2014) Definitions and History of Generational Epigeneti Inheritance. In: Tollefsbol TO, editor. Translational Epigenetics: Evidence and Debate. San Diego, CA: Academic Press, 11-24. |
| [44] | Hoile SP, Lillycrop KA, Grenfel LR, et al. (2014) Phenotypic and Epigenetic Inheritance Across Multiple Generations in Mammals through the Female Line. In: Tollefsbol TO, editor. Transgenerational Epigenetics: Evidence and Debate. San Diego, CA: Academic Press, 269-278. |
| [45] | Mashoodh R, Champagne FA (2014) Paternal Epigenetic Inheritance. In: Tollefsbol TO, editor. Transgenerational Epigenetics: Evidence and Debate, 231-236. |
| [46] | Toscano WA, Oehlke KP (2011) Nutrigenomics: A New Frontier in Environmental Health Sciences. In: Nriagu JO, editor. Encyclopedia of Environmental Health Burlington, VT: Elsevier, 202-205. |
| [47] | Toscano WA, Oehlke KP, Kafouri R (2010) An Environmental Systems Biology Approach to the Study of Asthma. In: Pawnkar R, Holgate S, Rosenwasser LJ, editors. Allergy Frontiers: Future Perspectives. Springer: NY, NY, 239-252. |
| [48] | Toscano WAJ, Lee P, Oehlke KP (2008) The Role of Environmental Health Research in Understanding Chronic Conditions. Minn Physician 21: 27-28. |
| [49] |
McLachlan JA (2001) Environmental Signaling: What Embryos and Evolution Teach us about Endocrine Disrupting Chemicals. Endocr Rev 22: 319-341. doi: 10.1210/edrv.22.3.0432
|
| [50] |
Vandenberg LN (2014) Non-monotonic Dose Responses in Studies of Endocrine Disrupting Chemicals: Bisphenol-A as a Case Study. Dose-Response 12: 259-276. doi: 10.2203/dose-response.13-020.Vandenberg
|
| [51] |
Vandenberg LN, Colborn T, Hayes TB, et al. (2012) Hormones and endocrine-disrupting chemicals: low-dose effects and nonmonotonic dose responses. Endocr Rev 33: 378-455. doi: 10.1210/er.2011-1050
|
| [52] | Vandenberg LN, Erlich S, Belcher SM, et al. (2013) Low-dose Effects of Bisphenol-A: An Integrated review of in vitro, Laboratory, Animal, and Epidemiology Studies. Endocrine Disruptors 1: 1-20. |
| [53] | CDC (2014) National Diabetes Statistics Report, 2014. Center for Disease Control, Atlanta, GA. |
| [54] |
Misra A, Bhardwai S (2014) Obesity and the Metbolic Syndrome in Developing Countries: Focus on South Asians. Nestle Nutr Inst Workshop Ser 78: 133-140. doi: 10.1159/000354952
|
| [55] |
Chase K, Sharma RP (2013) Epigenetic Developmental Programs and Adipogenesis - Implications for Psychotropic Induced Obesity. Epigenetics 8: 1133-1140. doi: 10.4161/epi.26027
|
| [56] | Janesick A, Blumberg B (2011) Minireview: PPARγ as the target of obesogens. J Ster Biochem & Mol Biol 127: 4-8. |
| [57] | Janesick A, Blumberg B (2011) The role of environmental obesogens in the obesity epidemic. Obesity before birth: Springer. 383-399. |
| [58] | Kishida K, Shimomura H, Nishizawa N, et al. (2001) Enhancement of the aquaporin adipose gene expression by a peroxisome proliferator-activated receptor gamma. J Biol Chem 276: 48572-48579. |
| [59] |
Grun F, Blumberg B (2009) Minireview: the case for obesogens. Mol Endocrinol 23: 1127-1134. doi: 10.1210/me.2008-0485
|
| [60] |
Kemper MF, Stirone C, N. KD, et al. (2014) Genomic and non-genomic regulation of PGC1 isoforms by estrogen to increase cerebral vascular mitochondrial biogenesis and reactive oxygen species. Eur J Pharmacol 723: 322-329. doi: 10.1016/j.ejphar.2013.11.009
|
| [61] |
Riu A, Grimaldi M, le Maire A, et al. (2011) Peroxisome Proliferator-Activated Receptor γ Is a Target for Halogenated Analogs of Bisphenol A. Environ Health Perspect 119: 1227-1232. doi: 10.1289/ehp.1003328
|
| [62] |
Swedenborg E, Ruegg J, Mäkelä S, et al. (2009) Endocrine disruptive chemicals: mechanisms of action and involvement in metabolic disorders. J Mol Endovrinol 43: 1-10. doi: 10.1677/JME-08-0132
|
| [63] |
Willy PJ, Murray IR, Busch BB, et al. (2004) Regulation of PPARγ coactivator 1α (PGC-1α) signaling by an estrogen-related receptor α (ERRα) ligand. Proc Natl Acad Sci (USA) 101: 8912-8917. doi: 10.1073/pnas.0401420101
|
| [64] | Alaynick WA (2006) Nuclear Receptors, Mitochondria and Lipid Metabolism. Mitochondrion 8: 329-337. |
| [65] |
Koo S-H, Satoh H, Herzig S, et al. (2004) PGC-1 Promotes Insulin Resistance in Liver through PPAR-alpha-dependent Induction of TRB-3. Nature Med 10: 530-534. doi: 10.1038/nm1044
|
| [66] |
Sargis RM, Johnson DN, Choudhury RA, et al. (2010) Environmental endocrine disruptors promote adipogenesis in the 3T3-L1 cell line through glucocorticoid receptor activation. Obesity 18: 1283-1288. doi: 10.1038/oby.2009.419
|
| [67] |
Unterberger C, Stapls KJ, Smallie T, et al. (2008) Role of STAT3 in Glucocorticoid-induced Expression of the Human IL-10 Gene. Mol Immunol 45: 3230-3237. doi: 10.1016/j.molimm.2008.02.020
|
| [68] |
Tanaka T, Nobuiaki Y, Kishimoto T, et al. (1997) Defective Adipocyte Differentiation in Mice Lacking the C/EBP or C/EBP Gene. EMBO J 16: 7432-7443. doi: 10.1093/emboj/16.24.7432
|
| [69] |
Berger J, Tanen M, Elbrecht A, et al. (2001) Peroxisome Proliferator-activated Receptor-γ Ligands Inhibit Adipocyte 11β-Hydroxysteroid Dehydrogenase Type 1 Expression and Activity. J Biol Chem 276: 12629-12635. doi: 10.1074/jbc.M003592200
|
| [70] |
Fukazawa H, Hoshino K, Shiozawa T, et al. (2001) Identification and Quantification of Chlorinated Bisphenol A in Wastewater from Wastepaper Recycling Plants. Chemosphere 44: 973-979. doi: 10.1016/S0045-6535(00)00507-5
|
| [71] |
Guerra P, Eljarrat E, Barcelo D (2010) Simultaneous Determination Hexabromocyclododecane, Tetrabromobisphenol A, and Related Compounds in Sewage Sludge and Sediment tsamples from the Ebro River Basin (Spain). Anal Bianal Chem 397: 2817-2824. doi: 10.1007/s00216-010-3670-3
|
| [72] | Kolpin DW, Furlong ET, Meyer MT, et al. (2002) Pharmaceuticals, Hormones, and Other Organic Wastewater Contaminants in U.S. Streams, 1999-2000: A National Reconnaissance. Environ Sci Technol 36: 1202-1211. |
| [73] |
Korshin G, Kim J, Gan L (2006) Comparative Study of Reactions of Endocrine Eisruptors, Bisphenol-A and Diethylstilbestrol in Electrochemical Treatment and Chlorination. Water Res 40: 1070-1078. doi: 10.1016/j.watres.2006.01.003
|
| [74] |
Gallard H, Leclercq A, Croué J-P (2004) Chlorination of Bisphenol A: Kinetics and By-products Formation. Chemosphere 56: 465-473. doi: 10.1016/j.chemosphere.2004.03.001
|
| [75] | Levian C, Ruiz E, Yang X (2014) The Pathogenesis of Obesity from a Genomic and Systems Biology Perspective. Yale J Biol Med 87: 113-126. |
| [76] |
Zhao S, Iyengar R (2012) Systems pharmacology: Network analysis to identify multiscale mechanisms of drug action. Ann Rev Pharmacol Toxicol 52: 505-521. doi: 10.1146/annurev-pharmtox-010611-134520
|
| [77] |
Toscano WA, Oehlke KP (2005) Systems Biology: New Approaches to Old Environmental Health Problems. Int J Environ Res Public Health 2: 4-9. doi: 10.3390/ijerph2005010004
|
| [78] | Leischow SJ, Best A, Trochim WM, et al. (2008) Systems Thinking to Improve the Public's Health. Amer J Prev Med Suppl 2 35: S-196 - S-203. |
| [79] |
Cox J, Mann M (2011) Quantitative, High-Resolution Proteomics for Data-Driven Systems Biology. Annu Rev Biochem 80: 273-299. doi: 10.1146/annurev-biochem-061308-093216
|
| [80] |
Li H, Deng H (2010) Systems Genetics, Bioinformatics and eQTL Mapping. Genetica 138: 915-924. doi: 10.1007/s10709-010-9480-x
|
| [81] |
Majithia AR, Flannick J, Shahanian P, et al. (2014) Rare variants in PPARG with decreased activity in adipocyte differentiation are associated with increased risk of type 2 diabetes. Proc Natl Acad Sci (USA) 111: 13127-13132. doi: 10.1073/pnas.1410428111
|
| [82] |
Hoffman S (2011) Computational Analysis of High Throughput Sequencing Data. Methods Mol Biol 719: 199-217. doi: 10.1007/978-1-61779-027-0_9
|
| [83] | Dondorp WJ, de Wert GM (2013) The Thousand-Dollar Genome: An Ethical Exploration. Eur J Hum Genet Suppl 1: S6-S-26. |
| [84] | Zhang L, Kim S (2014) Learning Gene Networks under SNP Perturbations using eQTL Datasets. PLoS Comput Biol 10: e1003420. |
| [85] |
Meaburn E, Schulz R (2012) Next generation sequencing in epigenetics: insights and challenges. Semin Cell Dev Biol 23: 192-199. doi: 10.1016/j.semcdb.2011.10.010
|
| [86] |
Pellegrini M, Ferrri R (2012) Epigenetic analysis: ChIP-chip and ChIP-seq. Methods Mol Biol 802: 377-387. doi: 10.1007/978-1-61779-400-1_25
|
| [87] |
Ku CS, Naidoo N, Wu M, et al. (2011) Studying the epigenome using next generation sequencing. J Med Genet 48: 721-730. doi: 10.1136/jmedgenet-2011-100242
|
| [88] |
Patel CJ, Battacharya J, Butte AJ (2010) An Environment-Wide Association Study (EWAS) on Type 2 Diabetes Mellitus. PLoS One 5: e10746. doi: 10.1371/journal.pone.0010746
|
| [89] |
Manolio TA, Brooks LD, Collins F (2008) A HapMap harvest of insights into the genetics of common disease. J Clin Invest 118: 1590-1605. doi: 10.1172/JCI34772
|
| [90] | Zeggini E, Weeden MN, Lindgren CM, et al. Replication of Genome-Wide Association Signals in UK Samples Reveals Risk Loci for Type 2 Diabetes. Science (Washington) 316: 1336-1341. |
| [91] | Ioannidis J, Loy EY, Poulton R, et al. (2009) Researching Nongenetic Determinants of Disease A Comparison and Proposed Unification. Sci Transl Med 1: 7ps8. |
| [92] |
Davis AP, Rosenstein MC, Wiegers TC, et al. (2011) DiseaseComps: a metric that discovers similar diseases based upon common toxicogenomic profiles at CTD. Bioinformation 7: 154-156. doi: 10.6026/97320630007154
|
Figures(3)
William A. Toscano, Jr., Hitakshi Sehgal, Emily Yang, Lindsey Spaude, A. Frank Bettmann. Environment-Gene interaction in common complex diseases: New approaches[J]. AIMS Molecular Science, 2014, 1(4): 126-140. doi: 10.3934/molsci.2014.4.126
DownLoad:
