Functional exploration and drug prediction on programmed cell death-related biomarkers in lung adenocarcinoma

Abstract

Background: Our study aims to perform functional exploration and drug prediction of programmed cell death (PCD)-related biomarkers in lung adenocarcinoma (LUAD).

Methods: UCSC-Xena obtained LUAD-related genes. DESeq2 screened PCD-specific differentially expressed genes (DEGs), and these DEGs were intersected with genes identified by weighted gene co-expression network analysis (WGCNA) to pinpoint the key genes. KOBAS-i was used for enrichment analysis. String and GeneMania were used to construct protein interaction networks and gene-gene interaction networks, respectively. Using two machine learning algorithms to screen for key genes, and taking the intersection as biomarkers, validating via receiver operating characteristic (ROC) and in vitro experiments. Building a diagnostic model with a nomogram. Construct transcription factor (TF) regulatory network. CIBERSORT was used for immune infiltration analysis. Enrichr predicts targeted drugs and AutodockTools simulates molecular docking.

Results: 120 hub genes related to PCD were identified, and an intersection of these genes with DEGs yielded 10 key genes, which were enriched in apoptosis-related pathways. Further machine learning screening of these genes led to the selection of 7 genes, among which 6 genes (FGR, LAPTM5, SIRPA, TLR4, ZEB2, and NLRC4) exhibited significant differences upon ROC validation, ultimately serving as biomarkers, in vitro experiments also confirmed. A nomogram demonstrated their excellent diagnostic performance. These six biomarkers are correlated with the infiltration status of most immune cells, suggesting that they affect LUAD through the immune system. TF regulation analysis identified the upstream miRNAs. Finally, drug prediction yielded three potential drugs: Lenvatinib, methadone, and trimethoprim.

Conclusion: PCD-related biomarkers in LUAD were explored, which may contribute to further understanding on PCD in LUAD.

Keywords: Lung adenocarcinoma; Machine learning; Programmed cell death; Transcriptome sequencing.

Fig. 1

Results of WGCNA analysis. (A) Quantified results of ssGSEA score on LUAD tissue (n = 510) and control tissue (n = 58) in LUAD. (B) The procedures sorting on the optimal soft threshold. (C) Gene dendrogram based on the topological overlap, together with the assigned modules colors. (D) Correlation between the traits and the modules (the numbers without the brackets were the correlation coefficient and those within the brackets were the P-values). (E) The gene-significance-module membership plot of MEgreen module.

Fig. 2

DEGs analysis. (A) Volcano plot displaying the DEGs (n = 333) based on the data from TCGA. Each dot annotates a gene, and the left and right parts represents the down-regulated and up-regulated DEGs, respectively. (B) Heatmap based on the top 20 up-regulated and down-regulated DEGs in tumor and normal tissue. (C) DEGs (n = 333) and Hubgene(n = 120) obtained by WGCNA were intersected to obtain 10 key genes. (D) Locations of the intersected PCD-related genes on the chromosome. The outermost circle in the figure is chromosome, the middle gene line represents the position of the gene on the chromosome, and the inner circle represents the DEGs, with red representing upregulation and blue representing downregulation.

Fig. 3

Gene enrichment analysis and gene-gene network construction. (A–B) GO and KEGG enrichment analyses were then carried out on the 10 common PCD-related genes using KOBAS-i. (C–D) Constructed PPI network (C) and GGI network (D) of these 10 PCD-related genes using String and GeneMania.

Fig. 4

Machine learning results. (A–B) Results on LASSO regression on the training set. (C) Relationship between generalization error and the characteristics in LUAD. (D) Venn diagram showing the common feature genes from LASSO regression test and SVM-RFE.

Fig. 5

Validation on the efficacy of feature genes in LUAD. (A–B) ROC curve on the feature genes in the training set (A) and testing set (B). (C–D) Box plot on the expression levels of feature genes in the training set (C) and testing set (D).

Fig. 6

Construction of the diagnostic model. (A) Nomogram on the 6 biomarkers (FGR, NLRC4, TLR4, LAPTM5, ZEB2 and SIRPA) in LUAD. (B) Plotted ROC curve of the nomogram. (C) Plotted calibration curve of the nomogram.

Fig. 7

Determination on the immune cell infiltration. (A) Differentially immune cells infiltration in both LUAD tumor tissue and normal tissue. (B) Correlation on the 6 biomarkers and the immune cells infiltration. *P ≤ 0.05, *P ≤ 0.01, *P ≤ 0.001, *P ≤ 0.0001.

Fig. 8

Immune checkpoint analysis results. (A) ssGSEA calculated scores for 79 immune checkpoint related genes. (B) The analysis of correlation between immune checkpoints and seven PCD genes (FGR, NLRC4, TLR4, LAPTM5, CYBB, ZEB2 and SIRPA) indicates that the higher the numerical value, the stronger the correlation.

Fig. 9

Plotted TF regulatory network on the biomarkers (NLRC4 was removed due to the zero number of TF predicted).

Fig. 10

Results on the molecular docking. (A) Molecular docking of FGR to Lenvatinib and the hydrogen binding sites located on the amino acids of MET-90 and LEU-95. (B) Molecular docking of TLR4 to Methadone and the hydrogen binding sites located on the amino acid of ALA-139. (C) Molecular docking of TLR4 to trimethoprim and the hydrogen binding sites located on the amino acids of PRO-202, LEU-204, MET-201 and LEU-198.

Fig. 11

In-vitro validation using LUAD cells. (A–F) Quantified expression levels of 6 biomarkers FGR (A), LAPTM5 (B), SIRPA (C) TLR4 (D), ZEB2 (E) and NLRC4 (F) in LUAD cells A549 and bronchial epithelial cells BEAS-2B via qPCR. (G) Effects of TLR4 silencing on the migration of LUAD cells A549 tested via wound healing assay. (H) Effects of TLR4 silencing on the invasion of LUAD cells A549 tested via Transwell assay. All experimental data of independent triplicates were expressed as mean ± standard deviation. ***P < 0.001; ****P < 0.0001.