Prognostic prediction of m6A and ferroptosis-associated lncRNAs in liver hepatocellular carcinoma

Introduction

Liver hepatocellular carcinoma (LIHC) is a prevalent malignancy, ranking sixth globally.^[1] Current treatments include resection, transplantation, radiotherapy, chemotherapy, and molecular targeted therapy. However, early diagnosis is challenging due to the absence of typical symptoms, leading to poor prognosis and a 5-year survival rate of 30%-40%.^[2] Early prediction models and therapeutic targets are urgently needed. N6-methyladenosine (m6A), abundant in messenger RNA (mRNA) and Long non-coding RNA (lncRNA), regulates various ribonucleic acid (RNA) processes and is crucial in cell fate, cycle, differentiation, and circadian rhythm. In LIHC, METTL3 promotes cell growth, migration, and tumorigenicity, while METTL14 downregulation indicates poor prognosis.^[3] YTHDF2 facilitates LIHC cell proliferation by recognizing m6A sites. Ferroptosis, an iron-dependent programmed cell death, is more sensitive in cancer cells due to their metabolic activity and reactive oxygen species (ROS) demand. m6A modification is strongly correlated with ferroptosis in LIHC, with roles in ATG5 expression and autophagy signaling.^[4-5] Studies suggest m6A’s pivotal role in regulating ferroptosis and liver cancer treatment. Prognostic models based on lncRNA expression have gained attention, with m6A-associated lncRNA signatures providing new rationale for LIHC diagnosis and therapy.^[6] However, additional prognostic markers are needed. This study identified ferroptosis- and m6A-associated lncRNAs, developed an improved prognostic model for LIHC, and investigated the relationship between tumor mutational burden (TMB) and prognosis. Functional enrichment analysis of differentially expressed genes was also conducted.

Materials and methods

Expression data and clinical information of 424 LIHC patients were obtained from the The Cancer Genome Atlas (TCGA) database. We focused on 22 m6A-related genes and 5 ferroptosis-related genes, and distinguished between mRNA and lncRNA, yielding 16,773 lncRNAs. After filtering, 370 clinical data were matched with expression data. Co-expression analysis was conducted using the “limma” package in R, with a correlation coefficient of 0.4 and P-value threshold of 0.001. A prognostic model was developed using m6A/ferroptosis-related lncRNAs, with samples randomly divided into train and test sets. Univariate and lasso-Cox regression analyses identified significant lncRNAs, and a risk score formula was constructed. Survival analysis showed the model’s effectiveness in stratifying high- and low-risk patients. The model was validated through univariate and multivariate Cox regression, with receiver operating characteristic curve (ROC) curves and calibration curves assessing accuracy. Mutation information was analyzed to determine gene mutation frequency, TMB, and tumor immune dysfunction and exclusion (TIDE) scores. Differential expression analysis was performed, followed by Kyoto Encyclopedia of Genes and Genomes (KEGG) and Gene Ontology (GO) enrichment analyses, visualized using various R packages.

Result

Screening of m6A and ferroptosis-related lncRNAs and construction of prognostic model

We used 5 ferroptosis-related genes and m6A genes for co-expression analysis, identifying 20 m6A/ferroptosis related genes (mFRGs). From TCGA-LIHC, 16, 773 lncRNAs were analyzed, revealing 593 m6A/ferroptosis related lncRNA (mFlncRNAs). Univariate cox analysis screened 80 mFlncRNAs, and lasso regression further reduced to 13, ultimately yielding 7 for the risk model. High-risk included ELFN₁-AS1, AL603839.2, etc., while low-risk were AC073573.1, AC069307.1 (Figure 1A-B) (Supplemental Figure 1A-B) (Supplemental Table 1).

Figure 1

Prognostic prediction of m6A and ferroptosis-associated lncRNAs in liver hepatocellular carcinoma. (A) Co-expression analysis of m6A/ferroptosis genes and lncRNAs. (B) Selection of the optimal penalty parameter for LASSO regression. (C) Prognostic curves in risk groups of the three data sets. (D) Univariate cox regression analysis. (E) Multivariate cox regression analysis. (F) ROC curves for all-sample group. (G) Nomogram for LIHC patient. (H) Survival curves for high-low TMB in all-sample group. (I) Survival curves at change in TMB and risk in all-sample group. (J) Comparison of TIDE in high-low risk groups.

Discussion

LIHC is a complex disease caused by multiple factors, including chronic hepatitis B virus (HBV), hepatitis C virus (HCV), aflatoxin-contaminated food, heavy alcohol consumption, obesity, smoking, and type 2 diabetes.^[7] Current early diagnosis relies on alPhafetoProtein (AFP), but its accuracy is limited.^[8] Thus, exploring additional molecules in LIHC’s molecular mechanism is crucial for better disease understanding and management. Studies have shown alterations in lncRNA levels are implicated in various cancers and could serve as diagnostic markers and therapeutic targets.^[9] Our study presents a novel risk-prognostic model for LIHC patients, focusing on lncRNAs, m6A modification, and ferroptosis. We identified ELFN₁-AS1, AC018690.1, AL603839.2, LINC02313, and AC004801.6 as high-risk factors, while AC069307.1 and AC073573.1 were low-risk factors. Our model, validated by PCA and ROC curve, showed higher predictive power than conventional clinical models. Additionally, we found that high TMB was associated with poor prognosis in LIHC patients, but our risk prognosis model demonstrated stronger predictive potential. KEGG and GO analysis identified gene enrichment in cell cycle regulation, bile secretion pathways, and xenobiotics metabolism. While our study identified mFlncRNAs’ potential in early LIHC diagnosis, further validation is needed. Leveraging mFlncRNAs can devise more rational therapeutic strategies, ultimately enhancing LIHC patients’ survival rates.