Pan-cancer prognostic model and immune microenvironment analysis of natural killer cell-related genes
Highlight box
Key findings
• The research developed a pan-cancer prognostic model using 63 genes related to natural killer (NK) cells. By screening the literature, DEPDC1 and ASPM were identified as potential areas for new tumor research.
What is known and what is new?
• NK cells play a significant role in antitumor immunity and are closely related to tumor prognosis and recurrence.
• A pan-cancer prognostic model was constructed with 63 NK cell-related genes and further identified two genes, DEPDC1 and ASPM, which were not reported to be associated with NK cell killing in previous studies.
What is the implication, and what should change now?
• In this paper, a pan-cancer prognostic model was constructed to analyze their role in the immune microenvironment, which may contribute new insights into tumor research.
• Two genes, DEPDC1 and ASPM, were further identified. By reviewing the literature, it is found that the possible associations of ASPM with thymoma and uveal melanoma have not yet been reported. Furthermore, there exists some disagreement about the relationship between DEPDC1 expression and the prognosis of stomach adenocarcinoma, which remains to be further explored in future studies.
Introduction
Natural killer (NK) cells are a type of cytotoxic lymphocytes that are essential for cancer surveillance and can act as effectors without prior sensitization (1). A study has demonstrated that infused allogeneic NK cells can safely cross the human leukocyte antigen barrier and avoid graft-versus-host disease reaction (2).
NK cell-based tumor immunotherapy, including NK cell-based immune checkpoint inhibition and CAR-engineered NK cells. Some researchers have proposed utilizing NK cells as novel targets for immune checkpoint inhibition, suggesting that the combination of anti-programmed death 1 (PD-1), anti-programmed death-ligand 1 (PD-L1) inhibitors with NK cell-specific checkpoint inhibitors may hold significant value for combination immune checkpoint therapy (3). Besides, drawing inspiration from the CAR-T immunotherapy, researchers have extended their focus to other immune cells, including CAR-NK, CAR-CIK, and CAR-MΦ. Among these, CAR-NK cells exhibit several advantages over CAR-T cells, including better safety, superior antitumor activity, and high efficiency for ‘off-the-shelf’ manufacturing (4-6).
A recent study demonstrates utilization of the nanotechnology in NK cell-based tumor immunotherapy. However, the full realization of engineered NK cells’ potential in clinical practice has been hindered by the absence of suitable models to comprehensively study human NK cell biology complexity (7). Additionally, it has been emphasized that to maximize patient benefits from immunotherapy, personalized analysis of cancer based on biomarkers are of paramount importance (3). Hence, there is an urgent need to develop NK cell-related models or identify associated biomarkers through diverse approaches. In this study, we construct a pan-cancer prognostic model based on 63 NK cell-related genes and screened two key genes, ASPM and DEPDC1, which may provide a new direction for future study to further analyze of the mechanisms underlying NK cell-mediated tumor immunity and lay the foundation for personalized drug development. We present this article in accordance with the TRIPOD reporting checklist (available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-434/rc).
Methods
Data download and preprocessing
The CRISPR/Cas9 library data related to NK cell killing were obtained from the literature (8-10), and the intersection was obtained via Venn diagram analysis (11). RNA sequencing (RNA-seq) [fragments per kilobase per million (FPKM) value] and the clinical data of 32 cancers from The Cancer Genome Atlas (TCGA) database were downloaded from the USCS Xena database (https://xena.ucsc.edu/), RNA-seq data (FPKM value) of normal tissues were downloaded from the GTEx database (https://gtexportal.org/home/), and pan-cancer prognostic data were downloaded from the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/gds). The study was conducted in accordance with the Declaration of Helsinki (as revised in 2013). The normalizeBetweenArrays algorithm in the “limma” R package (The R Foundation for Statistical Computing) (12) was used to correct the data.
Differentially expressed NK cells associated with tumor prognosis
We first extracted information on the messenger RNA (mRNA) expressions of NK cell-related genes from the RNA-seq data obtained from the GTEx and TCGA databases, analyzed the differences in these genes between the tumor and normal samples using the Wilcoxon test, and determined the common differentially expressed NK cell-related genes (DENKGs) in different tumors using the “RobustRankAggreg” R package (13). Univariate Cox analysis was employed to screen the prognosis-related DENKGs for solid tumors.
Prognostic model of NK cells in solid tumors
We used the “glmnet” R package (14) to screen the overfitted prognosis-related DENKGs via least absolute shrinkage and selection operator (LASSO) regression, and the prognostic model of NK cells in solid tumors was constructed via multivariate Cox analysis. Furthermore, the data were randomly categorized into training and test groups at a ratio of 8:2 for model verification. Additionally, six different cancer datasets from GEO were used to further validate our prognostic model.
Nomogram model
The “survival” R package (15) was used to conduct survival curve and receiver operating characteristic (ROC) analyses of the model, and independent prognostic analysis was employed to verify the effectiveness of the model. In addition, a nomogram model of solid tumor prognosis was constructed, and the nomogram calibration curve was used to test its prognostic effect. ROC and concordance index (C-index) were also used to analyze the accuracy of the nomogram model.
Analysis of the immune microenvironment
TCGA pan-cancer samples were classified into high- and low-risk groups according to the risk model. The “GSEA” R package (16) was used to analyze the enriched immune pathways in the high- and low-risk groups so as to verify the enrichment effect of the immune pathways on the model. The proportion of immune cells in pan-cancer samples was analyzed using the “CIBERSORT” R package (17), and the difference of immune cells between the high- and low-risk groups was analyzed.
Analysis of tumor mutational burden (TMB)
To characterize the differences in tumor mutation between the high- and low-risk groups, we downloaded pan-cancer mutation data from UCSC Xena, calculated the TMB using the Perl script, and categorized the patients into two groups. Furthermore, we separately analyzed the differences in the TMB value between the high- and low-risk groups.
Gene analysis of the NK cell prognostic model
First, the interactions between the genes in the model were analyzed using the Search Tool for the Retrieval of Interacting Genes/Proteins (STRING) (18) database (https://cn.string-db.org/), and those with binding scores >0.7 were selected as the core genes. Second, the relationship between core genes and immune microenvironment and stem cell scoring was further analyzed using the “CIBERSORT” R package (19). Prognostic survival analysis was also employed as an important method for core gene screening.
Statistical analysis
Data are presented as means ± standard error of the mean (SEM). Wilcoxon test was applied to analyze the differences of NK cell-related genes between the tumor and normal samples according to the TCGA and GTEx database. Statistical analyses were performed using R 4.1.2. P<0.05 was considered statistically significant. A Perl script was used to calculate the TMB value.
Results
Differential and prognostic NK cell-related genes
The design and process of our study are presented in Figure 1. Five CRISPR/Cas9 library results were obtained from the studies conducted by Kearney et al. (8), Freeman et al. (9), and Sheffer et al. (10). A total of 771 NK cell-related genes were obtained via Venn diagram analysis performed to determine the intersection of ≥ two datasets, which were found to be valid (Figure 2A,2B). After the combination of the TCGA and GTEx data, the differences in the 730 NK cell-related genes between tumor and control samples were tested. A total of 184 DENKGs were screened using the “RobustRankAggreg” R package (13) (Figure 2C). Univariate Cox analysis was employed to screen 136 DENKGs related to pan-cancer prognosis (Table 1).
Table 1
Gene | HR | HR 95L | HR 95H | P value |
---|---|---|---|---|
ASPM | 1.765335 | 1.64416 | 1.89544 | 2.63E−55 |
CRY2 | 0.300167 | 0.25784 | 0.34944 | 2.66E−54 |
DEPDC1 | 1.726787 | 1.60937 | 1.85277 | 3.33E−52 |
SNX1 | 0.176422 | 0.14077 | 0.2211 | 2.79E−51 |
TMEM158 | 1.614313 | 1.51535 | 1.71974 | 8.43E−50 |
GTSE1 | 1.747023 | 1.62188 | 1.88182 | 5.37E−49 |
LRRC27 | 0.391062 | 0.34329 | 0.44548 | 2.74E−45 |
MSANTD3 | 3.782003 | 3.13263 | 4.56599 | 1.46E−43 |
E2F8 | 1.629834 | 1.52008 | 1.74751 | 6.44E−43 |
ZIC2 | 1.396849 | 1.32864 | 1.46856 | 4.03E−39 |
B4GALT5 | 3.647815 | 3.00407 | 4.42951 | 5.31E−39 |
CFAP69 | 0.448298 | 0.39738 | 0.50574 | 6.93E−39 |
SRXN1 | 1.77672 | 1.62956 | 1.93717 | 8.34E−39 |
HERC1 | 0.334825 | 0.28332 | 0.39569 | 9.93E−38 |
SLC35C1 | 2.903049 | 2.46438 | 3.4198 | 3.08E−37 |
CRBN | 0.297112 | 0.24494 | 0.3604 | 7.18E−35 |
KANK3 | 0.558185 | 0.50789 | 0.61347 | 1.04E−33 |
TCF19 | 1.947962 | 1.73628 | 2.18545 | 6.58E−30 |
WDR20 | 0.264804 | 0.21015 | 0.33367 | 1.92E−29 |
NRBF2 | 4.169672 | 3.21529 | 5.40734 | 4.94E−27 |
CDC7 | 1.696141 | 1.54023 | 1.86783 | 6.61E−27 |
SYNE2 | 0.525487 | 0.46596 | 0.59263 | 9.76E−26 |
LMF1 | 0.561003 | 0.50312 | 0.62555 | 2.39E−25 |
PNO1 | 3.388568 | 2.65755 | 4.32067 | 7.32E−23 |
TMEM132A | 1.576354 | 1.43967 | 1.72601 | 7.96E−23 |
ANXA6 | 0.541267 | 0.47693 | 0.61429 | 1.97E−21 |
ERCC5 | 0.448042 | 0.37941 | 0.52908 | 2.94E−21 |
CCDC71 | 0.341456 | 0.27282 | 0.42736 | 6.29E−21 |
NMRK1 | 0.537575 | 0.47142 | 0.61301 | 1.96E−20 |
LAMB3 | 1.303835 | 1.22995 | 1.38216 | 4.94E−19 |
RAC1 | 6.708947 | 4.39132 | 10.2498 | 1.34E−18 |
ATOH8 | 0.759045 | 0.71341 | 0.8076 | 2.92E−18 |
DNTTIP1 | 2.82242 | 2.23025 | 3.57182 | 5.81E−18 |
BRIP1 | 1.439142 | 1.32391 | 1.56441 | 1.24E−17 |
ARF3 | 0.240691 | 0.17233 | 0.33618 | 6.58E−17 |
DLG2 | 0.63909 | 0.57383 | 0.71176 | 3.72E−16 |
LUC7L2 | 0.311179 | 0.2349 | 0.41223 | 4.07E−16 |
THBS3 | 1.748406 | 1.52642 | 2.00267 | 7.32E−16 |
PDK4 | 0.788456 | 0.74393 | 0.83565 | 1.11E−15 |
TRAF7 | 2.762126 | 2.15155 | 3.54598 | 1.57E−15 |
PPIL1 | 2.558938 | 2.03091 | 3.22425 | 1.61E−15 |
IKZF3 | 0.757816 | 0.70736 | 0.81187 | 3.04E−15 |
UACA | 0.617324 | 0.54721 | 0.69642 | 4.43E−15 |
PRR15L | 0.861562 | 0.82994 | 0.89438 | 5.67E−15 |
PEG3 | 0.787961 | 0.74086 | 0.83805 | 3.49E−14 |
CCAR2 | 0.399318 | 0.31483 | 0.50648 | 3.77E−14 |
TSEN15 | 2.341763 | 1.86901 | 2.9341 | 1.41E−13 |
CMYA5 | 0.741102 | 0.68439 | 0.80252 | 1.63E−13 |
ZNF331 | 0.674646 | 0.6074 | 0.74934 | 2.03E−13 |
RELB | 1.63665 | 1.43324 | 1.86893 | 3.45E−13 |
AP1M1 | 1.948639 | 1.62797 | 2.33247 | 3.53E−13 |
DMGDH | 0.779652 | 0.72862 | 0.83426 | 5.73E−13 |
C1QTNF1 | 1.412917 | 1.28567 | 1.55276 | 7.06E−13 |
GTF2H3 | 2.453881 | 1.92031 | 3.13571 | 7.19E−13 |
MRPL15 | 2.507044 | 1.94514 | 3.23128 | 1.26E−12 |
CDK5RAP3 | 0.480827 | 0.39191 | 0.58992 | 2.24E−12 |
PITPNM2 | 1.439583 | 1.29924 | 1.59509 | 3.36E−12 |
BCAN | 1.149548 | 1.10369 | 1.19731 | 1.95E−11 |
EDEM2 | 2.480786 | 1.89159 | 3.2535 | 5.12E−11 |
MOB1A | 2.455448 | 1.87707 | 3.21204 | 5.56E−11 |
EXD3 | 0.660047 | 0.58233 | 0.74813 | 8.02E−11 |
OXTR | 1.264809 | 1.17689 | 1.3593 | 1.65E−10 |
VEGFA | 1.339791 | 1.22426 | 1.46622 | 2.05E−10 |
DNAJB6 | 1.950501 | 1.57684 | 2.41271 | 7.40E−10 |
ITPKC | 1.751919 | 1.46413 | 2.09628 | 9.12E−10 |
SAC3D1 | 1.483356 | 1.30664 | 1.68398 | 1.11E−09 |
POLE3 | 2.436907 | 1.82458 | 3.25474 | 1.61E−09 |
SLC35A2 | 1.980095 | 1.58576 | 2.47248 | 1.65E−09 |
C2orf74 | 0.780065 | 0.719329 | 0.84593 | 1.90E−09 |
CYP51A1 | 1.345769 | 1.220756 | 1.48359 | 2.37E−09 |
FBLN1 | 1.230317 | 1.148977 | 1.31741 | 2.86E−09 |
BRI3BP | 1.506361 | 1.315343 | 1.72512 | 3.18E−09 |
STK4 | 1.891211 | 1.528935 | 2.33933 | 4.27E−09 |
PYGO2 | 0.451677 | 0.341712 | 0.59703 | 2.36E−08 |
EFEMP1 | 1.238168 | 1.148672 | 1.33464 | 2.39E−08 |
IKBKG | 1.491816 | 1.292539 | 1.72182 | 4.56E−08 |
CD34 | 0.752885 | 0.677572 | 0.83657 | 1.30E−07 |
PFKFB4 | 1.30387 | 1.180117 | 1.4406 | 1.84E−07 |
HFM1 | 0.671015 | 0.577273 | 0.77998 | 2.03E−07 |
TREM2 | 1.253687 | 1.151075 | 1.36545 | 2.11E−07 |
SHC2 | 0.836544 | 0.781848 | 0.89507 | 2.30E−07 |
SCAMP5 | 0.793343 | 0.726576 | 0.86624 | 2.45E−07 |
PID1 | 1.168711 | 1.100984 | 1.2406 | 3.08E−07 |
MET | 1.163761 | 1.097253 | 1.2343 | 4.39E−07 |
PDP1 | 0.718749 | 0.632071 | 0.81731 | 4.74E−07 |
DDR2 | 1.211411 | 1.123883 | 1.30576 | 5.38E−07 |
ME3 | 0.813547 | 0.749458 | 0.88312 | 8.26E−07 |
TAP1 | 1.43316 | 1.241208 | 1.6548 | 9.33E−07 |
NDEL1 | 1.685718 | 1.355619 | 2.0962 | 2.65E−06 |
PDCL | 0.53852 | 0.414899 | 0.69897 | 3.29E−06 |
ARHGAP10 | 0.769212 | 0.68814 | 0.85984 | 3.88E−06 |
CYB5A | 0.810169 | 0.738754 | 0.88849 | 7.78E−06 |
HERC2 | 0.667784 | 0.555699 | 0.80248 | 1.65E−05 |
SOX10 | 5.64E−36 | 2.26E−52 | 1.41E−19 | 2.52E−05 |
GP1BB | 5.64E−36 | 2.26E−52 | 1.41E−19 | 2.52E−05 |
DAZAP2 | 0.500619 | 0.358446 | 0.69918 | 4.92E−05 |
TM9SF2 | 0.554394 | 0.415449 | 0.73981 | 6.14E−05 |
MZF1 | 0.767031 | 0.67355 | 0.87349 | 6.34E−05 |
CC2D2A | 0.805539 | 0.723826 | 0.89648 | 7.42E−05 |
SPINT1 | 0.904722 | 0.860319 | 0.95142 | 9.64E−05 |
MYL3 | 0.869839 | 0.810814 | 0.93316 | 0.0001 |
GRB2 | 0.488561 | 0.338547 | 0.70505 | 0.000129 |
CDK10 | 0.752894 | 0.650412 | 0.87152 | 0.000144 |
POLR1B | 1.438897 | 1.19174 | 1.73731 | 0.000154 |
STAT1 | 1.357285 | 1.153264 | 1.5974 | 0.000237 |
TMEM81 | 1.304715 | 1.131917 | 1.50389 | 0.000243 |
CLIP4 | 1.161742 | 1.070908 | 1.26028 | 0.000307 |
MAP3K12 | 1.210209 | 1.090331 | 1.34327 | 0.000337 |
NIF3L1 | 1.671069 | 1.24791 | 2.23772 | 0.000568 |
YDJC | 1.239501 | 1.093703 | 1.40473 | 0.000772 |
LEO1 | 0.667489 | 0.525121 | 0.84845 | 0.000958 |
SON | 0.645861 | 0.491231 | 0.84917 | 0.001743 |
PLIN2 | 1.145417 | 1.049574 | 1.25001 | 0.002325 |
CNBD2 | 0.629512 | 0.463982 | 0.8541 | 0.002948 |
SLC35B3 | 0.692636 | 0.538824 | 0.89035 | 0.004152 |
ABCA8 | 0.901293 | 0.839446 | 0.9677 | 0.004166 |
TMEM209 | 1.355676 | 1.099138 | 1.67209 | 0.004467 |
ASPG | 1.100884 | 1.028925 | 1.17788 | 0.005325 |
NOP10 | 1.617456 | 1.135009 | 2.30497 | 0.007798 |
AKAP12 | 1.101045 | 1.025406 | 1.18226 | 0.008029 |
MBD2 | 0.749007 | 0.604766 | 0.92765 | 0.008095 |
TRAF2 | 0.785976 | 0.657473 | 0.93959 | 0.008193 |
PSMB3 | 1.516803 | 1.110832 | 2.07114 | 0.008759 |
ARHGAP24 | 0.900004 | 0.831823 | 0.97377 | 0.008763 |
IFNGR2 | 1.403803 | 1.082233 | 1.82092 | 0.010609 |
HSPA4 | 0.670607 | 0.48797 | 0.9216 | 0.013767 |
RNF31 | 1.255841 | 1.04685 | 1.50655 | 0.014168 |
UGP2 | 1.346234 | 1.058687 | 1.71188 | 0.015302 |
LTC4S | 0.701113 | 0.524051 | 0.938 | 0.016805 |
SPTBN1 | 1.227038 | 1.034371 | 1.45559 | 0.018888 |
MICAL3 | 1.126359 | 1.01696 | 1.24753 | 0.022455 |
UNG | 1.295687 | 1.036186 | 1.62018 | 0.023105 |
ZBTB12 | 1.123996 | 1.010022 | 1.25083 | 0.032132 |
TMEM87B | 0.860211 | 0.747686 | 0.98967 | 0.03528 |
NAALAD2 | 0.873619 | 0.766974 | 0.99509 | 0.041947 |
SAMD4A | 1.096398 | 1.000643 | 1.20132 | 0.04841 |
HR, hazard ratio; 95L, lower limit of 95% confidence interval; 95H, higher limit of 95% confidence interval.
DENKG prognostic risk model
After the overfitted DENKGs were screened via LASSO regression, multifactor Cox analysis was employed to construct the risk model. Furthermore, 63 prognostic genes (Table S1) were identified via multivariate Cox analysis and were used to construct a pan-cancer prognostic risk model according to the following risk formula: (expressing gene1 × β gene1) + (expressing gene2 × β gene2). In addition, the pan-cancer samples from TCGA were randomly categorized into groups at a ratio of 8:2 to test the effectiveness of the model.
Testing the prognostic risk model for pan-cancer
Survival and ROC analyses revealed that our model could well predict tumor prognosis (Figure 3A,3B), disease progression, and recurrence (Figure 3C,3D). This finding was also verified in the survival analysis of both the training and validation groups (Figure 3E,3F). Independent prognostic analysis was performed on age, sex, cancer stage, and risk score. Multivariate and univariate independent prognostic analyses revealed that our risk score could predict tumor prognosis independently of other clinical data (Figure 4A,4B).
Nomogram model
To further study the prognosis of pan-cancer, we constructed a nomogram model (Figure 4C) and used the nomogram (Figure 4D) and ROC curve (Figure 4E) to verify the effectiveness of the model. The line chart and ROC curve demonstrated that our model was effective, and the C-index (Figure 4F) was >0.7, which was ideal.
Analysis of the immune microenvironment
The enrichment of immune pathways showed that there were differences in most immune pathways between the high- and low-risk groups of our model (Figure 5A). These differences were analyzed using the CIBERSORT algorithm (19). As can be seen from the Figure 5B, our model was different in most immune cells. To further understand the relationship between NK cells and tumor prognosis, we analyzed the survival of resting and activated NK cells for pan-cancer prognosis, and activated NK cells predicted according to the algorithm cells were found to be associated with good tumor prognosis (Figure 5C,5D).
Analysis of TMB
The higher the TMB value is, the worse the tumor prognosis. We obtained pan-cancer mutation data from TCGA database and determined the correlation between pan-cancer TMB and the high- and low-risk groups of the model. Then we testified the accuracy of our model and classified each tumor to analyze the correlation between TMB and risk score and drew a radar map (Figure 5E,5F).
Prognostic core genes of pan-cancer NK cells
A total of 63 prognosis-related DENKGs from the STRING database were included to analyze the protein-protein relationship. A protein-protein interaction (PPI) network map was established and visualized using the Cytoscape software (version 3.9.1) (Figure 6A). Genes with a binding coefficient of ≥0.7 were selected as the core genes and were displayed in a correlation analysis map (Figure 6B).
Analysis of the immune microenvironment and prognosis of core genes
The expression values of the above-mentioned core genes in various tumors were correlated with immune, stromal, and stem cell scores, and the related heat map was created (Figure 6C-6F). The immune and stromal scores were negatively correlated with the expression levels of the core genes (Figure 6C,6D). Conversely, a significant positive correlation existed between the expression of the core genes and RNA stemness scores (RNAss) (Figure 6E). A total of 15 core genes were linked to tumor prognostic data, and the gene survival curves in each tumor were drawn to further screen the pan-cancer related genes that were linked to prognosis in NK cells. Among the 15 core genes, ASPM and DEPDC1 were found to play a key role in the prognosis of 12 and 14 tumor types, respectively (Table 2). It is noteworthy to highlight that high DEPDC1 gene expression was associated with better prognosis of colon adenocarcinoma (COAD), stomach adenocarcinoma (STAD), rectum adenocarcinoma (READ) and thymoma (THYM), whereas low DEPDC1 gene expression was associated with better prognosis in ten other types of cancer (Figure 7, Figure S1). In addition, five different tumor types were collected in GEO database to verify our findings (Figure 8). We also observed that the expression of DEPDC1 varied in its correlation with activated and resting NK cells across different cancer types (Figure S2).
Table 2
Gene | Cancer type | P value |
---|---|---|
HERC1 | KIRC | 1.46E−06 |
LGG | 9.58E−06 | |
LIHC | 0.045917 | |
UVM | 0.003986 | |
DNAJB6 | CESC | 0.018468 |
OV | 0.01384 | |
SKCM | 0.029439 | |
UCEC | 0.00546 | |
SHC2 | ACC | 0.010142 |
CESC | 0.048921 | |
KICH | 0.02189 | |
LGG | 0.000556 | |
OV | 0.005896 | |
SKCM | 0.035547 | |
SOX10 | KIRC | 0.003434 |
LIHC | 0.03873 | |
VEGFA | BLCA | 0.018589 |
CESC | 0.006546 | |
KIRP | 0.000228 | |
LGG | 0.001337 | |
PRAD | 0.026486 | |
SARC | 0.025041 | |
UCEC | 0.039948 | |
MET | ACC | 0.027873 |
LGG | 0.000381 | |
PAAD | 9.42E−05 | |
UVM | 0.016287 | |
POLE3 | ACC | 0.00025 |
ESCA | 0.02957 | |
THYM | 0.028938 | |
ASPM | ACC | 4.00E−09 |
KIRC | 0.000205 | |
KIRP | 0.000217 | |
LGG | 6.97E−08 | |
LIHC | 0.006928 | |
LUAD | 0.015197 | |
MESO | 0.002203 | |
PAAD | 0.043426 | |
PCPG | 0.007278 | |
THYM | 0.016468 | |
UCEC | 0.00157 | |
UVM | 0.047649 | |
E2F8 | BLCA | 0.048864 |
KIRC | 0.002008 | |
KIRP | 0.000346 | |
LGG | 0.000236 | |
LIHC | 1.64E−06 | |
MESO | 0.000717 | |
PAAD | 0.044548 | |
PCPG | 0.032873 | |
STAD | 0.003008 | |
THYM | 0.019958 | |
DEPDC1 | ACC | 2.32E−05 |
COAD | 0.049063 | |
KIRC | 0.015496 | |
KIRP | 0.000238 | |
LGG | 6.99E−08 | |
LIHC | 0.000543 | |
LUAD | 0.012019 | |
MESO | 2.40E−05 | |
PAAD | 0.010571 | |
PCPG | 0.006069 | |
READ | 0.027713 | |
STAD | 0.027779 | |
THYM | 0.007332 | |
UVM | 0.039939 | |
BRIP1 | ACC | 0.020404 |
COAD | 0.017831 | |
KIRP | 0.014144 | |
LGG | 5.95E−05 | |
LUAD | 0.014512 | |
MESO | 1.74E−06 | |
PAAD | 0.014325 | |
READ | 0.00895 | |
THYM | 0.019396 | |
SLC35C1 | COAD | 0.0358 |
PCPG | 0.030677 | |
STAD | 0.045321 | |
SLC35B3 | READ | 0.002272 |
SARC | 0.003276 | |
THCA | 0.027777 | |
THYM | 0.024145 | |
HSPA4 | LUAD | 0.046428 |
SARC | 0.030118 | |
SKCM | 0.011044 | |
STAD | 0.047031 | |
RAC1 | ACC | 0.008144 |
DLBC | 0.021017 | |
GBM | 0.000632 | |
KIRC | 0.002018 | |
LGG | 0.006966 | |
LIHC | 0.001216 | |
MESO | 0.002586 | |
PCPG | 0.02042 | |
SKCM | 0.041701 | |
UVM | 0.012914 |
ACC, adrenocortical carcinoma; BLCA, bladder urothelial carcinoma; BRCA, breast invasive carcinoma; CESC, cervical squamous cell carcinoma and endocervical adenocarcinoma; COAD, colon adenocarcinoma; DLBC, lymphoid neoplasm diffuse large B-cell lymphoma; ESCA, esophageal carcinoma; GBM, glioblastoma multiforme; KICH, kidney chromophobe; KIRC, kidney renal clear cell carcinoma; KIRP, kidney renal papillary cell carcinoma; LAML, acute myeloid leukemia; LGG, brain lower grade glioma; LIHC, liver hepatocellular carcinoma; LUAD, lung adenocarcinoma; MESO, mesothelioma; OV, ovarian serous cystadenocarcinoma; PAAD, pancreatic adenocarcinoma; PCPG, pheochromocytoma and paraganglioma; PRAD, prostate adenocarcinoma; READ, rectum adenocarcinoma; SARC sarcoma; SKCM, skin cutaneous melanoma; STAD, stomach adenocarcinoma; THCA, thyroid carcinoma; THYM, thymoma; UCEC, uterine corpus endometrial carcinoma; UVM, uveal melanoma.
Discussion
PD-L1 inhibitors have been approved by the US Food and Drug Administration (FDA) for the treatment of melanoma, lung cancer, and other diseases, and their application in immunotherapy has made a significant progress (20); however, immunotherapy is not effective in all tumors (21), and this limitation has spurred research into identifying the reasons underlying this lack of efficacy and into developing new therapeutic approaches.
NK cells are innate immune-related lymphocytes, which also play a particularly significant role in antitumor immunity (22). NK cell-based tumor immunotherapy, including immune checkpoint inhibition and CAR-engineered NK cells, is a promising area of research. However, there is a need for better NK cell-related models and associated biomarkers. Thus, we conducted a pan-cancer analysis focusing on the role of NK cell-related genes in pan-cancer. However, in comparison with other studies, the NK cell-related genes used in the present study were not obtained from the Kyoto Encyclopedia of Genes and Genomes (KEGG) database. Considering that there are still unknown NK cell-related genes to be mined, we selected the library of CRISPR/Cas9, a gene-editing technology which enables large-scale and in-depth sequencing (23), for analysis.
Because the library was sequenced in vitro, 14,148 samples were used from TCGA and GTEx pan-cancer data to analyze the expressions of the aforementioned genes in the tumor and control samples. After the DENKGs were analyzed, a pan-cancer prognostic model was constructed with 63 NK cell-related genes via univariate and multivariate Cox analyses. Based on the ROC results, the model had a 1-year area under the curve (AUC) value of 0.747 for predicting tumor prognosis, which was higher than that based on tumor staging (1-year AUC =0.673). The C-index of the model was >0.7, thus confirming its value in the study of pan-cancer prognosis. The high- and low-risk grouping of the prognostic model was verified through immune pathway and immune cell analysis.
Given that the range of 63 genes is still too large for researchers to practically examine, PPI analysis was employed to further screen the genes that play a key role in pan-cancer. Subsequently, 15 genes were selected with a binding coefficient of >0.7, and their correlations were analyzed in relation to the pan-cancer immune, matrix, RNA stem cell, and DNA stem cell scores. Among all these 15 genes, the one with the highest immune and matrix scores was BRCT repeats of breast cancer, type 1 (BRIP1), which was negatively correlated with the immune score of most tumors but positively correlated with the stem cell score. The protein encoded by this gene is a member of the RecQ DEAH helicase family and interacts with the BRCT repeats of breast cancer, type 1 (BRCA1). Previous studies have demonstrated that BRCA1 is a tumor-suppressor gene, and its mutations are known to increase susceptibility to many cancers, including breast, ovarian, pancreatic, and prostate cancer (24,25). The results of the present bioinformatics analysis indicated that BRCA1 is highly expressed in most tumors and may be related to tumor stem cells. Furthermore, survival analysis revealed that BRCA1 was associated with the poor prognosis of most tumors, such as pancreatic adenocarcinoma and adrenocortical carcinoma. This provides a new direction for the study of NK cell-related genes in pan-cancer.
Finally, the prognostic ability of the 15 genes was analyzed in pan-cancer, with the most prominent genes being ASPM and DEPDC1, as they were found to play a significant role in the prognosis of 12 and 14 tumor types, respectively.
ASPM (spindle microtubule assembly factor) is a protein-coding gene, which is mainly involved in cell mitosis, cell cycle progression, and DNA damage repair (26,27). Initially, research on ASPM focused on its mutations with autosomal recessive primary microcephaly [MicroCephaly Primary Hereditary (MCPH)], with mutations in ASPM accounting for over 40% of MCPH cases (28,29). In a recent study, Razuvaeva et al. hypothesized that mutations in ASPM inhibit the growth of neural progenitor cells, thereby impeding neurogenesis and leading to MCPH, thus providing a possible explanation for why ASPM mutations are the most commonly mutated genes in MCPH (30). Moreover, ASPM is also closely associated with the occurrence and development of various cancers (31-35). A study indicates that ASPM promotes the proliferation, migration, invasion, and stemness of malignant tumors via the WNT/β-catenin signaling pathway; for example, in the case of prostate cancer, ASPM maintains a subpopulation of prostate cancer stem cells by increasing the protein stability of disheveled-3 (Dvl-3), the cardinal upstream regulator of the canonical Wnt signaling pathway (36). Moreover, Tsai et al. suggest the clinical utility of APSM as a prognostic biomarker for cancer and propose viable molecular targeting and synthetic lethal approaches to leverage its therapeutic potential (27). The results of our bioinformatic analysis confirmed that ASPM is an oncogene that is upregulated in most tumors, and our study suggested, for the first time, that ASPM plays a significant role in the pan cancer immune microenvironment. Moreover, the possible associations of ASPM with THYM and uveal melanoma have not yet been reported. Further experiments should be conducted to confirm this result.
DEPDC1 is a DEP domain protein-coding gene containing 1, which is closely associated with poor prognosis in various malignant tumors, such as breast cancer, bladder cancer, osteosarcoma, and oral squamous cell carcinoma (37-40). Huang et al., through the analysis of glycolysis-related genes, confirmed that DEPDC1 promotes the malignant progression of oral squamous cell carcinoma through the WNT/β-catenin signaling pathway and suggest that DEPDC1 may be a novel biomarker and therapeutic target for oral squamous cell carcinoma (40). In our survival analysis, DEPDC1 was associated with the poor prognosis of most tumors, including adrenocortical carcinoma, kidney renal clear cell carcinoma, kidney renal papillary cell carcinoma, and liver hepatocellular carcinoma (HCC), but interestingly, the presence of DEPDC1 in COAD, STAD, and rectal adenocarcinoma was associated with a good prognosis. The good prognosis of STAD is consistent with another bioinformatics analysis study, in which a high level of DEPDC1 expression was associated with a good progression-free interval in cases of STAD (41). However, another study revealed that a higher expression of DEPDC1 was associated with poor prognosis in STAD, and further experimentation is needed to confirm whether the expression of DEPDC1 is correlated with tumor metastasis and differentiation (42). This inconsistency is likely due to differences in data processing and analytical tools; nevertheless, additional studies should be conducted to further clarify the relationship between DEPDC1 expression and the prognosis of STAD.
As identified the most prominent NK cell-related genes in this research, both ASPM and DEPDC1 can promote the malignant progression of cancers through the WNT/β-catenin signaling pathway, which plays a significant role in various physiological processes such as cell proliferation, differentiation, migration (43). Increasing research has revealed the correlation between dysregulation of the Wnt/β-catenin signaling pathway and the development and progression of tumors, such as colorectal cancer, melanoma, and leukemia (44-46). Through a comprehensive literature review, a significant correlation was unveiled between NK cells and the WNT/β-catenin signaling pathway. Emerging evidence highlights the participation of the Wnt/β-catenin signaling pathway in the development and differentiation of NK cells (47,48). For instance, one study demonstrated that the introduction of DKK1, a natural inhibitor of β-catenin-dependent Wnt signal, results in diminished NK cell counts (49). Nevertheless, there exists some inconsistency regarding the impact of Wnt/β-catenin signaling pathway inhibition on NK cell activation and cytotoxicity. In one study, Xiao et al. proposed that the suppression of NK cell activation mediated by DKK2, also a natural inhibitor of Wnt signal, may be independent of the Wnt/β-catenin signaling pathway (50). However, in gastrointestinal tumors, particularly HCC and gastric cancer (GC), ISG12a has been demonstrated to suppress Wnt/β-catenin signaling pathway, thereby downregulating PD-L1 expression and rendering cancer cells sensitive to NK cell-mediated death (51). Given the limited literature on the association between Wnt/β-catenin signaling and NK cells, this controversy deserves more attention and further exploration in the future.
Conclusions
This study investigated the role of NK cell-related genes in pan-cancer and constructed a prognostic model with 63 NK cell-related genes. Survival and ROC analyses employed prove the effectiveness of the model. In addition, the roles of the 63 NK cell-related genes in cancer were analyzed, and two significant genes were identified—DEPDC1 and ASPM—that may offer a potential direction in tumor immune research.
We also further discovered that DEPDC1 is variably related to the prognosis of 14 kinds of cancer; among these, the association between DEPDC1 expression and the prognosis of STAD remains to be further explored. The association between ASPM expression and the poor prognosis of THYM and uveal melanoma has been characterized, which have not been examined in previous bioinformatics analyses. Although our research still has certain limitations, including missing clinical cohort data and a lack of experimental verification to evaluate the analysis results, our findings potentially open new avenues of research in this field.
Acknowledgments
Funding: None.
Footnote
Reporting Checklist: The authors have completed the TRIPOD reporting checklist. Available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-434/rc
Peer Review File: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-434/prf
Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-434/coif). The authors have no conflicts of interest to declare.
Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. Besides, the study was conducted in accordance with the Declaration of Helsinki (as revised in 2013).
Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.
References
- Romee R, Foley B, Lenvik T, et al. NK cell CD16 surface expression and function is regulated by a disintegrin and metalloprotease-17 (ADAM17). Blood 2013;121:3599-608. [Crossref] [PubMed]
- Yu J, Venstrom JM, Liu XR, et al. Breaking tolerance to self, circulating natural killer cells expressing inhibitory KIR for non-self HLA exhibit effector function after T cell-depleted allogeneic hematopoietic cell transplantation. Blood 2009;113:3875-84. [Crossref] [PubMed]
- Khan M, Arooj S, Wang H NK. Cell-Based Immune Checkpoint Inhibition. Front Immunol 2020;11:167. [Crossref] [PubMed]
- Xie G, Dong H, Liang Y, et al. CAR-NK cells: A promising cellular immunotherapy for cancer. EBioMedicine 2020;59:102975. [Crossref] [PubMed]
- Marofi F, Saleh MM, Rahman HS, et al. CAR-engineered NK cells; a promising therapeutic option for treatment of hematological malignancies. Stem Cell Res Ther 2021;12:374. [Crossref] [PubMed]
- Chen Z, Hu Y, Mei H. Advances in CAR-Engineered Immune Cell Generation: Engineering Approaches and Sourcing Strategies. Adv Sci (Weinh) 2023;10:e2303215. [Crossref] [PubMed]
- Zhang H, Yang L, Wang T, et al. NK cell-based tumor immunotherapy. Bioact Mater 2023;31:63-86. [Crossref] [PubMed]
- Kearney CJ, Vervoort SJ, Hogg SJ, et al. Tumor immune evasion arises through loss of TNF sensitivity. Sci Immunol 2018;3:eaar3451. [Crossref] [PubMed]
- Freeman AJ, Vervoort SJ, Ramsbottom KM, et al. Natural Killer Cells Suppress T Cell-Associated Tumor Immune Evasion. Cell Rep 2019;28:2784-2794.e5. [Crossref] [PubMed]
- Sheffer M, Lowry E, Beelen N, et al. Genome-scale screens identify factors regulating tumor cell responses to natural killer cells. Nat Genet 2021;53:1196-206. [Crossref] [PubMed]
- Vyas J, Gryk MR, Schiller MR. VENN, a tool for titrating sequence conservation onto protein structures. Nucleic Acids Res 2009;37:e124. [Crossref] [PubMed]
- Smyth GK. limma: Linear Models for Microarray Data. In: Gentleman R, Carey VJ, Huber W, et al. editors. Bioinformatics and Computational Biology Solutions Using R and Bioconductor. Statistics for Biology and Health. New York, NY: Springer; 2005:397-420.
- Kolde R, Laur S. RobustRankAggreg: Methods for robust rank aggregation. Bioinformatics 2012;28:573-80. [Crossref] [PubMed]
- Hastie T, Tibshirani R, Tibshirani R. Best subset, forward stepwise or lasso? Analysis and recommendations based on extensive comparisons. Statist Sci 2020;35:579-92. [Crossref]
- Therneau T. A package for survival analysis in S. R package version 2015;2.
- Barbie DA, Tamayo P, Boehm JS, et al. Systematic RNA interference reveals that oncogenic KRAS-driven cancers require TBK1. Nature 2009;462:108-12. [Crossref] [PubMed]
- Newman AM, Steen CB, Liu CL, et al. Determining cell type abundance and expression from bulk tissues with digital cytometry. Nat Biotechnol 2019;37:773-82. [Crossref] [PubMed]
- von Mering C, Huynen M, Jaeggi D, et al. STRING: a database of predicted functional associations between proteins. Nucleic Acids Res 2003;31:258-61. [Crossref] [PubMed]
- Chen B, Khodadoust MS, Liu CL, et al. Profiling Tumor Infiltrating Immune Cells with CIBERSORT. Methods Mol Biol 2018;1711:243-59. [Crossref] [PubMed]
- Azad NS, Gray RJ, Overman MJ, et al. Nivolumab Is Effective in Mismatch Repair-Deficient Noncolorectal Cancers: Results From Arm Z1D-A Subprotocol of the NCI-MATCH (EAY131) Study. J Clin Oncol 2020;38:214-22. [Crossref] [PubMed]
- Marigo I, Zilio S, Desantis G, et al. T Cell Cancer Therapy Requires CD40-CD40L Activation of Tumor Necrosis Factor and Inducible Nitric-Oxide-Synthase-Producing Dendritic Cells. Cancer Cell 2016;30:377-90. [Crossref] [PubMed]
- Horowitz A, Strauss-Albee DM, Leipold M, et al. Genetic and environmental determinants of human NK cell diversity revealed by mass cytometry. Sci Transl Med 2013;5:208ra145. [Crossref] [PubMed]
- Shechner DM, Hacisuleyman E, Younger ST, et al. Multiplexable, locus-specific targeting of long RNAs with CRISPR-Display. Nat Methods 2015;12:664-70. [Crossref] [PubMed]
- Volcic M, Karl S, Baumann B, et al. NF-κB regulates DNA double-strand break repair in conjunction with BRCA1-CtIP complexes. Nucleic Acids Res 2012;40:181-95. [Crossref] [PubMed]
- Gupta R, Liu AY, Glazer PM, et al. LKB1 preserves genome integrity by stimulating BRCA1 expression. Nucleic Acids Res 2015;43:259-71. [Crossref] [PubMed]
- Jiang K, Rezabkova L, Hua S, et al. Microtubule minus-end regulation at spindle poles by an ASPM-katanin complex. Nat Cell Biol 2017;19:480-92. [Crossref] [PubMed]
- Tsai KK, Bae BI, Hsu CC, et al. Oncogenic ASPM Is a Regulatory Hub of Developmental and Stemness Signaling in Cancers. Cancer Res 2023;83:2993-3000. [Crossref] [PubMed]
- Kouprina N, Pavlicek A, Collins NK, et al. The microcephaly ASPM gene is expressed in proliferating tissues and encodes for a mitotic spindle protein. Hum Mol Genet 2005;14:2155-65. [Crossref] [PubMed]
- Létard P, Drunat S, Vial Y, et al. Autosomal recessive primary microcephaly due to ASPM mutations: An update. Hum Mutat 2018;39:319-32. [Crossref] [PubMed]
- Razuvaeva AV, Graziadio L, Palumbo V, et al. The Multiple Mitotic Roles of the ASPM Orthologous Proteins: Insight into the Etiology of ASPM-Dependent Microcephaly. Cells 2023;12:922. [Crossref] [PubMed]
- Li N, Chu J, Hu K, et al. ASPM overexpression enhances cellular proliferation and migration and predicts worse prognosis for papillary renal cell carcinoma. J Biosci 2023;48:17. [Crossref] [PubMed]
- Jiang L, Zhang S, An N, et al. ASPM Promotes the Progression of Anaplastic Thyroid Carcinomas by Regulating the Wnt/beta-Catenin Signaling Pathway. Int J Endocrinol 2022;2022:5316102. [Crossref] [PubMed]
- Yuan YJ, Sun Y, Gao R, et al. Abnormal spindle-like microcephaly-associated protein (ASPM) contributes to the progression of Lung Squamous Cell Carcinoma (LSCC) by regulating CDK4. J Cancer 2020;11:5413-23. [Crossref] [PubMed]
- Wang F, Li J, Liu J, et al. Controversial role of the possible oxyntic stem cell marker ASPM in gastric cancer. J Pathol 2017;241:559-61. [Crossref] [PubMed]
- Chen X, Huang L, Yang Y, et al. ASPM promotes glioblastoma growth by regulating G1 restriction point progression and Wnt-β-catenin signaling. Aging (Albany NY) 2020;12:224-41. [Crossref] [PubMed]
- Pai VC, Hsu CC, Chan TS, et al. ASPM promotes prostate cancer stemness and progression by augmenting Wnt-Dvl-3-β-catenin signaling. Oncogene 2019;38:1340-53. Erratum in: Oncogene 2019;38:1354. [Crossref] [PubMed]
- Harada Y, Kanehira M, Fujisawa Y, et al. Cell-permeable peptide DEPDC1-ZNF224 interferes with transcriptional repression and oncogenicity in bladder cancer cells. Cancer Res 2010;70:5829-39. [Crossref] [PubMed]
- Kretschmer C, Sterner-Kock A, Siedentopf F, et al. Identification of early molecular markers for breast cancer. Mol Cancer 2011;10:15. [Crossref] [PubMed]
- Shen L, Li H, Liu R, et al. DEPDC1 as a crucial factor in the progression of human osteosarcoma. Cancer Med 2023;12:5798-808. [Crossref] [PubMed]
- Huang G, Chen S, Washio J, et al. Glycolysis-Related Gene Analyses Indicate That DEPDC1 Promotes the Malignant Progression of Oral Squamous Cell Carcinoma via the WNT/β-Catenin Signaling Pathway. Int J Mol Sci 2023;24:1992. [Crossref] [PubMed]
- Jia B, Liu J, Hu X, et al. Pan-cancer analysis of DEPDC1 as a candidate prognostic biomarker and associated with immune infiltration. Ann Transl Med 2022;10:1355. [Crossref] [PubMed]
- Gong Z, Chu H, Chen J, et al. DEPDC1 upregulation promotes cell proliferation and predicts poor prognosis in patients with gastric cancer. Cancer Biomark 2021;30:299-307. [Crossref] [PubMed]
- Zhang Y, Wang X. Targeting the Wnt/β-catenin signaling pathway in cancer. J Hematol Oncol 2020;13:165. [Crossref] [PubMed]
- Yuan S, Tao F, Zhang X, et al. Role of Wnt/β-Catenin Signaling in the Chemoresistance Modulation of Colorectal Cancer. Biomed Res Int 2020;2020:9390878. [Crossref] [PubMed]
- Gajos-Michniewicz A, Czyz M. WNT Signaling in Melanoma. Int J Mol Sci 2020;21:4852. [Crossref] [PubMed]
- Wang H, Qiu Y, Zhang H, et al. Histone acetylation by HBO1 (KAT7) activates Wnt/β-catenin signaling to promote leukemogenesis in B-cell acute lymphoblastic leukemia. Cell Death Dis 2023;14:498. [Crossref] [PubMed]
- Bourayou E, Golub R. Signaling Pathways Tuning Innate Lymphoid Cell Response to Hepatocellular Carcinoma. Front Immunol 2022;13:846923. [Crossref] [PubMed]
- Haseeb M, Pirzada RH, Ain QU, et al. Wnt Signaling in the Regulation of Immune Cell and Cancer Therapeutics. Cells 2019;8:1380. [Crossref] [PubMed]
- Grzywacz B, Kataria N, Kataria N, et al. Natural killer-cell differentiation by myeloid progenitors. Blood 2011;117:3548-58. [Crossref] [PubMed]
- Xiao Q, Wu J, Wang WJ, et al. DKK2 imparts tumor immunity evasion through β-catenin-independent suppression of cytotoxic immune-cell activation. Nat Med 2018;24:262-70. [Crossref] [PubMed]
- Deng R, Zuo C, Li Y, et al. The innate immune effector ISG12a promotes cancer immunity by suppressing the canonical Wnt/β-catenin signaling pathway. Cell Mol Immunol 2020;17:1163-79. [Crossref] [PubMed]