ALDH2 promotes uterine corpus endometrial carcinoma proliferation and construction of clinical survival prognostic model

Data collection

Expression data and clinical annotation data of 33 types of cancer patients were downloaded from the TCGA data portal (https://gdc.cancer.gov). Perl software (Version 5.3.2) was used in the 33 types of cancer to sort the data to analyze the data for each patient. The expression profile data of normal human EC is obtained from the GTEx database (http://www.gtexportal.org). All database data were accessed in September 2020.

Analysis of overall survival

R with survival and survminer packages used the Kaplan-Meier method and Cox regression to analyze survival between different groups. P-value of less than 0.05 is considered statistically significant.

Analysis of different genes

R software with limma package was used to analyze the different genes between the data of each group, and log2 foldchange >0.5 and adjust P value <0.05 was used as the standard.

Analysis of gene and pathway enrichment

To determine the biological processes, cell components, molecular functions and biological pathways of differential gene enrichment, R software with clusterProfiler and org.Hs.eg.db packages were used for gene ontology (GO) enrichment and Kyoto Gene and Genome Encyclopedia The whole book (KGGE) approach analysis [17–19]. At the same time, Gene Set Enrichment Analysis (GSEA) was used to analyze related functions among different groups. The screening criterion is P < 0.05.

Construction of protein-protein interaction (PPI) network

The STRING database was used to construct a protein interaction network of differential genes [20]. The screening criterion was a comprehensive score of ≥0.9. Cytoscape software with cytoHubba and MCODE apps was used to visualize and screen highly critical genes and core sub-networks.

Model construction and verification

To construct the prognostic risk scoring system for UCEC, univariate Cox regression analysis was used to determine the prognostic gene of DEGs. P < 0.05 gene was considered as significant. Then, Lasso penalized Cox regression analysis was used to further select OS-related prognostic genes in UCEC. Finally, UCEC patients in the TCGA database were randomly divided into training sets and validation sets with equal numbers, and univariate Cox analysis was used to regress to establish a risk scoring model gradually. The C-index, the AUC value of the model, is calculated by the validation set and the verification of the Nomo diagram to test the accuracy of the risk score model.

Tumor immunity relevance of the ALDH2

Regarding the close association between the ALDH2 and immune cells of the UCEC microenvironment. Perl software is used to perform statistical analysis on the data downloaded by TCGA. The CIBERSORT deconvolution algorithm (http://cibersort.stanford.edu/) is used to estimate the proportion of different immune cells in cancer tissues, and the differences between immune cells in UCEC and normal adjacent tissues are analyzed. This tool was developed by Newman et al. [21] and has been verified to successfully quantify the abundance of specific cell types.

Cell culture

Three different cancer cell lines and 293T were purchased from ATCC (ATCC, USA). The highly differentiated EC cell line ISHIKAWA and the moderately differentiated EC cell line SPEC-2 were cultured with MEM medium (GIBICO, USA) containing 10% fetal bovine serum (BI, Australia), The poorly differentiated EC cell line KLE was cultured in DMEM/F12 (GIBICO, USA) mixed medium containing 10% fetal bovine serum. The 293T cell line and endometrial epithelial cell line hEEC were cultured in DMEM (GIBICO, USA) mixed medium containing 10% fetal bovine serum. All cells are cultured in a 37°C incubator containing 5% CO₂.

CCK-8 assay

The cells suspension was inoculated in a 96-well plate and continued to culture for 12 h, 24 h, 48 h, 72 h, and then measured cell proliferation ability with the Cell Counting Kit 8 (CCK-8) Kit (Dojindo, Japan). The absorbance values at OD 490 nm were measured using a plate reader (Biorad, USA).

Colony formation assay

The cells were seeded in a 6-well plate, and when they were cultured normally to form visible clones, the cells were fixed with 4% formaldehyde, stained with purple crystal (Meilunbio, China) and photographed for analysis.

Xenograft tumor in nude mice

Nude mice were randomly divided into the Control group and the ALDH2 overexpression group. Nude mice in each group were subcutaneously injected with 0.2 ml of UCEC cell suspension (1 × 10⁷ cells/ml) on the ventral side of the right hind limb. On the 21th day, the nude mice were sacrificed by cervical dislocation, the subcutaneous tumors were completely removed, and weighed with an electronic balance. The Wuhan University of Science and Technology Animal Ethics and Use Committee approved the tumor-forming experiment in nude mice.

Dual-luciferase reporter assay

3′-UTR of ALDH2 and mutant form were separately subcloned into a pmirGLO (Addgene, USA) vector to establish wt-ALDH2-luc and mut-ALDH2-luc respectively. The miRNA mimic, internal control and Luc plasmid were co-transfected into cells. Then, luciferase activities were detected 48 h after transfection by the Dual-Luciferase Reporter Assay System (Promega, USA).

RNA Pull-Down

RNA pull-down assay was carried out using Magnetic RNA-Protein Pull-Down Kit (Thermo, USA). The RNA-bound beads were added to the cell nuclear lysate. Then, the eluted proteins were detected by western blot analysis.

RNA Immunoprecipitation (RIP)

RNA-protein-antibody complexes were captured using Protein A/G (Thermo, USA). RNA was eluted by adding TRIzol directly to magnetic beads and isolated as per the manufacturer’s instructions. cDNA was synthesized using HiScript^® II 1st Strand cDNA Synthesis Kit (Vazyme, China) and analyzed by qRT-PCR.

RNA isolation and qRT-PCR

RNeasy Plus Universal Mini Kit (Qiagen, USA) was used to isolate total RNA from cell lines, and HiScript^® II 1st Strand cDNA Synthesis Kit (Vazyme, China) was used to reverse transcribe into cRNA. At the same time, the miRNeasy Micro Kit (Qiagen, USA) was used to isolate miRNA from cell lines, and the miRNA 1st Strand cDNA Synthesis Kit (Vazyme, China) was used to reverse the transcription of miRNA. For the extracted RNA and miRNA, qRT-PCR was performed on the Bio-Rad CFX-96 (Biorad, USA) system using the SYBR Green (Yisen, China) method to determine the relative expression level. β-actin and U6 were used as an endogenous control. The primer sequences of ALDH2 are as follows: ALDH2-F: 5-ATGGCAAGCCCTATGTCATCT-3, ALDH2-R: 5-CCGTGGTACTTATCAGCCCA-3. The synthesis of primers, plasmid sequencing, miRNA reverse transcription sequences, primers and probes are all designed by Ribobio (Ribobio, China).

Western blot

Western blot assay was performed according to the standard protocol. The ALDH2 antibody was purchased from CST (CST, USA), and the antibody was diluted according to the recommended ratio in the instructions.

Statistical analyses

Comparisons between groups and normality were performed using R 3.6.3 software. Comparisons were completed using a one-way analysis of variance (ANOVA), two-tailed Student’s t-test, non-parametric tests. Kaplan–Meier analyses with log-rank tests and the Cox proportional hazard model were used to analyze for survival. All data are presented as the mean ± SD. Statistical significance was set at P < 0.05.

Ethics approval and consent to participate

All mice experimental procedures and methods were evaluated and authorized in strict accordance with the guiding principles of the Animal Protection and Use Committee of Wuhan University of Science and Technology and accordance with the “Hubei Province Experimental Animal Management Regulations.”

Availability of data and material

The data generated during this study are included in this article and its supplementary, information files are available from the corresponding author on reasonable request.

The expression of the ALDH family in all cancers of TCGA

Through the TCGA database, the expression data of all ALDH families in all cancers were obtained. As shown in Figure 1, the expression genes of the ALDH family genes in cancer were very different (Figure 1A), and the expressions of the ALDH family genes were highly correlated (Figure 1B). In detail, ALDH9H1, ALDH18A1, ALDH3A2, ALDH1A1, ALDH1B1 and ALDH2 were expressed higher than other ALDH family genes in all cancers. The expression of ALDH2 and ALDH8A1 was highly positively correlated. The first six ALDH family genes were highly expressed in cancer was further analysis. According to the classification of normal tissues and tumor tissues, it was found that the expression of ALDH3B1, ALDH18A1, etc., increased in a variety of cancers, but ALDH9A1 and ALDH2 showed the opposite trend (Figure 1C) in UCEC. Except for ALDH3B1, which was statistically significant, the expression of the other 5 ALDH family members decreased in tumor tissues. Taken together, our data showed that the expression of ALDH was altered in a variety of cancers compared with normal tissues.

Figure 1. The expression of ALDH family in all cancers of TCGA. (A) The expression level of patients of ALDH family in 33 types of cancer in TCGA. (B) The correlations among the expressions of the ALDH family. (C) The relative expression levels of ALDH1A1, ALDH1B1, ALDH2, ALDH9H1, ALDH3B1, and ALDH18A1 in the TCGA database compared with normal tissues. (The number of normal tissue samples less than 3 was hidden). ^*P < 0.05, ^***P < 0.001.

The prognostic significance of the ALDH family

To further confirm the effect of ALDH family genes that were highly expressed in cancer on UCEC patients, we retrieved the survival information of each patient in the TCGA database, combined with the expression analysis of ALDH family genes. Changes in other ALDH family genes did not affect the overall survival of UCEC patients. Kaplan-Meier analysis results showed that only the expression level of ALDH2 and ALDH18A1 mRNA had a significant impact on the OS of UCEC patients (P = 0.003 to ALDH2, P = 0.017 to ALDH18A1). Compared with the low expression group, the overall survival time of UCEC patients in the high expression group was significantly prolonged in the lower expression group (Figure 2A). At the same time, the survival of DSS, DFI, and PFI of the low expression group of ALDH2 was the same as the overall survival result, all of which showed a worse survival rate (Figure 2B). For DSS in the ALDH18A1 low expression group, DFI and PFI were not different from the high expression group (Figure 2C). These results indicated that patients in the ALDH2 low expression group were associated with worse prognostic survival in different survival groups. In addition to survival time and survival status, TCGA data also contained complete clinical information of the patient. The results of the clinical analysis revealed that the expression of ALDH2 was not significantly correlated with the age and race of the patients (Figure 2D).

Figure 2. The prognostic value of the ALDH family. (A) Kaplan–Meier survival curves by OS. (B) Kaplan–Meier survival curves by PFI, DSS, and DFI of ALDH2. (C) Kaplan–Meier survival curves by PFI, DSS, and DFI of ALDH18A1. (D) The relationship between clinicopathological characteristics and the expression of ALDH2.

Analysis of GO, KGGE, GSEA function enrichment

In this study, we downloaded the expression data of all UCEC samples in the TCGA database. The deletion was divided into 2 groups according to the median value and the expression of ALDH2, and the different genes between the groups were analyzed. A total of 477 differential genes were identified. Among them, 322 genes were up-regulated, and 155 genes were down-regulated (Figure 3A, 3B). The enrichment of GO and KEGG helps us understand the potential regulatory functions of these differential genes. The six pathways with the highest correlation found in GO enrichment are, humoral immune response mediated by circulating immunoglobulin, complement activation, classical pathway, complement activation, immunoglobulin mediated immune response, B cell-mediated immunity, and protein activation cascade (Figure 3C). The 6 pathways with the highest correlation of KGGE enrichment are, Cell adhesion molecules, Type I diabetes mellitus, Intestinal immune network for IgA production, Hematopoietic cell lineage, Th1 and Th2 cell differentiation, and Phagosome (Figure 3D). GSEA enrichment showed that the GO processes of the main ALDH2 up-regulated group included Azurophil Granule Lumen, Response to Interferon Gamma, Myeloid Leukocyte Migration, Mmune Receptor Activity, Tetrapyrrole Metabolic Process (Figure 3E). The GO processes of the down-regulation group included Mitochondrial Genome Maintenance, Regulation of DNA Recombination, Cell Cycle Checkpoint, DNA Replication Checkpoint, Sulfur Amino Acid Biosynthetic Process (Figure 3F). The above enrichment results indicated that the alteration of ALDH2 expression might be related to cellular immunity. Therefore, the immune infiltration score of each UCEC patient was calculated by R software with CIBERSORT package. The results showed that patients in the ALDH2 low expression group had lower scores for CD8⁺ T cells and plasma cells (Figure 3G).

Figure 3. Analysis of GO, KEGG, GSEA function enrichment. (A) Heatmap showing differential gene expression (FDR < 0.05) between low expression group and high expression group of ALDH2. (B) Volcano plot of all differential gene expression analysis. (C) GO enrichment analysis. (D) KEGG enrichment analysis. (E) GSEA enrichment analysis of high expression group of ALDH2. (F) GSEA enrichment analysis of low expression group of ALDH2. (G) Enrichment scores for 22 immune cell subpopulations based on deconvolution by CIBERSORT between low expression group and high expression group of ALDH2. ^*P < 0.05, ^***P < 0.001.

PPI network construction

PPI of DEGs was constructed by using the String network tool. After hiding all the individual nodes, it was found that there were 392 nodes and 295 edges (Figure 4A). Statistics showed that the first six genes with the number of interaction relationships were ORM2, C3, CTSH, HLA-DPA1, HLA-DPB1 and HLA-DQA1 (Figure 4B). At the same time, the MCODE app was used to analyze the most significant interaction modules. The module network consisted of 24 nodes (Figure 4C). In addition, the hub node was calculated through cytoHubba app (Figure 4D). These genes had a potential relationship with ALDH2. It played the role of an oncogene in the occurrence and development of UCEC, which was related to the occurrence, development and prognosis of UCEC.

Figure 4. PPI network construction. (A) PPI network of al DEGs. The individual nodes were hidden. The interaction relationship prediction threshold is > 0.900. (B) Top 30 genes with the highest number of nodes. (C) The core subnet of the PPI network by using the MCODE app. (D) The core gene of the PPI network by using the cytoHubba app.

The risk score model of UCEC construction

Our data showed that the high expression ALDH2 group had a better OS, DSS, and PFI (Figure 2). Therefore, a prognostic model based on the overall survival score of ALDH2 expression was constructed. The expression profile of UCEC in the TCGA database was randomly divided into an equal number of training groups and test groups. COX univariate results revealed that a total of 62 genes were significantly related to the overall survival of UCEC patients. Further use lasso regression to screen out 15 genes to prevent the model from overfitting (Figure 5A, 5B). Finally, a 6-factor overall prognostic risk score model was established using COX stepwise regression. The coefficient of risk score model of overall survival was: CDKN2A × 0.18793 + WNT10A × 0.23536 + AQP5 × (−0.21242) + STX18 × (−0.30617) + GZMA × (−0.44271) + LCN2 × (−0.15148). COX multifactorial analysis found that they were all highly significantly related to prognosis (Figure 5C).

Figure 5. The risk score model of UCEC construction. (A) DEGs were identified by the LASSO logistic regression model with non-zero coefficients. (B) The graph of the relationship between likelihood deviation and log (Lambda). The vertical dashed lines indicate the λ value and the maximum λ value with the smallest error. (C) Forest plot of the hazard ratio for OS of parameters. (D) Predicted risk of overall survival by the risk score model. (E) Scatter plot of UCEC death predicted by risk score model. (F) The survival curve of the high expression group and the low expression progenitor in the risk score model. (G) ROC curve of risk score model. (H) Predicts the OS of patients with the Nomogram. (I) Calibration curve of the nomogram in the training dataset.

Next, use the test group to test the model. The risk scores of 268 UCEC patients in the verification group were shown in the figure (Figure 5D). The state chart showed that the number of deaths in the high-risk group was more (Figure 5E). Moreover, the survival prognosis of the high-risk group was worse (Figure 5F). Besides, the C-index of the model was calculated to be 0.798. The AUC value of 1 and 3 years was 0.738, and the AUC value of 5 years was 0.764 (Figure 5G). These results showed that the model had a better ability to predict. Therefore, a nomogram (Figure 5H) was gained according to the risk score model, and the verification result was shown in Figure 5I.

Effect of ALDH2 Overexpression on Tumor Progression in vitro and in vivo

To determine whether ALDH2 expression was reduced in EC cell lines, it was detected by Western Blot and qRT-PCR. Compared with normal human endometrial epithelial cells hEEC, the expression of ALDH2 was lower in ISHIKAWA, SPEC-2, and KLE (Figure 6A, 6B). In the meantime, a lentiviral plasmid pLKO.1-ALDH2 was constructed to overexpressed ALDH2. The overexpression efficiency was verified by Western Blot and qRT-PCR (Figure 6C, 6D). Next, the relationship between the expression of ALDH2 and the proliferation ability of EC was verified by CCK-8 and colony formation analysis (Figure 6E, 6F). The results showed that restoring the overexpression of ALDH2 reduced the proliferation ability of EC cell lines. Subsequently, a subcutaneous tumor model was used to assess the ability of tumor genesis and growth. As shown in Figure 6G, the tumors in overexpression of the ALDH2 group were significantly smaller than those of the control group (Figure 6G).

Figure 6. Effect of ALDH2 overexpression on tumor progression in vitro and in vivo. (A) Expression of ALDH2 in the UCEC cells by Western Blot. (B) Expression of ALDH2 in the UCEC cells by qRT-PCR. (C) Overexpression efficiency of ALDH2 by Western Blot. (D) Overexpression efficiency of ALDH2 by qRT-PCR. (E) CCK-8 assays OD 450 nm. (F) Colony formation ability was determined using colony formation assays. (G) Tumor xenografts from nude mice subsequent assays. ^*P < 0.05, ^**P < 0.01, ^***P < 0.001.

Hsa-miR-135b-3p binding the 3′UTR of ALDH2

In the previous analysis, the OS of UCEC patients in the ALDH2 low expression group was worse (Figure 2A). Therefore, whether ALDH2 genetic mutations were the cause of the decreased OS in UCEC patients was investigated. Analysis of gene mutation frequency was found that in most patients, it was not the mutation that caused the disappearance of ALDH2 effect (Figure 7A, 7B).

Figure 7. Hsa-miR-135b-3p binds to the 3′UTR of ALDH2. (A) Gene mutation distribution of UCEC patients detected in TCGA. (B) KM survival curve of the mutant group and non-mutation group. (C) The volcano plot of miRNAs in UCEC patients with TCGA: Red dots and green dots indicate differentially expressed miRNAs based on the fold change P < 0.05. (D) Venn diagram of up-regulated miRNA and predicted binding miRNA. (E) Correlation analysis between the expression of predicted bound miRNA and the expression of ALDH2. (F) Relative Luciferase activity was measured with a dual-luciferase reporter assay. (G) Dual-luciferase assay of the mutation group. (H) RNA pull-down analysis with ALDH2 antibody. (I) RIP assay was further verified for a direct association between hsa-miR-135-3p and ALDH2. (J) Overexpression efficiency of hsa-miR-135-3p mimic. (K) WB analysis the expression of ALDH2 with overexpression of hsa-miR-135-3p. (L) Relative expression of ALDH2 with overexpression of hsa-miR-135-3p by using Image J software. (M) Kaplan-Meier survival curve for patients with UCEC (According to the median value of aldh2 expression, it is divided into two groups of high and low expression). ^*P < 0.05, ^**P < 0.01, ^***P < 0.001.

Numerous studies had proved that miRNA could silence gene expression by targeting gene 3′UTR, consequently, it was speculated that ALDH2 might also be regulated by this mechanism. By screening the difference in miRNA between UCEC patients and the normal group, the volcano plot showed that a total of 122 up-regulated miRNAs were screened (Figure 7C). At the same time, the bioinformatics prediction website TargetScan was also visited. As a result, 312 miRNAs were predicted to bind to ALDH2’s 3′UTR. A total of 10 miRNAs were selected from the intersection, and the correlation between their expression and ALDH2 expression was shown in the figure (Figure 7D, 7E). The top three miRNAs with the highest negative correlation were hsa-miR-301b-5p, hsa-miR-3187-3p and hsa-miR-135b-3p.

To verify which miRNA regulated the expression of ALDH2, the Luciferase plasmid of ALDH2 3′UTR was constructed. The dual fluorescein report experiment results showed that only the hsa-miR-135b-3p transfection group had a decrease in relative fluorescence intensity. However, there was no significant difference between hsa-miR-301b-5p and hsa-miR-3187-3p after transfection (Figure 7F). The mutation group showed that the mutation could eliminate the inhibitory effect of hsa-miR-135b-3p on ALDH2 (Figure 7G). Simultaneously, RIP and RNA pull-down found that hsa-miR-135b-3p could enrich the expression of ALDH2, but the control group did not have this phenomenon (Figure 7H, 7I). QRT-PCR showed that after transfection of hsa-miR-135-3p mimic, the expression of hsa-miR-135b-3p was enhanced significantly (Figure 7J), while the results of WB and qRT-PCR showed that the expression of ALDH2 was decreased (Figure 7K, 7L). Besides, TCGA survival data showed that the high expression group of hsa-miR-135b-3p had a lower prognosis (Figure 7M). These results indicated that hsa-miR-135b-3p could down-regulate the expression of ALDH2.