Regulation of proline-directed kinases and the trans-histone code H3K9me3/H4K20me3 during human myogenesis

We present a system-level analysis of proteome, phosphoproteome, and chromatin state of precursors of muscle cells (myoblasts) differentiating into specialized myotubes. Using stable isotope labeling of amino acids in cell culture and nano-liqud chromatography-mass spectrometry/mass spectrometry, we found that phosphorylation motifs targeted by the kinases protein kinase C, cyclin-dependent kinase, and mitogen-activated protein kinase showed increased phosphorylation during myodifferentiation of LHCN-M2 human skeletal myoblast cell line. Drugs known to inhibit these kinases either promoted (PD0325901 and GW8510) or stalled (CHIR99021 and roscovitine) differentiation, resulting in myotube and myoblast phenotypes, respectively. The proteomes, especially the myogenic and chromatin-related proteins including histone methyltransferases, correlated with their phenotypes, leading us to quantify histone post-translational modifications and identify two gene-silencing marks, H3K9me3 and H4K20me3, with relative abundances changing in correlation with these phenotypes. ChIP–quantitative PCR demonstrated that H3K9me3 is erased from the gene loci of myogenic regulatory factors namely MYOD1, MYOG, and MYF5 in differentiating myotubes. Together, our work integrating histone post-translational modification, phosphoproteomics, and full proteome analysis gives a comprehensive understanding of the close connection between signaling pathways and epigenetics during myodifferentiation in vitro.

Human skeletal myogenesis is a complex process of differentiation of proliferating myoblasts into myotubes. As mononucleate myoblasts exit the cell cycle, they express differentiationspecific genes, elongate, and fuse into multinucleate myotubes. Transcriptome analysis in mouse myoblast cell lines and human myoblasts indicate extensive and highly coordinated changes in muscle-specific genes, homeodomain factors, and paired-box transcription factors during myogenesis (1). These transcriptional changes are accompanied by key epigenetic regulatory steps and complex temporal protein dynamics, the underlying mechanics of which are poorly understood. Alterations in these molecular control mechanisms lead to muscular diseases. Therefore, a more comprehensive characterization of the proteomic and epigenetic regulatory components is necessary and can be exploited to develop new therapies for muscular diseases.
During embryogenesis, the temporal and stage-specific expression of four closely related proteins, MyoD1, MyoG (myogenin), Myf5, and Myf6, collectively called the myogenic regulatory factors (MRFs), 3 commit cells to a myogenic fate (2). Transcriptome analyses have revealed new skeletal musclespecific genes (3,4), a few of which were validated using antibody-based protein profiling to have exclusive expression in skeletal muscles (5). Proteomic analysis using MS, while high throughput has challenging dynamic range, because tissue-specific proteins (e.g. myosin) are expressed in high abundance to other proteins. Previously, the SILAC/LC-MS/MS approach had demonstrated only ϳ1170 proteins of the mouse myotube proteome, of which 379 were shown to be differentially regulated during differentiation (6). Much lower sensitivity was reported for human myoblasts and myotubes (7,8). Here, by using the LHCN-M2 human skeletal myoblast cell line, we employed SILAC/nano-LC-MS/MS and quantified 5345 proteins by increasing sensitivity via peptide fractionation, i.e. strong-cation exchange chromatography.
In this work, we performed an initial screening of regulated motifs around phosphorylation sites to predict the major kinases involved in human myoblast differentiation; then we selected four inhibitors targeting different kinase families and investigated the proteome, phosphoproteome, and chromatin state (via histone PTM analysis) to characterize the proteomics changes during the canonical and the altered differentiation This work was supported by National Institutes of Health Grants R01GM110174 and R01AI118891 and Department of Defense Grant W81XWH-113-1-0426. The authors declare that they have no conflicts of interest with the contents of this article. The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health. This article contains Tables S1-S5 and Figs. S1-S3. The mass spectrometry proteomics data have been deposited to the Proteome Xchange Consortium via the PRIDE partner repository with the data set identifier PXD007589. 1 Present address: Regeneron Pharmaceuticals, Tarrytown, NY 10591.
process. Direct phosphorylation of transcription factors by mitogen-activated protein kinases (MAPKs) triggers recruitment of chromatin remodeling complexes, histone acetyl transferases, and histone deacetylases (HDACs) to target genes (9). In C2C12 cells, H3K4 methylations are known to increase across the MYOG gene locus during p38 MAPK activation (10). Similarly, phosphorylation of MyoD1 inhibits its acetylation in myoblasts, and MyoD1 dephosphorylation followed by acetylation is a critical determinant of myodifferentiation (11,12). Our results highlighted the abundance of two histone marks with known role in gene silencing, H3K9m3 and H4K20me3, regulated over the MRFs genes. This newly described role for these dual marks may provide additional insight into the mechanism of muscle differentiation.

Phosphoproteome profiling identifies kinases active during myodifferentiation
An in vitro model for studying human myogenesis has recently become available with the LHCN-M2 myoblast cell line that differentiates into myotubes in low-serum conditions (13). We set up forward and reverse label-switch experiments for controls in triplicate wherein light and heavy conditions are exchanged to identify variables attributable to labeling. In the forward experiment, we grew proliferating myoblasts in SILAC medium while maintaining matched cultures in light unlabeled medium to differentiate into myotubes (referred as forward controls) and performed phosphoproteomic analysis (Fig. 1). Myodifferentiation was confirmed by fused cell morphology and expression of myosin heavy chain (MHC) expression, both of which were absent in myoblasts (see Fig. 3A).
We identified and quantified 5362 phosphopeptide ion pairs, which allowed for the relative quantification of the phosphorylation changes between the two cell states (Table S1). Each experiment was performed in three biological replicates, and the obtained ratio was the average of these three replicates. We filtered out phosphopeptides with a phosphorylation site mapped with a nonconfident localization score (A score Ͻ 20) corrected by a false-discovery rate (confident if p Ͻ 0.01). We obtained a list of 3479 phosphopeptides (Table S2), which we then filtered further for fold change Ն1.4 or Յ0.71 between myoblast and myotube ( Fig. 2A). Approximately 88.5% of the phosphorylation events were mapped on serine residues, 10.6% on threonine, and 0.9% on tyrosine. Because kinases are the major players in cell signaling pathways, we used Motif-X (14) with significance threshold of 0.000001 to identify the enriched kinase motifs for different categories (up-and down-regulated in myoblasts and myotubes) using the human IPI database as background (full list of mapped motifs in Table S3). The proline-directed ([pS]-P) motif was the most over-represented in the regulated phosphopeptides in both myoblasts and myotubes (Fig. 2B). This is the consensus motif of the CMGC family of kinases, which includes cyclin-dependent kinases (CDKs), mitogen-activated protein kinases (MAPKs), glycogen synthase kinases (GSK), and CDK-like kinases. CDKs sustain cell-cycle progression, in association with cyclins in myoblasts; this activity requires down-regulation during cell-cycle exit and differentiation in myotubes. GSK3 is also known to repress myoblast differentiation and fusion (15). MAPKs typically stimulate growth-related signals during myoblast proliferation, implying that they should also be tamped down during differentiation. By performing Pathway Linker analysis using our data set, we con- Figure 1. Workflow for proteome-wide phosphorylation analysis. The "forward" experiment was performed in three biological replicates to predict kinases actively involved in myogenesis. Proliferating mononuclear myoblasts were SILAC-labeled. In the "reverse" experiment (drug treatment analysis), SILAC labeling was applied to differentiating myotubes, whereas unlabeled cell cultures were proliferating myoblasts (control) or cells treated with drugs. The kinase inhibitors were 1 M PD325901, 1 M GW8510, 10 nM roscovitine, and 3 M CHIR99021. Light and heavy amino acid-labeled cells were lysed and mixed in equal amounts prior sample processing. Peptides were purified and fractionated using SCX. The fractions collected were enriched for phosphopeptides using TiO 2 , and flow-through was analyzed for total proteome. The data shown are averages of at least three independent experiments. nLC, nano-LC.

Stage-specific protein and PTM dynamics in human myogenesis
firmed that protein kinases, and specifically cyclin-dependent protein kinases, were the most enriched protein categories (Fig.  2C). Based on this prediction, we chose four inhibitors of four kinases targeting regulated phosphopeptides during differentiation process aiming to alter the normal differentiation process: GW8510 and roscovitine targeting CDK/CDC2, PD325901 targeting MAPK, and CHIR99021 targeting GSK.

Kinase inhibition during proliferation alters differentiation phenotype
We then performed the reverse experiment to complement the results from forward experiment, as well as to examine the effects of kinase inhibition on myogenesis. For this, we cultured myoblasts in light culture medium with and without the four selected inhibitors. The myoblasts were harvested as control, and the four drug-treated myoblast cultures were differentiated in light medium. Control and drug-treated myoblasts proliferated at the same rate (observation based on confluence). A parallel myoblast culture maintained in drug-free, heavy SILAC proliferation medium and subsequently differentiated in heavy serum-poor SILAC medium served as the myotube control; the light myoblast and heavy myotube controls of the inhibitor study will be referred to as reverse controls. We observed growth and morphology using microscopy and MHC expression using fluorescent microscopy to confirm differentiation of myoblasts into myotubes (Fig. 3B). By day 6, myotubes show a clear phenotype as spindle-shaped multinucleate cells, expressing MHC (Fig. 3B, panel ii, 6 days of differentiation). Surpris-ingly, on day 2 of differentiation, whereas control myotubes had no MHC expression yet (Fig. 3B, panel i), cells treated with either PD325901 or GW8510 prematurely differentiated, evident by fused cell morphology and MHC. On the other hand, cells treated with roscovitine and CHIR99021 showed myoblast-like morphology and no MHC expression (Fig. 3B). Although both roscovitine and GW8510 inhibit CDK2, roscovitine additionally also inhibits CDK5, possibly contributing to its antidifferentiation effect. The involvement of CDKs in MyoD1 expression has been demonstrated extensively in previous studies. In fact, CDK2 can phosphorylate pRb and disrupt MyoD1 function (16), in addition to degrading MyoD1 directly by phosphorylating at serine 200 and initiating ubiquitinproteosomal pathway (17). Like CDKs, GSKs also have progrowth effects. When we plated differentiation-arrested cells treated either with roscovitine or CHIR99021, they did not become confluent, indicating that proliferation was indeed suspended (data not shown). To check whether differentiation was delayed, we maintained them in differentiation medium for a week longer than control, but it only led to cell death (results not shown). These observations prompted us to explore the functional significance of CMGC kinase inhibition during human myogenic differentiation.

Functional analysis validates phenotypic changes induced by kinase inhibition
From six experiments (four drugs plus forward and reverse controls; Fig. 1), we identified and quantified a total of 6976

Stage-specific protein and PTM dynamics in human myogenesis
proteins (Table S4). Of these, ϳ1166 showed regulated phosphosites, and ϳ623 proteins were regulated during normal versus perturbed myogenesis (Table S5). To assess the sensitivity of the analysis in detecting muscle-related proteins (acquired using the GeneCite tool (18)), we plotted the quantified proteins as fold change of their relative abundance in differentiated myotubes versus proliferating myoblasts (Fig. S1A). We then overlaid proteins mapped by GeneCite to be related to muscle formation (PubMed repository). A third of these proteins were not known to be associated with myogenesis by this search. Proteins such as ARSB, BCAM, SPB9, and TPPP3 have no reported protein expression in any muscle tissue so far (www. proteinatlas.org) (5). 4 One of the "novel" myogenesis-associated proteins from our study, AK1C1 (20-␣-hydroxysteroid dehydrogenase) has only been demonstrated at the mRNA level in many nonmuscle tissues (19), whereas another, AK1C2 (Aldo-keto reductase family 1 member C2) has only recently been demonstrated at protein level in dermis by MALDI-MS, despite its known mRNA overexpression in skeletal muscle (20). Another interesting new myogenesis-associated protein, ANS1A, has reported function only in neuron remodeling, mediated through EphA8 signaling (21), also important for myogenesis. More work will be needed to confirm the functional significance of these proteins in myogenesis.
Further, when we filtered the proteins from this study against Uniprot proteins, we found ϳ2801 proteins annotated with the term "muscle." Plotting their fold changes across all the quantified proteins, many proteins associated with muscle formation were regulated in our data set (Fig. S1, A and B). Of the regulated proteins, those with up-regulated levels in myotubes enriched for steroid pathway signaling (AK1C1 and DHI2), cytochrome reductase (NB5R1), ER transmembrane protein (DJB14), probable methyltransferases (MET7A, MET7B), and nucleolar protein (NPA1P). Such regulation in protein abun-dance were also noted during kinase inhibition (to be discussed in Fig. S2).
In particular, plotting the myosin proteins that were specific for muscle lineage, we found that treatment with PD325901 and GW8510 resulted in similar patterns of protein expression, as did treatment with roscovitine and CHIR9902; this similarity in their effect on protein abundance was consistent with phenotypes they induce (Fig. S1C). Next, we used iGPS (Groupbased Prediction System with the interaction filter) (22, 64) (http://igps.biocuckoo.org/) 4 to predict the potential protein kinases corresponding to identified phosphosites. We plotted the log2 fold change of the phosphopeptides (myotube/myoblast) averaged depending on their predicted kinase, and the graph was sorted from smallest to largest fold change (Fig. S1D). We found that phosphoproteins with phosphorylated sites relatively more abundant in proliferative conditions mapped to the MAPK pathway, whereas those more abundant in differentiating myotubes were enriched in the GSK pathway. This denotes that, when the kinase activities by PD325901 and CHIR99021 are disrupted, differentiation would accelerate and stall, respectively. This was completely consistent with our morphological observations (Fig. 3) and an independent confirmation of our hypothesis generated by motif analysis.
Next, we performed Gene Ontology (GO) enrichment analysis of the proteins with significantly regulated relative abundance in control and kinase inhibition conditions (Table S4) using ClueGO (23) from Cytoscape (24) (Fig. S2, A-E). GO term enrichment revealed actin filament-based process, cellular component assembly, and movement of cell or subcellular component to be overrepresented in myoblasts (Fig. S2A), whereas proteins important for muscle-related functions enriched in myotubes (Fig. S2B). However, the rest of the inhibition conditions did not reflect this trend, especially accelerated differentiation not including enrichment of muscle proteins (Fig. S2, D and E). Because of ambiguous results obtained from the observation of the proteome, we decided to analyze

Stage-specific protein and PTM dynamics in human myogenesis
the phosphorylations to verify whether their regulation might be independent from protein abundance.
By monitoring phosphorylation events on myo-related proteins, we found that the phosphorylation of specific sites of myogenic proteins and their relative abundance conferred differentiation specificity (Fig. S2F). Myo9B associated with actinfilament and ATP (movement) binding was differentially phosphorylated in the responses to drug treatment. On this protein, we identified two erstwhile unknown phosphorylation sites, Ser-1290 and Thr-1271, which were significantly down-regulated in control and GW8510-treated myotubes (promoting differentiation). Both these sites were prominently phosphorylated in the presence of CHIR99021 (arresting differentiation). Another myogenic protein, Myo18A (actin-filament binding) was phosphorylated at a novel phosphosite, Ser-1639 in cells treated with PD325901, whereas its paralog, Myo18B (intracellular trafficking), was phosphorylated at an unreported site: Ser-729 in control myotubes. Interestingly, CHIR99021 was the only condition that has significant phosphorylation on the major myogenic determinant protein, MyoD1 (at Ser-278); this correlated with muscle structure development inhibited significantly by CHIR99021 (Fig. S2B). Among proteins regulated in CHIR99021 treatment (Fig. S2F), we show that MyoD1 with phosphorylation at Ser-278. MyoD1 is known to interact with some important muscle-related proteins such as MYH9 (cellular myosin) and MAPK14 (myoblast differentiation involved in skeletal muscle regeneration). One of the associated proteins, RRAS2 (plasma membrane GTP-binding protein) is highly expressed in skeletal muscles, whereas HSP90AA1 (cardiac muscle cell apoptotic process) has no documented skeletal muscle protein expression or function. Transcription factors (STAT3 and SNW1) and epigenetic factors (HDAC5) seemingly were co-regulated along with MyoD1. Although STAT3 is a known mediator of differentiation interacting directly with MyoD1 (25) and HDAC5 is known to repress MyoD1 (26), SNW1 seemingly has no reported functional association with MyoD1 so far.
Because straightforward GO analysis of the conditions teased out muscle-related terms only in CHIR9902 (as being down-regulated), we clustered the different protein regulations without any paired comparisons using k-means clustering of the SILAC fold changes of the proteome. From trends seen across all the analyzed conditions, we obtained an ideal grouping of 19 classes that displayed unique trends (Fig. 4). A ClueGO-Cytoscape representation of the proteins relatively more abundant in both PD325901 and GW8510 (Fig. 4A) revealed enrichment for muscle development. In contrast, another prominent trend (Fig. 4B, left panel) with PD325901, GW8510, and CHIR99021 showed down-regulation of nonmuscle proteins, further validating our model system.

Kinase modulation of myogenesis through chromatinassociated proteins
We next checked the known protein-protein interactions for highly regulated proteins in our data sets using STRING (Fig. S3) (27). Most of these proteins were reported muscle cells according to the Human Protein Atlas (www.proteinatlas.org) 4 (data not shown). Further, the implicated targets of kinase inhibitors we used were represented in these networks. For example, PD325901-treated cells featured MAPK1 that was strongly associated with MAPK3 and interacted with SFPQ, FBL, PPP1CC, and CDK4. Most of these interactors have been annotated in the context of muscle development: SFPQ is a binding partner of DUX4, a gene linked to facioscapulohumeral muscular dystrophy (28); nucleolar localization of normal proteins to fibrillarin (FBL) has been known to be affected in myositis (29); PPP1CC is possibly attributed to muscle strength (30); and CDK4 maintained myogenic population (13). All of these proteins were up-regulated in accelerated differentiation by PD325901 and GW8510. STRING revealed that there was a close functional network of transcription factors and histone variants differentially regulated among the kinase-perturbed conditions. One protein that was down-regulated in PD325901 treatment is the heterochromatin protein CBX3, known to bind to H3K9me3 to regulate gene expression (31). The histone H2A variant H2AFV was remarkably up-regulated in myotubes, as reported earlier (32), forming a hub with several prominent histones and chromatin-interacting proteins. Similarly, H3F3B (H3.3B), H2B variants (HIST1H2BD and HIST1H2BH), and those of H1 (HIST1H1E, HIST1H1B, and H1FX) were up-regulated in myotubes and all drug-treated conditions. Although STRING analysis of our data set revealed that proteins known to interact from other studies are potentially co-regulated, functional associations of other correlating, new interactor proteins mapped in this study need to be validated further by protein-protein interaction studies.

Kinase inhibition of myogenesis regulates histone PTMs
Because histones were highly regulated and featured in the STRING interactome generated for different condition (Fig.  S3), we quantified PTMs from histone H3 and H4 peptides to find the extent to which the chromatin state (histone marks) was modified to effect the phenotypic change (differentiation of myoblast to myotube). From a total of 76 post-translationally modified forms of histone H3 and H4 peptides, single PTM abundances were calculated for each peptide, and 40 were found to be significantly regulated between the control and kinase-inhibited conditions during myogenic differentiation ( Fig. 5A; analysis of variance p value Ͻ 0.05). H3K27me3, H4K16ac, and H4K20me1 appeared to be uniquely enriched in myoblasts, whereas H3K4me1, H3K18me1, H3K18ac, H3K23me1, H3K23ac, H3K79me3, and H3K122ac up-regulation was typical for myotubes. The drugs that accelerated differentiation were specifically characterized by certain marks: the repressive H3K9me1 for PD325901 and active acetylation marks such as H3.3K27, H3K79, H4K5, H4K8, and H4K12 for GW8510. This is consistent with our earlier finding that very different histone PTM abundance patterns may interplay to result in similar phenotypes (33). Next, we observed that CHIR99021 and roscovitine treatment did not have a signature of their own separately or together for drug effects but functionally co-trended with myoblasts with respect to H3K9me2/3, H3K14ac, H3K27me1, and H3K36me1. Overall, we observed an increased methylation in proliferating myoblasts, with total global methylations trending uniquely for both proliferation and proliferation arrest (roscovitine and Stage-specific protein and PTM dynamics in human myogenesis CHIR99021) (Fig. 5B). In contrast, acetylation predominated in normal (control myotubes) and accelerated differentiation (PD325901 and GW8510). These observations clearly reiterate that proliferation/myoblast-like states have a poised chromatin configuration as opposed to active chromatin of the differentiation/myotube-like conditions, as reported previously (34). Uniquely, a repressive mark, H3K9me3 showed a higher abundance in myoblast/CHIR99021/roscovitine, compared with myotube/PD325901/GW8510; other marks such as H3.3K27me3 and H4K20me3 were enriched in myoblast and CHIR99021.
Apart from this, we correlated histone marks and respective histone writers especially ϳ57 methyltransferases from our proteome data and found 10 that were regulated in abundance at least in one condition (Fig. 5C). SETD7 (monomethylation of H3K4) was more abundant in differentiated myotubes, whereas SETD3 (H3K4 and 36 methyltransferase) was up-regulated in GW8510. Both these methyltransferases were reported to promote muscle differentiation when transfected in mouse fibro-blasts and myoblasts (35,36). In contrast, another H3K4 methyltransferase, SMYD1 was increased in CHIR99021-treated cells and decreased drastically in PD325901-treated cells. Thus, abundance of catalyzing enzymes and marks did not correlate, because marks were catalyzed by redundant enzymes, as in the case of H3K4 methylations (by KMT2A, MLL3, KMT2C, SMYD1, SMYD2, SMYD3, SETMAR, SETD3, and SETD7) and one methyltransferase catalyzing more than one mark (Fig. 5C). Because ultimately we wanted to understand the role of these kinases on transcriptional regulation of muscle cell differentiation, we chose two marks, H3K9me3 and H4K20me3, that showed stage-specific expression and were regulated consistently by kinases and examined their abundance at muscle-specific genes.

Genetic mapping of the trans-histone code H3K9me3/H4K20me3
We observed that the dual repressive marks, H3K9me3/ H4K20me3 co-trended in proliferating myoblasts/myoblast- . k means clustering analysis of protein trends across analyzed conditions. A, k means analysis grouped data into 18 clusters with different trends estimated using z score normalization of protein abundance across conditions. The first represented cluster (left panel) shows proteins more abundant in GW8510 (GW) and PD0325901 (PD) treatment, because the y axis represents the normalized z score for ratio of normal differentiation/drug treatment. A negative value indicates that the function is promoted by the drug. On the right, the biological process was enriched using the Gene Ontology database. p values are illustrated on the right. Connector lines indicate similar biological functions. Ball size represents relative number of genes that belong to that function. B, left panel, cluster including proteins down-regulated in GW8510, PD0325901, and CHIR99021 (CHIR) treatment. Right panel, most of these proteins are unrelated to muscle development.

Stage-specific protein and PTM dynamics in human myogenesis
like cells (CHIR99021 and roscovitine) versus normal (myotube)/accelerated differentiation (PD325901 and GW8510) (Fig. 6A). This prompted us to look for their abundance on the major myogenic genes. Using qPCR, we showed that MYOD1, MYOG, and MYF5 not only showed stage-specific expression but were also perturbed during kinase inhibition (Fig. 6B). Although all the three genes were expressed in PD325901treated cells and myotubes, in GW8510 the expression of the master regulator of differentiation, MYOD1, alone was enough to sustain accelerated differentiation. Compared with these, the expression in CHIR99021-and roscovitine-treated cells was lower. This consistency between morphology and gene expression confirmed the direct effect of drugs on phenotype. Based on this, we hypothesized that depletion of repressive marks in differentiation-promoting conditions may allow transcription of the myogenic genes in myotubes.
So far, the role of H3K9me3 in the transcriptional regulation of muscle cell differentiation is poorly understood, although it has been reported to be crucial for lineage restriction in other cell fates (37)(38)(39). Using ChIP, H3K9 methylation was shown to be enriched in promoters of cell cycle gene, Rb/E2F as myoblasts differentiate into myotubes (40). On the other hand, others showed that H3K9me3 was neither enriched in cell cycle genes nor differentially regulated myogene expression (41,42). Recently, using Suv39h1 knockdown experiments, Jin et al. Similarly, very few studies are available about H4K20 methylation in myogenesis. Recently, it has been shown that a related H4K20 dimethyltransferase, Suv4 -20H1, promoted muscle stem cell quiescence and inhibited differentiation by suppressing MyoD1 (44). By using ChIP-qPCR, we investigated the abundance of H3K9me3 and H4K20me3 in three regulatory regions (promoter, transcription start site (TSS), and 3Ј end) of the MRF genes, MYOD1, MYOG, and MYF5 (Fig. 6C). We observed that both of these marks were depleted in control myotubes and generally less abundant in PD325901. On the other hand, both CHIR99021 and roscovitine had almost identical enrichment patterns for H3K9me3, abundant at promoter and thereby suppressing myogenic gene expression and differentiation. CHIR99021 had some depletion of H4K20me3 across myoG but still did not express mRNA possibly because of other repressive marks and combinatorial effects on the myogenic genes. Exceptionally, GW8510 was an outlier, because it was enriched for both these repressive marks across the MRFs and paradoxically still expressing at least MyoD1 and inducing differentiation. It is possible that combinatorial effects by H3K4ac (activating mark) and H3K27ac (enhancer-associated) that are abundant exclusively in myotubes and GW8510, as well as other dynamic and transient marks, may have a role in opening

Stage-specific protein and PTM dynamics in human myogenesis
chromatin. H3K9me3 is known to co-localize with H4K20me3 in ES cells, the redundancy possibly to ensure heterochromatin formation and interact with a wider range of other repressors (45,46). In addition to histone PTMs, other epigenetic factors such as DNA methylation over the genes themselves may regulate differentiation (47). It has been shown that other proteins such as pRb associate with and regulate SUV39H1 (H3K9 methylation) and Suv4 -20h (H4K20me3) (48 -50). Collectively, our results linking H3K9me3 and H4K20me3 to myodifferention provide new insights into transcriptional regulation of skeletal muscle development.

Conclusions
We combined MS and ChIP-qPCR experiments to characterize the modulation of proteins and histone PTMs during myogenesis. We identified correlations between phenotypes and proteomic profiles during treatment with kinase inhibitors, which we used to boost or delay myotube formation. We uncovered that two repressive histone marks, H3K9me3 and H4K20me3, co-regulated on MRF genes, adding crucial insight to our understanding of myodifferentiation.
Although reported as specific for particular kinases, we cannot exclude some off-target effects with the inhibitors we used. We used these inhibitors at a much higher nonlethal dose than the reported IC 50 values, so as to see their effects on myogenesis (1 M PD325901, 250 M GW8510, 10 M roscovitine, and 3 M CHIR99021). Other MEK inhibitors, U0126, MIIC, and PD98059 have been known to mediate responses later attributed to their effects on proliferation or calcium homeostasis rather than inhibition of ERK1/2 activity (51,52). However, the MAPK inhibitor we used, PD325901, was found by others to be specific without such off-target effects (52). Both GW8510 and roscovitine inhibit cell-cycle progression by G 2 /M arrest; how-ever, we observed two different phenotypic effects in this study. This may be due to the downstream pathways involved, conferring some specificity in action. In pancreatic ␣ cells, GW8510 may have JNK-and p38-dependent stimulatory effects on transcriptional activity (53); at the 1-40 M range, roscovitine inhibits a wide variety of kinases including ERK1 and ERK2 (42,54,55). Similarly, GSK3␣/␤ inhibitor, CHIR99021 is known to have off-target effects on kinases within the CDK2-cyclin A2/E cell-cycle complex (56). Despite this, we believe that the effects we observe have some specificity because they are replicated in GO enrichment in the K-clustering analysis, kinase interaction partner module, and abundance of dual repressive mark, H3K9me3/H4K20me3, on MRF genes and clearly correlate with myogenesis. Because genetic approaches with RNAi and CRISPR/Cas9 also suffer from this disadvantage, the exact mechanism underlying kinase-inhibited modulations need to be interpreted with caution.
In this study, inhibition of the MAPK pathway seems to be accelerating differentiation by acting on unique suite of proteome and histone marks. Although PD325901 is not useful as such for differentiation therapy because of toxicity observed in clinical trials (ClinicalTrials.gov identifiers NCT00147550 and NCT00174369), analogs producing similar phenotypic effects without inducing toxicity or even combination therapy at low dose with GW8510 may hold promise. So far, the effectiveness of HDAC inhibitors has led us to believe that histone acetylation is key to the reversal of muscular dystrophy such as Duchenne (57,58). Here, we propose that depletion of repressive histone methylation marks locally on the muscle-specific genes MYOD, MYOG, and MYF5 might be a new path to regulate myodifferentiation. Together, we speculate that a combination of different strategies including the kinase-directed global and

Stage-specific protein and PTM dynamics in human myogenesis
gene-specific histone PTM abundances may enhance the outcome of therapy for dystrophy and muscle injury.

Cell culture
The human myogenic cell line LHCN-M2 was a kind gift from Dr. Woodring Wright (University of Texas Southwestern Medical Center at Dallas, Dallas, TX). The immortalized myoblasts were cultured in proliferation medium and, once confluent, were differentiated in low-serum-containing medium as described previously (13). For heavy labeling of myoblasts and myotubes, Dulbecco's modified Eagle's medium deficient in arginine and lysine (Thermo Scientific) was supplemented with [ 13 C, 15 N] L-arginine (98% purity) (Cambridge Isotope Laboratories) and light lysine. Culture medium was prepared using 10% dialyzed fetal bovine serum as previously described and then filter-sterilized through a 0.22-m filter. Myoblasts were cultured on gelatinized 15-cm culture dishes in proliferation medium for 3-5 days to 80% confluence and then switched to low-serum medium to differentiate into myotubes for 6 days. Both forward and reverse experiments with heavy isotopes in either myoblast or myotube culture and vice versa were performed. To harvest, the medium was removed; then the cells quickly washed twice with sterile PBS at 37°C, trypsinized, collected off the dishes, pelleted by centrifugation, snap-frozen in liquid nitrogen, and stored at Ϫ80°C.

Treatment with kinase inhibitors
All drugs except roscovitine (Adipogen, San Diego, CA) were purchased from Calbiochem-Millipore. Drugs were dissolved in DMSO and titrated for dose response. LHCN-M2 cells were plated at 1 ϫ 10 6 cells in 15-cm culture dishes, and the drugs were used in the following concentrations: 1 M PD325901, 1 M GW8510, 10 nM roscovitine, and 3 M CHIR99021. The final DMSO concentration was Ͻ0.1%.

Proteome and phosphoproteome studies
Cell membranes were disrupted in five volumes of lysis buffer (50 mM Tris-HCl, pH 8.2, 8 M urea, 75 mM NaCl, and 1ϫ HALT complete protease inhibitors mixture; Roche) and sonicated at 15 W output with three bursts of 30 s each, followed by cooling on ice for 10 s between each sonication. The lysates were cleared by centrifugation at 16,000 ϫ g for 20 min. Protein concentrations were determined with Bradford assay. The proteins were reduced with 5 mM dithiothreitol at 55°C, cooled to room temperature, and then alkylated with 14 mM iodoacetamide at room temperature in the dark. Equal amounts of protein (ϳ5 mg each) from heavy and light extracts were mixed. Reduced and alkylated proteins were diluted 5-fold with 0.25 M Tris-HCl, pH 8.0, and digested overnight at 37°C with 1:100 enzyme: protein ratio using sequencing grade modified trypsin (Promega, Madison, WI). The resulting peptide solution was desalted using C 18 reverse phase cartridges (Sep-Pak, Waters) and lyophilized. The dried peptides were dissolved in 5 mM aqueous ammonium acetate containing 30% acetonitrile, 7 mM KH 2 PO 4 , pH 2.9 (solvent A) for fractionation with strong cation-exchange chromatography (SCX). The SCX separation was performed on a 9.4-mm ϫ 100-mm column packed with 5-m particle size, and 200 Å pore size SCX resin (PolySulfoethylA, PolyLC). The peptides were eluted using the following gradient program: 25% to 100% B (30% acetonitrile, 7 mM KH 2 PO 4 , 350 mM KCl, pH 2.65) for 33 min followed by re-equilibration of column at 100% A. The flow rate was 2.5 ml/min, and 10 fractions were collected. Each fraction was lyophilized and redissolved in 0.1% acetic acid; peptides were extracted with C 18 reversed phase cartridges (Oasis, Waters) prior to phosphopeptide enrichment. Phosphopeptide enrichment was performed using a TiO 2 protocol described previously (59). In short, ϳ0.5-2 mg of lyophilized peptides were mixed with 500 l of loading buffer (65% acetonitrile and 2% TFA saturated with glutamic acid). Approximately 2 mg of TiO 2 beads (Titansphere TiO 2 5 m, GL Sciences Inc., Japan) were equilibrated in loading buffer for 15 min (slurry), and slurry was mixed with 0.5 mg of peptide (60) at room temperature. TiO 2 beads were then packed by centrifugation in equilibrated C 18 spin columns (PepClean C 18 Spin Columns; Thermo Scientific). The beads were sequentially washed with 800 l of 65% acetonitrile and 0.5% TFA followed by 65% acetonitrile and 0.1% TFA. Unbound flow-through comprised the total protein. For phosphopeptide elution, the beads were incubated for 1 min in 300 mM ammonium hydroxide, 50% acetonitrile at room temperature and centrifuged. This sample was immediately acidified by the addition of glacial acetic acid to a final concentration of 3% as phosphosites were unstable at basic pH). Two extra elution steps in 500 mM ammonium hydroxide in 60% acetonitrile followed. The three eluates of each fraction were pooled and dried using a SpeedVac, and the pellets were stored at Ϫ80°C. Concentrated peptides were diluted with 0.1% acetic acid for LC-MS-/MS analysis.
Dried samples were resuspended in buffer A (0.1% formic acid) and loaded onto an Easy-nLC system (Thermo Fisher Scientific), coupled online with an Orbitrap Fusion Tribrid mass spectrometer (Thermo Scientific). Peptides were loaded into a PicoFrit 18-cm-long fused silica capillary column (75-m inner diameter) packed in-house with reversed-phase Repro-Sil Pur C 18 -AQ 3 m resin. A gradient of 165 min was set for peptide elution from 2 to 26% buffer B (100% acetonitrile, 0.1% formic acid) for phosphopeptides and from 4 to 28% B for unmodified peptides at a flow rate of 300 nl/min. The MS method was set up in a data-dependent acquisition mode. For full MS scan, the mass range of 350 -1200 m/z was analyzed in the Orbitrap at 120,000 full width at half maximum (FWHM) (200 m/z) resolution and 5 ϫ 10e5 AGC target value. MS/MS was performed in the ion trap in the rapid mode using the TopSpeed mode (3 s). HCD collision energy was set to 27, AGC target was set to 10e4, and maximum injection time was set to 200 ms.

Histone PTM analysis
Histone purification and analysis was performed as previously described (61). Briefly, histones were acid-extracted from nuclei with 0.2 M H 2 SO 4 for 2 h and precipitated with 33% TCA for 1 h. Purified histones were dissolved in 30 l of 50 mM NH 4 HCO 3 , pH 8.0. Derivatization reagent was prepared by mixing propionic anhydride with acetonitrile in a ratio of 1:4 (v/v) and was added to the histone sample in 1:2 (v/v) ratio for 15 min at 37°C. This reaction was performed twice. Histones were then digested with trypsin (enzyme:sample ratio 1:20, overnight, room temperature) in 50 mM NH 4 HCO 3 . After digestion, the derivatization reaction was performed again twice to cap peptide N termini. The samples were desalted prior LC-MS analysis by using C 18 Stage-tips. Samples were analyzed by using a nano-LC-MS-/MS setup. Chromatography was configured with the same type of column and HPLC as for the proteomics analysis. The HPLC gradient was as follows: 2% to 28% solvent B (A ϭ 0.1% formic acid; B ϭ 95% MeCN, 0.1% formic acid) over 45 min, from 28% to 80% solvent B in 5 min, 80% B for 10 min at a flow-rate of 300 nl/min. nano-LC was coupled to an LTQ-Orbitrap Elite mass spectrometer (Thermo Scientific). Full scan MS spectrum (m/z 290 -1400) was performed in the Orbitrap with a resolution of 60,000 (at 400 m/z) with an AGC target of 10e6. The acquisition method contained both data-dependent and targeted scans. The targeted signals were the isobaric peptides of histone H3 and H4, performed to discriminate their relative abundance with the fragment ion profiles. MS/MS was performed with collision-induced dissociation with normalized collision energy of 35, an AGC target of 10e4, and a maximum injection time of 100 ms. MS/MS data were collected in centroid mode. Precursor ion charge state screening was enabled, and all unassigned charge states, as well as singly charged species, were rejected.

Immunofluorescence microscopy
The cells were washed in PBS, fixed in 4% formaldehyde, stained for myosin heavy chain (MF-20; Developmental Studies Hybridoma Bank, University of Iowa), and visualized using Alexa Fluor 488 secondary antibody (Invitrogen). Fluorescent microscopy was performed using Nikon Eclipse TE2000-U.

Proteomics data analysis
MS data were analyzed using standard proteomics pipelines for protein identification/quantification. Specifically, the spectra obtained from the forward experiment ( Fig. 1) were identified by using Mascot 2.3 (Matrix Science). Searches were performed against human entries in the IPI database to assign peptide sequences and identify their proteins of origin. Precursor ion mass tolerance was set to 10 ppm; fragment ion mass tolerance was set to 0.6 Da. Enzyme was set to trypsin, and as many as two missed cleavages were allowed. Dynamic modifications included carbamidomethylation of cysteine, and oxidation of methionine, as well as phosphorylation of serine, threonine, and tyrosine. SILAC labeling was set for the search, considering heavy lysines (ϩ8 Da) and heavy arginines (ϩ10 Da), allowing a maximum of three labeled amino acids (Lys and Arg) per peptide. The data were filtered for 1% false discovery rate, and only Mascot scores Ͼ20 were considered (ϪLog10 match probability). SILAC ratios were only calculated for peptides that could be uniquely mapped to a given protein group.
Phosphopeptides were determined based on Mascot delta scores, and the mass spectra were peak-picked and the centroid data were read from raw tandem mass spectra. These were searched in Mascot 2.3 using the same sequence database, enzyme specificity, modifications, and mass tolerance as described above. The ratio heavy/light H/L was reported to discriminate regulated phosphorylation events. The A score (PTM localization score) was used to define unambiguously localized phosphorylations, which was essential for our motif analysis. The H/L ratio of the phosphopeptides was divided by the quan-

Stage-specific protein and PTM dynamics in human myogenesis
tified H/L of the respective protein identified from the flowthrough of the TiO 2 chromatography, also identified and quantified using Mascot tool. We filtered for phosphorylation events with at least 1.5-fold deviation between heavy and light sample after correcting for protein regulation (no statistical p value adopted). The reverse experiment (Fig. 1) including the four drug treatments was performed by using MaxQuant (62). Same settings from Mascot 2.3 were used for database searching, except that UniProt human database was used instead of IPI. The results were filtered for 1% false discovery rate. Phosphopeptide and protein quantification was performed by using the algorithm of the software, i.e. extracted ion chromatography of the precursor masses and definition of the intensity of the chromatogram (62). The protein H/L ratio was assessed by averaging the H/L ratio of all unmodified peptides quantified from the flow-through of the TiO 2 chromatographic enrichment. Statistics was performed by using the two tails heteroscedastic t test (significance when p value Ͻ0.05); once proteins or phosphopeptide ratios were log2 transformed, the t test was used to assess whether the replicates were significantly different from 0. We used the t test because of the limited number of replicates (n ϭ 3 or 4), where nonparametric tests would have not provided sufficient strength to define statistically regulated proteins and peptides. All calculations were made in Microsoft Excel. The MS proteomics data have been deposited to the ProteomeXchange Consortium via the PRIDE partner repository with the data set identifier PXD007589 and 10.6019/PXD007589.
For histone analysis, EpiProfile examined extracted ion chromatography of a list of known modified and unmodified histone peptides. The list of the masses corresponding to the selected histone peptides was defined previously (63). The peptide relative ratio was calculated using the total area under the extracted ion chromatograms of all peptides with the same amino acid sequence (including all of its modified forms) as 100%. For isobaric peptides, the relative ratio of two isobaric forms was estimated by averaging the ratio for each fragment ion with different mass between the two species. Statistics was assessed by using the two-tailed homoscedastic t test (significant when the p value is Ͻ0.05).