Unraveling lung cancer dynamics: a new metabolic signature improving the prediction of recurrence in resected lung adenocarcinoma.

Background
Lung cancer remains a leading cause of cancer-related death worldwide, with tobacco smoking being the main risk factor. The high mortality rate is mainly due to locally advanced or metastatic stage at diagnosis and high relapse rate following therapies with radical or palliative intent [1]. Therefore, efforts to improve prevention, early detection, and treatment strategies remain essential in addressing this major health challenge [2].
Non-Small Cell Lung Cancer (NSCLC) is the most common histological type of lung cancer and includes several distinct subtypes, the most common being adenocarcinoma (LUAD) and squamous cell carcinoma (LUSC). These histotypes are characterized by different molecular and clinical features, even when classified into the same pathological stage [3–5]. The current approach for patients with upfront radically resected early stage LUAD includes the use of adjuvant systemic therapies based on pathological stage and several predictive biomarkers, including molecular alterations such as Epidermal Growth Factor Receptor (EGFR) activating mutations and Anaplastic Lymphoma Kinase (ALK) gene rearrangements, as well as Programmed Death Ligand protein 1 (PD-L1) protein expression [6]. Nonetheless, adjuvant treatment strategies are still far from being really tailored to the individual patient, with many patients still being over- or under-treated, experiencing relapse even after adjuvant therapies. Therefore, there is an urgent need for novel predictive and prognostic biomarkers that can be integrated with those already available [7, 8].
Metabolic reprogramming is necessary to support the relentless proliferation of tumor cells, their migration and growth in metastatic sites [9, 10]. Moreover, it influences the tumor microenvironment, immune milieu, extracellular matrix deposition/remodeling and prepares a permissive soil for colonization of metastatic cells or drug resistance [11]. High expression of key metabolic enzymes is coupled with the upregulation of genetic and epigenetic factors, oncogenes, and proliferative pathways [12].
Several investigators have developed and proposed different metabolic signatures for LUAD, detecting a significant correlation with patients’ prognosis in some series [13–20]. Although these studies were encouraging, to date no metabolic biomarker for LUAD has been introduced in clinical practice, even if tumor metabolism is tightly connected with therapeutic interventions in lung cancer [13–20]. The only cancer metabolic drugs approved for LUAD treatment are several cytotoxic agents classified as “antimetabolites”, including those targeting pyrimidine synthesis [21].
Considering the above background, this study was designed to identify a novel metabolic gene signature using clinical and transcriptomic data together with bioinformatics tools. The potential predictive role of the signature was then evaluated.

Methods

Study design and population
This is an observational, retrospective, multicenter study and male and female patients were considered. Enrolment was conducted at San Luigi Gonzaga and Città della Salute e della Scienza University Hospitals in Orbassano and Turin (Italy), respectively. DEFLeCT patients were part of the prospective observational clinical trial PROMOLE approved by the Ethical Committee of the San Luigi Gonzaga University Hospital (protocol n.14057 approved on September 28, 2018, n.73/2018, version-3 January 30, 2023). MAGA patients were included in the PROFILING protocol (No. 001-IRCC-00IIS-10) of Candiolo Cancer Institute FPO IRCCS, Candiolo, Turin Italy. The informed consent was collected for each patient at the time of the surgical resection. Demographic characteristics, smoking habit, clinical and radiological data were collected, together with data on pathologic staging (according to the VIII edition of the AJCC/UICC TNM staging system) and molecular analysis.
Some analyses were performed on LUAD integrated cohort (LUAD-IC), a dataset reported by Gyorffy et al. [22]. Briefly, the dataset gathered together 2852 lung tumour specimens from 17 independent patient cohorts obtained from NCBI Gene Expression Omnibus (https://www.ncbi.nlm.nih.gov/geo/) and the Genomic Data Commons Data Portal (https://portal.gdc.cancer.gov/). The cohorts included both LUAD (55%) and LUSC (40%) tumors. Most tumors were predominantly stage I (60%), with a significant proportion of stage II samples (32%), while stage III samples were fewer (8%) [22]. For more detailed and comprehensive information, the original study provides a full description of the dataset and its characteristics [22]. To perform our investigation, we restricted the analysis on LUAD histology, which consists of 1308 LUAD samples.

Transcriptomic analysis and exome sequencing
As part of DEFLeCT project, we reanalyzed mRNA-sequencing of the DEFLeCT cohort collected as previously described [23]. RNA samples and mRNA-sequencing data of MAGA cohort were instead both collected and analyzed ex novo. To ensure data consistency and uniformity across all results, we performed the analysis on the new samples using the same methodology as previously applied [23]. Briefly, tissues were isolated at the time of surgical resection and RNA extraction was performed by Maxwell-RSC simply RNA tissue kits (Promega Corporation). RNA samples were purified from DNA by treatment with RNase-Free DNase Set (QIAGEN). The quality of samples was verified by 2100 Bioanalyzer (Agilent Technologies) and Qubit assays (Thermo-Fisher Scientific).
Libraries for RNA sequencing were generated using TruSeq Stranded mRNA Library Prep (Illumina Inc.) following manufacturer’s instructions, using 1 µg of total RNA as input material and 15 PCR cycles for DNA amplification. Each library was analyzed with the DNA High Sensitivity chip using Agilent 2100 Bioanalyzer (Agilent Technologies) and quantified by Qubit 2.0 Fluorometer using ds DNA High Sensitivity Qubit Assay (Thermo-Fisher Scientific). Libraries were pooled together in equimolar amounts and run at the concentration of 1000 pM on the NextSeq1000 sequencer (Illumina Inc.) in 75 nts paired end sequencing mode following manufacturer instruction. Only genes with a value of Transcripts Per Million (TPM) > 2 in at least 80% of patients were considered for the analysis. TPM normalisation has been computed by Python script calculating the reads per kilobase (RPK) ratio for each gene and then scaling the results to the sum of the RPK per sample divided by 10^6. The sum of TPM for each sample, by construction, resulted equal to 10^6”. As described by Gyorffy et al. [22] that performed the analysis, LUAD-IC samples were examined using the in-situ oligonucleotide array platforms GPL96 (Affymetrix Human Genome U133A Array), GPL3921 (Affymetrix HT Human Genome U133A Array) and GPL570 (Affymetrix Human Genome U133 Plus 2.0 Array). Data obtained from these analyses were subjected to two rounds of normalization [22]. For further details, please refer to Gyorffy et al. [22].
Exome sequencing was carried out using the SureSelectXT Human All Exon Kit, following the manufacturer’s guidelines (Agilent Technologies). The sequencing was conducted on the Illumina NovaSeq 6000 (Illumina Inc.) platform, with 35 million and 70 million paired-end reads for normal and tumor samples, respectively. To identify tumor-specific mutations, each tumor sample was compared to the corresponding normal DNA from the normal tissue. After the initial quality check, NGS sequencing of exomes was processed through a bioinformatics pipeline to assess the mutational landscape of LUAD patients. For each sample, TMB was determined by counting the total number of single nucleotide variants (SNVs) per megabase within the sequenced genomic region. In this study, we employed the IDEA pipeline workflow, which enables the identification of SNVs, insertions or deletions (INDELs), and gene copy-number alterations (CNAs). The results of somatic variant calling were compared with those from the established Genome Analysis Toolkit (GATK) Mutect2 and found to be consistent. SNV and INDEL annotations were based on the Catalogue of Somatic Mutations In Cancer (COSMIC) (https://cancer.sanger.ac.uk/cosmic).

Metabolic genes screening and analysis
For metabolic gene enrichment analysis two databases were used: “Metabolomics Workbench tool” of Enrichr software and the Molecular signatures database (MSigDB) of the gene set enrichment analysis (GSEA) software [24–26]. The Gene Ontology biological processes (GO:BP) with a cut-off of p < 0.05 were considered in our analysis.

Methodology for sample stratification
To define the cutpoints for the different metabolic subsets, we applied an unsupervised hierarchical clustering approach using Morpheus software (https://software.broadinstitute.org/morpheus) [27], allowing the data to naturally segregate according to the expression levels of the genes composing the signature. This unbiased clustering procedure identified three stable and biologically meaningful groups that reflected low, intermediate, and high expression patterns. The resulting dendrogram and cluster structure provided clear separation points, which were then used to assign samples to the low-, medium-, and high-expression subsets. This data-driven strategy ensured that the cutpoints were determined objectively based on intrinsic transcriptional variability rather than arbitrary thresholds, thereby preserving the biological relevance of the metabolic subgroups.
The transcriptomic data were also analysed and visualized as Principal Component Analysis (PCA) plots via the web tool ClustVis [28].
To identify DEGs between the subset high and low of patients, we considered a log2 TPM ratio ≥ 1 and performed a 2-sample t test with a p value threshold of 0.05.

Survival analysis
For all cohorts considered, according to the definitions adopted in the AJCC Cancer Staging Manual (Amin MB, Edge SB, Greene FL, Byrd DR, Brookland RK, Washington MK, et al., editors. AJCC Cancer Staging System, Version 9. Chicago, IL: American College of Surgeons; 2025) DFS represents the interval from curative-intent treatment to the first occurrence of disease recurrence (either local or distant).There was no hierarchy among the events, the first that occurred (either local or distant) was considered as an event for the Kaplan-Meier estimation. Comparisons were made with log-rank test and hazard ratio (HR) estimates calculated by means of Cox regression. Multivariate analysis was performed using the Cox proportional hazards model including in the multivariate analysis only covariates with p-value ≤0.1 at the univariate analysis. Whenever indicated, Tests were two-sided and results were reported with 95% confidence intervals (95%CI) or interquartile ranges (IQR).
To perform our investigations on the LUAD-IC cohort, we used the Kaplan Meier plotter platform (https://www.kmplot.com) [22], restricting the analysis on LUAD histology. Depending on the selected clinical parameters, the software analysed only the samples for which the requested clinical information was available, excluding cases with missing data. Therefore, the reported data for the different analyses referred to the subset of samples with available information. The software indicated this number as”Using the selected parameters, the analysis will run on n patients”.
GraphPadPrism 8.0.1 and IBM SPSS v29.0 were used to identify the correlation between the expression of metabolic genes and the clinical outcomes of DEFLeCT and MAGA patients.

Immunohistochemistry (IHC)
All cases were analysed using IHC (n = 11 low and n = 30 high representative samples for LUAD metabolic “subtype”). Briefly, 5 µm-thick serial paraffin sections from representative paraffin blocks were processed using an automated platform (Ventana BenchMark Ultra, Roche) with a primary rabbit anti-TK1 (JF0970) antibody (Thermo-Fisher Scientific). TK1 tumor positivity was measured by two pathologists (L.R. and M.V.) using a semiquantitative histological score (H-score) and calculated as previously described [29].

Flow chart
To provide a clearer and more comprehensive illustration of the study design, analytical workflow and data integration steps, we have included a detailed flow chart summarising the entire methodology (see Fig. 1). Created in BioRender. https://BioRender.com/h37y798

Statistical analysis
The statistical analyses were performed using IBM SPSS v29.0 and GraphPadPrism 8.0.1 software. Descriptive statistics were used for patients’ characteristics. Qualitative variables were compared using the χ2 and Fisher’s exact tests. Appropriate statistical tests were performed as indicated in the Results section, and p < 0.05 was considered to indicate statistical significance in all analysis.

Results

Study population
In the present study we included 87 patients from a larger cohort of early-stage resectable LUAD patients who were enrolled in a prospective observational clinical study, called PROMOLE, within the research project DEFLeCT (Digital Technology For Lung Cancer Treatment) [23, 30]. All patients were ≥18 years of age and had pathological stage I-III LUAD (according to the VIII edition of the AJCC/UICC TNM staging system) [31]. No patients received neo-adjuvant therapies. A summary of patient characteristics is provided in Fig. 2, panel A.
Thirty-one patients (35.6%) were females and the median age at diagnosis was 66.7 years (range 37–81). All patients underwent radical surgery: 78 (89.7%) received lobectomy, 3 (3.5%) bi-lobectomy, and 3 (3.5%) pneumonectomy, while for 3 (3.5%) patients details on the type of surgical procedure were missing. Twenty-two patients (25.3%) were active smokers, 47 (54%) former smokers, and 18 (20.7%) never smokers. Pathologic stage was I in 47 patients (54.1%), II in 20 patients (22.9%), and III in 17 patients (19.5%). Data were missing for 3 (3.5%) patients who were excluded from survival analysis. Maximum standardized uptake (SUVmax) values at the [(18)F] fluorodeoxyglucose positron emission tomography/computed tomography (PET/CT) scan performed at the time of diagnosis were also collected, with a median SUVmax of 7.7 (range 2.1–16.1). Seventeen patients (19.5%) received adjuvant chemotherapy; of these, 4 (23.5%) received cisplatin plus vinorelbine, 5 (29.4%) cisplatin plus gemcitabine, while the remaining 8 patients (47.1%) received different adjuvant treatments within clinical trials. A total of 7 patients (8.1%) received postoperative thoracic radiotherapy.

Identification of metabolism related gene clusters in LUAD
We performed a bulk mRNA-seq analysis of the 87 samples to identify metabolic genes that characterized LUAD cohort. After gene expression quantification and normalization to exclude genes with irrelevant expression and background noise, a total of 5612 genes were considered for the further analysis. To identify metabolic genes from this gene set we interrogated the specific “Metabolomics Workbench tool” of Enrichr software [24, 25] and the Gene Set Enrichment Analysis (GSEA) database [26]. We found 461 genes involved in the metabolism of lipids, carbohydrates, amino acid, nucleotides as highlighted by the analysis of the 100 most significantly enriched Gene Ontology categories related to “Biological Process” (GO:BP) (FDR-adjusted p < 0.05) (Supplementary figure 1, Supplementary data 1).
We performed a prompt investigation for the expression of the 461 metabolic genes in the 87 tumor samples. By using Morpheus clustering web application [27], genes were categorized in seven different metabolic clusters (A-G) through an unsupervised hierarchical clustering (Supplementary Figure 2). A subsequent functional analysis of these 7 clusters using the “Metabolomics Workbench” metadata [25] showed no overlap of metabolic function among the clusters (Fig. 2 panel B).
The most significant GO:BP for each cluster were Amino Sugar Catabolic Process (GO:0046348) for Cluster A, Branched-Chain Amino Acid catabolic Process (GO:0009083) for Cluster B, Steroid Metabolic Process (GO:0008202) for Cluster C, Purine-Containing Compound Biosynthetic Process (GO:0072522) for Cluster D, Phosphatidylinositol Biosynthetic Process (GO:0006661) for Cluster E, Glycerophospholipid Biosynthetic Process (GO:0046474) for Cluster F and Glycolytic Process (GO:0006096) for Cluster G (Fig. 2 panel B, Supplementary data 5).

Identification of a putative metabolic cluster of genes for LUAD stratification
We explored the potential association of the genes within the A-G clusters for Disease-Free Survival (DFS). Considering the mean expression of selected genes, we documented that only the high expression levels of “Cluster G” of genes was significantly associated with worse DFS of DEFLeCT LUAD patients (log-rank p-value = 0.04)(Fig. 3 panel A). We also analyzed a larger LUAD cohort, the “Integrated cohort” (LUAD-IC), encompassing 1308 tumor specimens derived from 12 independent LUAD cohorts as previously described [22]. In LUAD-IC cohort the high expression level of “Cluster G” of genes was significantly associated with worse DFS (log-rank p-value < 0.0001) while the high expression levels of “Cluster B” (log-rank p-value < 0.0001) or “Cluster F” of genes (log-rank p-value = 0.01) were associated with an improved DFS (Fig. 3 panel B).

“Cluster G” of genes identifies different metabolic subsets of LUAD
By means of an unsupervised hierarchical clustering using Morpheus software [27], we categorized the 87 DEFLeCT LUAD tumor samples according to expression levels of “Cluster G” of genes into 3 different subsets: low- (subset low, n = 15 patients), medium- (subset medium, n = 36 patients), and high-expression (subset high, n = 36 patients), as shown in Fig. 4 by the Heat map (panel A) and Principal Component Analysis (PCA) plot (panel B).
We then deepened into the clinical parameters of patients in each metabolic subset (Supplementary data 6). Although the “subset high” was characterized by a higher percentage of smokers, stage-III or high-grade tumors, the differences were not statistically significant among the three metabolic subsets for smoking habits (p-value = 0.283), gender (p-value = 0.554), or TNM stage (stage I vs. stages II-III, p-value = 0.227).
The SUVmax values were categorized into tertiles and we observed that the majority of patients in the metabolic subset-high showed SUVmax values in the higher tertile (SUVmax > 8) (Supplementary data 6). Nonetheless, SUVmax values were not significantly associated with the three metabolic subsets, regardless of whether they were categorized into tertiles (SUVmax < 4, 4–8, >8) or using the median value of 7.7 to differentiate the population into “high” and “low” SUVmax groups (p-values of 0.779 and 0.698, respectively). Furthermore, baseline SUVmax values did not show a significant correlation with the metabolic subsets, even when considering high and low categories only, yielding two-sided p-values of 0.423 and 0.541 for the SUVmax tertiles and median thresholds, respectively.

Mutational profile of LUAD
For 84 out of 87 tumors (96.5%) included in DEFLeCT LUAD cohort, the mutational status was available. We analyzed the distribution of Tumor Mutation Burden (TMB) in the different metabolic subsets and we evidenced that the high expression of “Cluster G” of genes was positively associated with TMB (p-value < 0.0001)(Fig. 4 panel C). Focusing on the distribution of the most important pathological mutated oncogenes characterizing LUAD (i.e. EGFR, KRAS, BRAF, MET) we observed an significant enrichment of KRAS mutant tumors in metabolic subset-high (% of KRAS mutant tumors vs all tumors in subset medium 5%, high 14%, p-value = 0.022)(Fig. 4 panel D).

Correlation between metabolic subsets and DFS in LUAD patients
We investigated the putative correlation between the metabolic subsets and the DFS of DEFLeCT LUAD cohort. After exclusion of 7 patients due to missing relevant data or inadequate follow-up, the survival analysis was performed on the remaining 80 patients. Overall, with a median follow-up of 56.5 months (95% CI: 56.6 -62.3), the median DFS for the entire cohort was 49 months (95% CI: not estimable), with an estimated DFS at 3 years of 70.6%. Among these 80 patients, DFS showed a statistically significant difference according to metabolic subset (p-value = 0.025) (Fig. 5 panel A).
Namely, patients in “subset low” exhibited a significantly better DFS compared to those in “subset high”, with a median DFS that was not reached (NR) versus 32 months (95% CI: 5–67, log-rank test p-value = 0.040) with a HR of 0.156 (95% CI: 0.020–1.201, p-value = 0.074) (Fig. 5 panel A).
As expected, TNM staging demonstrated a significant association with DFS (p-value = 0.017). It is noteworthy that, although with the limitations due to the small sample size, the metabolic signature had a greater impact in stages I-II (Fig. 5 panel B). In contrast, the correlation between metabolic signature and outcomes in stage III patients was not statistically significant (p-value = 0.799). Specifically, in the early stages (stage I-II), none of the patients in “subset low” experienced relapse, resulting in a significantly better DFS compared to patients in “subset high” (median DFS NR vs 36.1 months; p-value = 0.022) (Fig. 5 panel B). The correlation with metabolic subset was particularly significant in stage II patients, where all but two patients with metabolic “subset high” experienced an early relapse, while none of the patients in metabolic “subset low” relapsed. The median DFS for “subset high” was only 8 months vs NR for “subset low” (p-value = 0.002) (Fig. 5 panel B). Patients with grade 3 tumors exhibited worse DFS compared to those with grade 1–2 tumors (p-value = 0.005). Notably, a higher tumor grade was not associated with “subset high” (p-value = 0.224).
In line with these findings across the cohort, gender did not correlate with DFS or metabolic subsets (p-value = 0.293). However, metabolic subsets were significantly associated with DFS in both males and females (p-value = 0.049) (Fig. 5 panel C).
Smoking history was not associated with DFS (p-value = 0.471) (Fig. 5 panel D). Nevertheless, patients with a history of smoking had a higher likelihood of belonging to “subset high”, although this difference did not reach statistical significance (OR 3.100, 95% CI: 0.742–12.953, p-value = 0.135). Notably, high metabolic subset was significantly associated with DFS also in patients with a smoking history (p-value = 0.035) (Fig. 5 panel D). Finally, SUVmax values were not associated with DFS (p-value = 0.956).
We investigated the putative correlation between the oncogene mutational profile and DFS of DEFLeCT cohort (Fig. 5 panel E). The comparison between non-mutated and mutated tumors indicated a not significant difference in terms of DFS (p-value = 0.97). Probably due to the small number of data, the significant difference in DFS between “subset low” and “subset high” is lost by separating patients into non-mutated and mutated ones (p-value = 0.056 and 0.29 respectively). It is noteworthy that DFS in non-mutated and mutated patients within both the low subset and the high subset was not statistically significant (p-value = 0.73 and 0.08 respectively).
In the multivariate analysis, we included variables that had a p-value < 0.1 at the univariate analysis: pathological staging, grading, and metabolic subset. Pathological staging and metabolic subset maintained a statistically significant association with DFS (p-value = 0.036 and 0.033 respectively), whereas grading did not maintain a significant correlation with the survival outcome (p-value = 0.44)(Fig. 5 panel F).

Identification of Lung Metabolic Signature “LMetSig” related to poor prognosis in LUAD
Based on exploratory analyses in DEFLeCT LUAD patients, we compared the gene expression level between “subset high” vs “subset low” (5612 genes) and we identified 429 differentially expressed genes (DEGs) (Supplementary figure 3, Supplementary data 7). The intersection of DEGs with “Cluster G” of genes revealed a list of 23 differentially expressed genes (gDEGs) up-regulated in tumors of “subset high” versus “subset low” (Fig. 6 panel A and B, Supplementary data 7). The high expression of 10 out of the 23 gDEGs was significantly associated with a reduced DFS (p-value = 0.002) (Fig. 6 panel C). This set of 10 genes, that was named “Lung Metabolic Signature” (LMetSig) is composed by Aldolase-A (ALDOA), 7-Dehydrocholesterol Reductase (DHCR7), Deoxythymidylate Kinase (DTYMK), Glyceraldehyde-3-Phosphate Dehydrogenase (GAPDH), NAD (P) Dependent Steroid Dehydrogenase-Like Protein (NSDHL), 6-Phosphofructokinase Type C-F (PFKP), Phosphoglycerate Kinase 1 (PGK1), Glucose Phosphate Isomerase (GPI), Triosephosphate Isomerase (TPI1) and Thymidine Kinase (TK1) (Fig. 6 panel C).
We analyzed the expression level of LMetSig genes in tumor and normal tissues of the same patient, if available. We did not detect any significant difference between the expression of LMetSig genes in the normal tissues of the “subset high” and “subset low” (Fig. 7) while, as expected, the transcription of these genes increased significantly in tumor tissues compared to normal tissues in both subsets (except for DHCR7 and NSDHL in “subset low”) (Fig. 7). Notably, this increase was higher in tumors of patients belonging to “subset high” when compared to “subset low” (Fig. 7).
The predictive ability of LMetSig was verified in the LUAD-IC cohort (n = 528 patients): considering the mean expression of LMetSig genes patients with low expression of LMetSig showed a significantly better DFS compared to high expressors (median DFS of 20.07 vs 50.56 months, log-rank p-value < 0.0001)(Supplementary figure 4 panel A). Analyzing LUAD-IC cohort according to TNM stage (n = 274 patients), LMetSig showed a significant impact on DFS in stage-I patients (log-rank p-value = 0.0017) but not in stage-II patients (n = 98 patients)(log-rank p-value = 0.76)(Supplementary figure 4 panel B). Sample sizes of other stages were too low for meaningful analysis (n = 8 patients for stage-III and n = 0 for stage-IV).
LMetSig expression level was significantly associated with DFS in both sexes (n = 240 female patients, log-rank p-value = 0.004; n = 288 male patients, log-rank p-value = 0.0021) (Supplementary figure 4 panel B). The high metabolic subset was associated with a worse DFS for both never-smoker and smoker (former and current) patients although not significantly for the latter ones (n = 140 never-smoker, log-rank p-value = 0.024); n = 229 smoker, log-rank p-value = 0.135) (Supplementary figure 4 panel B).
We performed a multivariate analysis, including stage and LMetSig as covariates showing that both variables maintained a statistically significant correlation with DFS also in the LUAD-IC cohort (n = 380 patients, stage log-rank p-value < 0.0001, LMetSig log-rank p-value = 0.014)(Supplementary figure 4 panel C).

Exploration of LMetSig expression in a new early-stage radically resected LUAD cohort
LMetSig expression distribution was investigated in a different independent prospective cohort of early-stage radically resected LUAD (MAGA cohort n = 48 patients, Supplementary data 6). A summary of patient characteristics data is provided in Figure 8 panel A.
After mRNA-seq analysis and unsupervised clustering according on LMetSig expression level, we profiled tumors in “subset low”, “subset medium”, and “subset high” (Fig. 8 panels B and C).
As previously defined for the DEFLeCT LUAD cohort, the MAGA cohort did not exhibit a statistically significant association between metabolic subsets and smoking habit (p-value = 0.409), gender (p-value = 0.305), or TNM stage (stage I vs. stages II-III, p-value = 0.173). A borderline significant association was observed with tumor grading (p-value = 0.061). Indeed, in this cohort, high-grade tumors accounted for 64, 53, and 21% of the samples classified as metabolic subsets high, intermediate, and low, respectively (Supplementary data 6).

TK1 as a putative new prognostic biomarker for LUAD
Comparing subsets “low” and “high” in DEFLeCT and MAGA LUAD cohorts, TK1 emerged as the most significant up-regulated DEG both among the metabolic genes considered and among LMetSig genes (Fig. 9 panel A, Supplementary data 7)
The high expression of TK1 was associated with a shorter DFS in both DEFLeCT (median DFS 42 vs NR months, log-rank p-value = 0.001)(Fig. 6 panel C) and LUAD-IC patients (median DFS 11.2 vs 32 months, log-rank p-value < 0.0001) (Supplementary figure 4 panel D).
Immunohistochemistry analyses on DEFLeCT samples highlighted the specific TK1 protein expression in tumor cells only and also validated the sample classification in “subset low” and “subset high” according to gene expression level (Fig. 9 panel B).

Discussion
Cancer cell metabolism deregulation leads to uncontrolled proliferation of tumor cells and it is also associated with other malignant “phenotypes”, such as migration and invasion [11]. This deregulation is generated through a complex set of genetic alterations, aberrant signaling pathways, crosstalk with stromal cells, and the tumor microenvironment [10–12]. Because of the relevant role of metabolic networks and the superior performance of multi-gene models compared to single gene tests, we could reasonably expect that metabolism-related prognostic models might play a quite significant role in this setting. Some attempts have been made to propose genomic signatures with prognostic roles based on metabolic genes in LUAD, although the results have been highly variable and poorly standardized [13–20].
In this study we identified “LMetSig”, a new metabolic signature, able to stratify early-stage LUAD patients and to predict prognosis. Due to the limited number of events at the time of analysis, OS data were not mature. Furthermore, OS might be greatly influenced by the type of treatment received in case of relapse that was driven by the tumor profile (e.g., oncogene addiction, PDL1 expression, etc). Therefore, our main focus was on the correlation between the metabolic signature and DFS. The high expression of LMetSig was associated with shorter DFS in two different independent early-stage LUAD cohorts: we detected a statistically significant correlation of LMetSig with DFS in DEFLeCT LUAD cohort composed by 87 patients but also in the LUAD-IC composed by more than 1300 patients, further highlighting its potential prognostic value.
Notably, metabolic subsets had no significant association with clinical parameters, but had an independently significant correlation with DFS. With the limitation of the sample size, LMetSig up-regulation showed a greater impact on DFS in patients with stage I-II than in those with stage III disease. Our data suggest that the metabolic signature might play a pivotal role in identifying high-risk patients within stage I-II. The lack of a significant correlation between metabolic signature and outcomes in stage III patients deserves further investigation. Indeed, this result may be due to the small sample size in our study or it might suggest that, in more advanced stages, other tumor hallmarks overshadow the predictive value of metabolic signatures. Indeed, at a more advanced stage of the disease, the substantial tumour burden, increased genomic instability, and multiple concurrent pathophysiological processes can accelerate tumour progression. These factors may introduce confounding variables that interfere with the reliability of analytical measurements.
The absence of significant associations between LMetSig and standard clinical parameters such as smoking status, gender, TNM stage or tumor grade suggests that the signature may act as an independent biological determinant rather than reflecting conventional clinical features. This independence might reinforce its potential value as a complementary prognostic marker, capable of capturing punctual molecular aspects of tumor behavior not explained by routinely assessed clinical variables.
“Subtype high” tumors exhibited a trend towards higher uptake values (quantified as SUVmax) of the primary tumor at PET/CT scan. To date, several studies have indicated that SUVmax may represent an independent risk factor for recurrence but available evidence is not conclusive [32, 33].
The study also suggested a putative relationship between the metabolic profile and the mutational landscape of tumors. The significant correlation between high expression of LMetSig and TMB is worthy of further investigation as TMB is a potential biomarker for response to immune checkpoint inhibitors, overall tumor aggressiveness and often associated with the occurrence of neoantigens. In addition, the analysis of the distribution of the key LUAD mutated oncogenes (such as EGFR, KRAS, BRAF, and MET) showed that KRAS mutations are enriched in the metabolic “subset high” tumors. Only 30% of all patients in DEFLeCT LUAD cohort who relapsed harbored KRAS mutations, and all had a medium or high expression of LMetSig. The lack of a significant difference in DFS between mutated and non-mutated tumours, both overall and within metabolic subsets, emphasizes the complexity of the genetic and epigenetic networks in tumor progression.
LMetSig is composed by ten genes: ALDOA, DHCR7, DTYMK, GAPDH, NSDHL, PFKP, PGK1, GPI, TPI1 and TK1. We examined the degree of overlap between LMetSig and previously published metabolic signatures. Indeed, GPI, ALDOA, TPI1 and PFKP were indicated as metabolic markers related to the poor prognosis of LUAD patients in the studies of Xue C. et al. [15] and Zhao Z. et al. [16], as well as in the study of Wang Z. et al. [17]. In particular, these studies found that these genes were associated with a reduction in OS of LUAD patients. Although the overall overlap was limited, these data emphasised the prognostic value of these metabolic genes.
ALDOA gene, encodes the enzyme that facilitates the reversible transformation of fructose-1,6-bisphosphate into glyceraldehyde-3-phosphate and dihydroxyacetone phosphate and has a clear role in the progression of several cancer types, including lung cancer [34]. Recently, ALDOA has also emerged as a potential therapeutic target [35]. The PFKP gene is responsible for encoding the enzyme that facilitates the conversion of fructose-6-phosphate into fructose-1,6-bisphosphate. PFKP may be involved in metabolic reprogramming of lung cancer and it is often found to be overexpressed [36]. Specifically, Shen et al. have shown that higher PFKP expression in lung tumors correlates with worse OS while the reduction of PFKP level leads to decreased glucose uptake, which has potential antitumor effects [36]. PGK1 is a key gene in glycolysis, known for its overexpression in LUAD, which has been associated with worse prognosis [37]. Its function extends beyond metabolism, influencing hypoxia mechanisms, DNA synthesis, the recruitment of immune cells, including M2 macrophages and exhausted T cells, contributing to the immunosuppressive tumor microenvironment [37]. TPI1 has shown a significant correlation between overexpression and unfavorable prognosis in LUAD [38]. Studies in mouse models indicate that inhibiting TPI1 can significantly reduce cell migration, colony formation, and tumor growth. Finally, two genes of the LMetSig are involved in nucleotide synthesis: DTYMK and TK1. Regarding DTYMK, its prognostic role has already been demonstrated in patients with NSCLC [39]. In LUAD, high expression levels of DTYMK had been significantly associated with reduced DFS and OS, with an increase of immune-infiltration, tumor cell resistance to chemotherapy, and cell proliferation and invasion [40, 41].
Among LMetSig, TK1 gene was identified as a candidate prognostic biomarker for LUAD due to three peculiar features: i) TK1 was undetectable in normal tissue; ii) TK1 emerged as one of the most significant DEGs comparing different tumor tissues; iii) TK1 showed a significant correlation with worse DFS.
TK1 could serve as a putative therapeutic target, a predictive biomarker, and a marker of treatment response in LUAD. TK1 encodes an enzyme involved in DNA synthesis, particularly of pyrimidines, during the cell replication phase. Its expression significantly increases in active proliferating tumor cells where the “salvage pathway” of synthesis of pyrimidines supports the normal synthesis by “de novo pathway” [42, 43]. TK1 is particularly overexpressed during the S phase of the cell cycle, making it an ideal biomarker for tumor proliferation. In vitro evidence suggests that TK1 is not only expressed in the cytoplasm but may also be present on the cell membrane of lung cancer cells [44]. Membrane TK1 is rarely detected in normal tissue, making it a potential tumor-specific target for antibody-based or cell-based therapies [44], particularly in patients with immune checkpoint inhibitors-resistant LUAD [45].
Inhibition of pyrimidine synthesis, particularly the de novo pathway, is an established therapeutic strategy in oncology. The most important group of antitumor agents acting at this level includes nucleoside analogous and antimetabolites, such as 5-fluorouracil, gemcitabine, and pemetrexed. Unfortunately, the use of these drugs is associated with the rapid onset of resistance, which is also due to the parallel activation of the “salvage pathway” [42, 43]. Currently, there are no specific TK1 inhibitors approved or under clinical investigation. Inhibition of pyrimidine recovery mediated by TK1 is currently possible only through siRNA. Interestingly, silencing the cytosolic TK1 via siRNA enhances the effect of thymidylate synthase inhibition involved in “de novo” synthesis, conferring greater sensitivity to pemetrexed in cell lines [46]. Silencing of TK1 gene reduces tumour cell proliferation, invasion and immune evasion. Concurrent knockdown of both KPNA2 and TK1 restores STAT3 regulation and enhances antigen presentation [47]. After pemetrexed administration, the transient increase in TK1 (as shown in salvage-pathway “flare” PET studies) suggests that TK1 activity is dynamically modulated [48].
In the context of identifying new predictors or tracers of lung cancer, molecules that reflect the course of the disease and can be detected in a patient’s bloodstream are of fundamental importance. From this perspective, the levels of serum TK1 activity and serum TK1 protein concentration are related to tumor progression and might be used for monitoring outcomes both in early and advanced disease [49–56].
Epigenetic variables, including those linked to tumor cell metabolism, are increasingly recognized as key determinants of cancer progression. Alterations in metabolic pathways can reshape the epigenetic landscape by modulating the availability of essential cofactors, thereby influencing gene expression programs that promote tumor growth, adaptation, and therapeutic resistance. As a result, metabolic interactions are emerging as critical drivers of malignant evolution and represent a growing area of interest for both biomarker development and targeted therapeutic strategies [57].
Despite their limited numbers, the DEFLeCT, MAGA and LUAD-IC cohorts still include patients with a variety of demographic characteristics, smoking habits, clinical and radiological parameters, pathological staging, and molecular features. These characteristics collectively reflect the spectrum observed in most LUAD patients. Interestingly, LMetSig consistently clustered early-stage resected LUAD of DEFLeCT and MAGA cohorts: the comparison of the “subset low” and of “subset high” of the two cohorts showed that the two populations are similar in terms of smoking habit, stage and pathological grading. As previously described, lung tumors, especially LUAD, are heterogeneous, and prognosis varies even within the same TNM stage. The current staging system does not fully capture the molecular complexity of the disease, indicating the need for additional prognostic/predictive factors. A better understanding of tumor metabolism could improve predictions of recurrence and guide more tailored therapeutic strategies, enabling better patient selection for personalized treatments and optimizing follow-up care.

Putative clinical implementation strategy
In summary, LMetSig is a promising tool for improving the prognostic assessment of early-stage lung adenocarcinoma (LUAD) beyond traditional clinical and pathological parameters. It could help to guide risk-adapted management by identifying individuals who could benefit from intensified surveillance or adjuvant treatment despite having early-stage disease.
Although the present study is exploratory, it aims to generate preliminary evidence and lay the groundwork for future investigations based on larger, prospective cohorts.
LMetSig ability to predict disease-free survival independently and to stratify patients within the same TNM stage suggests that metabolic profiling could address a major unmet clinical need: it might distinguish between early-stage patients who are truly at low risk and those who might potentially benefit from more intensive surveillance or adjuvant therapy. Notably, the lack of association between LMetSig and most clinical variables, including tumor grade, highlights its potential as a complementary, biology-driven biomarker.
If prospectively validated, LMetSig could be incorporated as an additional molecular test performed on surgical or biopsy specimens, alongside established clinicopathological parameters. In particular, TK1 assessment might represent a surrogate and reliable biomarker of LMetSig.
Baseline TK1 expression (or its dynamic induction) could influence responsiveness to pemetrexed chemotherapy. Moreover, concurrent targeting or inhibition of TK1 may enhance the anti-tumor efficacy of pemetrexed. Thus, TK1 warrants investigation as both a predictive biomarker and a therapeutic target in combination with pemetrexed-based regimens in LUAD. Since TK1 levels reflect proliferation and salvage-pathway activity and given the fact that TK1 is overexpressed in many LUAD tumors, measuring TK1 expression prior to chemotherapy might represent a novel predictive biomarker for pemetrexed activity (alone or in combination). For example, tumors with high TK1 might show different degrees of sensitivity than tumors with low TK1, depending on how salvage pathways compensate for Thymidine Synthase inhibition. As shown in salvage-pathway “flare” PET studies, the transient increase in TK1 after pemetrexed administration suggested that TK1 activity is dynamically modulated [48]. Given the presence of TK1 on the cell membrane in lung cancer but not in healthy lung tissue, targeted therapies (e.g., antibodies, ADCs, CAR-T cells) against TK1 might selectively kill LUAD tumor cells either alone, or in synergy with pemetrexed chemotherapy. This dual approach could both reduce tumor burden and overcome adaptive resistance. Finally, combining TK1 evaluation (as DFS predictor marker) with pemetrexed might limit metastatic potential, since TK1 appears to support not only proliferation but also invasion and metastasis in LUAD.
From a health economic perspective, the signature relies on a limited number of genes and could be investigated using widely available proteomic (e.g. immunohistochemistry), or transcriptomic platforms (e.g. RT-qPCR), which could support its feasible and cost-effective integration into existing testing pipelines. TK1 is particularly amenable to investigation by IHC, as its expression can easily be assessed using routinely archived FFPE.
Furthermore, given the evidence that serum levels of TK1 can serve as a non-invasive biomarker in different tumors [49–56], integrating serum TK1 assessment alongside LMetSig may represent an additionally, minimally invasive method for monitoring tumour dynamics, improving risk stratification and potentially informing longitudinal clinical decision-making. A prospective evaluation of serum TK1 levels may be carried out using a high-sensitivity commercial Enzyme-Linked Immunosorbent Assay (ELISA).
Overall, when combined with a serum-based evaluation of TK1 where appropriate, LMetSig has the potential to be incorporated into clinical workflows as a practical, scalable and cost-effective molecular classifier that strengthens evidence-based decision-making in LUAD, particularly where current prognostic tools are insufficient. This is particularly relevant given the heterogeneity of LUAD. Therefore, incorporating metabolic signatures such as LMetSig into clinical workflows could help personalise post-surgical management, reduce overtreatment, and ultimately improve patient outcomes.

Electronic supplementary material
Below is the link to the electronic supplementary material.