Alternative Polyadenylation Utilization Results in Ribosome Assembly and mRNA Translation Deficiencies in a Model for Muscle Aging.

Aging-associated muscle wasting is regulated by multiple molecular processes, whereby aberrant mRNA processing regulation induces muscle wasting. The poly(A)-binding protein nuclear 1 (PABPN1) regulates polyadenylation site (PAS) utilization, in the absence of PABPN1 the alternative polyadenylation (APA) is utilized. Reduced PABPN1 levels induce muscle wasting where the expression of cellular processes regulating protein homeostasis, the ubiquitin-proteasome system, and translation, are robustly dysregulated. Translation is affected by mRNA levels, but PABPN1 impact on translation is not fully understood. Here we show that a persistent reduction in PABPN1 levels led to a significant loss of translation efficiency. RNA-sequencing of rRNA-depleted libraries from polysome traces revealed reduced mRNA abundance across ribosomal fractions, as well as reduced levels of small RNAs. We show that the abundance of translated mRNAs in the polysomes correlated with PAS switches at the 3'-UTR. Those mRNAs are enriched in cellular processes that are essential for proper muscle function. This study suggests that the effect of PABPN1 on translation efficiency impacts protein homeostasis in aging-associated muscle atrophy.

Aging is marked by reduced regulation of the fundamental cellular machinery including transcription and translation that together control protein production. The level of translation can be studied by monitoring the extent to which mRNAs are associated with polysomes (1). RNA sequencing (RNA-seq) analysis of mRNAs in the heavy polysome fractions is often employed to gain global insight into the cell translatome (2). Translational efficiency measures the ratio of translated mRNAs relative to the steady-state levels of the mRNAs (3). Changes in response to a wide range of stress conditions can be associated with reduced translation efficiency (3–5).

Generally, mRNA transcription rates positively correlate with translation efficiency (6). Translation efficiency is affected by the transcription start site at the 5′-untranslated region (UTR) of mRNA (4,7). The cytoplasmic mRNA pool is strongly influenced by the 3′-UTR lengths, which is determined by the polyadenylation site cleavage site (7,8). The canonical polyadenylation site (cPAS) motif is often found at a distal region of the 3′-UTR, resulting in a long 3′-UTR, whereas the APA motifs can be found at both the 3′-UTR or intragenic regions (introns or exons). APA utilization at the 3′-UTR often results in shorter 3′-UTR, affecting mRNA stability, nuclear export, and overall cytoplasmic mRNA levels (9). Together, this determines the pool of mRNA that could bind to the ribosome and be translated. The effect of cPAS and APA utilization on translation efficiency is not sufficiently understood (3,7).

The poly-(A) nuclear RNA binding protein 1 (PABPN1) is a key regulator of APA utilization. PABPN1-mediated APA utilization at the 3′-UTR affects ~70% of mRNA, resulting in higher mRNA stability (10,11). PABPN1 also affects APA utilization in intragenic regions (exons and introns) leading to alternative splicing and the production of transcript variants (12). PABPN1 is ubiquitously expressed in all eukaryotic cells (13), and its abundance declines with age, predominantly in human skeletal muscle (14). Hallmarks of muscle aging, including atrophy and impaired regeneration are caused by reduced PABPN1 levels (14–17), indicating PABPN1 regulatory role in muscle aging. A short alanine-expansion mutation in PABPN1 is the genetic cause for oculopharyngeal muscular dystrophy (OPMD), an adult-onset myopathy (18). The pathogenic PABPN1 protein forms nuclear aggregates (19). Levels of nuclear PABPN1 are depleted by PABPN1 aggregation limiting levels of functional nuclear PABPN1 (20,21). Below a critical level, PABPN1 expression could induce muscle wasting (22). The molecular mechanisms underlying PABPN1-mediated muscle wasting have been explored using the PABPN1-dependent proteome in mouse muscles (17). Although reduced PABPN1 levels have been proposed to negatively affect nuclear export and translational efficiency of target mRNAs (23,24), the effect of PABPN1 on translation efficiency has not been studied in-depth. Here we investigated PABPN1 role on translation efficiency and translation profiles in muscle cell cultures.

Materials and Methods

Cell Culture

Human 7304.1 immortalized myoblasts and human primary myoblasts from Rf1034.5 (50,51) were propagated in F12 medium (Invitrogen) supplemented with 20% fetal bovine serum, 5 ng/ml hEGF (R&D systems, Minneapolis, MN, USA), and 1 µm dexamethasone. Myoblasts were cultured at 90% confluence. The parental lines were not transduced with a lentivirus. Stable cell lines shPAB and shScram were transduced with shRNA to PABPN1 or scramble shRNA as described in (14,24). A summary of experimental procedure for stable or transient PABPN1 down-regulation is found in Supplementary Figure S1. PABPN1 downregulation was confirmed by Western blot analysis and immunofluorescence.

Translation Efficiency

Protein translation was determined in living cells using the Protein synthesis assay kit (Cayman chemicals, Ann Arbor, MI). The assay was performed in a 96-well plate following the manufacturer’s protocol, with the following modifications: the OPP Alexa-flour-488 Azide (named here OPP-488) staining solution was diluted twice compared to the company protocol. During the washing steps, Hoechst 33342 (Thermo Fisher Scientific, Waltham, MA) was added to stain the nucleus. Cells were pretreated with 100 µg/ml cycloheximide (CHX) were used as a negative control. Living cells were directly imaged with the High-Content Screening (HTC VTI platform, Thermo Fisher Scientific). Cell-based image quantification was carried out with the Compartmental Analysis BioApplication. The MFI values were corrected to the background (cultures without Azide-488 labeling). Image quantification was carried out per cell, using nuclear segmentation to mark an individual cell. After live-cell imaging, cells were fixed with 2% formaldehyde and immunolabelled with antibodies to PABPN1 as described in (14). Plates were reimaged with the HTC VTI platform and PABPN1 MFI was quantified within the nuclei in all objects and OPP-488 positive objects. Experiments were carried out in triplicate and were repeated at least 4 times.

Ribosome Fractionation and RNA Isolation

Fractionation was performed as described in (1). In brief, before cell collection, cultures were treated with 100 µg/ml CHX for 10 minutes. The cytoplasmic fraction of the collected cells were lysed in PEB (20 mM Tris-HCl pH 7.5, 100 mM KCl, 5 mM MgCl2, 0.3 % Igepal CA-630) supplemented with protease and RNase inhibitors + CHX. Cell lysis was carried out on ice for 10 minutes. Subsequently, nuclei and cell debris were removed by centrifugation at 10 000 rpm for 10 min. Fractionation was carried out through a 10%–50% sucrose gradient, centrifuged for 90 min at 39 000 rpm at 4°C. Twelve fractions, 1 ml each, were collected per gradient; real-time data acquisition, traces recording, and quantification was carried out with Peak Chart (Brandel, Inc.), polysome traces including controls are shown in Supplementary Figure S2. The polysome fractions were stored at −80 before RNA extraction. RNA and protein extractions were performed with Trizol Reagent, according to the manufacturer’s protocol (Thermo Fisher Scientific). RNA quality was assessed with the Lab-on-a-Chip technology, and RNA quantity was determined with Qubit RNA HS assay kit (Thermo Fisher Scientific). Fraction 12 had too little RNA and was omitted from further analysis. The protein pellet was dissolved in a protein loading buffer.

Ribosomal RNA Depletion, RNA Library Preparation, and Sequencing

Library preparation was carried out using 2 protocols for all ribosomal fractions. Ribosomal RNA depletion was done using the NEBNext rRNA Depletion Kit (NEB). After ribosomal RNA depletion, first-strand synthesis was performed with random hexamers and Superscript III (Thermo Fisher Scientific). Second-strand synthesis was performed in the presence of uracil to generate stranded libraries. Sequencing libraries were generated using NEBNext chemistry (NEB) and library quality was assessed with the Agilent 2100 Bioanalyzer. Poly(A) RNA was amplified using an oligo(dT) primer or the TeloPrime Full-Length cDNA Amplification Kit V2 (Lexogen) according to the manufacturer’s protocol. Libraries were quantified using a High Sensitivity assay on the Qubit (Invitrogen) and library quality was assessed on a Bioanalyzer 2100 (Agilent). Sequencing of the oligo(dT) libraries was not performed due to extremely low cDNA output in shPAB samples. Libraries were sequenced on the Illumina HiSeq4000 2 × 150bp PE to a depth of 30 million reads per sample. Before sequencing cDNA samples were pooled: F3―3–4; F4―5–7; F5―8–11, fractions 1 and 2 remained separately as depicted in Supplementary Figure S3.

Sequencing and RNA-seq Analyses

All RNA-seq files were processed using the open-source BIOWDL RNA-seq pipeline developed at the LUMC (https://github.com/biowdl/RNA-seq/). This pipeline performs read preprocessing (including quality control, quality trimming, and adapter clipping), alignment, expression quantification, and transcript assembly. FastQC was used for checking raw read QC. Adapter clipping was performed using Cutadapt (v2.4) with default settings. Alignment of RNA-seq reads was performed using STAR (v2.6.0c) with the GRCh38 reference genome without alternative contigs. Ensembl gene annotation v.94 was used for both transcript expression quantification.

RNA-seq analysis was carried out on transcripts using TPM normalization, correcting read counts for transcript length. Transcript quantification was performed using StringTie (1.3.4d) with the “-e” mode to only quantify the expression of known transcripts. Expression values were normalized using the TPM procedure (52). The PCA of CPM, the variance, the hierarchical clustering and heatmap, and TPM MDS were calculated and plotted in R.

The assignment of biotype to a category is found in Supplementary Table S1. In total, 206 049 transcripts were included in further downstream analysis. Subsequently, 2 filtering steps were added: (a) transcripts that have a “zero” TPM in at least 9 out of the 10 samples were removed, reducing the data set to 92 925 transcripts; and (b) transcripts with an average TPM <1 (fractions 3–5) were removed and the data set was reduced to 36 883 transcripts. Both filtering criteria minimized the impact of transcripts whose abundance between fractions was inconsistent. All subsequent analyses were carried out on this data set.

Fold changes between shPAB and parental samples were calculated per fraction using a log 2 (TPM + 1) transformation. The average fold change was calculated on a weighted fold change using the following weighting factor: 1, 1.1, 1.8, 1.9, and 2 for fractions 1–5, respectively. Transcripts with an average fold-change >2, absolute value, were then selected (<−0.98568 or >0.99645 in log 2) as differential abundance. Transcripts after fold-change selection are found in Supplementary Table S2.

3′-UTR Analysis

First, adapter sequences were removed from the raw RNA-seq reads using cutadapt v.1.16, and the read quality was checked using FastQC (v0.11.7) and MultiQC (v1.7). The processed RNA-seq reads were aligned to the human reference transcriptome based on GRCh38 and Ensembl gene annotation v.103 using STAR (version 2.7.5a) with the following parameters: ―outReadsUnmapped Fastx―twopassMode Basic―quantMode TranscriptomeSAM―chimJunctionOverhangMin 15―chimSegmentMin 15. We then collected 81 764 transcripts, for which the 3′-UTR is annotated in Ensembl v.103. Using a customized Python script, we created a GTF file annotating 1 canonical polyadenylation signal sequence (cPAS; AATAAA) and 2 alternative polyadenylation signal sequences (APA; CCYTCY and CWGGYC) in the 3′ UTR of these 81 764 transcripts. These 3 polyadenylation signal sequences are selected because they are the most significant signals in muscles based on our earlier study (11). Next, we used featurecounts from the subread package (v2.0.1) to quantify reads in all samples with the identified 3′ UTR polyadenylation signal sites. Multiple alignments are included in the quantification as a fraction using the featurecounts parameters of “-M –fraction”. The normalized reads aligned to cPAS or APA were clustered together and the average counts were calculated per PAS motif, per transcript, and sample. A zero was assigned if one of these polyadenylation signals was not found at the 3′-UTR, and if reads were not mapped to these sequences. Subsequently, TPM was applied to normalize this raw average expression table of polyadenylation signals where library size is calculated as all aligned reads in the STAR alignment files. To select transcripts where a clear signal can be observed, only transcripts having nonzero average expression at both distal and proximal polyadenylation signal sites in at least 3 parental samples are kept. This resulted in 15 298 transcripts out of the original 81 764 transcripts. We further divided them into 2 groups based on their biotype: the translated CDS consisted of 12 104 transcripts and the untranslated CDS consisted of 3 194 transcripts. These filtered transcripts were used for analysis. A ratio was also calculated between the average expression of distal and proximal polyadenylation signal sites per sample for these filtered transcripts. The 3′-UTR ratios are found in Supplementary Table S4.

RT-qPCR Analysis

RT-qPCR analysis was performed on RNA from the pooled ribosomal fractions in triplicate. Fractions 1 and 2 were pooled together. RNA (500 ng) was converted to cDNA with the RevertAid First Strand cDNA Synthesis Kit (Thermo Fisher), according to the manufacturer’s protocol using the oligo dT primer. Three nanogram RNA were used for qPCR amplification using technical duplicates using a standard amplification protocol (TM 60 degrees). The average of technical duplicate CT values were normalized to the values in parental fractions 1 + 2. Per gene, the primer set covered the last 2 exons, ensuring that only mRNA was amplified. Primer sets are found in Supplementary Table S5. Samples with CT values >35 were excluded.

RNAscope

Spatial localization of a single RNA molecule was carried out in adherent muscle cells using the RNAscope Fluorescent Multiplex Assays kit according to the manufacturer’s protocol (ACD biotech) with the following modifications: the protease solution was diluted 1:30, and the amplifier solutions were diluted 1:2. Ten genes were included in the RNAscope: The positive control probe mix, provided by the company (POLR2A, UBC, and PPIB). The second probe mix included LMNA, CAV1, and DES; the third probe mix was included PABPN1 and DMD, and the fourth probe mix included TP53 and MYF5. All probes (Supplementary Table S5) are from the Human probe set (ACD biotech); The PABPN1 probe is from the mouse set, the homology between human and mouse is >95% for the probe regions. Experiments were carried out in a 96-well plate, and the labeled plate was imaged and quantitatively analyzed with the Cellomics HCS, 7-filter with LED light source (Thermo Fisher Scientific). Cell-based image quantification was carried out with the SPOT toolbox. Every probe was hybridized in 3 biological replicates. Results are from a representative experiment, mean and variation are from >1 000 cells.

Protein Extraction and Western Blot Analysis

Lysis and separation by SDS–PAGE were carried out as described in (17). Protein aliquots from the ribosome fractionation were separated by SDS–PAGE on 10% acrylamide gels. Proteins were transferred into PVDF membranes and incubated with primary antibodies, followed by antirabbit or antimouse secondary IRDye 800CW or IRDye 680RD conjugated (LI-COR) antibodies and were subsequently detected using the Odyssey CLx Infrared imaging system (LI-COR, NE). Antibodies are listed in Supplementary Table S5.

Ribosome Structural and Enrichment Network Analyses

The structure of the human 80S ribosome (protein database code [PDB] 4UG0) was annotated using UCSF Chimera (53). The mRNA was placed in the model in the same position and orientation as found in PDB code 6HCF. The ribosomal proteins in the shPAB proteome (17) were included in the model.

GSEA was carried out as detailed in https://www.gsea-msigdb.org/gsea/doc/GSEAUserGuideFrame.html. We selected for GSEA gene networks that were found as significant in DAVID https://david.ncifcrf.gov/ and showed a strong protein–protein network connectivity in STRING (https://string-db.org), V.9.1 (32). Per gene set, the protein–protein interaction of the genes showing significance enrichment with shPAB was demonstrated in STRNG using the high confidence (0.7) setting.

Statistical Analyses

Statistical analyses were performed in R and GraphPad Prism 8.4.2. Statistical analyses for TPM differences between fractions were assessed with Brown-Forsythe and Welch Anova, assuming no equal standard deviations (SDs), corrected for multiple testing for all biotypes n > 50 transcripts (CDS-T; CDS-UT; pseudogene transcripts and long ncRNAs). For small ncRNAs and misc-RNA, a parametric t test was used to assess statistical differences between parental and shPAB. Statistical assessment of differences in the number of transcripts per biotype between parental and shPAB was carried out using the chi-square test. Differences in dCT and the number of RNA foci (RNAscope) were estimated using the student’s t test.

Results

Changes in Translation Efficiency Upon Reduction in PABPN1 Levels

Global translation efficiency (25) was measured with the fluorescently tagged O-propargyl-puromycin (OPP) (26), and the effect of PABPN1 was investigated in human myoblast cultures with a constitutive or transient PABPN1 downregulation (stable cell-line expressing PABPN1 shRNA [shPAB]). PABPN1 levels in constitutive shPAB cells were reduced in both nuclear and cytosolic fractions by close to 60% or 50%, respectively (Figure 1A). In transient down-regulated cell cultures nuclear PABPN1 levels were reduced in 50% of controls (Supplementary Figure S1). In cells with reduced PABPN1 levels we found a significant reduction in translation efficiency (Figure 1B and C and Supplementary Figure S1). High differences in translation activity were found between cells, in both parental and shPAB (Figure 1B). Therefore, our analysis was cell-based, fluorescence intensity was measured only from the green-segmented cells (Figure 1C and Supplementary Figure S1). In parental cells, 40% of cells were gated as OPP-positive, but in the shPAB cell culture, only 15% were gated as OPP-positive. We then asked if PABPN1 levels correlate with the percentage of OPP-positive cells. We measured PABPN1 levels from the same cultures stained with an antibody to PABPN1 (Figure 1D). Confirming the Western blot analysis, PABPN1 mean fluorescence intensity (MFI) was nearly 50% lower in shPAB compared with parental cells (Figure 1E). Yet, in the OPP-488 positive cells, PABPN1 levels were nearly doubled (Figure 1F). This indicates that high PABPN1 levels are associated with translation activity.

Figure 1.

Analysis of translation efficiency in stable muscle cells with reduced PABPN1 levels. (A) Western blot analysis of PABPN1 levels in nuclear and cytoplasmic fractions from parental to shPAB cells in 2 biological replicates. glyceraldehyde-3-phosphate dehydrogenase and Histone 3 are loading controls for the cytoplasmic or nuclear fractions, respectively. (B and C) Nascent protein analysis. Representative images (B) after incubation with OPP-488 (green) and nuclear counterstaining (blue). Cell cultures pretreated with cycloheximide (CHX) are negative controls. Bar charts (C) show the average of OPP-488 mean fluorescence intensity (MFI) per cell. Average and standard deviations are from N = 3. The average number of cells that were included in that analysis is depicted above/in the bars. (D–F) Images and analysis of OPP-488 and PABPN1 in fixed cells. (D–F) images and analysis in fixed and PABPN1 immunofluorescence. Representative images (D) of OPP-488 (green), PABPN1 (red), and nuclear counterstain (blue). The lower panel shows examples of segmented OPP-488 positive cells (green counter). Bar charts show the average nuclear PABPN1 MFI per cell (E), or PABPN1 MFI ratio between shPAB and parental cultures in all cells or the OPP-488 positive cells (F). Statistical analyses are from 3 biological replicates, the average number of cells is denoted within the bars in panel C. p < .05 is denoted with *, not significant (NS). Scale bar is 30 µm. OPP-488 = OPP Alexa-flour-488 Azide; PABPN1 = poly(A)-binding protein nuclear 1.

PABPN1 Differentially Affects the Presence of RNA Biotypes on Polysome Profiles

To systematically determine how PABPN1 downregulation reduces translation efficiency, we analyzed RNA from polysomes prepared from parental and shPAB cell cultures. Polysomes were fractionated on sucrose gradients and the RNAs present in different fractions were investigated. Polysome fractions were collected from 5 independent experiments. In general, fraction 1 (F1) contained RNA not bound to ribosomes, fractions 2–4 contained the 40S and the 60S and the monosomes (80S), and fractions 5–11 contained polysomes of increasing molecular weights (Figure 2A). The polysome traces were abolished in the presence of ethylenediamine tetraacetic acid (EDTA), whereas cycloheximide or puromycin abolished polysomes (Supplementary Figure S2). Overall, the polysome profiles of parental and shPAB cells indicate that the translation machinery is generally functional and intact compared to other conditions like EDTA or puromycin. The fractions #1, #2–4, and #5–11 were distinguished based on RNA size distribution (Supplementary Figure S3A) and differed from cytoplasmic or nuclear fractions (Supplementary Figure S3B). Based on bulk RNA profiles, fractions 1 and 2 (F1 and F2) were kept separate, and we pooled fractions 3 and 4 (F3), fractions 5–7 (F4) containing the light polysomes, and fractions 8–11 (F5; Figure 2A). The bulk RNA, as detected at 260 nm, in F4 and F5 was higher in shPAB compared with parental, warranting further investigation (Figure 2B). The polysome traces profiles were comparable within an experiment but diverge between technical replicates (Supplementary Figure S2). Therefore, we did not combine the polysome traces profiles from different experiments, and sets were studied separately as is customary. We used 1 polysome traces experiment for RNA-seq and another for validation.

Figure 2.

RNA biotype distribution in polysome profiling fractions. (A) Absorbance at 254 nm in polysome profiling fractions from parental and shPAB cell cultures as collected from sucrose gradient in 12 fractions. Fraction 12 contained no RNA and was omitted. Pooled fractions are denoted under the sucrose gradient fractions. Fractions 1 and 2 were kept separate, fractions 3 and 4 contained the monosomes, fractions 5–7 contained the light polysomes, and fractions 7–11 the heavy polysomes. See also Supplementary Figures S2 and S3. (B) Bar chart shows the percentage of the area under the curve (AUC) across pooled polysome fractions in parental and shPAB cell culture. Average and standard variations are from N = 4 experiments. A statistical significance between parental and shPAB is denoted with *. (C) A pipeline to filter transcripts with low signal levels. (D) Hierarchical clustering and heatmap of correlation RNA data set from parental (P) and shPAB (sh) polysome fractions 1–5. See also Supplementary Figure S3. (E) Multidimensional scaling (MDS) plot shows a pairwise distance similarity between the polysome fractions for MDS1 and MDS2. (F) Box plots show the median of the total TPM of transcripts in the 7 biotype groups: translated coding sequences (CDS), untranslated (UT) CDS, long noncoding RNAs (nc), small nc, pseudogene transcripts, rRNA, and misc-RNAs. Fraction number and parental or shPAB are annotated with colors or shapes (circle or triangle), respectively. The percentage of transcripts is depicted between brackets. (G) The distribution of RNA biotype across ribosomal fractions is demonstrated by TPM proportions. TPM = transcript per million.

Considering the limitation of 1 experimental RNA-seq data set, we applied 3 steps to filter our transcripts that are not consistent across fractions: (a) transcripts with <50 counts after normalization were excluded; (b) transcripts with inconsistent counts across the samples were removed, only transcripts with counts in 9 or 10 fractions were included; and (c) transcripts with low counts in the polysome fractions were excluded as well. The principal component analysis (PCA) showed that variations between the ribosomal fractions were larger than the variations between parental and shPAB (Supplementary Figure S3D). The hierarchical clustering of correlations showed high similarity between fractions F2 and F5, whereas F1 was distant (Figure 2D). The transcript per million (TPM) multidimensional analysis (MDS) further indicated that F1 and F2 were more distant from F3, F4, and F5. Altogether, most differences were found between fractions rather than between parental and shPAB (Figure 2E).

Since we used rRNA-depleted libraries for sequencing, we investigated the distribution of different RNA biotypes across the ribosomal fraction, and then asked which biotypes’ distribution are affected by PABPN1 (Figure 2F; Supplementary Table S1). The lowest abundance was found for rRNA, supporting the successful depletion of rRNA (Figure 2F, Supplementary Table S1). RNA biotype groups are characterized by transcript length (27), we confirmed that the long RNAs were enriched in the heavy polysome fractions, whereas the short RNAs were enriched in F1 and F2 (Supplementary Figure S4). Therefore, the enrichment analysis was corrected for transcript length. Protein coding sequences (CDS) transcripts were the most abundant and were mostly enriched in the light and heavy polysome fractions (Figure 2G). This observation is consistent with previous studies (1,28,29). We discriminated between the translated (T)-CDS and the untranslated (UT)-CDS. UT-CDS lacks elements that are required for translation. Because PABPN1 affects APA at internal gene regions, which can result in untranslated transcripts, we assessed this group separately from T-CDS. The proportions of UT-CDS were 3-folds lower than T-CDS, but the relative abundance across fractions was similar (Figure 2G). High read counts were also found for small miscellaneous RNAs (misc; Figure 2G), represented by only a small number of transcripts (Supplementary Table S1). The abundance misc-RNAs across the fractions was high in the unbound and monosome fractions but low levels in the polysome fractions (Figure 2G). The function of misc-RNAs in the ribosome is predominantly unknown (30). Other biotypes, pseudogenes (transcripts that have homology to CDS but have a disrupted open reading frame), small and long noncoding RNAs, had relatively small read count across the ribosomal fractions (Figure 2G and Supplementary Table S1). In contrast to the bulk RNA measurements suggesting accumulation of RNA in the heavy fractions in the shPAB polysomes (Figure 2B), the bulk TPM reads did not show similar differences between parental and shPAB (Figure 2G). This difference could result from rRNA depletion and tRNA exclusion from library preparation and biotype analysis, of the RNAseq-based approach, respectively.

PABPN1 Primarily Affects the Abundance of Translated CDS and misc-RNAs

The effect of PABPN1 on RNA abundance in the ribosomal fraction was investigated with a difference of the weighted fold change across the fractions. A weighted value was calculated to correct for the differences in abundance between fractions. We applied an average weighted fold-change ≥2 to identify RNAs whose abundance in the ribosome is affected by PABPN1. From those, 53.6% (4 269) transcripts had low abundance (fold change ≤0.5) and 46.4% (3 699) had high abundance (fold change ≥2) in shPAB (Figure 3A). Furthermore, we focused on differences that were significant with a chi-square test (Figure 3B). Overall, T-CDS and misc-RNAs were the most significantly affected by PABPN1 levels across polysomal fractions (Figure 3C and Supplementary Figure S5). Both biotypes had lower abundance in shPAB samples in both monosomes and polysomes fractions (Figure 3C). The abundance of UT-CDS in the ribosomal fractions was unaffected by PABPN1 levels (Figure 3C). A higher abundance for both T-CDS and UT-CDS was found in F1, however, the number of transcripts is relatively low compared with the ribosomal fractions (Figure 3C).

Figure 3.

Altered RNA abundance in biotype groups across ribosomal fractions is affected by PABPN1. (A) Volcano plot shows average TPM (log 2) versus average weighted fold-change (FC) in shPAB versus parental combined ribosomal fractions. The PABPN1-affected transcripts are considered with an FC ≥ |2| cutoff and are denoted in red (FC ≥2) or blue (FC ≤2). The unchanged transcripts are marked in gray. (B) RNA biotypes of the PABPN1-affected transcripts. For each biotype, the p value was calculated using the chi-square test for significant enrichment in shPAB cells. The biotype groups with reduced or enriched transcripts in shPAB are denoted in blue or red, respectively. See also Supplementary Figure S5. (C) Dot plots show TPM (log 2) in parental versus shPAB in each ribosomal fraction and biotype group. A statistical difference between parental and shPAB is denoted with *p < .05; **p < .005; ***p < .001 the mean difference is between brackets. The plus or minus values of the mean difference indicate lower or higher levels in shPAB, respectively. NS denotes not significant. (D and E) Examples of TPM in the pooled polysome fractions in parental (black) and shPAB (red) for misc-RNAs (D) Y RNA (sum of 14 transcripts), 7SK small nuclear RNAs (sum of 4 transcripts), or T-CDS (CDS-T; E) PABPN1, RPL3, Tubulin (TUBB), and NCL. For the translated CDS protein levels were assayed with a Western blot. See also Supplementary Figure S6. (F) Plots show RT-qPCR (left y-axis, dCT normalized to F1 in parental) and RNA-seq (right y-axis, sum TPM for all transcripts, the transcript number in between brackets) for 5 genes across ribosomal fractions. Fractions 1 + 2 were averaged due to low expression levels as determined by RT-qPCR analysis. (G) Bar chart shows the mean single-molecule RNA foci per cell in parental (black) and shPAB (red) cell cultures. Genes with a lower or higher fold change in shPAB RNA-seq are highlighted in blue or red, respectively. Unchanged genes are denoted in gray. RNA foci were quantified from RNA scope experiments, averages, and standard error are from N > 1 000 cells from 3 independent experiments. Statistical significance (p < .05, Student’s t test) is denoted with an asterisk. Negative controls (nc, probes to bacterial DapB) are shown per fluorophore. Representative images are shown in Supplementary Figure S7. (H) Bar charts show the mean weighted fold change in RNA-seq (gray) and RNAscope (purple). Between brackets is the number of transcripts that were annotated with TPM. CDS = coding sequences; PABPN1 = poly(A)-binding protein nuclear 1; TPM = transcript per million.

Examples of TPM quantification across ribosomal fractions are shown in Figure 3D and the reads mapped to genomic regions in 4 T-CDS are shown in Supplementary Figure S6. Protein levels across the ribosomal fractions were found to be consistent with TPM levels (Figure 3E). Reduced PABPN1 levels in shPAB F3–F5 are consistent with the shRNA expression to PABPN1 (Figure 3E). Higher ribosomal protein L3 (RPL3) levels in shPAB polysome profiling (Figure 3E) are consistent with higher levels in the shPAB proteome. In addition, TUBB mRNA and nucleolin (NCL) mRNA levels were unchanged between parental and shPAB polysomes, and protein levels were also unchanged (Figure 3E).

Validation of mRNA distribution using reverse transcriptase quantitative PCR (RT-qPCR) was carried out in an independent polysome profiling set. Expression values were normalized to F1 + F2 in parental. F1 + F2 were combined due to high CT values. Overall, the differences between parental and shPAB that were observed with the RNA-seq were recapitulated by the RT-qPCR results (Figure 3F). Nevertheless, changes in expression levels across the fractions were not always consistent between RNA-seq and RT-qPCR (Figure 3F). These differences are likely due to differences between the polysome traces.

To gain additional confidence on the effect of PABPN1 on mRNAs, that potentially bind to the ribosome and affect the PABPN1 translatome, we employed the RNAscope procedure (31) and analyzed only cytosolic RNA puncta (representative images are in Supplementary Figure S7). PABPN1 effect on the mean number of puncta per cell in cytosolic regions was found to be consistent with the RNA-seq fold-change direction (Figure 3G): we confirmed PABPN1, as well as reduced levels of DES, MYF5, and UBC mRNAs or CAV1 and TP53 elevated levels the levels of POLR2A, LMNA, and PPIB mRNAs, were consistently unaffected by PABPN1 (Figure 3G). Overall, the mean weighted fold change was consistent with the fold-changes measured by the RNAscope for 9 out of the 10 tested genes (Figure 3G). Only for DMD, the difference in TPM counts between parental and shPAB was not confirmed (Figure 3H). Cell imaging showed nuclear accumulation of DMD mRNA (Supplementary Figure S7), and the number of DMD puncta in the cytoplasm was too low to assess differences between parental and shPAB (Figure 3G). The nuclear accumulation of DMD mRNA is consistent with previous studies (32,33). To verify differential gene expression by reduced PABPN1 levels, we performed RNAscope in cell cultures with a transient PABPN1 down-regulation (Supplementary Figure S8). Overall similar changes in mRNA foci were found in shPAB constitutive and transient cells.

To assess the functional impact of the PABPN1-translatome we performed an enrichment of translated transcripts in the gene network. We first performed a global pathway enrichment of T-CDS whose abundance in F3–F5 is FC > |2| using GO terms, and after removal of redundancies we investigated enrichment of the most significant GO terms using Gene Set Enrichment Analysis (GSEA). The DNA repair and nuclear chromosome gene sets showed a significant positive enrichment score in shPAB (Figure 4A). Within the nuclear chromosome gene set, the histone deacetylates and the DNA damage formed a distinguished protein–protein network (Figure 4A). Both gene sets were not significantly affected in muscles with reduced PABPN1 levels or in OPMD. In culture, PABPN1 affects cell division (14,34), which is consistent with the dysregulation of cell cycle gene networks. The significant negatively enriched gene sets in shPAB included muscle contraction and ribosome (Figure 4B). Within muscle contraction, the actin–cytoskeleton gene network was prominent, and PABPN1 was also noted (Figure 4B-ii). The cytoskeleton and actin-binding proteins were found to be affected by PABPN1 levels (7,14,17). The protein–protein network of the ribosome showed a remarkable connectivity (Figure 4B-iv). A prominent enrichment of the ribosome protein network was also found in the PABPN1-dependent proteome in mouse muscles (17). Modeling of the commonly PABPN1-regulated ribosomal small subunit proteins in both mouse and cell model showed that ribosomal proteins that bind to mRNA, specifically S5 that binds to the 3′-UTR, are affected by PABPN1 levels (Figure 4B-v). We suspected that the abundance of transcripts from APA at the 3′-UTR in polysomes might be affected by PABPN1.

Figure 4.

Enrichment of protein networks in the PABPN1 translatome. (A) The positively enriched gene networks. (B) The negatively enriched gene networks. (i) and (iii) show a graphic presentation of GSEA enrichment results, and the main enrichment features: p value (corrected), enrichment score (ES) and the number of genes (N), and the percentage (%) that are specifically enriched in each gene network are depicted. (ii) and (iv) show the protein–protein network of the specifically enriched genes. The most abundant daughter networks are colored. (V) A model of the ribosome small subunit, the PABPN1-affected proteins in both shPAB cells and shPAB mice are highlighted. The mRNA is in green. GSEA = Gene Set Enrichment Analysis; PABPN1 = poly(A)-binding protein nuclear 1.

A PABPN1-mediated Polyadenylation Signal Utilization at the 3′-UTR Impacts the Abundance of Translated CDS Across Polysome Fractions

We then investigated if a change in polyadenylation site (PAS) utilization at the 3′-UTR influenced mRNA enrichment in polysome fractions after silencing PABPN1. We focused on canonical PAS (cPAS; AATAAA) and alternative PAS (APA; CCYTCY/CWGGYC) sequences at the 3′-UTR of transcripts and the ratio between the 2 parameters, which were found to be the most significantly enriched in muscles (11). The mean count per million (CPM) in cPAS or APA motifs at the 3′-UTR (examples are shown in Figure 5A and B) was used for analyses across polysome fractions. Examples of the mean CPM in cPAS and APA at the 3′-UTR are shown for unchanged or down-regulated transcripts in Figure 5A and B, respectively. The mean normalized counts in T-CDS with cPAS/APA motifs at the 3′-UTR significantly differed between parental and shPAB in the monosome fractions (F2–F3, Figure 5C and D). In contrast, a difference was not found for UT-CDS with cPAS/APA motifs at the 3′-UTR (Figure 5C and D). As expected, the number of transcripts with reads at cPAS or APA at the 3′-UTR was highest in the heavy polysomes, and nearly 3-folds higher in T-CDS compared with UT-CDS (Supplementary Figure S9 and Supplementary Table S1). The lower mean read counts for PAS or APA in UT-CDS compared with T-CDS in the polysome fractions suggested less translation for UT-CDS and further suggested that PABPN1 has a more pronounced effect in translated transcripts than on untranslated transcripts.

Figure 5.

Analysis of the 3′-UTR in shPAB coding sequences. (A and B) Screenshots of normalized read counts that mapped to cPAS (AATAAA) or APA at the 3′-UTR of ENST00000377482 (ERRFI1; A) and ENST00000228644 (MYF5; B) in parental or shPAB across polysome fractions. ERRFI1 transcript levels are unaffected by PABPN1 (fold-change 0.92), MYF5 transcript levels are 4.2 folds lower in shPAB. The location of cPAS or APA (CCYTCY/CWGGYC) is depicted in pink or blue lines, respectively. The 3′-UTR is depicted with a black bar and the length is in bp, orientation 5′–3′. The bar chart shows the mean counts per gene for cPAS (pink) or APA (blue). (C and D) Dot plot shows log mean CPM to cPAS (C) or APA (D) across the polysome traces fractions in parental (black) or shPAB (red). (E) Dot plot shows the mean ratio cPAS to APA across the polysome traces fractions. (C–E) T-CDS are denoted with a continuous line and UT-CDS with a dashed line. Error bars show SEM. Parental and shPAB samples are denoted in black or red, respectively. Statistical significance was determined with a one-way ANOVA test (Sidak’s multiple comparison test), ***q < .005; ****q < .001. (F) Dot plot shows Pearson correlation (R) between the mean ratio (cPAS to APA) and the mean weighted average fold change of T-CDS across the polysome traces fractions. ****p < .001.

A switch between canonical and alternative PAS utilization is a robust readout for the function of PABPN1 (11). Calculating the mean cPAS to APA per ribosomal fractions revealed a significant difference between parental and shPAB for the T-CDS but not UT-CDS (Figure 5E). Notably, the ratio differed only in F2–F5, polysome fractions, but not in the unbound fraction (Figure 5E). Together, these results suggest that the cPAS-to-APA ratios may affect T-CDS abundance in the polysomes and may affect translation efficiency. To test this, we correlated cPAS-to-APA ratios with the average weighted fold changes of T-CDS in F2–F5 and found a significant correlation in F3–F5 from parental F3–F5 (Figure 5F). In contrast, in shPAB samples, a correlation was not found (Figure 5F). This result suggests that a low mean cPAS-to-APA ratio could positively affect T-CDS abundance in the polysome fractions and translation efficiency. In contrast, a high mean cPAS-to-APA ratio, driven by reduced PABPN1 levels, seems to negatively affect translation efficiency.

Discussion

The disbalance of protein homeostasis plays a central role in a wide range of adult-onset pathologies. Whereas the involvement of protein breakdown pathways has been extensively studied, resulting in potential therapeutical approaches, the involvement of protein synthesis and translation efficiency have been far less studied. Translational control changes have been reported understress and cell senescence (35,36). Here we investigated translation efficiency in muscle cells with reduced PABPN1, as a model for adult-onset OPMD and sarcopenia. The ribosome content was found to be aging-associated in rat skeletal muscles (37) and aging human muscles (38). Also in muscles with reduced PABPN1 levels, the translation network was significantly dysregulated (17). Here we show that reduced translation efficiency in a human muscle cell with reduced PABPN1 levels translation efficiency is highly impaired (Figure 1). In both mouse and cell models, which exhibit features of muscle aging (14,15,17,39), we found the ribosome protein network to be highly affected (Figure 4B and (17)). Although expression of orthologous ribosomal proteins allow for a high degree of plasticity to ensure translation (40), insufficiency or overexpression of only a few ribosomal proteins can lead to reduced translation (41,42). Aberrant expression of ribosomal proteins in skeletal muscles affects muscle mass (43,44). For example, an increase in RPL3 and its paralog RPL3l led to a global reduction in protein synthesis (42). Higher RPL3 levels were also found in both the PABPN1-dependent proteome and the PABPN1-dependent translatome. An association between ribosome biogenesis and muscle hypotrophy or atrophy, via the mTOR and Wnt/β-catenin/c-myc signaling pathways, has been reported (45). Yet, how translation efficiency contributes to muscle atrophy and muscular pathologies requires further investigations.

We present evidence that PABPN1-mediated APA utilization at the 3′-UTR of coding transcripts is linked to reduced translated CDS in polysomes and a global reduced translation (Figure 5). In human skeletal muscles, PABPN1 levels are reduced from midlife onwards (14). Reduced PABPN1 levels below a critical expression level cause muscle atrophy and wasting, mimicking muscle aging (15,17). Since the PABPN1 effect on translation efficiency is highly prominent, it suggests a causal role for translation in muscle aging. Indeed, a global translation deceleration was reported in a wide variety of aging tissues and organisms (46). In short-lived organisms (eg, Caenorhabditis elegans, Drosophila melanogaster), alterations in ribosome stoichiometry and reduced translation efficiency are associated with longevity (46). Yet, in-depth studies are warranted to fully understand the molecular and structural changes of the ribosome and translation efficiency, and the role of cytoplasmic PABPN1 in aging-associated pathologies.

In addition to the prominent role of mRNA in translation efficiency, we found that the levels of misc-RNAs were significantly lower in shPAB polysomes (Figure 3B–D). Misc-RNAs are ~120 nucleotides-long ncRNAs whose functions are mostly unknown (30). The most enriched misc-RNAs were Y RNAs, a small RNA family bearing a step-loop structure with conserved binding sites for proteins such as the Ro60 ribonuclear protein (47). The Y RNA family has been implicated in multiple cellular functions, including RNA stability, but their precise cellular functions remain unknown (47). Interestingly, RNY5 binds to the 5S (48), which is among the ribosomal proteins that are significantly upregulated in shPAB cells. The effects of small, unknown RNAs, including Y RNAs on translation efficiency, are waiting for future studies.

The PABPN1-dependent translatome that we report here is limited to a single polysome profile experiment. Due to this limitation, the translatome includes only transcripts whose abundance is high and the average fold-change > |2|. The major enriched gene networks were found to be highly similar to the protein networks that were found in muscles with PABPN1 down-regulation (15,17,39) or in muscles expressing the expanded PABPN1 (14,49). The reduced abundance of muscle proteins, due to PABPN1 silencing (Figure 4B), agrees with the myogenic defects (34). In both muscle cells, culture and muscles reduced PABPN1 levels led to dysregulation of the ubiquitin-proteasome system and muscle atrophy (15). Moreover, upregulation of the nucleosome protein network in cell culture was also consistent with histone upregulation in the PABPN1-dependent proteome (17). The DNA damage protein network in the PABPN1-translatome (Figure 4A) agrees with an aberrant cell proliferation caused by reduced PABPN1 levels (14,34). Together, the PABPN1-translatome agrees with the PABPN1 effect on cell proliferation and muscle cell biology.

In summary, we show that translation efficiency is highly affected by PABPN1 levels. PABPN1 is predominantly nuclear-localized, affecting mRNA processing at the 3′-UTR. Translation efficiency is known to be affected by levels of cytosolic mRNA, and sequence elements at the 5′-UTR. Our study suggests the APA utilization also affects translation efficiency. Moreover, we suggest a role for small unknown RNAs in ribosome biogenesis and or translation.

Supplementary Material

Supplementary data are available at The Journals of Gerontology, Series A: Biological Sciences and Medical Sciences online.

Multimedia

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.