Retrieval of data from the Genotype-Tissue Expression (GTEx) project database and The Cancer Genome Atlas (TCGA)

Retrieval of data from TCGA (16) and GTEx (17) was performed using the University of California, Santa Cruz (UCSC) Xena platform (https://xenabrowser.net/datapages/?cohort=TCGA%20TARGET%20GTEx&removeHub=https%3A%2F%2Fxena.treehouse.gi.ucsc.edu%3A443) (18). RNA sequencing (RNA-seq) data presented as log₂(transcript per million +0.001) values were used and compared. Survival data, including progression-free interval (PFI; restricted to 56 months, the longest time without late-stage cross-over) and overall survival (OS) data, were obtained for primary LUAD cases in TCGA (n=472 with PFI data and n=510 with OS data) for Kaplan-Meier survival analysis. Median gene expression was used for the low and high MRPL51 expression cut-off.

MRPL51 protein expression data in the Human Protein Atlas (HPA)

MRPL51 expression at the protein level in normal lung and LUAD tissues was examined using immunohistochemistry (IHC) images deposited in HPA (19,20).

Cell culture and treatment

The A549 and Calu-3 human LUAD cell lines were obtained from the Cell Resource Center, Peking Union Medical College [which is the headquarters of the National Science & Technology Infrastructure-National BioMedical Cell Line Resource (NSTI-BMCR)]. 293T cells were obtained from Procell Life Science & Technology Co., Ltd. In brief, all these cell lines were cultured in Dulbecco's Modified Eagle Medium (DMEM; Gibco; Thermo Fisher Scientific, Inc.), supplemented with 10% fetal bovine serum (FBS; Gibco; Thermo Fisher Scientific, Inc.) and 1% penicillin-streptomycin solution (100 U/ml penicillin and 100 µg/ml streptomycin). Cells were cultured in a humidified environment with 5% CO₂ at 37°C.

Lentiviral MRPL51 and FOXM1 short hairpin (sh)RNAs were produced using the pLKO.1-puro vector (Addgene, Inc.). The following shRNA sequences were used: shMRPL51#1, 5′-CGCATCCGCTATCTCTACAAA-3′; shMRPL51#2, 5′-GTTCGGAGTGTATGACAACAT-3′; shMRPL51#3, 5′-CTTAATAAACGCATCCGCTAT-3′; shFOXM1#1, 5′-GCCCAACAGGAGTCTAATCAA-3′; shFOXM1#2, 5′-GCACTATCAACAATAGCCTAT-3′; and shFOXM1#3, 5′-CGCCGGAACATGACCATCAAA-3′. A scramble sequence was used as the negative control (NC): shNC, 5′-CCTAAGGTTAAGTCGCCCTCG-3′.

The lentiviruses used for infection were generated using the 2nd generation system, which included the psPAX2 packaging plasmid, the pMD2.G envelope plasmid and 293T cells. In total, 1 µg pLKO.1 shRNA plasmid, 750 ng psPAX2 packaging plasmid and 250 ng pMD2.G envelope plasmid were resuspended in 20 µl serum-free OPTI-MEM (Gibco; Thermo Fisher Scientific, Inc.). Transfection into 293T cells was performed using FuGENE 6 transfection reagent (Promega Corporation) at room temperature for 30 min, according to the manufacturer's instruction. Media from the cell culture was harvested and centrifuged at 245 × g for 5 min at 4°C. Then, the media was filtered through a 0.45 µm filter to remove the cells.

Reverse transcription-quantitative PCR (RT-qPCR)

The knockdown efficiency of lentiviral MRPL51 and FOXM1 shRNAs in A549 and Calu-3 cells was examined 48 h after infection (MOI=10) by RT-qPCR. Total RNAs were extracted from cells using TRIzol Reagent (Thermo Fisher Scientific, Inc.). cDNA was reversely transcribed using the High-Capacity cDNA Reverse Transcription Kit (Thermo Fisher Scientific, Inc.) according to the manufacturer's instruction. qPCR was performed using FastStart Universal SYBR Green Master Mix (Roche Diagnostics), with a 7900HT Fast Real Time PCR System (Thermo Fisher Scientific, Inc.). The following primers were used: MRPL51 forward, 5′-GATCGTTGGAACGAGAAAAGGGC-3′ and reverse, 5′-CCCTTTCCAACCTCGAAGCCAT-3′; FOXM1 forward, 5′-TCTGCCAATGGCAAGGTCTCCT-3′ and reverse, 5′-CTGGATTCGGTCGTTTCTGCTG-3′; and β-actin (ACTB) forward, 5′-CACCATTGGCAATGAGCGGTTC-3′ and reverse, 5′-AGGTCTTTGCGGATGTCCACGT-3′. The following thermocycler conditions were used: Stage 1, 50°C for 2 min and 95°C for 10 min; Stage 2, 95°C for 15 sec and 60°C for 1 min, 40 cycles; Stage 3 (dissociation curve), 95°C for 15 sec, 60°C for 15 sec and 95°C for 15 sec. Relative gene expression was calculated using the 2^−ΔΔCq method (21), and ACTB was used as the control gene.

Western blotting

In brief, A549 and Calu-3 cell samples were collected and lysed using RIPA lysis buffer following the manufacturer's instruction (Beyotime Institute of Biotechnology). The concentration of proteins in the supernatant was measured using a BCA kit (Beyotime Institute of Biotechnology). A total of 25 µg per lane of each sample was separated by 10% SDS-PAGE and transferred onto nitrocellulose membranes. The membranes were washed with 1X TBST (0.1% Tween-20 in 1X TBS), blocked with 1X TBST containing 5% non-fat dry milk for 1 h at room temperature, and incubated with primary antibodies at 4°C overnight. The following primary antibodies and dilutions were used: Anti-MRPL51 (1:1,000; cat. no. PA5-58988; Thermo Fisher Scientific, Inc.), anti-E-cadherin (1:1,000; cat. no. 20874-1-AP; Proteintech Group, Inc.), anti-N-cadherin (1:1,000; cat. no. 22018-1-AP; Proteintech Group, Inc.), anti-Vimentin (1:1,000; cat. no. 10366-1-AP; Proteintech Group, Inc.), anti-FOXM1 (1:1,000; cat. no. 13147-1-AP; Proteintech Group, Inc.) and anti-β-tubulin (1:1,000; cat. no. 10094-1-AP; Proteintech Group, Inc.). Subsequently, the membranes were washed with 1X TBST and incubated with HRP-conjugated AffiniPure goat anti-rabbit IgG (H+L) (1:5,000; cat. no. SA00001-2; Proteintech Group, Inc.) for 1 h at room temperature. The protein bands were developed using BeyoECL Star reagent (Beyotime Institute of Biotechnology) with a ChemiDoc Imaging System (Bio-Rad Laboratories, Inc.). β-tubulin served as the endogenous protein expression control.

Cell Counting Kit-8 (CCK-8) assays

CCK-8 assays were conducted using a commercial CCK-8 (Beyotime Institute of Biotechnology). In brief, A549 and Calu-3 cells with or without lentiviral-mediated MRPL51 knockdown were seeded into 96-well plates (2,000 cells/well) and allowed to attach for 2 h; 2 h after seeding was considered the 0 h time point. Cell viability was measured at 0, 24, 48 and 72 h after seeding according to the manufacturer's protocol. For this, 10 µl CCK-8 solution was added to each well and incubated in a cell incubator at 37°C for 2 h. Absorption at 450 nm was measured using a microplate reader. All CCK-8 experiments were repeated three times with three technical repeats.

Immunofluorescent staining

In brief, A549 and Calu-3 cells (1×10⁵ cells per well/24-well plate) grown on coverslips were fixed with 4% paraformaldehyde for 15 min at room temperature, permeabilized using 0.1% Triton X-100 in PBS for 15 min at room temperature, and blocked using 2% BSA (Beyotime Institute of Biotechnology) in PBST (0.1% Tween-20) for 30 min at room temperature. For MRPL51 and mitochondria double labeling, cells were incubated with MRPL51 rabbit polyclonal antibody (1:400; cat. no. PA5-58988; Thermo Fisher Scientific, Inc.) overnight at 4°C. After washing, cells were incubated with Alexa Fluor® 488-conjugated anti-rabbit IgG (H+L), F(ab')₂ (1:2,000, cat. no. 4412, Cell Signaling Technology, Inc.) for 1 h at room temperature in the dark. Mitochondria were labeled using MitoTracker Red CMXRos (Beyotime Institute of Biotechnology) according to the manufacturer's instructions. For E-cadherin and Vimentin double labeling, cells were incubated with mouse anti-Vimentin (1:500; cat. no. 60330-1-Ig; Proteintech Group, Inc.) and rabbit anti-E-cadherin (1:250; cat. no. 20874-1-AP; Proteintech Group, Inc.) overnight at 4°C. After washing, cells were incubated with Alexa Fluor® 488-conjugated anti-mouse and Alexa Fluor® 594-conjugated anti-rabbit (1:2,000; cat. nos. 4408 and 8889; Cell Signaling Technology, Inc.) secondary antibodies for 1 h at room temperature in the dark. The nucleus was stained with DAPI (1 µg/ml) at room temperature for 1 min. Fluorescence images were captured using an Olympus IX81 widefield fluorescence microscope (Olympus Corporation).

Gene Set Enrichment Analysis (GSEA)

GSEA was performed using GSEA Desktop version 4.1.0 (https://www.gsea-msigdb.org/gsea/index.jsp). In brief, the RNA-seq data of primary LUAD cases from TCGA were separated into high and low MRPL51 expression groups based on the median MRPL51 expression. The difference in gene set enrichment was explored in the Human (H): hallmark gene sets within the Molecular Signatures Database (MSigDB_v2022.1; http://www.gsea-msigdb.org/gsea/msigdb/index.jsp). The following parameters were set to run GSEA: Number of permutations, 1,000; permutation type, phenotype. Only the gene sets with |normalized enrichment score (NES)|>1, nominal (NOM) P<0.05 and false discovery rate (FDR) q<0.05 were included.

Characterization of the correlation between MRPL51 expression and cellular functional states in LUAD cells at the single-cell level

MRPL51 expression at the single-cell level was characterized using the single-cell RNA-seq dataset from Kim et al (22) based on patient-derived xenograft (PDX) tumor cells from patients with LUAD (n=126). Subsequently, its correlation with 14 single-cell functional states of each tumor cell in this dataset, including ‘stemness’, ‘invasion’, ‘metastasis’, ‘proliferation’, ‘EMT’, ‘angiogenesis’, ‘apoptosis’, ‘cell cycle’, ‘differentiation’, ‘DNA damage’, ‘DNA repair’, ‘hypoxia’, ‘inflammation’ and ‘quiescence’, was calculated. single-cell functional state scores were obtained from CancerSEA (23). In this database, the 14 single-cell functional states were scored by individual cells. Subsequently, the correlations between gene expression and the functional states were assessed by calculating the Pearson's correlation coefficients.

Transwell assay of cell invasion

A Transwell assay of A549 and Calu-3 cells (with or without MRPL51 knockdown) invasion was performed using a 6.5-mm Transwell chamber (8-µm pore size; Corning, Inc.) coated with Matrigel (1 mg/ml), according to a previously described protocol (24). In brief, 2×10⁵ cells in 100 µl serum-free DMEM were seeded into the top chamber. The bottom chamber was filled with 0.6 ml DMEM supplemented with 10% FBS as the attractant. The plate was further incubated at 37°C for 24 h. Then, the invaded cells on the bottom side were visualized using a light microscope (Shenzhen Aowei Optical Instrument Co., Ltd.) and counted manually.

Cell cycle distribution

In brief, 48 h after lentiviral infections, A549 and Calu-3 cells (1×10⁶) were collected, fixed using pre-chilled 70% ethanol (−20°C) for 3 h, washed with FACS buffer [Multisciences (Lianke) Biotech Co., Ltd.] twice and resuspended in 500 µl DAPI (1 µg/ml), according to a previously described protocol (25). Cell cycle distribution was examined using a BD FACSCalibur flow cytometer (BD Biosciences). Cell cycle distribution was analyzed using the Watson Pragmatic algorithm with NovoExpress software (v.1.5.4; Agilent Technologies, Inc.).

Dual luciferase assays

The wild-type (WT) MRPL51 promoter sequence (Fig. S1) and the mutant (MT) sequence were chemically synthesized using an automated DNA synthesizer [General Biology (Anhui) Co., Ltd.] and inserted into the pGL3-basic plasmid (Promega Corporation). A549 and Calu-3 cells with or without lentiviral-mediated FOXM1 knockdown were seeded into 24-well plates (5.0×10⁴ cells/well). After 24 h, cells were transfected with 0.5 µg of either WT or MT luciferase reporter plasmid together with the Renilla luciferase reporter plasmid pRL-TK (0.02 µg; Promega Corporation, using FuGENE 6 transfection reagent (Promega Corporation). After 24 h, cells were washed, lysed using the passive lysis buffer (Promega Corporation) and then subjected to the measurement of firefly and Renilla luciferase activities using the Dual-Luciferase Reporter Assay System (Promega Corporation). Firefly luciferase activity was normalized to Renilla luciferase activity. The experiments were performed in triplicate.

Kaplan-Meier survival analysis in the Kaplan-Meier plotter

To validate the prognostic significance of MRPL51 expression in patients with LUAD, it was meaningful to perform survival data from other databases. Therefore, validation analysis was performed using survival data collected from Kaplan-Meier plotter (http://kmplot.com/analysis/) (26), which integrated data from both TCGA and multiple datasets in the Gene Expression Omnibus (GEO). Accession numbers of the specific databases extracted from the GEO can be obtained from the supplemental materials section of the Kaplan-Meier plotter website (http://kmplot.com/analysis/index.php?p=download). To avoid analysis duplication, survival data collected from TCGA were excluded. First progression-free survival (FPS) and OS were compared between patients with the top 50% and bottom 50% of MRPL51 expression.

Chromatin immunoprecipitation (ChIP)-qPCR assay

To understand the mechanisms of aberrant upregulation of MRPL51 in LUAD tumor cells, the transcription factors with potential binding to the MRPL51 promoter were investigated using available ChIP sequencing (ChIP-seq) data in the Cistrome Data Browser (CistroDB; http://cistrome.org/db/; accession no. 53260) (27). In CistroDB, no ChIP-seq data were available for LUAD cells. Therefore, 293 cells exhibiting epithelial morphology were selected for bioinformatics prediction. A ChIP assay was performed using the SimpleChIP Plus Enzymatic Chromatin IP Kit (Agarose Beads; cat. no. 9004; Cell Signaling Technology, Inc.) according to the manufacturer's protocol. In short, 4×10⁶ cells with or without FOXM1 knockdown were collected, fixed 1% with formaldehyde for 10 min at room temperature and lysed using RIPA lysis buffer following the manufacturer's instruction (Beyotime Institute of Biotechnology). Subsequently, fragmented chromatin was prepared via partial digestion with micrococcal nuclease for 20 min at 37°C. ChIP-validated rabbit anti-FOXM1 antibody (3 µg/10⁶ cells; cat. no. 13147-1-AP; Proteintech Group, Inc.) and ChIP-Grade Protein G Agarose Beads were used. Rabbit IgG (3 µg/10⁶ cells, cat. no. 30000-0-AP; Proteintech Group, Inc.) served as the negative control. DNA samples were collected, purified and subsequently used as template for qPCR as described above. The primers covering the predicted FOXM1 binding site were designed as follows: Forward, 5′-TTTCAAGTCGGCGGAGAAT-3′ and reverse, 5′-GAACTCGGTATTTGTGGCTTTG-3′.

Statistical analysis

Statistical analysis was performed using GraphPad Prism (version 8.10; Dotmatics). For quantitative data, mean ± SD values were calculated based on three biological repeats of three technical repeats. For samples with assumed normal distribution by Shapiro-Wilk test, one-way ANOVA with Tukey's post hoc test was used for multiple group comparisons. For data without normality assumption, Kruskal-Wallis test with Dunn's post hoc test was performed. Pearson's correlation coefficient was calculated to estimate correlations. The log-rank test was conducted to estimate the difference between survival curves. P<0.05 was considered to indicate a statistically significant difference.

MRPL51 is aberrantly upregulated in LUAD tissues compared with normal lung tissues

RNA-seq data were retrieved for normal lung tissues in the GTEx database and for the LUAD dataset in TCGA. MRPL51 expression was compared among normal lung tissues (n=287), LUAD tumor-adjacent normal tissues (n=59) and LUAD tissues (n=513; Fig. 1A). The results demonstrated that MRPL51 expression was significantly higher in LUAD tissues compared with that in either normal lung or LUAD tumor-adjacent normal tissues (Fig. 1A). IHC staining data of MRPL51 expression in normal lung and LUAD tumor sections were retrieved from HPA (28), an open and accessible database providing protein expression data in human cells, tissues and organs. MRPL51 expression in alveolar epithelial cells was nearly undetectable (Fig. 1D; left). However, medium to high MRPL51 expression was observed in LUAD tissue sections (Fig. 1D; middle and right). To visualize the distribution of MRPL51 in LUAD cells, mitochondria were double-labeled using MitoTracker and MRPL51 in A549 and Calu-3 cells (Fig. 1B and C). Florescence images confirmed the localization of MRPL51 in the mitochondria found in the cytoplasm of these two cell lines (Fig. 1B and C).

GSEA results

To explore the differences in the behaviors of LUAD tumors with high and low MRPL51 expression, patients with primary tumors in TCGA-LUAD (n=513) were separated into two groups according to the median MRPL51 expression. The enrichment in the hallmark gene sets was compared between the two groups. Only the gene sets with |NES|>1, NOM P<0.05 and FDR q<0.05 were considered significantly enriched (Fig. 2A). The detailed GSEA plots indicating the significantly enriched gene sets are shown in Fig. 2B, and the peak point of the green plot is the enrichment score (ES). The results demonstrated that the high MRPL51 expression group exhibited significantly higher expression of genes (reflected by the positive ES values in Fig. 2A) enriched in multiple gene sets compared with the low MRPL51 expression group, including ‘DNA_ REPAIR’, ‘UNFOLDED_ PROTEIN_RESPONSE’, ‘MYC_ TARGETS_V1’, ‘OXIDATIVE_ PHOSPHORYLATION’, ‘MTORC1_SIGNALING’, ‘REACTIVE_OXYGEN_SPECIES_PATHWAY’, ‘MYC_ TARGETS_ V2’, ‘E2F_ TARGETS’ and ‘G2M_CHECKPOINT’ (Fig. 2A and B). In the plots of Fig. 2B, the location of the peaks suggested that the genes in these gene sets were upregulated in the high MRPL51 expression group compared to the low MRPL51 expression group.

Characterization of the association between MRPL51 expression and LUAD cellular behaviors at the single-cell level

To explore the association between MRPL51 expression and LUAD cellular behaviors at the single-cell level, MRPL51 expression was investigated in the single-cell RNA-seq dataset from Kim et al (22) based on PDX tumor cells (n=126) from patients with LUAD (Fig. 3A). A total of 14 single-cell functional state scores of every tumor cell in this dataset, including ‘stemness’, ‘invasion’, ‘metastasis’, ‘proliferation’, ‘EMT’, ‘angiogenesis’, ‘apoptosis’, ‘cell cycle’, ‘differentiation’, ‘DNA damage’, ‘DNA repair’, ‘hypoxia’, ‘inflammation’ and ‘quiescence’, were obtained from CancerSEA (23) (Fig. 3B). The correlation between MRPL51 expression and the functional state scores was assessed by calculating Pearson's correlation coefficient. Since the cell number was relatively small (n=126), and the bioinformatics analysis only provided directions for validation studies, weak correlation (Pearson's r≥0.2) was set as the cut-off. In total, 6 out of the 14 analyzed functional states showed weak but significant correlation with MRPL51 expression (Pearson's r≥0.2; P<0.05). In detail, MRPL51 expression was positively correlated with ‘cell cycle’, ‘DNA damage’, ‘DNA repair’, ‘EMT’, ‘invasion’ and ‘proliferation’ of LUAD cells at the single-cell level (Fig. 3D).

Knockdown of MRPL51 increases epithelial properties, induces G1 phase arrest and weakens the invasion of LUAD tumor cells

To validate the regulatory effect of MRPL51 on cellular behaviors, lentiviral shRNA was used for MRPL51 knockdown in A549 and Calu-3 cells (Fig. 3C and E). shMRPL51#1 and shMRPL51#2 were used for further studies due to a more optimal suppression of MRPL51 expression. MRPL51 knockdown was associated with decreased N-cadherin and Vimentin expression, but increased E-cadherin expression in A549 and Calu-3 cells (Fig. 3E and F). A Transwell invasion assay confirmed that A549 and Calu-3 cells with MRPL51 knockdown exhibited significantly decreased invasion (Fig. 3J-K). Flow cytometry and CCK-8 assays demonstrated that MRPL51 knockdown induced G1 phase arrest and slowed cell proliferation (Fig. 3G-I and L-M).

MRPL51 upregulation is associated with unfavorable OS of patients with LUAD

Kaplan-Meier survival analysis was performed using PFI and OS data from TCGA-LUAD. Patients were grouped by median MRPL51 expression. The log-rank test showed that patients with higher MRPL51 expression (top 50%) had a similar rate of PFI [hazard ratio (HR), 1.108; 95% CI, 0.840-1.462; P=0.47] compared with the patients with lower MRPL51 expression (bottom 50%; Fig. 4A). Patients with higher MRPL51 expression had a significantly worse OS rate than the patients with lower MRPL51 expression (HR, 1.462; 95% CI, 1.094-1.953; P=0.011; Fig. 4B). Using the same cut-off of MRPL51 expression, the prognostic significance was validated using data from LUAD cases in the Kaplan-Meier plotter. The log-rank test indicated that patients with higher MRPL51 expression had a significantly worse FP probability (HR, 1.55; 95% CI, 1.12-2.14; P=0.0077; Fig. 4C) and OS rate (HR, 2.58; 95% CI, 1.97-3.38; P<0.001; Fig. 4D).

MRPL51 transcription is activated by FOXM1 in LUAD tumor cells

The transcription factors with potential binding to the MRPL51 promoter were investigated using available ChIP-seq data in Cistrome Data Browser (CistroDB; http://cistrome.org/db/) (27). Only transcription factor genes with a regulatory potential >0.85 in epithelial tissues were identified as candidates, including PATZ1, POLR2A, MYC, FOXM1, STAT1, ERG and MBD2 (Fig. 5A). The expression correlation between the transcription factor genes and MRPL51 was assessed using RNA-seq data from TCGA-LUAD (Fig. 5A). Among the transcription factor genes, FOXM1 had the strongest correlation with MRPL51 expression (Pearson's r=0.52; Fig. 5B, data not shown for the other genes since their correlation with MRPL51 were all <0.2). FOXM1 expression was also significantly higher in LUAD tissues than in normal lung and LUAD tumor-adjacent normal tissues (Fig. 5D). In CistroDB, no ChIP-seq data were available for LUAD cells. Therefore, 293 cells exhibiting epithelial morphology were selected for bioinformatics prediction. Using the UCSC Xena platform and the anti-FOXM1 based ChIP-seq data for 293 cells (11), it was demonstrated that the binding sites of FOXM1 in the MRPL51 promoter (5′-ACAAATA-3′) match the FOXM1 binding motif (5′-RYAAAYA-3′; Fig. 5C and E) (11). To verify the regulatory effects of FOXM1 on MRPL51 transcription, lentiviral FOXM1 shRNAs were generated. FOXM1 knockdown decreased the mRNA and protein expression levels of MRPL51 (Fig. 5F and G). shFOXM1#1 and shFOXM1#3 were used for further studies due to a more optimal suppression of FOXM1 expression. To validate the activating effect of FOXM1 on the MRPL51 promoter, a recombinant pGL3-basic reporter plasmid carrying either wild-type (WT) [MRPL51-promoter(p)-WT] or mutant (MT) (MRPL51-p-MT, 5′-ACCCCTA-3′) MRPL51 promoter fragments was generated (Fig. 5E and I). Dual-luciferase assays confirmed that the knockdown of endogenous FOXM1 in A549 and Calu-3 cells significantly impaired the luciferase activity of MRPL51-p-WT but had limited influence on MRPL51-p-MT (Fig. 5I). Subsequently, ChIP-qPCR was conducted using anti-FOXM1 in A549 and Calu-3 cells with or without FOXM1 knockdown. Primers covering the predicted FOXM1 binding site were designed (blue arrows; Fig. 5E). qPCR data revealed that the MRPL51 promoter fragments with the predicted FOXM1 binding site were enriched in the samples precipitated using anti-FOXM1 compared with IgG controls (Fig. 5H). However, FOXM1 knockdown significantly impaired the enrichment (Fig. 5H).