Genome-wide identification of FHL1 as a powerful prognostic candidate and potential therapeutic target in acute myeloid leukaemia.

BACKGROUND: Acute myeloid leukaemia (AML) is a malignant haematological tumour with high heterogeneity and mortality. A reliable prognostic assessment is critical for treatment strategies. However, the current prognostic evaluation system of AML is insufficient.
METHODS: Genome-wide univariate Cox regression analysis was performed on three independent AML datasets to screen for the prognostic-related genes. Kaplan-Meier survival analysis was employed to verify the efficacy of FHL1 in evaluating overall survival in 1298 de novo AML patients, 648 non-acute promyelocytic leukaemia AML patients and 407 cytogenetically normal AML patients; the data for some of these patients were also used for EFS and RFS validation. Multivariate Cox regression was performed to validate FHL1 as an independent prognostic indicator. WGCNA, GSEA, and gene correlation analysis were applied to explore the mechanism of FHL1 in AML. The synergistic cytocidal effect of FHL1 knockdown was verified in in vitro experiments.
FINDINGS: Comprehensive genome-wide analyses and large-sample validation showed that FHL1 is a powerful prognostic candidate for overall survival, event-free survival, and relapse-free survival in AML and is independent of prognosis-related clinical factors and genetic abnormalities. The molecular mechanism may occur through regulation of FHL1 in leukaemia stem cells, tumour-associated signalling pathways, and transmembrane transport of chemotherapeutic drugs. FHL1-targeted intervention enhances the sensitivity of AML cells to cytarabine.
INTERPRETATION: FHL1 may serve as an evaluation factor for clinical strategy selection, and its targeted intervention may be beneficial for chemotherapy in AML patients.

AML is a haematological malignancy derived from haematopoietic stem cells that has high heterogeneity and mortality. Risk stratification and prognosis assessment are critical for the clinical strategy selection of AML patients. The current evaluation system mainly depends on patient-related factors, such as increasing age, and prognosis-related cytogenetics and gene mutations, which are very inadequate. However, to some extent, research based on abnormalities in gene expression profiles may be able to compensate for the deficiency. Although studies have shown that gene expression signature-derived scoring systems have clinical significance for the prognosis of AML, they are not widely used in clinical practice due to the large number of genes involved. Thus, confirmation of effective single-gene candidates may constitute a breakthrough in AML prognosis assessment and targeted intervention. Therefore, we screened single-gene AML prognostic candidates at the genome level. In addition, the value of

Our research creatively integrated genome-wide gene expression data and clinical information for relatively large samples from public databases to search for potential prognostic-related genes. Bioinformatics analysis and

The findings of our research suggest that

Alt-text: Unlabelled box

Introduction

Acute myeloid leukaemia (AML) is a complex and heterogeneous disease with different genetic backgrounds and responses to treatment. Risk stratification and prognosis assessment are of great significance in the choice of treatment for AML patients [1], [2], [3]. Pretreatment prognostic factors for AML mainly include patient-related factors, such as increasing age, and AML-related genetic factors, such as cytogenetic abnormalities and gene mutations. For example, monosomal karyotype and RUNX1, ASXL1, and TP53 mutations have been added as characteristics of the adverse-risk population due to their independent associations with risk [4], [5], [6]. However, about 30% of AML patients are classified as having intermediate risk according to the 2017 European LeukemiaNet (ELN) recommendations, many of which do not carry karyotype abnormalities or gene mutations with prognostic value, and the choice of ideal treatment remains unclear [4,7,8]. Therefore, the current risk stratification and prognosis assessment of AML patients need to be further improved.

Changes in the expression profile of patients with AML based on high-throughput sequencing and gene microarray technology have prognostic value [9]. Indeed, studies have shown that gene expression signature-derived scoring systems have clinical significance for the prognosis of AML, such as the 86-probe-set gene expression signature [10], the 24-gene prognostic signature [11] and the 11-gene risk scoring system [12]. However, due to technical problems in implementing a large number of gene signatures, the wide use of these scoring systems in clinical practice is difficult and not conducive to in-depth mechanistic research. Although some studies suggest that individual genes have prognostic value in AML, the findings are far from adequate [13], [14], [15]. To avoid differences between cohorts, we selected 1298 de novo adult AML patients from multiple independent datasets, including those with normal and aberrant karyotypes, to screen and validate candidate genes and compared them with the reported prognostic scoring systems to obtain the most powerful and generally effective prognostic indicators of AML, thus laying the foundation for subsequent research into the relevant cellular functions and molecular mechanisms.

In this study, we used bioinformatic tools to screen and validate prognosis-related markers at the genome-wide level in a large-scale AML patient cohort and found that high expression of a single candidate gene, four-and-a-half LIM domain 1 (FHL1), which is generally considered to be a tumour suppressor gene, can be used as a powerful and independent indicator of prognosis for AML patients, including non-acute promyelocytic leukaemia (non-APL) AML patients and cytogenetically normal AML (CN-AML) patients. High expression of FHL1 is associated with multiple clinical features of AML and may be involved in drug resistance and relapse in AML patients by regulating leukaemia stem cell (LSC) function, affecting multiple tumour-associated signalling pathways and the transmembrane transport of chemotherapeutic drugs. In addition, targeted intervention of FHL1 enhances the sensitivity of AML cells to cytarabine, suggesting that it may be a new strategy for AML treatment.

Materials and methods

Patients and treatment

The first cohort was derived from The Cancer Genome Atlas (TCGA), which provided 151 de novo AML patients representing the major morphological and cytogenetic subtypes of AML, with RNA expression profiles based on high-throughput sequencing (RNA-Seq) and detailed clinical information. All gene expression data were available through the data portal (https://portal.gdc.cancer.gov/). The clinical characteristics, cytogenetic and molecular information, treatments, and survival status of AML patients were summarized by TCGA Research Network [16].

The other two independent cohorts were provided by the German AML Cooperative Group (AMLCG) and included 162 and 78 untreated CN-AML patients (MDS-RAEB was excluded). All patients received intensive chemotherapy. Clinical and survival information as well as microarray data were obtained from Gene Expression Omnibus (GEO, accession number: GSE12417).

Several other published AML data sets were also used for further analyses, as follows: GSE37642-GPL96 (n = 417), GSE37642-GPL570 (n = 136) and GSE106291 (n = 250) for survival analyses of all de novo or non-APL AML patients; GSE71014 (n = 104) for survival analyses of CN-AML patients; GSE83533 (n = 19) for relapsed AML case analysis; and GSE30029 (n = 121) to investigate gene expression patterns in stem cells. All the expression information and clinical data mentioned above are publicly available in the GEO database; cases with incomplete clinical data were omitted.

For further verification, bone marrow samples from 28 newly diagnosed AML patients were collected with informed consent from the Department of haematology at Qilu Hospital of Shandong University in China; of these cases, 15 were sensitive to the subsequent standard chemotherapy regimen, and 13 were resistant. This study was approved by the Medical Ethics Committee of the School of Medicine at Shandong University. Informed consent of the patients was obtained.

Univariate Cox regression analyses

To identify prognostic genes, genome-wide univariate Cox regression analyses for overall survival (OS) were performed in the cohort from TCGA as the training set. With FDR<0.05 as the statistical boundary, prognostic candidate genes were obtained. The candidate genes were validated by univariate Cox regression in two AMLCG cohorts representing CN-AML patients to identify potential prognostic indicators. Gene expression levels were dichotomized based on the median expression level of each gene in the individual cohorts as the cutoff value.

Weighted gene co-expression network analysis (WGCNA)

WGCNA is a systematic biology method used to describe the pattern of gene correlation between different samples. The association between gene sets and phenotypes identifies candidate biomarker genes or therapeutic targets. In this study, a gene expression matrix for TCGA data was constructed, and 25% of the genes with the largest variance were selected as input data. WGCNA was performed using the “WGCNA” package in R software with construction of an adjacency matrix and a topological overlap matrix (TOM) and calculation of the corresponding dissimilarity (1-TOM). Gene dendrogram construction and module identification were performed with a dynamic tree cut, and correlations between the module eigengenes and survival conditions were calculated [17]. Further analyses of the module containing FHL1 were performed for gene expression and functional enrichment analyses.

Gene signature analyses

The biological function related to FHL1 was analysed by Gene Set Enrichment Analysis (GSEA, http://software.broadinstitute.org/gsea/index.jsp). GSEA was performed on samples with high (top quartile) and low (bottom quartile) FHL1 expression in TCGA. The cutoff values for GSEA were nominal P <0·05 and false discovery rate (FDR)<0·25.

Cell culture and lentivirus transduction

The AML cell lines Kasumi-3, U937 and Kasumi-1 were obtained from American Type Culture Collection (ATCC) and cultured in RPMI 1640 medium (Life Technologies, Carlsbad, CA) supplemented with 10% foetal bovine serum (Gibco, Carlsbad, CA) in a humidified atmosphere with 5% CO2 at 37 °C [18]. Kasumi-3 and U937 cells were transduced with lentiviruses expressing short hairpin RNAs (shRNAs) against FHL1 or with scramble sequences obtained from GenePharma (Shanghai, China) to compare the biological effects of FHL1 knockdown. Puromycin was used to select cells positive for viral infection. The interference sequences targeting FHL1 were sh-FHL1-1 5′-GGACTTCTACTGCGTGACTTG-3′ and sh-FHL1-2 5′-GCTGTGGAGGACCAGTATTAC-3′.

Plasmid transfection

The 3Flag-tagged FHL1 eukaryotic expression plasmid and its vector (Genechem, Shanghai, China) were transfected into Kasumi-1 and U937 cells with Roche Transfection Reagent (Roche, Switzerland) according to the manufacturer's protocol.

Cytarabine treatment and cell viability

Cytarabine was purchased from Sigma-Aldrich (Sigma-Aldrich, St Louis, MO), dissolved in dimethyl sulfoxide (DMSO) and kept frozen at −20 °C. AML cells were treated with cytarabine at a concentration of 100 nM or 500 nM. Cell Counting Kit-8 (CCK-8) assays were performed to detect cell viability. Cells were dispensed into 96-well plates at a density of 8 × 103 cells in 100 µL of complete medium with different concentrations of cytarabine for 72 h, after which 10 µL of CCK-8 reagent was added to each well and incubated for 3 h in an incubator. Absorbance at 450 nm was measured using a microplate reader, and cell viability was calculated.

RNA extraction and quantitative real-time PCR (qRT-PCR)

Total RNA in cells and clinical samples was extracted using the Trizol method as previously reported (Invitrogen, Carlsbad, CA, USA). Complementary DNA (cDNA) was synthesized using PrimeScript™ RT Reagent Kit with gDNA Eraser (TaKaRa, Japan). FHL1 gene expression was verified by PCR using TB Green™ Premix Ex Taq™ (TaKaRa, Japan) with β-actin as a control. The primer sequences were as follows: FHL1 forward 5′-CCAACACCTGTGTGGAATG-3′ and reverse 5′-GAGTCCTCCCGAGTGGTG-3′; β-actin forward 5′-AGTTGCGTTACACCCTTTCTTG-3′ and reverse 5′-CACCTTCACCGTTCCAGTTTT-3′.

Statistical analysis

FHL1 expression was dichotomized based on either a higher (FHL1high) or lower (FHL1low) expression value compared to the median of the cohort to which the patient belongs as a cutoff. Overall survival (OS) was defined as the time from the initial diagnosis to death for any reason or the end of observation. Event-free survival (EFS) comprised the interval between the initial diagnosis and relapse or death. The probabilities of OS, EFS and relapse-free survival (RFS) were calculated with the log-rank test.

Multivariate Cox regression models containing multiple prognostic variables, including age, white blood cell (WBC) count, adverse cytogenetic risk and gene mutations that were statistically significant in univariate Cox regression, were used in analyses of OS, EFS or RFS in total de novo AML patients, non-APL AML or CN-AML patients (adverse cytogenetic risk was removed). C-index is calculated by the “survcomp” package of R software, and C-index is compared using the “cindex.comp” package. The time-dependent receiver operating characteristic (ROC) analysis and the calculation of the area under the curve (AUC) are performed by the “survcomp” and “survival” packages of R software. The comparison of integrated area under the curves (IAUC) uses “iauc.comp” package.

Analyses between the two groups were performed using Student's t-tests, Welch's t-tests, paired t-tests or Mann–Whitney tests, and the relationship between gene expression and clinical features was determined using Mann–Whitney tests, Chi-square tests or Fisher's exact tests. Gene expression correlation analyses and all statistics were completed using R 3.5.3, Stata/IC 15.0 and GraphPad Prism 8.0.2 software. Each experiment was repeated three times, and all data are presented as the mean ± standard error. For all tests, the level of significance was p-value <0·05.

Results

Integrated genomic screening for new prognostic markers in AML.

To identify robust AML prognosis-related genes for prognosis prediction and subsequent targeted intervention, we synthesized gene expression data and survival status in patients with AML and performed genome-wide univariate Cox regression analysis to screen for prognosis-related genes in three independent cohorts, including TCGA-LAML and two CN-AML datasets: AMLCG (1999-2003) and AMLCG (2004). First, we identified 12,311 genes that exist in all the three independent data sets. Then using TCGA data set as the training set, the above-mentioned genes were subjected to the univariate Cox regression after dividing into two groups based on the median gene expression level. With FDR<0.05 as the statistical boundary, 394 candidate genes were obtained. The candidate genes were validated by univariate Cox regression in AMLCG (1999-2003) and AMLCG (2004). Finally, three genes, FHL1, HOPX and FAM124B were verified in both validation data sets and evaluated further in prognosis assessment (Fig. 1a).

Fig. 1

Identification of new prognostic markers in AML based on a genome-wide approach. a. The schematic diagram to identify potential prognosis-related genes in AML by univariate Cox regression based on TCGA data and two other cytogenetically normal AML datasets, AMLCG (1999-2003) and AMLCG (2004) (univariate Cox regression). b. Multivariate Cox regression analyses of the three potential prognosis-related genes FHL1, HOPX and FAM124B screened above in TCGA data set (multivariate Cox regression; FHL1 p-value = 0·0001, HOPX p-value = 0·005, FAM124B p-value = 0·048). c. The C-index of FHL1, HOPX and FAM124B in TCGA data set (Student's t-tests; compared with FHL1, HOPX p-value = 0·089, FAM124B p-value = 0·020). d. The time-dependent receiver operating characteristic (ROC) analysis of the three potential prognosis-related gene. The area under the curve (AUC) was calculated and compared (Mann–Whitney tests; compared with FHL1, HOPX p-value = 0·011, FAM124B p-value<0·0001). HR, hazard ratio; CI, confidence interval.

To compare the prognostic evaluation capacity of these three genes, we performed the multivariate Cox regression analysis of FHL1, HOPX and FAM124B on the survival of AML patients. The results showed that all the three genes were statistically significant in the Cox model, in which FHL1 showed a larger Hazard Ratio (HR), might play an important role in the model (Fig. 1b). Then we calculated the C-index for the survival assessment of these three genes. It is shown that the C-index of FHL1 was 0·762 (0·683–0·842), HOPX was 0·690 (0·601–0·779), and FAM124B was 0.652 (0·558–0·746). The C-index of FHL1 was higher than that of FAM124B (Student's t-test; P= = 0·02), and slightly higher than that of HOPX although there was no statistical significance (Student's t-test; P= = 0·089) (Fig. 1c). The time-dependent receiver operating characteristic (ROC) analysis on the three genes was also performed, and the results showed that FHL1 had a larger area under the curve (AUC), especially for 10-year survival. Statistical analysis of the integrated area under the curves (IAUC) showed that the AUC of FHL1 is higher than that of HOPX (Wilcoxon rank-sum test; P= = 0·011) and FAM124B (Wilcoxon rank-sum test; P <0·0001) (Fig. 1d). Considering the specificity of APL in AML, we also performed survival assessment analysis of FHL1, HOPX and FAM124B in non-APL AML patients. Combining the results of univariate Cox regression, multivariate Cox regression results and time-dependent ROC, we found that FHL1 also has a better survival evaluation performance in non-APL AML patients (Fig. S1a–c). Based on the results above, FHL1 has a better prognostic assessment efficacy of AML in three screened genes, but its role in leukaemia is unclear, so we chose FHL1 for further research.

Higher FHL1 expression is powerful for predicting poor clinical outcomes in de novo AML patients

A total of 1298 adult patients with de novo AML from multiple trials were used to assess the prognostic value of FHL1. Among them, total de novo AML patients with higher FHL1 expression had a shorter OS than did those with lower FHL1 expression in four independent datasets (GSE37642-GPL96 (median OS 8·7 vs. 16·4 months; log-rank test, P = 0·011), GSE37642-GPL570 (median OS 11·4 vs. 24·1 months; log-rank test, P = 0·016), TCGA (median OS 8·1 vs. 46·5 months; log-rank test, P<0·0001) and GSE106291 (median OS 15·1 vs. 28·8 months; log-rank test, P = 0·04)) (Fig. 2a, Fig. S2a). Similarly, in the non-APL group, high FHL1 expression also suggests shorter survival and worse prognosis for AML patients (Fig. 2b, Fig. S2b). The prognostic value of FHL1 expression was also evaluated in the intermediate cytogenetic risk category. We found unfavourable prognostic effects of higher FHL1 expression in 407 de novo AML patients with normal karyotypes, including those in TCGA CN-AML (median OS 10·7 vs. 25·3 months; log-rank test, P = 0·021), GSE71014 (median OS 25·0 months vs. undefined; log-rank test, P = 0·029) datasets, AMLCG (1999-2003) (median OS 7·9 vs. 33·3 months; log-rank test, P <0·0001) and AMLCG (2004) (median OS 11·4 months vs. undefined; log-rank test, P = 0·0024) (Fig. 2c, Fig. S2c). In addition, compared with FHL1-low patients, FHL1-high patients had worse EFS (median OS 6·4 vs. 17·0 months; log-rank test, P <0·0001) and RFS (median OS 13·4 vs. 34·1 months; log-rank test, P = 0·005), which was also observed in CN-AML patients (EFS median OS 7·5 vs. 14·8 months, log-rank test, P = 0·0045; RFS median OS 8·3 vs. 17·0 months, log-rank test, P = 0·025) and the EFS of non-APL AML patients, but not RFS (EFS median OS 7·2 vs. 12·95 months, log-rank test, P = 0·0018; RFS median OS 14·9 vs. 17·0 months, log-rank test, P = 0·23) (Fig. 2d–f). Moreover, consistent with prognosis predictions, patients with high FHL1 expression had worse cytogenetic and molecular risk classifications (Fig. S3a–c). Therefore, FHL1 is an effective prognostic factor for survival in patients with de novo AML and is also useful in intermediate-genetic risk patients.

Fig. 2

High FHL1 expression is associated with poor survival in AML. a. Overall survival (OS) curves according to high expression levels of FHL1 compared with low expression levels of FHL1 divided by the median in 2 independent de novo AML cohorts GSE37642-GPL96 (log-rank test, p-value = 0·011) and GSE37642-GPL570 (log-rank test, p-value = 0·016). b. OS curves according to dichotomized expression of FHL1 as mentioned above in non-APL patients in GSE37642-GPL96 (log-rank test, p-value = 0·019) and GSE37642-GPL570 (log-rank test, p-value = 0·0045). c. OS curves grouped by median of FHL1 in 2 independent cytogenetically normal (CN)-AML cohorts TCGA CN-AML (log-rank test, p-value = 0·021) and GSE71014 (log-rank test, p-value = 0·029). d-f. Analyses of event-free survival (EFS) and relapse-free survival (RFS) with dichotomized expression of FHL1 in the cohorts from total de novo AML (log-rank tests; EFS p-value<0·0001, RFS p-value = 0·005), non-APL AML (log-rank tests; EFS p-value = 0·0018, RFS p-value = 0·23) and CN-AML patients (log-rank tests; EFS p-value = 0·0045, RFS p-value = 0·025).

FHL1 is an independent predictor of poor outcomes in AML

Multiple clinical factors have a significant impact on the prognosis of AML, represented by age, white blood cell (WBC) count, and cytogenetic risk. In addition, many gene mutations also have important indications for the prognosis of AML. Univariate Cox regression was used for genes with higher mutation frequencies to select prognosis-related mutant genes that were included in the multivariate Cox regression (Table S1–3). For total de novo AML patients, multivariable analyses were performed with adjustment for age, WBC count, adverse cytogenetic risk and several parameters that exhibited significant impacts on OS, including MLL-PTD, FLT3, DNMT3A, TP53 and RUNX1 mutations in total de novo AML patients in TCGA. When combined with all of these prognostic factors in multivariate analysis, higher FHL1 expression was still an independent poor prognostic factor for OS (HR = 2·194; 95% confidence interval (CI), 1·409-3·416; P  =  0·001), EFS (HR = 2·025; 95% CI, 1·326-3·093; P  =  0·001) and RFS (HR = 2·240; 95% CI, 1·300-3·859; P  =  0·004) in AML (Table S1, Fig. 3a). High FHL1 expression was still an independent poor prognostic factor for OS in non-APL AML patients (Table S2, Fig. S4). In the subtype of CN-AML, higher expression of FHL1 was also associated with OS (HR = 2·143; 95% CI, 1·134–4·050; P  =  0·019), EFS (HR = 2·498; 95% CI, 1·355–4·604; P  =  0·003) and RFS (HR = 2·504; 95% CI, 1·235–5·078; P  =  0·011), independent of age, WBC count, MLL-PTD, and FLT3 and DNMT3A mutation status (Table S3, Fig. 3b).

Fig. 3

FHL1 is an independent prognostic factor of poor outcomes in AML. Multivariate Cox regression analyses of OS, EFS and RFS in TCGA (a) and TCGA CN-AML (b). The p-values are obtained by multivariate Cox regression and shown in the figures above. No Cox regression assumptions were violated assessed using the Stata software. HR, hazard ratio; CI, confidence interval; WBC, white blood cells; FHL1 high, the expression value of FHL1 is greater than the median.

Comparison of FHL1 with other published predictive gene sets for prognostic assessment

Several genome-wide AML prognostic gene sets based on gene expression have been reported. To compare the prognostic value of FHL1 with other prognostic models, we performed pairwise multivariate Cox regression on FHL1 and 3-gene, 7-gene and 24-gene models [11,19,20]. According to multivariate Cox regression, FHL1 was still independent and statistically significant (P <0·05) for OS (8/9), EFS (3/3) and RFS (2/3) in most comparisons (Table 1) and for OS (3/3) and EFS (3/3) in non-APL group (Table S4). These results show that FHL1 is an effective, independent and simple AML prognostic indicator.

Table 1

Comparison of FHL1 with other published predictive gene sets in AML prognostic evaluation.

Predictors	TCGA-OS (n = 151)	TCGA-EFS (n = 151)	TCGA-RFS (n = 151)	GSE12417-GPL96 (n = 162)	GSE12417-GPL570 (n = 78)
FHL1	<0.0001 (2.718;1.821- 4.056)	<0.0001 (2.302;1.573- 3.369)	0.015 (1.828;1.125- 2.970)	0.0003 (2.216;1.444- 3.401)	0.068 (1.840;0.956- 3.544)
3-gene score (Wilop et al.)	<0.0001 (1.740;1.387- 2.182)	<0.0001 (1.662;1.338- 2.064)	0.002 (1.569;1.181- 2.083)	0.475 (1.135;0.801- 1.608)	0.025 (1.895;1.082- 3.320)
FHL1	<0.0001 (3.044;1.893- 4.894)	<0.0001 (2.899;1.832- 4.589)	0.008 (2.251; 1.238- 4.094)	0.045 (1.642;1.011- 2.668)	0.018 (2.124;1.139- 3.961)
7-gene score (Marcucci et al.)	0.882 (0.966;0.611- 1.528)	0.32 (0.798;0.511- 1.246)	0.457 (0.798;0.441- 1.445)	0.012 (1.870;1.146- 3.053)	0.098 (1.680;0.908- 3.106)
FHL1	<0.0001 (2.411;1.584- 3.670)	0.0003 (2.092;1.406- 3.112)	0.056 (1.636; 0.987- 2.712)	0.001 (2.116;1.338- 3.349)	0.018 (2.293;1.155- 4.552)
24-gene score (Li et al.)	0.002 (1.944;1.277- 2.958)	0.002 (1.896;1.247- 2.823)	0.013 (1.903;1.144- 3.168)	0.376 (1.224;0.782- 1.917)	0.697 (1.141;0.588- 2.214)

Statistical test: multivariate Cox regression.

FHL1 expression is associated with multiple clinical and molecular characteristics of AML

Patients with AML in TCGA were divided into two groups according to the median expression level of FHL1. Higher expression of FHL1 was associated with older age (Mann-Whitney test, P  =  0·0054) and complex cytogenetics (Chi-square test, P  =  0·022). Furthermore, the FHL1-high and -low groups revealed different distribution characteristics of the French-American-British (FAB) classification (Fisher's exact test, P  =  0·001). Patients with M0 more frequently had higher FHL1 expression, whereas those with M3 mostly showed lower FHL1 expression. We analysed genes with higher mutation frequencies in TCGA and found that patients with higher FHL1 expression had higher incidences of mutations in RUNX1 (Chi-square test, P <0·0001) but that those with lower FHL1 expression had mutations in genes with a better clinical prognosis, such as NPM1 (Chi-square test, P  =  0·025) and CEBPAdouble mutations (Fisher's exact test, P  =  0·028) (Table 2).

Table 2

FHL1 expression is associated with multiple clinical and molecular characteristics of AML.

	Expression of FHL1		P-value
	Low (n = 75)	High (n = 76)
Age			0.0054⁎
Median (range)	51 (21–81)	61 (21–88)
Gender, no. (%)			0.745⁎⁎
Male	42 (56.0)	40 (52.6)
Female	33 (44.0)	36 (47.4)
AML-FAB subtype, no. (%)			0.001†
M0	1 (1.3)	14 (18.4)	0.001⁎⁎
M1	16 (21.3)	20 (26.3)	0.567⁎⁎
M2	20 (26.7)	17 (22.3)	0.575⁎⁎
M3	13 (17.3)	2 (2.6)	0.003⁎⁎
M4	16 (21.3)	13 (17.1)	0.541⁎⁎
M5	8 (10.7)	7 (9.2)	0.792⁎⁎
M6	1 (1.3)	1 (1.3)	1.000†
M7	0 (0.0)	1 (1.3)	1.000†
N.D.	0 (0.0)	1 (1.3)	1.000†
WBC (× 109/L)			0.562⁎
Median (range)	15.2 (0.4–223.8)	18.75(0.6–171.9)
BM blast (%)			0.278⁎
Median (range)	74 (30–100)	71 (32–97)
PB blast (%)			0.313⁎
Median (range)	34 (0–97)	48.5 (0–91)
Cytogenetics, no. (%)
PML-RARA	12 (16.0)	2 (2.6)	0.0049⁎⁎
BCR-ABL	0 (0.00)	3 (3.9)	0.245†
Complex cytogenetics	4 (5.3)	14 (18.4)	0.022⁎⁎
Mutation, no. (%)
MLL-PTD	2 (2.7)	6 (7.9)	0.276†
FLT3-ITD	14 (18.7)	19 (25.0)	0.432⁎⁎
NPM1	25 (33.3)	13 (17.1)	0.025⁎⁎
DNMT3A	20 (26.7)	16 (21.1)	0.450⁎⁎
RUNX1	0 (0.0)	14 (18.4)	<0.0001⁎⁎
CEBPAdouble	5 (6.7)	0 (0.0)	0.028†
TP53	2 (2.7)	9 (11.8)	0.056⁎⁎
NRAS	4 (5.3)	3 (3.9)	0.719†
KRAS	3 (4.0)	4 (5.3)	1.000†

Mann–Whitney test.

Chi-square test.

Fisher's exact test.

WGCNA and identification of AML prognosis-related modules

WGCNA, a systematic biological method that divides genes into different clusters or modules based on the similarity of gene expression patterns, can be used for analysis of the potential functions of genes [21]. To identify AML prognosis-related modules, we performed WGCNA based on gene expression profiles and survival information, including OS time, OS status, EFS time, and EFS status, for AML patients in TCGA. Samples with large differences in gene expression patterns were removed (Fig. 4a). We chose the scale-free R2 = 0·8 as a soft threshold to implement a scale-free network (Fig. 4b, c). As a result, 32 gene co-expression modules were identified after dynamic tree cut merging (Fig. 4d). The heat map plotted the TOM among 1000 selected genes, indicating that each module was verified independently of each other (Fig. 4e), and correlation between each module and the survival of AML patients was calculated. We identified two modules with statistically significant associations with OS time, OS status, EFS time, and EFS status (P < 0·05): one for tumour suppression and the other for tumour promotion. The green module correlated positively with the survival time of AML patients (Pearson correlation coefficient, OS P  =  2 × 10−4; EFS P  =  8 × 10−5) but negatively with the survival state of patients (defining death as the positive event; Pearson correlation coefficient, OS P  =  7 × 10−6, EFS P  =  6 × 10−6); the red module correlated negatively with survival time (Pearson correlation coefficient, OS P  =  0·03; EFS P  =  0·04) but positively with survival state (Pearson correlation coefficient, OS P  =  0·01, EFS P  =  0·02), in which FHL1 was located (Fig. 4f).

Fig. 4

Weighted gene co-expression network analysis (WGCNA) in de novo AML patients. a. Hierarchical clustering dendrogram of TCGA samples, which is a clustering based on the similarity of gene expression profiles and shows the similarity of gene expression pattern between samples. The OS time, OS status, EFS time and EFS status of each patient are displayed at the bottom. b, c. Analysis of the scale-free fit index and the mean connectivity for the determination of soft-thresholding powers. d. Hierarchical clustering dendrogram of genes in different modules. e. The heat map describes the TOM among 1000 selected genes in WGCNA. f. Correlation analysis between the module eigengenes and survival traits. Each row corresponds to a module eigengene, and columns represent different survival traits. The p-values of each module's correlation with the corresponding survival condition are shown in parentheses (Pearson correlation coefficient); FHL1 is included in the red module, which is not conducive to the survival of AML patients. g. Heat map shows gene expression of the red module in which FHL1 is located in 151 AML samples in TCGA. h. KEGG functional enrichment analysis was used to enrich the biological function of the red-module genes.

As genes in the same module usually have similar biological functions, we analysed the genes in the same module as FHL1. This red module contained some genes that have been reported as prognostic markers for AML, such as CALCRL, DOCK1 and HOPX [8,13,22]. It also included some stem cell-related genes, such as CD34 and FAM30A [23,24], and some transcription factor-related genes, such as JAK1, SMAD1 and Notch1 (Fig. 4g) [25], [26], [27]. Kyoto encyclopedia of Genes and Genomes (KEGG) enrichment analysis also revealed that the red-module genes are involved in transcriptional misregulation in cancer (Fig. 4h). WGCNA of TCGA data revealed the module containing FHL1 to be significantly associated with the prognosis of AML, and functional analysis of the module-related genes suggested that FHL1 might be involved in the regulation of gene transcription and stem cells and may act as a candidate prognostic marker.

FHL1-related genes are enriched in LSC signatures, pathways in cancer and the transcellular transport of chemotherapeutic drugs

To explore the relevance of FHL1 expression to tumour cell hallmarks, we performed GSEA for patients with AML based on the level of FHL1 expression. AML, LSC signatures [24], pathways in cancer (including the Wnt signalling pathway, MAPK signalling pathway and JAK/STAT signalling pathway), the cellular response to drugs, and the transcellular transport of chemotherapeutic drugs were significantly enriched (Fig. 5a–h). To further validate the results, we examined the correlation between FHL1 and well-known LSC-related genes in 151 de novo AML patients and found that expression of FHL1 correlated significantly and positively with ADGRG1, FAM30A, CD34, ZBTB46 and NYNRIN, which are reported to be important components of LSC signatures and also upregulated in relapsed patients than in newly diagnosed AML patients (Fig. 5i, Fig. S5a–i) [24]. This positive correlation of FHL1 with ADGRG1, FAM30A, CD34 persisted in drug-resistant AML patients (Fig. S5j). In addition, expression of FHL1 and that of SMAD3 and TCF4 in the Wnt signalling pathway [28,29], MEF2C and TAB2 in the MAPK signalling pathway [30,31], and IL2RA and JAK1 in the JAK/STAT signalling pathway showed an obvious positive correlation (Fig. 5j–l) [25,32]. The GSEA results also showed that FHL1 expression-related genes were enriched in the transmembrane transport of chemotherapeutic drugs; moreover, expression of ABCC1 and ABCC4, which are responsible for transporting chemotherapeutic drugs such as cytarabine out of the cell, was elevated in the FHL1-high expression group (Fig. 5m–n) [33], [34], [35]. Additionally, in patients with chemotherapeutic drug-resistant disease, expression of FHL1 correlated negatively with the SLC29A1 gene, which is responsible for transporting cytarabine into cells (Fig. 5o) [36,37]. Therefore, FHL1 may be involved in the resistance of AML cells to chemotherapeutic drugs by reducing the uptake of chemotherapeutic drugs into cells and increasing their excretion out of cells.

Fig. 5

FHL1 is associated with LSCs, vital pathways in cancer and the transport of chemotherapeutic drugs in AML. a–h. GSEA of the molecular characteristics of AML, LSCs, multiple cancer-related signalling pathways, and the response to drugs in FHL1-highvs -low samples from TCGA (permutation tests). NES, normalized enrichment score. i. Correlation between FHL1 and LSC-related genes in AML samples from TCGA (Pearson correlation coefficient). j-l. Correlation between FHL1 and key genes in the Wnt, MAPK and JAK/STAT signalling pathways (Pearson correlation coefficient). m, n. Analysis of expression of drug outward transport genes, ABCC1 and ABCC4, in the FHL1-high and -low groups in TCGA (Student's t-tests; ABCC1 p-value=0·00,037, ABCC4 p-value=0·015). o. Correlation between FHL1 and SLC29A1, which transfers chemotherapeutic drugs into cells, in relapsed AML cases in GSE106291 (Pearson correlation coefficient; p-value=0·0069).

FHL1 is highly expressed in LSCs, and its targeted intervention enhances the cytotoxic effect of cytarabine

To further clarify the function of FHL1 in LSCs and chemotherapy resistance of AML, we compared expression of FHL1 in FAB-M0, also known as minimally differentiated acute myeloblastic leukaemia, with that in FAB-M4 and -M5, which have higher degrees of differentiation (Fig. 6a). A previous study in normal haematopoiesis reported that FHL1 is overexpressed in umbilical cord blood (CB)-derived primitive haematopoietic stem progenitor cells (HSPCs) compared with their progeny cells [38]. Therefore, we explored whether FHL1 is related to LSCs and found that FHL1 expression was increased in cells with a high 17-gene LSC score (LSC17), a prognostic biomarker related to stemness (Fig. 6b) [24]. Moreover, expression of FHL1 in AML CD34+ cells was higher than that in AML CD34- cells (Fig. 6c). FHL1 is highly expressed in LSCs and may participate in the regulation of LSCs, and the high recurrence rate of AML is attributed to the persistence of LSCs, which have many stem cell characteristics associated with therapeutic resistance, such as quiescence [24,39]. Our study and another study based on high-throughput screening show that FHL1 expression is elevated at initial diagnosis in chemotherapy-resistant AML patients compared with chemotherapy-sensitive patients (Fig. 6d) [40]. We also performed qRT-PCR of FHL1 in samples from AML-sensitive and drug-resistant patients we collected at the time of initial diagnosis and obtained the same results (Fig. 6e). In addition, we found that expression of FHL1 in the same AML patient was higher at the time of relapse than at the first diagnosis (Fig. 6f). Therefore, FHL1 may be involved in the chemotherapy resistance and relapse of AML.

Fig. 6

FHL1 is highly expressed in LSCs and may be involved in the chemotherapy resistance and relapse of AML. a. Comparison of FHL1 expression for FAB-M0 (n  =  15) and FAB-M4/M5 (n  =  44) subtypes in TCGA (Mann–Whitney test; p-value<0·0001). b. Difference in FHL1 expression between patients with a higher (n  =  76) and a lower (n  =  75) LSC17 score (Mann–Whitney test; p-value<0·0001). c. Analysis of FHL1 expression in normal CD34+ cells (n  =  31), CD34+ cells (n  =  46) and CD34- cells (n  =  44) in AML patients (Mann-Whitney tests; normal CD34+ cells vs. AML CD34+ cells p-value = 0·075, AML CD34+ cells vs. AML CD34- cells p-value<0·0001). d. Differential expression of FHL1 at initial diagnosis in patients with chemotherapy-sensitive (n  =  164) or drug-resistant (n  =  71) disease in the GSE106291 dataset (Student's t-test; p-value = 0·003). e. Expression of FHL1 in chemotherapy-sensitive (n  =  15) or drug-resistant AML cases (n  =  13) at the initial diagnosis was further verified using qRT-PCR (Mann–Whitney test; p-value = 0·0022). f. Comparison of FHL1 expression in paired de novo and relapsed AML samples (n  =  19) in GSE83533 (paired t-test; p-value = 0·017). g, h. OS and EFS curves of patients who underwent induction chemotherapy in TCGA grouped by expression of FHL1 (log-rank tests; OS p-value = 0·0071, EFS p-value = 0·0083). i. The OS curve of AML patients with bone marrow transplantation according to dichotomized expression of FHL1 in TCGA (log-rank test; p-value = 0·00,081). j. FHL1 mRNA levels in Kasumi-3 and U937 cells transduced with FHL1 shRNAs (Welch's t-tests). k, l. Percentage of surviving cells after transduction with control or FHL1 shRNAs and treatment with different concentrations of cytarabine (100 nM and 500 nM) or without cytarabine (DMSO) in Kasumi-3 and U937cells (Welch's t-tests). m. FHL1 mRNA levels in Kasumi-1 and U937cells transduced with overexpression plasmid of FHL1 (Welch's t-tests). n, o. Percentage of surviving cells after transduction with vector or overexpression plasmid of FHL1 and treatment with different concentrations of cytarabine (100 nM and 500 nM) or without cytarabine (DMSO) in Kasumi-1 and U937cells (Welch's t-tests). P-values were shown in the figures above.

We divided AML patients receiving cytarabine-based standard induction chemotherapy into two groups according to FHL1 expression at the first diagnosis, and survival analysis showed that the FHL1-high expression group had shorter OS and EFS (Fig. 6g, h). Interestingly, de novo AML patients with high FHL1 expression also had worse predicted outcomes after transplantation (Fig. 6i). We knocked down expression of FHL1 in the AML cell line Kasumi-3 and U937 (Fig. 6j), and detected the effects on the sensitivity of these cells to cytarabine. According to the results, knocking down FHL1 with two different shRNAs increased the cytocidal effect of cytarabine on AML cells compared with the control group (Fig. 6k, l). We overexpressed FHL1 in AML cells and found that FHL1 can enhance the survival of AML cells with cytarabine treatment (Fig. 6m-o). The results suggested that targeted intervention of FHL1 might enhance the sensitivity of AML cells to cytarabine. Therefore, the expression level of FHL1 at initial diagnosis may be useful for prognostic evaluation and targeted interventions for AML patients with chemotherapy.

Discussion

In this study, we used genome-wide expression data and clinical data from three independent AML and CN-AML datasets to screen for genes with significant and general prognostic value and identified several potential prognostic indicators for AML. The prognostic value of HOPX has been confirmed in previous studies. For example, Lin et al. showed that higher HOPX expression is associated with distinct features and predicts poor prognosis in de novo AML [13]. HOPX also acts as a member of a 4-gene expression prognostic signature which might guide post-remission therapy in patients with intermediate-risk cytogenetic acute myeloid leukaemia and participates primitive haematopoiesis [41,42]. A genome-wide association study of anorexia nervosa showed that FAM124B related to anorexia nervosa, but its function in tumours, especially in leukaemia, requires more research [43]. Combining the results of different cohorts, we found that FHL1 is most closely related to the prognosis of AML, though the role of FHL1 in AML remains unclear. Therefore, we further explored the functions and mechanisms of FHL1 in AML.

FHL1 belongs to the FHL protein family, which is characterized by a combination of four and a half highly conserved LIM domains and acts as an important mediator of protein-protein interactions. FHL1 is highly expressed in cardiomyocytes and skeletal muscle and participates in various cardiovascular diseases and skeletal myopathies [44], [45], [46]. FHL1 has long been recognized as a tumour suppressor gene and is downregulated in various tumours, such as liver cancer, lung cancer and breast cancer, inhibiting cell proliferation, invasion and metastasis [47], [48], [49]. However, a study by Wang et al. showed that cytosolic tyrosine kinase Src-phosphorylated FHL1 promotes cell proliferation, converting FHL1 from a tumour suppressor to a tumour promoter [50]. In addition, Xu et al. confirmed that FHL1 expression correlates positively with radioresistance in cancer patients [51]. Therefore, the function of FHL1 in tumours is complex and controversial.

We evaluated the prognostic value of FHL1 in 1298 de novo AML patients from four independent datasets and 407 CN-AML patients from four independent CN-AML datasets. The results showed that high expression of FHL1 in all of the above cohorts indicated a worse OS. Similar results also obtained in non-APL group. Studies of EFS and RFS also suggested that high expression of FHL1 was associated with worse EFS and RFS. To investigate whether the effects of FHL1 on prognosis are independent of other disease-related factors, we performed multivariate Cox regression analyses for patients with total de novo AML, non-APL AML and CN-AML, and the results confirmed that FHL1 has prognostic value independent of increasing age, higher WBC count, adverse cytogenetic risk and prognosis-related gene mutations. Several reports on the prognostic evaluation of AML based on gene expression profiling have shown good predictive ability. In this study, we compared the efficacy of FHL1 with the reported 3-gene, 7-gene, and 24-gene prognostic scoring systems for AML survival assessment and found that FHL1 has independent and effective assessment capabilities. Overall, the above results indicate the strong positive prognostic impact of high FHL1 expression in AML, which might be used as an independent prognostic indicator for non-APL AML and CN-AML patients.

Based on comprehensive analysis of FHL1 expression levels and clinical characteristics of AML patients, we found that expression of FHL1 in minimally differentiated AML (AML-M0) was higher than that in highly differentiated AML (AML-M4/5). FHL1 expression in HSPCs has been reported to be increased compared with that in their progeny cells [38], but the relationship between FHL1 and LSCs remains unclear. We compared the levels of FHL1 in LSCs and non-LSCs from AML and found the FHL1 level in AML LSCs was higher than that in non-LSCs of AML. Through a comprehensive analysis of the correlation between WGCNA, GSEA, and the relatedness of FHL1 with well-known stem cell genes, we identified a significant positive correlation between FHL1 and LSCs.

We and Heuser et al. observed that FHL1 is elevated in patients with AML-resistant disease [40]. Our results also show that high expression of FHL1 predicts worse clinical outcomes in patients with AML who receive cytarabine-based induction chemotherapy; although the mechanism remains unclear, regulation of LSCs by FHL1 may a reason. In addition, the results of WGCNA and GSEA suggested that FHL1 may be involved in the transcellular transport of chemotherapeutic drugs, as represented by cytarabine. The FHL1-high expression group showed upregulation of ABCC1 and ABCC4 genes, which transfer drugs out of cells, whereas expression of FHL1 in the chemotherapy-resistant AML group correlated negatively with the SLC29A1 gene, which is responsible for the transfer of cytarabine into AML cells [35,37]. Therefore, FHL1 may partly mediate AML resistance through regulation of the transmembrane transport of chemotherapeutic drugs, but its more in-depth mechanism needs further exploration. Our cellular experiments also demonstrated that knockdown of FHL1 enhances the sensitivity of AML cells to cytarabine, suggesting that targeted intervention of FHL1 may be used as a potential auxiliary to traditional chemotherapy in AML patients.

In conclusion, this study identified FHL1 as a powerful prognostic indicator independent of and complementary to existing clinical or genetic factors for prognostication in AML, and its high expression suggests a worse survival and chemotherapy response. Furthermore, knockdown of FHL1 enhances the sensitivity of AML cells to cytarabine. Hence, high expression of FHL1 may serve as an evaluation factor for clinical strategy selection, and its targeted intervention may be beneficial for the treatment of AML patients, especially those with chemotherapy-resistant disease.

Declaration of Competing Interest

The authors declare no conflicts of interest.