Work

Author response: Proteome-wide systems genetics identifies UFMylation as a regulator of skeletal muscle function

Year: 2022

Type: peer-review

Abstract: Article Figures and data Abstract Editor's evaluation Introduction Results Discussion Materials and methods Appendix 1 Data availability References Decision letter Author response Article and author information Metrics Abstract Improving muscle function has great potential to improve the quality of life. To identify novel regulators of skeletal muscle metabolism and function, we performed a proteomic analysis of gastrocnemius muscle from 73 genetically distinct inbred mouse strains, and integrated the data with previously acquired genomics and >300 molecular/phenotypic traits via quantitative trait loci mapping and correlation network analysis. These data identified thousands of associations between protein abundance and phenotypes and can be accessed online (https://muscle.coffeeprot.com/) to identify regulators of muscle function. We used this resource to prioritize targets for a functional genomic screen in human bioengineered skeletal muscle. This identified several negative regulators of muscle function including UFC1, an E2 ligase for protein UFMylation. We show UFMylation is up-regulated in a mouse model of amyotrophic lateral sclerosis, a disease that involves muscle atrophy. Furthermore, in vivo knockdown of UFMylation increased contraction force, implicating its role as a negative regulator of skeletal muscle function. Editor's evaluation This manuscript will be of broad interest to those working in the genetics of complex diseases, with the results strongly supporting the author's primary claims. Overall, this is an important study that demonstrates the power of proteomics-based systems genetics studies in the mouse. https://doi.org/10.7554/eLife.82951.sa0 Decision letter Reviews on Sciety eLife's review process Introduction Building and maintaining healthy skeletal muscles is crucial for all stages of life. Skeletal muscle makes up 30–40% of an adult human’s body mass and is vital not just for breathing and movement, but also for metabolism and longevity. The maintenance of muscle function is one of the best predictors for overall health, with sarcopenia being the major contributor of age-associated frailty (McGregor et al., 2014). Identifying factors that regulate skeletal muscle function has great potential to improve the quality of life for humans and animals. Systems genetics is a population-based approach that links genetic variation to complex traits (Baliga et al., 2017). These forward genetics approaches integrate genomics and other multi-omic data to phenotypic traits using various modelling approaches including genetic mapping and correlation analysis. Advances in high-throughput proteomic techniques have enabled the analysis of genetically diverse populations, and integration of these data with systems genetics has emerged as a powerful approach to identify novel associations between protein abundance and complex phenotypes (for a review, see Molendijk and Parker, 2021a). Recent large-scale plasma proteomic studies in human populations have begun to unravel complex genetic variations and their contribution to proteome diversity and disease-relevant phenotypes (Suhre et al., 2017; Benson et al., 2018; Sun et al., 2018; Emilsson et al., 2018). However, unlike studies in humans, the use of genetic reference panels (GRPs) enables accurate control of the environment, breeding patterns, and access to a range of tissues for molecular analysis. Systems genetic analyses incorporating transcriptomics, lipidomics, and/or metabolomics in GRPs have led to the discovery of a range of novel regulators of complex phenotypes ranging from insulin resistance (Parks et al., 2015), insulin secretion (Keller et al., 2019), atherosclerosis (Bennett et al., 2015), lipid metabolism (Jha et al., 2018a; Jha et al., 2018b; Linke et al., 2020), cardiac hypertrophy (Rau et al., 2015), cardiac diastolic dysfunction (Cao et al., 2022), and many more. The inclusion of proteomics into systems genetic analysis provides information on an important biological layer and has been performed in a range of GRPs including yeast (Foss et al., 2007; Picotti et al., 2013; Parts et al., 2014), worms (Singh et al., 2016), fruit fly (Okada et al., 2016), plants such as maize (Hu et al., 2017; Jiang et al., 2019), and several livestock such as cattle (Boudon et al., 2020) and pig (Bovo et al., 2018). The use of proteomics to analyse the liver proteome of mouse GRPs has also gained popularity and been used to analyse the BxD panel (Wu et al., 2014; Williams et al., 2016), the Hybrid Mouse Diversity Panel (HMDP) (Ghazalpour et al., 2011; Parker et al., 2019), and cohorts of the Collaborative Cross/Diversity Outbred (CC/DO) (Chick et al., 2016). More recently, several studies have performed proteomic analysis of additional tissues from cohorts of the BXD (Williams et al., 2018) and CC/DO (Xiao et al., 2022) and include further phenotypic associations. Functional screening of genetic perturbations with high-throughput phenotypic measurements have identified novel regulators of muscle biology and disease. These include forward genetic mutagenesis screens to identify regulators of skeletal muscle development and locomotion in zebrafish (Birely et al., 2005; Horstick et al., 2013; Johnson et al., 2013; Bennett et al., 2018) and worms (Beron et al., 2015), RNAi screening of muscle size and function in fruit fly (Kao et al., 2021; Graca et al., 2021), and CRISPR/Cas9 screening of muscle cells to identify regulators of myogenesis and cell survival in dystrophy models (Bi et al., 2017; Lek et al., 2020; Ashoti et al., 2022). Here, we combined forward genetics via proteomic analysis of a diverse mouse panel with a targeted reverse genetics screen via AAV6 vector-mediated expression of shRNAs to knock down specific genes in bioengineered skeletal muscle to identify candidate regulators of skeletal muscle function. Our approach identified UFMylation as a regulator of skeletal muscle function that was validated in vivo. Results Genetic regulation of the mouse skeletal muscle proteome To begin to understand how variations in the genome drive changes in the skeletal muscle proteome, we performed a proteomic analysis of gastrocnemius muscle from 73 inbred mouse strains of the HMDP that were fed a chow diet and housed under identical environmental conditions (n=2–4; 161 mice). The proteomic data were integrated with previously acquired genomic and various molecular/phenotypic data via systems genetics analysis (Figure 1A). Proteomics was performed with eighteen 10-plex tandem mass tag (TMT) experiments each consisting of nine strains plus a pooled common internal reference. Peptides were analysed by 2D-liquid chromatography coupled to tandem mass spectrometry (2D-LC-MS/MS) resulting in the quantification of 5350 proteins with 4027 quantified in >50 mice and 2069 proteins quantified in all 161 animals (Figure 1B and Supplementary file 1). Biological replicates showed similar proteomes, as evidenced by the hierarchical clustering dendrogram, where mice of the same strains are typically neighbouring or located closely (Figure 1—figure supplement 1). Median intra-strain coefficient of variation was 9.2% while inter-strain coefficient was significantly larger, suggesting reproducible quantification and genetically driven variation in the proteome was captured in the data (Figure 1C). Figure 1 with 2 supplements see all Download asset Open asset Proteome-wide systems genetics analysis of the mouse skeletal muscle proteome. (A) Overview of the experimental design. (B) Number of proteins identified. (C) Intra- and inter-strain coefficient of variation. (D) Protein-quantitative trait loci (pQTL) Manhattan plot. (E) pQTL variant and gene location density. (F) Ribosomal proteins correlation and variant network (upper), and scatterplots expressed as Log2(ratio to control) showing correlation coefficient calculated using biweight midcorrelation (n=161) (lower). (G) Genetic associations of variant hotspot on chromosome 13 associated with mitochondrial complex V subunits in trans. (H) Intragenic variants associated to ACADL abundance. (I) Boxplot showing variant allele associated to EPHX1 abundance (Student’s t-test). (J) EPHX1 Arg338Cys mutation DynaMut protein flexibility analysis. To identify possible genetic factors regulating protein abundance, we next associated single nucleotide polymorphisms (SNPs) to the abundance of 4027 skeletal muscle proteins quantified in >50 mice via a protein-quantitative trait loci (pQTL) mapping. We identified significant cis-regulation of the proteome with local SNPs associated to the abundance of 527 unique proteins (±10 Mb of gene; local adjusted p<1 × 10−4) and trans-regulation with distant SNPs (>10 Mb of the gene or on a different chromosome) associated to the abundance of 170 unique proteins (global adjusted p<5 × 10−8) (Figure 1D and Supplementary file 2). Visualizing SNP densities identified regions with a high number of variants (Figure 1E). Of particular interest are the genetic ‘hotspots’ (observed as vertical ‘streaks’ on chromosomes 7, 10, 12, and 14) abundant in trans-associations affecting genes in various locations. These regions agree with those observed in eQTLs of the HMDP (Lusis et al., 2016). We further performed variant effect predictions and annotated the position of variants relative to gene location (within a gene, i.e., intragenic or between genes, i.e., intergenic SNPs) to investigate possible mechanisms regulating protein abundance (McLaren et al., 2016). We observed a complex array of pQTLs on chromosome 1 with 46 and 21 proteins regulated in -cis and -trans, respectively. The majority of these pQTLs are in the distal region of the chromosome previously described as the QTL-rich region on chromosome 1 (Qrr1) with separate linkage disequilibrium (LD) blocks within 1qH2.1, 1qH5, and 1qH3 (Mozhui et al., 2008). The various loci in Qrr1 containing these pQTLs have been associated with a range of neural, behavioural, and cardiometabolic phenotypes, and form complex co-regulatory networks modulating a range of pathways such as RNA metabolism and translation. Here, we show a cis-pQTL on chromosome 1 (lead SNP rs31934459) is associated to the abundance of RPL7 and is also a trans-pQTL for RPL19 and RPL23. These proteins are correlated in abundance and form a co-regulated network with other 60S ribosomal proteins suggesting genetic variants regulating RPL7 abundance subsequently regulate protein complex stability/assembly (Figure 1F and Supplementary file 3). The distal region of chromosome 13 also contains a complex series of cis- and trans-pQTLs giving rise to ‘hotspot’ co-regulation. A cis-pQTL associated with MOCS2 is also a trans-pQTL for eight proteins all members of mitochondrial complex V (Figure 1G). MOCS2 is involved in molybdopterin biosynthesis and it is unclear how MOCS2 might regulate mitochondrial complex V abundance. Patients with mutations in MOCS1 that result in mild molybdenum cofactor deficiency also display reduced mitochondrial respiration suggesting a link between molybdopterin biosynthesis and ATP production (Grings et al., 2019). Among the 527 proteins with a cis-pQTL association, 212 had an intragenic association. For example, non-coding SNPs in the fifth and sixth intron of Acadl were cis-pQTLs for ACADL (Figure 1H). We also identified 14 missense variants as pQTLs such as the rs32746574 variant (GG >AA; R338C), which was associated with significantly lower abundance of EPHX1 (Figure 1I). The R338C mutation was found to be deleterious (PROVEAN score: –5.6) (Choi and Chan, 2015), destabilizing (FoldX:+1.1 ΔG) (Delgado et al., 2019) and possibly damaging (PolyPhen2 HumDiv: 0.62) (Adzhubei et al., 2013; Figure 1—figure supplement 2). R338C increases the flexibility of EPHX1 due to the missing interaction (hydrogen bond) between R338 and Y291 (Figure 1J). Aligning human and mouse EPHX1 protein structures revealed that the missense mutation causing R338C (yellow) is in close proximity to the catalytic (blue) and disease-related sites (red) identified by Gautheron et al., 2021; Figure 1—figure supplement 2. Taken together, our proteomic analysis of the HMDP has helped define the genetic factors potentially regulating the abundance of hundreds of skeletal muscle proteins. Systems genetics integration of the skeletal muscle proteome with molecular and phenotypic traits The renewable nature of extensively characterized inbred mouse strains from genetic references panels such as the HMDP and the BxD allow for the integration of data across multiple cohorts (Lusis et al., 2016; Ashbrook et al., 2021). We focused our analysis on a subset of 300 molecular or phenotypic traits incorporating various plasma metabolites, lipids, and cytokines; whole body measurements such as glucose/insulin sensitivity and body composition/organ weights; and muscle phenotypes such as cardiac and skeletal muscle function previously quantified in the same strains of mice in the HMDP (Lusis et al., 2016). Note that data integration is performed at the strain level, since the proteomic data was not generated from the same mice, as those used in previous studies. Supplementary file 4 summarizes all the phenotypic data integrated in the current study and includes data sources. Cumulatively, these previous systems-level mouse studies identified hundreds of loci associated to molecular or phenotypic traits (mol/pheQTLs) (Figure 2A and Supplementary file 5). We next devised a three-step systems genetic analyses to associate these molecular or phenotypic traits to the skeletal muscle proteome that included: (1) identification of SNPs shared between skeletal muscle cis-pQTLs and mol/pheQTLs, (2) correlation and supervised multivariate associations between the abundance of skeletal muscle proteins and each molecular or phenotypic trait, and (3) comparison of protein abundance and molecular or phenotypic differences between allelic variations (Figure 2B). All data can be browsed at muscle.coffeeprot.com and includes querying at the protein- or phenotype-level followed by several interactive visualizations. Analysis of this new resource identified hundreds of associations allowing for the prioritization of potential causal proteins regulating molecular or phenotypic traits. We present a few examples of phenotype-protein associations below. For example, a locus on chromosome 7 contained cis-pQTLs associated to MCEE that were shared with fasting glucose in mice subject to high-fat/high-sugar feeding (Figure 2C). Homozygous allelic variation of SNPs were associated to MCEE protein abundance and also a trend for greater visceral adiposity and fasting glucose on a chow diet (Figure 2D–E). MCEE functions as a methylmalonyl-CoA epimerase important for propionyl-CoA metabolism. Rare autosomal recessive missense mutations in MCEE have been identified in patients with methylmalonic aciduria (Bikker et al., 2006). We also identified cis-pQTLs associated to the abundance of both OCIAD1 and OCIAD2; two neighbouring genes located on chromosome 5. Intragenic SNPs in both genes were associated with several cardiometabolic parameters of adiposity including percentage body fat assessed by nuclear magnetic resonance, retroperitoneal fat mass, circulating free fatty acids, HDL, and cholesterol. Furthermore, the abundance of OCIAD1 and OCIAD2 were positively correlated to several of these traits, and homozygous allelic variation of SNPs in both the Ociad1 and Ociad2 loci were associated to protein abundance and cholesterol (Figure 2F–G). Very little is known about the functions of OCIAD1/2 but recent data have revealed a role in mitochondrial complex III assembly (Le Vasseur et al., 2021) and human GWAS analysis has identified variants in both the OCIAD1 and OCIAD2 loci are associated with susceptibility to type 2 diabetes (Vujkovic et al., 2020). We also identified genetic variants associated to the abundance of the ER resident glucose-6-phosphate transporter SLC37A4 that co-localize to fasting glucose and fat pad mass in mice fed a chow diet (Figure 2H–I). SLC37A4 has enhanced expression in the liver, gut, and kidney, and plays a role in the regulation of glycogenolysis and gluconeogenesis, however, its role in skeletal muscle metabolism has been comparatively less studied. Patients with mutations in SLC37A4 present with several metabolic complications, particularly hepatomegaly and enlarged kidneys due to the accumulation of glycogen but also often have signs of dyslipidaemia and hypoglycemia. Figure 2 Download asset Open asset Protein and phenotype quantitative trait locus (QTL) analysis. (A) Manhattan plot and genomic location distribution of mol/pheQTLs. (B) Overview of the three-step integrative analysis approach. (C) Mirrored Manhattan plots of MCEE and glucose QTLs. (D) Allelic variant boxplots of rs31160203 for MCEE and visceral fat. (E) Allelic variant boxplots of rs50173258 for MCEE and glucose. (F) Correlation scatterplot of OCIAD2 abundance expressed as Log2(ratio to control) and plasma cholesterol concentrations. (G) Allelic variant boxplots of rs33256997 for OCIAD2 and plasma cholesterol. (H) Mirrored Manhattan plots of SLC37A4 and glucose QTLs. (I) Mirrored Manhattan plots of SLC37A4 and fat pas mass QTLs. (J) Average distribution of lean mass per mouse strain. (K) Orthogonal partial least-squares (OPLS) loading plot of proteins explaining the variance related to strain lean mass. Separation on the x-axis shows variation related to the predictive component (p1), whilst the y-axis shows the orthogonal component (o1). Highlighted points reflect Student’s correlation p-values for multiple biweight midcorrelations of proteins correlated with lean mass (–0.3< r > 0.6, p<0.05). Correlation of lean mass and the protein abundance expressed as Log2(ratio to control) of INMT (L), MOCS2, (M) and TACC2 (N). Allelic variant boxplots of selected single nucleotide polymorphisms (SNPs) with lean mass and INMT (rs49460035) (O), MOCS2 (rs28163611) (P), and TACC2 (rs32292483) (Q). *p<0.05, **p<0.01, ***p<0.001, ****p<0.0001 by Student’s t-test. We next analysed associations of the skeletal muscle proteome to lean mass which showed a range of genetic variation in female mice fed a chow diet (Figure 2J). Using supervised multivariate and pairwise correlation analysis, we identified 300 proteins positively or negatively associated to lean mass (–0.3 < r > 0.6, p<0.05) (Figure 2K). Proteins involved in translation including ribosomal subunits and elongation initiation factors were positively correlated while we observed several extracellular matrix proteins including collagens and cathepsins to be negatively correlated with lean mass. Allelic variations were integrated to highlight several trends such as negative correlations between INMT, MOCS2, and TACC2 versus lean mass (Figure 2L–N), where alleles which reduced protein abundance (rs49460035, rs28162611, and rs32292483) generally displayed lower lean mass (Figure 2O–Q). We also investigated the relationships between mol/pheQTL’s and skeletal muscle cis-pQTLs in or around the Qrr1 region, located distal on chromosome 1. Figure 3A displays regional association plots of pQTLs (arrows indicating gene location) shared with mol/pheQTLs including HOMA-IR, lean mass, plasma concentrations of keratinocyte-derived growth factor (KC; CXCL1), and plasma cholesterol. Two genomic regions are associated with HOMA-IR which contains cis-pQTLs for MPZ, NIT1, PCP4L1, and UFC1 in the first region while EPHX1 is located in the second more distal region. It is important to note that these regions were also associated with plasma insulin concentrations and as such, variations in the expression of these proteins may play greater roles in other tissues such as the pancreas. The QTL for plasma cholesterol and KC concentrations were only shared to the first region of the HOMA-IR blocks while the QTL for lean mass was distributed across both regions. To further investigate these associations, we correlated the abundance of the proteins in skeletal muscle to the molecular or phenotypic traits (Figure 3B). Our data identified a negative correlation between UFC1, lean mass, and plasma cholesterol. UFC1 is the E2 ligase for the post-translational modification of ubiquitin-fold modifier 1 (UFM1) to target proteins through UFMylation. We also identified negative correlations between MPZ and plasma cholesterol/triacylglycerols (TAGs). It is unclear if MPZ plays a causative role in the regulation of lipid metabolism, however diabetic peripheral neuropathy is associated with reduced expression of Mpz, myelin abnormalities, and several defects in lipid metabolism (Cermenati et al., 2012). Positive correlations were found between BPNT1 and plasma insulin, HOMA-IR, and visceral fat while also negatively correlated with plasma cholesterol, TAGs, and lean mass. BPNT1 is involved in phosphatidylinositol phosphate and adenosine phosphate metabolism. Mice lacking Bpnt1 display liver defects but the role of this in whole body metabolism and insulin sensitivity is et al., we observed an overall negative correlation between the abundance of EPHX1 and HOMA-IR, visceral plasma insulin, and cholesterol. Our correlation results are further by associations accessed through the et al., Figure Of particular interest are the associations between 2 and The associations of UFC1 with and lean mass a potential role of this protein in skeletal muscle and/or functional Figure Download asset Open asset associations of Qrr1 region on chromosome 1. (A) Manhattan plot of selected genes located the Qrr1 region, with traits. (B) correlation (C) associations from the and and functional genomic screen in bioengineered skeletal To potential causal regulators of muscle function, we targeted genes novel skeletal muscle pQTLs and molecular/phenotypic associations and performed a targeted functional genomic screen in human skeletal et al., 2019). We focused on proteins with negative associations to lean mass, or other metabolic traits, and generated a of to knock down the expression of these proteins in an were around to and and to on the (Figure was to a contraction for of or with lower for of this were analysed by proteomics which quantified targets with 13 targets significantly reduced in abundance by (Figure of UFC1, and was at the protein and in significant increases in production the (Figure These of knockdown were with the observed negative correlation between lean mass and the abundance of UFC1, and in skeletal muscle of the HMDP for a causal regulation of muscle function. The in in muscle function in control (Figure knockdown of UFC1 increased production the of also but we were to identify in the proteomics analyses to Taken together, these data identify potential causal negative regulators of muscle function. Figure 4 Download asset Open asset Functional screening of skeletal muscle function. (A) Overview of experimental design. (B) of target proteins (C) force, and (D) of and target proteins. (Student’s relative to with UFMylation skeletal muscle function Our data that the regulation of UFC1 and changes in UFMylation may play a role in muscle function. We first investigated the regulation of UFMylation in a mouse model of amyotrophic lateral a disease that involves muscle et al., 2014). The model involves muscle from of due to expression of an mutation in 1 et al., of we observed a significant in and free UFC1, and in gastrocnemius skeletal muscles of mice, which was of changes in the abundance of as a of ER (Figure These data the overall of increases and is an in UFMylation We the first of changes in UFMylation muscle in vivo. We next in vivo UFMylation by into and skeletal muscles of mice in a experimental where muscles of the and muscles of the (Figure of we observed significant in the abundance of UFC1 and (Figure was a but significant in whole muscle but in muscle mass (Figure muscles were to vivo of muscle function using of UFC1 increased specific production in response to a single contraction and generated a of to (Figure Furthermore, both the and specific increased across all knockdown of UFC1 (Figure Figure Download asset Open asset UFMylation is regulated in and skeletal muscle function. (A) and (B) of UFMylation and in a gastrocnemius muscle of a mouse model of amyotrophic lateral (C) Overview of the experimental design. (D) of muscles with and (E) of (F) (G) mass and (H) mass vivo analysis of contraction in muscles showing (I) single contraction to (J) contraction to and (K) and specific to or Student’s Student’s Figure data 1 file containing as and a the used to Figure and all results for each The of each is 1. Download Figure data 2 file containing as and a the used to Figure and all results for each The of each is 1. Download We next performed a more and molecular analysis. we analysed of muscles by of which revealed differences in the abundance of but was a trend for a in the number of muscle knockdown of UFC1 (Figure we performed a proteomic analysis of muscles which quantified proteins of which were regulated in abundance knockdown of UFC1 (Figure and Supplementary file The up-regulated pathways proteins associated with muscle and to the while was observed in pathways and cell (Figure The of proteins knockdown of UFC1 includes which and the protein with the (Figure The up-regulated proteins in the

Source: eLife

DOI: 10.7554/elife.82951.sa2