First profiling of lysine crotonylation of myofilament proteins and ribosomal proteins in zebrafish embryos

Zebrafish embryos are translucent and develop rapidly in individual eggs ex utero; they are widely used as models for embryogenesis and organ development for human diseases and drug discovery. Lysine crotonylation (Kcr) is a type of histone post-translational modifications discovered in 2011. Kcr dynamics are involved in gene expression regulation and acute kidney injury; however, little is known about the effects of Kcr on non-histone proteins. In the present study, we conducted the first proteome-wide profiling of Kcr in zebrafish larvae and identified 557 Kcr sites on 218 proteins, representing the Kcr event in zebrafish. We identified two types of Kcr motifs containing hydrophobic (Leu, Ile, Val) and acidic (Asp and Glu) amino acids near the modified lysine residues. Our results show that both crotonylated proteins and sites of crotonylation were evolutionarily conserved between zebrafish embryos and humans. Specifically, Kcr on ribosomal proteins and myofilament proteins, including myosin, tropomyosin and troponin, were widely enriched. Interestingly, 55 lysine crotonylation sites on myosin were distributed throughout coiled coil regions. Therefore, Kcr may regulate muscle contraction and protein synthesis. Our results provide a foundation for future studies on the effects of lysine crotonylation on aging and heart failure.


Results
Profiling lysine crotonylation in zebrafish embryos. We investigated lysine crotonylation (Kcr) modification using larvae at 72-120 h post-fertilization (hpf). This developmental stage was examined because all larvae organs are well-developed at this point. Purified proteins were examined by immunoblot assay with a specific pan-Kcr antibody (Fig. 1a). We detected multiple major protein bands with molecular weights greater than those expected for histones, indicating Kcr modifications on non-histone proteins. To obtain the global crotonylome in zebrafish larvae, proteins were prepared from 72 and 120 hpf larvae. Lysine-crotonylated peptides were immune-enriched with anti-crotonyl lysine antibody-conjugated agarose beads and identified by nano-LC-MS/ MS (Fig. 1b). The obtained MS raw data were analyzed using MaxQuant software with the zebrafish database from UniProt (41,001 sequence). MaxQuant results were filtered by MaxQuant scores of more than 40, false discovery rate of less than 1% for both protein and peptide and site localization probability of greater than 0.75. For quality control validation of the MS data, we evaluated the mass error of all identified peptides. The distribution of mass error for precursor ions was close to zero and most values were less than 0.03 Da, indicating acceptable mass accuracy of the MS data (Fig. S1a). All identified Kcr peptides exhibited different abundances depending on their lengths (Fig. S1b). In this study, 557 Kcr sites in 218 proteins were identified in pooled larvae among 508 crotonylated peptides (Table S1). In all detected peptides, 154 Kcr sites, 194 Kcr peptides and 97 Kcr proteins were identified in individual triplicate experiments (Fig. 1c). Among our Kcr results from zebrafish embryos, Kcr proteins and sites converted to human were compared with recent studies to profile non-histone protein crotonylation in HeLa and H1299 cell lines, respectively 14,15 (Fig. S1c). To compare crotonylation and acetylation in zebrafish, we used a previously acquired Kac dataset in zebrafish 26 (Fig. S1d). Among the detected Kcr, only 67 (30.7%) Kcr proteins and 52 (9.3%) Kcr sites overlapped with Kac sites.
Our data set, including 484 surrounding sequences, was evaluated to identify site-specific sequence motifs from the −7 to the +7 positions surrounding the crotonylated lysine using the Motif-X program 27 . Of all surrounding sequences, 324 sequences were matched to a total of six definitively conserved motifs (Fig. 2a). The six motifs can be divided into two types: the first type includes hydrophobic residues at the +2 position relative to Kcr (Kcr-X-L, Kcr-X-V and Kcr-X-I), while the second type contains acidic residues at the −5, −1 and +2 positions relative to Kcr (E-X-X-X-Kcr, DKcr and Kcr-X-E). Approximately 56.6% of all motif peptides showed hydrophobic amino acid motifs and 43.5% showed acidic amino acid motifs (Fig. 2b). Kcr-X-L was the most common combination, accounting for 26.9% (87) of the motifs in zebrafish larvae.

Functional enrichment of Kcr by bioinformatics.
To comprehensively analyze the distribution and function of Kcr non-histone proteins, we performed GO functional annotation and enrichment analysis. Proteins were classified by GO annotation into three categories: biological process, cellular component and molecular function derived from the UniProt-GOA database (http://www.ebi.ac.uk/GOA/). Analyses of Kcr proteins and their subcellular localization are indicated in Fig. S2. Kcr proteins were in cells (37%), organelles (27%) and macromolecular complexes (27%) in cellular components. Molecular function analysis revealed protein functions in binding (47%), catalytic activity (27%), structural molecule activity (13%) and transporter activity (10%). The analysis of biological processes showed that Kcr proteins are primarily involved in cellular processes (28%), metabolic processes (22%) and single organism processes (19%). The identified Kcr proteins were localized in the cytosol (58%), mitochondria (11%), extracellular matrix (8%) and plasma membrane (6%). Figure 3a shows the preferred targets of Kcr proteins by GO enrichment analysis. Ribosomes, non-membrane-bound organelles and macromolecular complexes were significantly enriched among cellular components and lipid transporter activity, structural constituents of ribosomes and structural molecular activity were enriched in molecular function. Translation, metabolic process and diverse regulation of skeletal muscle contraction were strongly enriched in biological processes.
Protein domain enrichment analysis was performed such that domains related to myosin were given top ranks, which refer to GO analyses results (Fig. 3b). Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis was also performed to identify the metabolic pathways involving crotonylated proteins (Fig. 3c). Kcr occurs on many proteins involved in protein synthesis such as ribosome and protein processing. The cytoskeleton was involved in Kcr modification of tight junctions. Moreover, crotonylated proteins were identified in energy metabolism such as oxidative phosphorylation, carbon metabolism, glycolysis and glycogenesis.

Evolutionarily conserved Kcr in zebrafish and humans. To analyze the conservation of homologous
Kcr between zebrafish and humans, crotonylated proteins were analyzed with BLASTP 28 . We identified 213 (97.7%) Kcr proteins in zebrafish that significantly overlapped with 189 human proteins (Table S2). Moreover, sequences of amino acids upstream and downstream of Kcr sites in zebrafish larvae were analyzed with PSI-BlastP to identify conserved surrounding Kcr sites. We found that 428 zebrafish Kcr sites (76.8%) overlapped with 339 human surrounding sequences. Furthermore, 323 human Kcr sites contained the same lysine residues as 401 Kcr sites (72.0%) in zebrafish. These results indicate that Kcr sites and modifications are highly orthologous between human and zebrafish at the protein level.
Next, to compare potential cross-talk/interplay between Kcr and other PTMs such as Kac and lysine ubiquitination (Kubi), Kcr results were compared with the results in humans using the Kac and Kubi database from PhosphoSitePlus 29 . Of the total 323 Kcr sites converted to humans except for redundant sites, 95 Kac (29.4%) and 42 Kubi (13.0%) sites were detected in the same lysine residue. Kcr is expected to be involved in PTM cross-talk/ interplay with competitive Kac rather than with Kubi.
According to the results of functional enrichment and BLASTP, Kcr in zebrafish is abundant in ribosomal proteins and myofilament proteins, including myosin, TM and troponin. A summary of crotonylated myofilament and ribosomal proteins among zebrafish and humans is presented in Tables 1 and 2, respectively. We identified 194 crotonylation sites in myofilament proteins in zebrafish larvae, including 156 Kcr sites on myosin and myosin   (Table 1). Remarkably, myosin1 was high conserved from zebrafish to human and contained 63 crotonylation sites in zebrafish larvae, accounting for approximately 11% of the total crotonylation sites without overlap between myhz1.1 (B8A568) and 1.2 (B8A561). Interestingly, 55 crotonylation sites in myosin were distributed intensively in the coiled coil motif (amino acids 841-1937) compared to the myosin motor (amino acids 87-780) (Fig. 4). Furthermore, myosin light chain, TM and troponin were largely conserved from zebrafish to humans and contained several crotonylation sites that are also highly similar to site-specific lysine residues in humans ( Fig. S3).
Furthermore, 34 lysine sites on 21 ribosomal proteins were crotonylated in zebrafish larvae (Table 2). Crotonylated ribosomal proteins in zebrafish were found to share at least 85% homology and highly conserved site-specific lysine residues with humans. These results demonstrate that Kcr is evolutionarily well-conserved between zebrafish larvae and humans, at both the protein and amino acid levels.

Discussion
Lysine crotonylation (Kcr) was reported as a new type of PTM in 2011 and serves as a powerful indicator of active cellular genes near histones 11 . Recent studies showed that Kcr levels can be regulated by crotonylating P300 and decrotonylating sirtuin 3 12,13 . Although increased histone crotonylation is related to acute kidney injury due to cell stress and crotonate availability 30 , previous studies on Kcr have focused on the histone regulation of epigenetics. Non-histone Kcr substrates have not been identified. Here, we evaluated Kcr modification of non-histone proteins in zebrafish larvae. We identified and validated 557 Kcr sites on 218 proteins by immunoprecipitation with MS-based proteomics in zebrafish. This is the first large-scale dataset for crotonylation of non-histone proteins in zebrafish larvae. Next, we compared the Kcr results in zebrafish larvae with those of recent Kcr studies in human cell lines. However, Kcr sites did not overlap significantly because 1) our Kcr results from zebrafish embryos and reference results from human cell lines are of different biological status and 2) the Kcr-antibodies used for enrichment may have been different. Thus, lysine crotonylation may function in a variety of species, from zebrafish to human. To investigate PTM crosstalk, we compared Kcr modifications with previously reported acetylation in zebrafish. Only 67 Kcr proteins (30.7%) and 52 Kcr sites (23.9%) overlapped with Kac in zebrafish. These results indicate that Kcr is non-competitive with Kac and engages with biological pathways and interactions that differ from those of Kac sites in zebrafish. However, when analyzing total 323 Kcr sites converted to the human and human Kac, Kubi dataset from PhosphoSitePlus, 95 of the 323 Kac and 42 Kubi were consistent. Based on these results, there may be a mutual relationship between Kcr and other PTMs.
We found six crotonylated motifs in zebrafish larvae that include specific hydrophobic (leucine, valine and isoleucine) or acidic (glutamate and aspartate) amino acid residues flanking the modified lysine residue. The motif results suggest an interaction with specific enzymes, such as kinases and acetylases. For example, acetylation of the KxGS motif can regulate tau assembly by HDAC6 and these sites are hypoacetylated in patients with Alzheimer's disease 31 . Basophilic motifs such as RxxS are phosphorylated by protein kinase A and proline-directed motifs such as PxS and PxTPP are target sequences for mitogen-activated protein kinase 32 . Therefore, crotonylation-related proteins such as crotonylases and decrotonylases likely prefer hydrophobic or acidic flanking sequences. Because Kcr-X-E was also confirmed in Kac, Kcr-X-E may recognize an enzyme similar to Kac 33 . A recent study of Kcr on non-histone proteins showed that many acidic amino acid motifs such as Kcr-E, E-Kcr and Kcr-D were detected in H1299 and HeLa cells 14 . Furthermore, we identified only 5 crotonylated sites in histones: H2A K241, H2B K6, H3 K123, H4 K60 and H4 K78. Our bioinformatics results for subcellular localization also suggest that 89% of proteins were crotonylated in the cytosol, mitochondria and extracellular matrix, but not in the nucleus. Therefore, our results indicate that Kcr is distributed in various subcellular locations. To further investigate the biological regulatory effects of Kcr in zebrafish embryos, we carried out GO and KEGG pathway analysis. The results suggest that Kcr serves as a diverse regulatory factor in cellular and metabolic processes.
Additionally, Kcr sites and proteins were evolutionarily conserved between humans and zebrafish. A total of 97.7% of Kcr proteins and 76.8% of Kcr sites in zebrafish significantly overlapped with humans. Our previous study showed that 69% of zebrafish phosphoproteins were conserved in humans 9 . In addition, 51.7% of zebrafish Kac sites overlapped with humans and 34.5% of Kac sites were identified as human Kac sites 26 . Particularly, our dataset revealed that crotonylation of ribosomal proteins and myofilament proteins was highly enriched and evolutionarily conserved. Thus, we focused on myofilament proteins and ribosomal proteins for Kcr.
Numerous studies have examined the correlation between myofilament proteins and PTMs, such as short-term phosphorylation at multiple sites in myosin light chain (MLC), troponin, TM and myosin binding protein-C, which is associated with modulation of contractility 34 . Previous studies on PTMs with myofilaments showed that increased phosphorylation of MLC2 is well-known to increase Ca 2+ sensitivity 35 . Moreover, Tyr nitration and Cys S-nitrosylation of MLC1 is induced by oxidative stress or hypoxia-reoxygenation. As a result, nitrated and S-nitrosylated MLC1 may be prone to degradation by matrix metalloprotease-2 19 . Phosphorylation at Thr64 and Ser194 or 195 of human MLC1 is closely related to the stability of the myosin head 36 . Recently, Meishan et al. studied the relationship between myosin and PTMs in old age and found that modifications to myosin heavy chain type I and II (MYH1 and 2) in old age are associated with significant slowing of motility speed. They detected eight age-specific myosin PTMs: carbonylation of Pro79, Asn81, Asp900, Asp904 and Arg908; methylation of Glu1166; and deamidation of Gln1164 and Asn1168. Thus, these PTMs may be involved in disordered myosin organization and the slowing of motility 37 .
Other studies on TM and troponin showed that N-terminal acetylation of TM increased protein stability and strongly enhanced affinity to actin 38 . Acetylation enhances TM function, thereby regulating myosin activity 39 . In addition, TM can be phosphorylated by phosphoinositide 3-kinase, which activates myosin Mg 2+ ATPase and remodeling of the actin cytoskeleton [40][41][42] . Phosphorylation of troponin I on Ser23 and 24 by protein kinase A has been shown to reduce myofilament Ca 2+ sensitivity and is associated with heart failure 20 . Troponin T can be phosphorylated by several kinases such as protein kinase C, Ca 2+ /calmodulin-dependent protein kinase II and apoptosis signal-regulating kinase 1 43,44 . Consequently, phosphorylated troponin T at Ser209, 285 and Thr213, 294 by protein kinase C-α reduces tension, ATPase activity and Ca 2+ sensitivity 45 .
Thus, Kcr of myofilament proteins may play an important role in regulating Ca 2+ sensitivity, remodeling the actin cytoskeleton and modulating contractility. Furthermore, Kcr is significantly associated with heart failure, myocardial infarction and aging. Notably, Kcr of myosin is concentrated in myosin coiled coils, which are expected to contribute to intracellular transport.
Finally, we identified 34 crotonylated sites on 21 ribosomal proteins that are highly conserved between zebrafish and humans. Ribosomal proteins are among the major sources for protein synthesis and are responsible for translation. Since the 1970s, researchers have predicted that PTMs (such as acetylation) of ribosomal proteins are important for biological functions 46,47 . Recently, some studies revealed that the large ribosomal subunit L28 is substantially ubiquitinated during S phase in yeast and shows active ribosomal function during translation without targeting the protein for degradation 48 . In addition, protein N-terminal acetylation of ribosomal proteins by N-acetyltransferase is necessary to maintain protein synthesis in yeast 49 . Therefore, our data also indicate that Kcr of ribosomal proteins is important for the regulation of protein synthesis and ribosome assembly.
In conclusion, we determined the first large-scale crotonylome of zebrafish embryos. These crotonylated proteins and sites are widely distributed in non-histone proteins. Notably, our study revealed that Kcr is evolutionarily conserved between zebrafish and humans and is particularly enriched in ribosomal proteins and myofilament  proteins such as myosin, TM and troponin. Therefore, our results provide a foundation for future studies of the effects of crotonylation on aging and heart failure.

Danio rerio cultures. Danio rerio (wild-type) were acquired from the Korean Zebrafish Organogenesis
Mutant Bank (Daegu, South Korea). Zebrafish were maintained in a 14-h light/10-h dark cycle at 28.5 °C, with a recirculating filtration system using mechanical and biological filtration and fed with baby brine shrimp (Advanced Hatchery Technology, Inc., Salt Lake City, UT, USA) twice daily. The eggs were obtained by pair mating and kept at 28.5 °C in embryo media (0.3 mg/mL sea salt and 1 µg/mL methylene blue). Developmental stages are described as hours post fertilization (hpf) based on morphological features in standard embryogenesis 50 .
Protein extraction. The eggs were grown until the early larval period after hatching (72-120 hpf) in embryo media (0.3 mg/mL of sea salt and 1 µg/mL of methylene blue) and individually collected 3 times. Manually dechorionated embryos were deyolked by pipetting as previously described and pooled 51 . Deyolked embryos were mixed in lysis buffer containing complete RIPA buffer, protease inhibitor cocktails and histone deacetylase inhibitors and the mixtures were sonicated on ice. The supernatants were separated after centrifugation at 14,000 × g for 10 min at 4 °C. For protein purification, embryonic proteins were precipitated in 10% trichloroacetic acid overnight at 4 °C and then centrifuged at 12,000 × g for 7 min at 4 °C. Precipitated pellets were washed with −20 °C acetone twice and then dissolved in 50 mM ammonium bicarbonate buffer. Re-suspended proteins were quantified using BCA assay kits.
In-solution tryptic digestion. Before protein digestion, 10 mM DTT was added to reduce the protein lysates UniProt D. rerio database (41,001 sequences) concatenated with a reverse decoy database. Protein and peptides were acquired using the following parameters: trypsin/P for cleavage enzyme allowing up to 4 missing cleavages; 10 ppm for precursor ions and 0.02 Da for fragment ions of mass error; carbamidomethylation on Cys for fixed modification and oxidation on Met, crotonylation on lysine and acetylation on the protein N-terminus for variable modifications. False discovery rate for protein, peptide and Kcr site were specified at 1%. The minimum peptide length was set to 7. For selected specific Kcr sites, site localization probability was set to >0.75. All other parameters in MaxQuant were used as default.
Bioinformatics analysis for Gene Ontology annotation. Gene Ontology (GO) is a major bioinformatics initiative to unify the representation of gene and gene product attributes across all species. The GO annotation proteome was determined using the UniProt-GOA database (http://www.ebi.ac.uk/GOA/). Identified protein IDs were converted to UniProt ID and then mapped to GO IDs. If proteins were not annotated in the UniProt-GOA database, InterProScan software was used to annotate the protein's GO function based on protein sequence alignment. Next, identified proteins were categorized using GO annotation based on three classification: biological process, cellular component and molecular function.
For subcellular localization, we used Wolfpsort, a subcellular localization predication software that predicts subcellular localization. Wolfpsort is an updated version of PSORT/PSORT II for predicting eukaryotic sequences.
To investigate the KEGG pathway, identified proteins annotated by the KEGG database. First, we used the KEGG online service tool KAAS to annotate the protein's KEGG database. Next, we mapped the annotation results on the KEGG pathway database using the KEGG online service tool KEGG mapper.
Bioinformatics analysis for enrichment of GO and KEGG pathway analysis. For three GO annotation categories, biological process, cellular component and molecular function, we used the Functional Annotation Tool of DAVID Bioinformatics Resources 6.7 to identify GO enrichments against the background of zebrafish. Additionally, to identify enriched pathways, the KEGG database was used with the Functional Annotation Tool of DAVID against the background of zebrafish. To test the enrichment of protein-containing UniProt entries against all UniProt proteins, we used a two-tailed Fisher's exact test. Corrections for multiple hypothesis testing were performed using standard false discovery rate control methods. GO terms with a Scientific REPORts | (2018) 8:3652 | DOI:10.1038/s41598-018-22069-3 corrected p-value less than 0.05 were regarded as significant. Identified pathways were classified into hierarchical categories according to the KEGG website.
Motif and homologous analysis. Motif-X software was used to analyze the model of sequences with amino acids in specific positions of modifier-15-mers (7 amino acids upstream and downstream of the site) in all protein sequences. All database protein sequences were used as background database parameters and other parameters were used as default.
To analyze the conservation of Kcr, homologous proteins and sites between zebrafish and humans were examined using BLASTP 28 . The detailed procedure for examining the conservation of proteins and modification sites was previously described 52 . BLASTP parameters for humans were obtained from the UniprotKB database and p-value < 0.001 was considered as high conservation. To analyze potential cross-talk among Kcr, Kac and lysine ubiquitination, Kcro results converted to human were compared using database sets downloaded from PhosphoSitePlus 29 .
Immunoblot analysis. For SDS polyacrylamide gel electrophoresis, 20 μg of protein from each developmental phase was loaded in 10% gels. Separated gels were transferred onto polyvinylidene difluoride membranes on wetting blot systems and blocked with 5% bovine serum albumin with TBST buffer (20 mM Tris, 137 mM NaCl and 0.5% Tween-20 pH 7.4) for 5 h at RT. Membranes were incubated with anti-Kcr primary antibodies (PTM Biolabs, #PTM-501, 1:1000) overnight at 4 °C. After washing the membranes with TBST five times, they were incubated with anti-rabbit IgG horseradish peroxidase-linked secondary antibody (#7074, 1:2000; Cell Signaling Technology, Danvers, MA, USA) for 2 h at RT. Finally, the membranes were washed with TBST three times and target proteins were detected using ECL reagent with Image Quant LAS-4000 mini (GE Healthcare, Little Chalfont, UK).