Lysine Acetylation is an Important Post-Translational Modification that Modulates Heat Shock Response in the Sea Cucumber Apostichopus japonicus

Heat stress (HS) is an important factor for the survival of the marine organism Apostichopus japonicus. Lysine acetylation is a pivotal post-translational modification that modulates diverse physiological processes including heat shock response (HSR). In this study, 4028 lysine acetylation sites in 1439 proteins were identified in A. japonicus by acetylproteome sequencing. A total of 13 motifs were characterized around the acetylated lysine sites. Gene Ontology analysis showed that major acetylated protein groups were involved in “oxidation–reduction process”, “ribosome”, and “protein binding” terms. Compared to the control group, the acetylation quantitation of 25 and 41 lysine sites changed after 6 and 48 h HS. Notably, lysine acetyltransferase CREB-binding protein (CBP) was identified to have differential acetylation quantitation at multiple lysine sites under HS. Various chaperones, such as caseinolytic peptidase B protein homolog (CLBP), T-complex protein 1 (TCP1), and cyclophilin A (CYP1), showed differential acetylation quantitation after 48 h HS. Additionally, many translation-associated proteins, such as ribosomal proteins, translation initiation factor (IF), and elongation factors (EFs), had differential acetylation quantitation under HS. These proteins represented specific interaction networks. Collectively, our results offer novel insight into the complex HSR in A. japonicus and provide a resource for further mechanistic studies examining the regulation of protein function by lysine acetylation.


Introduction
Post-translation modifications (PTMs) are as important as the transcriptional ways in regulating protein activity [1]. PTMs can regulate allosteric effects of proteins and can affect binding sites and domains for protein-protein interactions [2]. PTMs consist of a variety of protein modifications, including phosphorylation, acetylation, ubiquitylation, methylation, and hydroxylation. Lysine acetylation, as a critical type of PTM, was first discovered in histone proteins [3]. Histone proteins interact with DNA and are the fundamental building blocks of eukaryotic chromatin. Lysine acetylation may direct the assembly of histones and help regulate the transcription of genes [4,5]. Outside of histones, there has been a rapid proliferation in the attention of nonhistone targets of lysine acetylation. For example, lysine acetylation is already extensively investigated in a wide range of transcription factors, molecular chaperones, cytoskeletal proteins, and nuclear import factors [6,7]. The acetylation of nonhistone proteins can regulate enzymatic activity, protein interactions, protein stability, and metabolic pathways [6,8].
Living organisms are frequently faced with abiotic stress conditions such as extreme temperatures, low nutrient availability, salinity, drought, and high ultraviolet irradiation. The physiology, behavior,

A. japonicus Had a Large Number of Acetylated Proteins and Sites by Proteome-Wide Analysis
To map lysine acetylation sites in A. japonicus, proteins were extracted and digested into peptides with trypsin as described in Materials and Methods ( Figure 1A). The length of most acetylation peptides obtained was between 6 and 30 amino acids, and most peptides consisted of 8-16 amino acids ( Figure 1B). Altogether, 4028 lysine acetylation sites in 1439 proteins were identified. We also analyzed the distribution of these acetylated sites. We identified 585 proteins that had one acetylated site, 202 proteins having two acetylated sites, and 106 proteins possessing three sites ( Figure 1C).

Acetylation Motifs Were Characterized in A. japonicus
We identified acetylation motifs using motif-X software (Boston, MA, USA). A substantial bias in the amino acid distribution was observed from positions −7 to +7 around the acetylated lysine residues. A total of 13 motifs were characterized in A. japonicus (Figure 2A). Among them, K ac XXK was the most common motif, followed by K ac L, KXXXXXXK ac , K ac F, and K ac Y ("K ac " represents the acetylated lysine, and "X" represents a random amino acid residue) ( Figure 2B). The sequence logos showed that small hydrophobic residues were frequently around K ac , including leucine (L), valine (V), alanine (A), and glycine (G) ( Figure 2C). Besides, lysine appears at a very high probability around K ac , especially from +3 to +7 and from −7 to −4 positions.

13 Acetylation Motifs Were Characterized in A. japonicus
We identified acetylation motifs using motif-X software (Boston, MA, USA). A substantial bias in the amino acid distribution was observed from positions −7 to +7 around the acetylated lysine residues. A total of 13 motifs were characterized in A. japonicus (Figure 2A). Among them, K ac XXK was the most common motif, followed by K ac L, KXXXXXXK ac , K ac F, and K ac Y ("K ac " represents the acetylated lysine, and "X" represents a random amino acid residue) ( Figure 2B). The sequence logos showed that small hydrophobic residues were frequently around K ac , including leucine (L), valine (V), alanine (A), and glycine (G) ( Figure 2C). Besides, lysine appears at a very high probability around K ac , especially from +3 to +7 and from −7 to −4 positions.

Functional Annotation of Acetylated Proteins
The Gene Ontology (GO) analysis showed that a total of 951 acetylated proteins were annotated to GO terms. Three ontologies, including molecular function (MF), cellular component (CC), and biological process (BP), were further analyzed ( Figure 3A). In the BP category, the major subcategories were "oxidation-reduction process" (140, 14.7%) and "metabolic process" (94, 9.9%). In the CC subcategory, a significant proportion of the clusters were classified as "ribosome" (51, 5.4%) and "intracellular" (40,4.2%). In the MF subcategory, the largest subcategory was "protein binding" (164, 17.2%), followed by "ATP binding" (84, 8.8%) and oxidoreductase activity (74, 7.8%). Clusters of Orthologous Groups (COG) annotation was also applied for the identified proteins. The results showed that the most abundant process was "translation, ribosomal structure, and biogenesis", which was enriched with 120 proteins ( Figure 3B). Both "lipid transport and metabolism" (with 108 proteins) and "post-translational modification, protein turnover, chaperones" with (77 proteins) were also of great abundance. To better describe the identified proteins, we also searched them in the Kyoto Encyclopedia of Genes and Genomes (KEGG) and the Interpro (IPR) database. The detailed information is listed in Table S1.

The Number of Acetylation Sites in Many Functional Proteins Changed Significantly under Heat Stress (HS)
By the quantification of the acetylated sites, a total of 25, 41, and 41 differential quantified acetylated sites were found in these three comparisons respectively (Figure 4). The results also showed that there were some overlapping sites with differential acetylated quantitation among the three comparisons. Notable, 3 acetylated sites (K178 of translation initiation factor 2 subunit 1, K29 of N-terminal kinase-like protein, and K120 of spectrin beta) had different amounts in both HS6h versus C and HS48h versus C comparisons.

HS6h versus C Comparison
A total of 25 differential modified acetylated sites were identified in the HS6h versus C comparison. The results showed that 20 acetylated sites were significantly up-regulated, while only 5 acetylated sites were significantly down-regulated. Three proteins with acetyltransferase activity, including histone acetyltransferase KAT5, KAT8 regulatory NSL complex subunit 1 (KANSL1), and CREB-binding protein (CBP), showed increased acetylation quantitation in some lysine sites in the 6 h group (Table 1). Two ribosomal proteins displayed up-regulated acetylation quantitation, including large subunit ribosomal protein L3e (Ribosomal_L3e) and small subunit ribosomal protein S15e (Ribosomal_S15e). Moreover, the acetylation quantitation of some translation-related proteins, such as translation initiation factor 2 (IF2), elongation factor 1α (EF1α), and elongation factor 2 (EF2), were up-regulated in the HS6h group compared the control group. Regarding the metabolic genes, the acetylation quantitation of K387 site in serine hydroxymethyltransferase (SHMT) decreased in the 6

HS6h versus C Comparison
A total of 25 differential modified acetylated sites were identified in the HS6h versus C comparison. The results showed that 20 acetylated sites were significantly up-regulated, while only 5 acetylated sites were significantly down-regulated. Three proteins with acetyltransferase activity, including histone acetyltransferase KAT5, KAT8 regulatory NSL complex subunit 1 (KANSL1), and CREB-binding protein (CBP), showed increased acetylation quantitation in some lysine sites in the 6 h group (Table 1). Two ribosomal proteins displayed up-regulated acetylation quantitation, including large subunit ribosomal protein L3e (Ribosomal_L3e) and small subunit ribosomal protein S15e (Ribosomal_S15e). Moreover, the acetylation quantitation of some translation-related proteins, such as translation initiation factor 2 (IF2), elongation factor 1α (EF1α), and elongation factor 2 (EF2), were up-regulated in the HS6h group compared the control group. Regarding the metabolic genes, the acetylation quantitation of K387 site in serine hydroxymethyltransferase (SHMT) decreased in the 6 h group. The acetylation amount of cholinesterase 1 (CHLE) was also down-regulated in the 6 h group.

HS48h versus C Comparison
The changes in acetylation amount in HS48h was much more obvious, with 41 differential quantified acetylated sites compared with the control group (Table 2). Of these sites, 20 sites increased, while the acetylation amount of the other 21 sites declined.  Three lysine sites of CBP were involved in the changes: both K1583 and K1590 sites showed up-regulated acetylation quantitation, while the acetylation amount of K1766 was down-regulated. Besides, the K60 site of chromatin modification-related protein EAF6 exhibited a decreased acetylation amount.
Remarkably, various lysine sites of some chaperones were also involved in the different acetylation amounts. A total of 9 lysine sites (K147, K225, K316, K395, K440, K469, K476, K590, and K596) were identified in heat shock protein 90 (HSP90). Of these sites, K440 of HSP90 showed an increased acetylation amount in the HS48h group. T-complex protein 1 (TCP1), belonging to distant homologues of the HSP60 family, also showed an up-regulated acetylation amount in K206. On the contrary, the acetylation amounts of caseinolytic peptidase B protein (CLPB) and cyclophilin A (CYP1) were down-regulated in the HS48 group. The acetylation quantitation of the K87 site in glutathione S-transferase (GST) was up-regulated.
Ribosomal proteins showed different acetylation amounts in the HS48 group, but with opposite regulation of the HS6h group. Three lysine sites, respectively, in 60S ribosomal protein L23a (Ribosomal_L23a), large subunit ribosomal protein L3e (Ribosomal_L3e), and ribosomal protein S5 (Ribosomal_S5) showed declined acetylation amounts in the HS48 group. Besides, the acetylation amount of EF1α still increased, while the acetylation quantitation of EF2 decreased in the HS48 group.
Solute carrier family 39 member 10 (SLC39A10), also known as zinc transporter ZIP10, was found to show an increasing acetylation amount at the K43 site. Similarly, sodium-dependent organic anion transporter (NPT) had a higher acetylation amount at K114 in the HS48h group.

HS48h versus HS6h Comparison
There were 41 lysine sites (21 up-regulated and 20 down-regulated) that had different acetylation amounts in the HS48h versus HS6h comparison. The results in this comparison shared some acetylated sites with those in the above two comparisons. For example, HSP60, FANCJ, EF1α, and some ribosome proteins were still included (Table 3). It is worth noting that various proteins involved in ion transport regulated the acetylation amount significantly in this comparison. For example, diverse proteins in solute carrier families (SLCs), such as NPT, SLC39A10, SLC25 and SLC28, increased the acetylation quantitation, while the acetylation amount of proteins in the transferrin family declined.

Differentially Acetylated Proteins Represented Specific Interaction Networks
For the purpose of understanding how these acetylated proteins interacted with each other, we used the STRING database (https://string-db.org) and Cytoscape software to construct the protein-protein interaction (PPI) networks for the proteins with differentially acetylated sites. The results showed that 37 proteins were mapped to the protein interaction database. These proteins were clustered into one big cluster and 5 protein pairs (Figure 4).

Lysine Acetylation Profiles in A. japonicus
In this study, we used global mass spectrometry based acetylproteome profiling to identify 4028 lysine acetylation sites arising from 1439 proteins in the intestine tissue of A. japonicus. A total of 13 motifs were characterized. K ac XXK was the most common motif, followed by K ac L, KXXXXXXK ac , K ac F, and K ac Y. Of these motifs, K ac XXK, K ac F, and K ac Y were reported in mammals and bacteria, which indicated extensively conserved sequences of these motifs [22][23][24]. Sequence logos showed that lysine was over-represented immediately around acetylated lysine in A. japonicus. The preference was also observed in many other species [25,26]. These motifs may be specially recognized by diverse lysine acetyltransferases (KATs). Further studies are needed to identify the corresponding KATs that recognize these acetylation motifs.
By GO functional analysis of the acetylated proteins, we found "protein binding" was one of the most abundant terms in the MF category in A. japonicus, which was consistent with the findings in bacteria [23], plants [25], and mammals [24]. The terms associated with "ribosome", "nucleus", "cytoplasm", and "membrane" were enriched in the CC category, which suggested the wide distribution of lysine acetylated proteins [25]. The acetylated proteins we identified also represented a wide range of BP terms. The result of COG analysis also implied that ribosomal proteins exhibited wide acetylation in A. japonicus, as reported in other species [25,27]. Besides, acetylation targeted proteins involved in "post-translational modification, protein turnover, chaperones" and "lipid transport and metabolism". These proteins had broad, important regulatory roles in HSR [28]. Overall, the functional analysis implied that lysine acetylation existed in proteins of A. japonicus widely, including the stress response proteins.

Differential Number of Acetylated Sites under HS
We studied the expression of protein in our previous proteomic study [28]. In that study, a total of 3432 proteins were identified, and 127 proteins showed significant difference in the C versus HS48h comparison. In the present study, 4028 lysine acetylation sites in 1439 proteins were identified by acetylproteome sequencing, and 41 acetylation sites had differential amounts in the C versus HS48h comparison. We summarized these data in the two studies, and found that only five proteins were significantly changed in both protein and acetylation amounts. The five proteins included heat shock protein 90 (HSP90), ELL-associated factor 2-like (EAF2), extracellular signal-regulated kinase 2 (ERK2), gamma-butyrobetaine dioxygenase (BBOX), and phosphodiesterase family member 7-like (ENPP7), and the changes in protein and acetylation amounts are listed in Table S2. Thus, we suppose that the quantitative changes of most acetylated sites resulted from protein modification rather than from protein abundance. The differential amounts of acetylated sites generally reflected the acetylation levels of these sites [24,27].

Lysine Acetyltransferases (KATs) and Deacetylases (KDACs)
KATs and KDACs are often referred to as histone acetyltransferases (HATs) and histone deacetylases (HDACs), respectively, because lysine acetylation was initially identified in histones [7]. Genome-wide sequence analyses in A. japonicus have indicated the presence of diverse KATs, such as MYST proteins and CREB-binding protein (CBP), as well as diverse kinds of KDACs [29].
Many KATs and KDACs themselves are autoacetylated [7]. Prominent examples are CBP and its homologue p300. CBP and p300 proteins are transcriptional co-activators of various transcription factors that are involved in a wide array of cellular activities, such as DNA repair, cell growth, differentiation, and apoptosis [30][31][32]. A reduced catalytic activity of p300 was achieved by autoacetylation of a cluster of acetylated sites within an apparent activation loop motif [33]. In the present study, we found that CBP (NCBI sequence ID: PIK57779.1) possessed conserved modular domains with the protein in other species. Notably, four lysine sites with close locations (K1583, K1590, K1766, and K1804) had different acetylation amounts under HS. These four lysine sites were around the zinc-binding sites as well as hydrophobic binding surface ( Figure S1). The three-dimensional domain structures of partial CBP predicted by the SWISS-MODEL server also proved that these lysine residues were very close to zinc ion binding sites, especially K1804. We, therefore, proposed that these acetylated lysine sites could alter the activities and functions of CBP and then influence a wide variety of signaling and transcriptional events under HS [33].
HSP90 is critical for the maturation of many crucial transcription factors, and acetylation regulation of HSP90 has drawn the most attention among chaperons. HSP90 interacts with transcription factor glucocorticoid receptor (GR), which helps GR to achieve a conformation suitable for binding its ligand and subsequent translocation into the nucleus. Acetylation of HSP90 prevents this interaction, ultimately leading to a loss of transcriptional activity of GR [35]. HDAC6, which binds and deacetylates HSP90, has been identified as a regulator of HSP90 acetylation and its chaperone activity [36][37][38]. A conserved HSP90 residue (K294), which is located in the middle region and C-terminal from the charged linker, is regarded as a key regulation site for HSP90 [7]. In yeast, K294 acetylation weakens the association with co-chaperones and client proteins [39]. In A. japonicus HSP90, K327 is equivalent to K294 in the human protein. Surprisingly, K327 was not an acetylation site in our data. In our study, we identified a total of nine lysine sites (K147, K225, K316, K395, K440, K469, K476, K590, and K596) in HSP90, which proved the presence of multiple acetylation sites in HSP90. Under HS48h, only the acetylation amount of K440 was up-regulated, to 1.94-fold. From our previous result, the protein expression of HSP90 rose dramatically in the HS48h group, by a factor of 6.11 [28]. Hence, we suppose the acetylation level of HSP90 (K440), which was represented with a ratio of the number of acetylated sites and amount of protein expression, was also reduced under HS in A. japonicus [35,39]. The results also showed that the regulation of HSP90 in A. japonicus under HS occurred at both translation and PTM levels. The specific characteristics and regulation mechanism in the acetylation of HSP90 still needs exploring.
CYP1, a peptidyl-prolyl isomerase, binds the immunosuppressive drug cyclosporin A and modulates T cell activation by inhibiting the phosphatase calcineurin [40]. Disruption of CYP1 decreases survival of cells following exposure to high temperatures, indicating that CYP1 plays a role in the stress response [41]. CYP1 is acetylated at K125, a residue located in the loop surrounding the drug pocket and a potential site of ionic interaction with calcineurin [34]. Acetylation of K125 could affect the interaction of cyclophilin with calcineurin. In our present study, CYP1 had decreasing acetylation amounts at both K47 and K80 sites in the HS48h group compared with the control group. No significant changes in the protein expression of CYP1 were examined in HSR of A. japonicus. We speculated that the deacetylation of CYP1 contributed to its isomerase activity under HS.

Translation-Associated Factors and Ribosome Proteins
Our results showed that the acetylation amounts of three translation-related factors (IF2, EF1, and EF2) increased significantly under HS in A. japonicus. Previous reports proved that lysine acetylation occurred in diverse EFs [42,43]. However, rare evidence is obtained about the impact of acetylation on the activity of these translation factors. Besides, some ribosome proteins had differential acetylation amounts in this study. Generally, the acetylation amounts of Ribosomal_L3e and Ribosomal_S15e increased at HS6h, while those of Ribosomal_L23e and Ribosomal_S5 declined at HS48h. Modifications of ribosomal proteins are important for protein synthesis [44]. For example, phosphorylation of ribosomal protein S6 has key regulatory roles in cell size and glucose homeostasis [45]. The acetylation level of mitochondrial ribosomal protein L10 influences mitochondrial protein synthesis [46]. We supposed that the dynamic acetylation changes of translation-related factors and ribosome proteins acted as the complex regulation mode on the protein synthesis under HS.

Samples
A. japonicus (100-120 g) was collected from an aquaculture farm in Weihai (Shandong, China) in April 2018. The sea cucumbers were transported to the tanks in our lab of Qingdao. The temperature and the salinity of sea water were around 16 • C and 30% , respectively, before the treatment experiment started. Acclimation in tanks lasted for 2 weeks. During both the acclimation and treatment experiment, the sea cucumbers were fed with regular fodder, and half of the water was changed daily [47].
When acclimation finished, three individuals taken from the tank were regarded as the control group. The methods of heating and sampling were followed by our previous procedures. Briefly, a heating rod with 2 kW of power was placed in the tank. The rate of heating was about 2 • C/h. The heating temperature was set to 26 • C, and then water temperature was maintained at 26 ± 0.5 • C during the subsequent experiment. The initial time was considered when the temperature just reached 26 • C. Three individuals were taken randomly after 6 and 48 h stimulations. The intestine tissue of these individuals was tweezed and washed with ddH 2 O, which helped us clear the intestine contents as much as possible. All samples were frozen quickly in liquid nitrogen and then stored at -80 • C.

Protein Extraction
The sample was ground with pestles in liquid nitrogen, respectively, and the sample powder was mixed with extraction buffer containing 8 M urea, 1% Triton-100, 65 mM dithiothreitol (DTT), 0.1% protease inhibitor cocktail, 3 µM trichostatin A, and 50 mM nicotinamide. Then, the sample was sonicated three times on ice using a high-intensity ultrasonic processor. After centrifugation at 20,000× g at 4 • C for 10 min, the protein in the supernatant was precipitated with 15% cold trichloroacetic acid (TCA) at −20 • C for 2 h. The protein was washed with cold acetone three times. Finally, the protein was re-dissolved in buffer (8 M urea, 100 mM NH 4 HCO 3 , pH 8.0), and the protein concentration was determined by a bicinchoninic acid (BCA) quantitative kit (Beyotime, Shanghai, China).

Trypsin Digestion and ENRICHMENT of Lysine-Acetylated Peptides
A total of 10 mg protein was reduced with 10 mM DTT for 1 h at 37 • C and alkylated with 20 mM iodoacetamide for 30 min at room temperature in darkness. For trypsin digestion, the protein sample was diluted in 0.1 M triethylammonium bicarbonate (TEAB). Finally, the diluted protein samples were digested with trypsin at 1:50 trypsin-to-protein mass ratio for 12 h first, and then 1:100 trypsin-to-protein mass ratio for 4 h.
To enrich acetylated peptides, the tryptic peptides were dissolved in NETN buffer (100 mM NaCl, 1 mM EDTA, 50 mM Tris-HCl, 0.5% NP-40, pH 8.0) and incubated with anti-acetyl lysine antibody beads (PTM Biolabs, Hangzhou, China) at 4 • C for 12 h with gentle shaking. The beads were washed four times with NETN buffer and twice with ddH 2 O. The bound peptides were eluted with 0.1% trifluoroacetic acid (TFA). In preparation for analysis, the peptides were desalted using peptide desalting spin columns (ThermoFisher, Waltham, MA, USA).

Liquid Chromatography-Tandem Mass Spectrometry (LC-MS/MS) Detection and Data Analysis
The eluted peptides were cleaned with C18 ZipTips (Millipore, Bedford, MA, USA) according to the manufacturer's instructions. LC-MS/MS analyses were performed according to the previous report [48]. Briefly, tryptic peptides were dissolved in solvent A (0.1% formic acid) and loaded onto a homemade reversed-phase analytical column (15 cm length, 75 µm inside diameter). Solvent B was 0.1% formic acid in 80% acetonitrile. The gradient elutions were carried out at a constant flow rate of 400 nL/min on an EASY-nLC 1000 UPLC system. The peptides were subjected to an NSI source followed by MS/MS analysis, which was performed in EASY-nLC 1200 UHPLC (Thermo Fisher Scientific, Madison, WI, USA) coupled to an Orbitrap Q Exactive HF-X mass spectrometer (Thermo Fisher Scientific, Madison, WI, USA) operating in the data-dependent acquisition (DDA) mode. For data-dependent analysis (DDA), a Q-Exactive HF-X mass spectrometer was operated in positive polarity mode with a spray voltage of 2.3 kV and a capillary temperature of 320 • C. Full MS scans from 350 to 1500 m/z were acquired at a resolution of 60,000 (at 200 m/z) with an automatic gain control (AGC) target value of 3 × 10 6 and a maximum ion injection time of 20 ms. From the full MS scan, a maximum number of 40 of the most abundant precursor ions were selected for higher energy collisional dissociation (HCD) fragment analysis at a resolution of 15,000 (at 200 m/z) with an AGC target value of 1 × 10 5 , a maximum ion injection time of 45 ms, a normalized collision energy of 27%, and an intensity threshold of 8.3 × 10 3 , and the dynamic exclusion parameter was set at 60 s. The MS spectra were collected in the 350-1500 m/z range.
The MS/MS data were processed using Proteome Discoverer 2.2 software (PD2.2, Waltham, MA, USA), and the processed data were searched against the A. japonicus genome database concatenated with a reverse decoy database [29]. The identification parameters were set as follows: precursor mass tolerance of 10 ppm and fragment mass tolerance of 0.02 Da. Trypsin/P was specified as the cleavage enzyme, allowing up to two missed cleavages, five modifications per peptide, and five charges. Carbamidomethylation on Cys was specified as a fixed modification, and oxidation on Met, phosphorylation on Ser, Thr, Tyr, and acetylation on the protein N-terminus were specified as variable modifications [49]. False discovery rate thresholds for peptide, protein, and modification site were specified at 5%. Minimum peptide length was set at seven residues. The site localization probability was set as >0.5. All other parameters in PD 2.2 software were set to default values. We checked the mass error of all identified peptides; the distribution was close to zero, and most were less than 0.02 Da, confirming that the mass accuracy of the MS data was acceptable.
The raw quantification values were calculated and validated using the Minora Feature Detector node in the Processing workflow of PD2.2. A quantification value was the intensity, and the area was detected for a given quantification channel depending on the setting in the precursor ions quantifier in the consensus method. A t-test with a p-value < 0.05 was applied to distinguish the differential quantitation of acetylated sites [26,49,50]. Fold changes (FCs) >1.2 or <0.83 were also regarded as the criteria of differentially acetylated sites [51]. Data are available via ProteomeXchange with the dataset identifier PXD013929.

Bioinformatics Analysis
Software motif-x (Boston, MA, USA) was used to analyze the model of sequences constituted with amino acids in specific positions (7 amino acids upstream and downstream of the acetylation site) in all acetylated proteins identified. Function annotations, such as GO and KEGG analyses, were performed as previously described [46]. COG annotation and IPR annotation were also performed to better classify the functions of proteins [52]. PPI for the identified differentially acetylated proteins were performed using the STRING database and Cytoscape software.
The amino acid sequence of CREB-binding protein (CBP) was analyzed for similarity with the BLAST programs at the National Center for Biotechnology Information (http://www.ncbi.nlm.nih.gov/ BLAST/). The three-dimensional domain structure of partial CBP was predicted by the SWISS-MODEL server (https://swissmodel.expasy.org/interactive).

Conclusions
In summary, this study presents a comprehensive analysis of acetylproteome in the sea cucumber A. japonicus. The profiling of acetylation proteins involved in HSR was identified in the species. Moreover, key proteins, such as lysine acetyltransferases, chaperones, translation-associated factors, and ribosome proteins, showed differential acetylation amounts under HS, suggesting complex post-translational modification was involved in the process. Collectively, our results offer novel insight into HSR in A. japonicus and provide a resource for further mechanistic studies examining the regulation of protein functions by lysine acetylation.