The blue mussels Mytilus galloprovincialis and M. trossulus are competing species with biogeographical ranges set in part by environmental exposure to heat and hyposalinity. The underlying cellular mechanisms influencing interspecific differences in stress tolerance are unknown, but are believed to be under regulation by sirtuins, nicotinamide adenine dinucleotide (NAD+)-dependent deacylases that play a critical role in the cellular stress response. A comparison of the proteomic responses of M. galloprovincialis and M. trossulus to an acute heat shock in the presence and absence of the sirtuin inhibitor suramin (SIRT1, 2 and 5) showed that sirtuins affected molecular chaperones, oxidative stress proteins, metabolic enzymes, cytoskeletal and signaling proteins more in the heat-sensitive M. trossulus than in the heat-tolerant M. galloprovincialis. Interactions between sirtuin inhibition and changes in the abundance of proteins of β-oxidation and oxidative stress in M. trossulus suggest a greater role of sirtuins in shifting metabolism to reduce the production of reactive oxygen species near thermal limits. Furthermore, RNA-binding proteins initiating and inhibiting translation were affected by suramin in M. galloprovincialis and M. trossulus, respectively. Western blot analysis showed that the levels of mitochondrial sirtuin 5 (SIRT5) were generally three times higher and increased with acute heat stress in response to sirtuin inhibition in M. trossulus but not in M. galloprovincialis, suggesting a possible feedback response in the former species and a greater reliance on SIRT5 for its stress response. Our findings suggest that SIRT5 plays an important role in setting interspecific differences in stress tolerance in Mytilus by affecting the stress proteome.
Climate change is altering sea surface and air temperatures and seawater pH, and is predicted to increase precipitation, which will alter seawater salinity (IPCC, 2014). Most vulnerable to changing climate conditions will be species at the edge of their tolerance limits and ranges (Tomanek, 2010, 2012). In blue mussels (genus Mytilus) along the Pacific coast, species range shifts are already being observed in response to a changing climate and according to species-specific stress tolerances (Braby and Somero, 2006a,b; Lockwood et al., 2015). Mytilus occur worldwide from subtropical to cold-temperate latitudes on intertidal rocky substrates. As abundant filter feeders, they are keystone species in coastal ecosystems (Gosling, 1992). The Mediterranean Mytilus galloprovincialis invaded southern California during the last century and replaced the native Mytilus trossulus from its southern range to San Francisco Bay, likely because of temperature increases due to climate change (Sarver and Foltz, 1993). Previous studies and current biogeographic distribution indicate that M. galloprovincialis is the more heat-tolerant species but is unable to withstand fluctuations in seawater salinity, while M. trossulus is more cold tolerant and is able to withstand hypo-saline conditions (Braby and Somero, 2006a; Lockwood et al., 2015; Tomanek et al., 2012). Mytilus physiological responses and adaptations have been characterized for various environmental conditions; however, understanding the underlying cellular mechanisms regulating an organism's ability to adapt to and withstand environmental stress is lacking.
Recent studies have shown widespread occurrence of acetylation and succinylation (or more generally acylation) of internal lysine residues of non-histone proteins involved in a wide range of cellular processes (Choudhary et al., 2014; Kim et al., 2006; Park et al., 2013; Rardin et al., 2013; Zhao et al., 2010). Deacylation of non-histone proteins is catalyzed by sirtuins (SIRT1–7), or class III lysine deacetylases (KDAC), using NAD+ as a cofactor (Bheda et al., 2016). While sirtuins are implicated in the regulation of cellular stress, genomic stability and energy metabolism (Houtkooper et al., 2012), their potential role in coping with environmental stress is largely unexplored.
To date, only histones and heat shock proteins (HSPs) have shown a response to deacylation by sirtuins during heat shock. For example, KDAC inhibitors increased the expression of heat shock genes upon heat treatment during Xenopus development (Ovakim and Heikkila, 2003). Additionally, inhibition of deacylases induces HSP90 acetylation and decreases its activity (Bali et al., 2005). Finally, deacetylation of heat shock factor 1 (HSF1) by SIRT1 extends the time HSF1 is bound to the HSP70 promoter, thereby increasing HSP70 transcription (Westerheide et al., 2009).
Our prior proteomic study found significant protein abundance changes in a number of non-HSPs involved in proteolysis, energy metabolism, scavenging of reactive oxygen species (ROS), cytoskeletal maintenance and cell signaling in Mytilus mussels during heat shock (Tomanek and Zuzow, 2010). These findings suggested that heat shock may induce shifts in metabolic pathways in Mytilus to increase nicotinamide adenine dinucleotide phosphate (NADPH) levels to scavenge ROS (Tomanek, 2015; Tomanek and Zuzow, 2010). In addition, mitochondrial SIRT5 and two electron transport system (ETS) proteins known to interact with sirtuins decreased in abundance simultaneously with heat shock, suggesting a role for sirtuins in setting thermal tolerance limits in Mytilus, possibly by regulating the response to the increased production of ROS during heat stress. Similarly, in the intertidal limpet Cellana toreuma exposed to a heat ramp up to 40°C, gene expression of SIRT1 increased at 30°C, indicating an influence of heat shock on cellular energy homeostasis and a shift towards catabolic metabolism (Han et al., 2013).
two-dimensional gel electrophoresis
D-3 hydroxybutyrate dehydrogenase
cytosolic/mitochondrial malate dehydrogenase
enoyl CoA hydratase
extracellular signal-regulated kinase
electron transport flavoprotein
electron transport system
GDP dissociation inhibitor
generalized linear model
heat shock cognate
heat shock factor
heat shock protein
immobilized pH gradient
nicotinamide adenine dinucleotide
nicotinamide adenine dinucleotide phosphate
Na+/H+ exchange regulatory factor
principal component analysis
protein disulfide isomerase
pentose phosphate pathway
reactive oxygen species
small heat shock protein
Tris-buffered saline with Tween
T-complex protein 1
To assess the role of sirtuin-dependent deacylation during heat shock on a proteomic scale, the current study compared the effects of suramin, an inhibitor of SIRT1, 2 and 5 (Lawson et al., 2010), on the heat shock responses of two marine mussel species of the genus Mytilus with different thermal tolerances. We acclimated M. galloprovincialis and M. trossulus mussels to 13°C under constant submersion, after which we excised gill tissue and exposed it to a range of stressful temperatures for 1 h. Tissues were subsequently placed into 13°C seawater for a 24 h recovery (chosen for its comparability to natural in situ recovery; Hofmann and Somero, 1996), and protein abundance changes were analyzed using 2D gel electrophoresis. Suramin was applied to half of the tissues under each temperature treatment for the duration of the experiment. Based on our previous work, we hypothesized that sirtuin inhibition by suramin would affect changes in the abundance of key cellular stress-response proteins (molecular chaperones and antioxidants), and that the effect would be different between the two Mytilus congeners. Our study is the first of its kind to use a known sirtuin inhibitor to investigate the role sirtuins may play in regulating the cellular stress response in organisms with different stress tolerances. Our findings reveal significant relationships between sirtuins and molecular chaperones, antioxidants, proteins involved in metabolic regulation and cell signaling, which vary in Mytilus congeners with different thermal tolerances. These data suggest a potential key regulatory function for sirtuins in the adaptation to thermal stress in Mytilus mussels.
MATERIALS AND METHODS
Animal collection, maintenance and experimental design
Specimens of M. galloprovincialis Lamarck 1819 were collected subtidally from Santa Barbara, CA, USA (34°24′15″N, 119°41′30″W) and those of M. trossulus Gould 1850 from Newport, OR, USA (44°38′25″N, 124°03′10″W). In a separate study, PCR was used to confirm that each site was occupied by only a single species (i.e. there were no hybrids present; Lockwood et al., 2010). Following acclimation to 13°C and full-strength (32 psu) seawater for 4 weeks, mussel gill tissue was dissected and individually placed in a 15 ml conical tube (constantly aerated) and heated at 6°C h−1 to the target temperatures of 28, 32 or 35°C [henceforward, heated treatments are referred to as HS (heat shock)], with 35°C being close to the thermal limit of Mytilus mussels (Braby and Somero, 2006b). Control mussel gill was left at the control temperature of 13°C and not heated. Exposure temperature duration was 1 h. Afterwards, gill tissues within tubes were brought back to the control temperature and left at 13°C for 24 h (under constant aeration and without a loss of ciliary activity), mimicking a recovery response from heat exposure in the intertidal zone. Half of the samples were kept in seawater with 100 μmol l−1 of suramin (inhibitor of SIRT1, 2 and 5) from the beginning of the experiment (henceforward, suramin treatments are referred to as HS+suramin) (Lawson et al., 2010). At the end of the experiment, gill tissue was immediately frozen using liquid nitrogen (N=6 per species for each treatment). Tissues were subsequently kept at −80°C until homogenization.
Sample preparation and analysis for proteomics
Frozen gill tissue was lysed by homogenization using an ice-cold ground-glass homogenizer and protein conformation was disrupted using a urea (7 mol l−1)-based buffer with 0.5% immobilized pH 4–7 gradient (IPG) buffer (GE Healthcare, Piscataway, NJ, USA). Proteins were subsequently precipitated using 10% trichloroacetic acid in acetone and protein concentration was determined with a 2D Quant kit (GE Healthcare), following the manufacturer's instructions.
Two-dimensional gel electrophoresis (2DGE)
Proteins were separated by isoelectric point using immobilized pH gradient (IPG) strips (pH 4–7, 11 cm; GE Healthcare). For second-dimension SDS-PAGE electrophoresis, IPG strips were placed on top of an 11.8% SDS-polyacrylamide gel, which was run (Criterion Dodeca; BioRad, Hercules, CA, USA) for 55 min. Gels were subsequently stained with colloidal Coomassie Blue (G-250) overnight and destained for 48 h. The resulting gel images were scanned for image analysis of protein abundance changes.
2D-gel image analysis of protein abundance
Digitized images of 2D gels were analyzed using Delta2D (version 3.6; Decodon, Greifswald, Germany) (Berth et al., 2007). Spots were detected on gels by fusing all images from one species into a composite image (Fig. 1). The relative amount of protein in each spot (i.e. spot volume) was quantified by normalizing against total spot volume of all proteins in the gel image.
All proteins detected were excised from gels and digested with trypsin (Promega, Madison, WI, USA). Digested proteins were extracted, concentrated and spotted on an Anchorchip™ target plate (Bruker Daltonics Inc., Billerica, MA, USA). We obtained peptide mass fingerprints (PMFs) using a matrix-assisted laser desorption ionization tandem time-of-flight (MALDI-ToF-ToF) mass spectrometer (Ultraflex II; Bruker Daltonics Inc., Billerica, MA, USA). To identify proteins, we combined PMFs and tandem mass spectra in a search against two databases, using Mascot (v2.2; Matrix Science Inc., Boston, MA, USA): (1) an NCBI-EST library limited to Mytilus that contained approximately 407,000 entries (May 2011), and (2) an NCBI dataset of non-redundant metazoan sequences (1,041,948 sequences; April 2011) (Table S1).
To confirm that changes in protein acetylation occurred during heat shock and differed between species, we conducted 2DGE and western analysis with acetyl-lysine antibodies.
Acetyl-lysine protein visualization
Following 2DGE, we transferred proteins to an 8.5×13.5 cm nitrocellulose membrane at 100 V for 30 min, using a Criterion Blotter (BioRad, Hercules, CA, USA). Transfer buffer (pH 8.8) consisted of 25% methanol, 0.2 mol l−1 glycine and 0.1% SDS. Membranes were blocked by rotating at 30 rpm at 18°C in 5% (w/v) bovine serum albumin (BSA) in Tris-buffered saline with Tween (TBST; 0.88% (w/v) sodium chloride, 0.25% Tris (w/v) and 0.6% Tween-20) for 1 h immediately following transfer. Then, the membrane was rinsed in TBST 4 times for 5 min. Following the rinses, membranes were incubated with a 1:500 solution of anti-acetylated lysine antibody (SPC-157F, Stressmarq Bioscience Inc., Victoria, BC, Canada) in 2.5% TBST with rotation at 30 rpm at 4°C for 17 h. The membranes were rinsed in TBST 4 times for 5 min. A secondary antibody was not necessary as the antibody was conjugated with horseradish peroxidase. Finally, we incubated membranes with ECL Plus (GE Healthcare) for 5 min at low light before they were scanned with a blue fluorescent laser at 430 and 200 μm resolution.
Western analysis of SIRT5
Protein samples (5 μg) were separated by one-dimensional SDS-PAGE with an 11.8% separating and a 4% stacking gel (0.5 mol l−1 Tris-HCl, pH 6.8) for 45 min at 200 V (Criterion, BioRad) and subsequently transferred to nitrocellulose membrane as described above for acetyl-lysine protein visualization. After blocking with 5% BSA in TBST, membranes were incubated for 1 h first with a primary antibody against SIRT5 (BML-SA464, Enzo Life Sciences, Inc., Farmingdale, NY, USA) diluted to 1:2500 with TBST and, subsequently, with a secondary antibody (goat anti-rabbit IgG; ADI-SAB-300, Enzo Life Sciences, Inc.) diluted at 1:5000 with TBST, with four 5 min washes with TBST following each incubation. The internal standard used was a recombinant human SIRT5 (5 ng; BML-SE555, Enzo Life Sciences, Inc.). The protocol followed all other steps and used reagents as described in the protocol for acetyl-lysine protein visualization.
To determine which proteins changed in volume significantly in response to HS+suramin treatment, we used an analysis of variance (two-way permutation ANOVA; P≤0.02) within each species and with HS and suramin treatment as the two main effects. We generated a null distribution for the two-way ANOVA (1000 permutations) to account for the unequal variance and non-normal distribution of some of the response variables. We chose a P-value of 0.02 to limit type I error. Post hoc testing to compare treatments was conducted using Tukey's analysis (P≤0.05), using MiniTab (version 15; State College, PA, USA), to support conclusions about differences in single protein abundance profiles.
To associate proteins with similar abundance patterns across samples, we employed hierarchical clustering with average linking (Delta2D), using a Pearson correlation metric. The resulting clusters are specific to each heat map only and were used to characterize associated or complementary changes in protein abundance within a protein functional category. To further assess the importance of specific proteins in differentiating the proteomes of mussels exposed to different temperatures and affected by the treatment with suramin, we also employed principal component analyses (PCA; Delta2D). Component loadings, which quantify the contribution of each protein in the separation of samples along a given component, were used to compare changes in protein abundance (Tables 1 and 2).
Western analysis of SIRT5
One-dimensional gels stained with a SIRT5-specific antibody were analyzed using TotalLab TL120 (version 2006f; Nonlinear Dynamics, Inc., Durham, NC, USA). The background noise was subtracted (‘rolling ball’ method). Bands were detected, identified and labeled relative to the control band. Corresponding bands in each lane were matched and calibrated based on the concentration of the known standard in the reference lane (average of three reference lanes). Band 2 was identified as the most comparable to the reference band while band 1 was located slightly above and band 3 slightly below the reference lane. Both bands 1 and 3 may represent post-translational modification (PTM) of SIRT5 but were not present in all samples. We focused our quantitative analysis on band 2, which was present in all samples. We used a three-factor generalized linear model (GLM) with species, heat shock and suramin set as fixed effects to determine differences in band 2 concentration (P≤0.05). Comparisons between treatments were conducted with a post hoc Tukey multiple comparisons test (P≤0.05).
We conducted PCA to assess the variation in protein abundance patterns based on statistical significance for HS, sirtuin inhibition by suramin or interaction effects (HS+suramin) using a two-way permutation ANOVA (P≤0.02), and identified protein candidates involved in the treatment response (Fig. 2). PCA of these abundance changes for the first two principal components (PCs) for the interaction effects (PCAs for main effects are not shown) in M. galloprovincialis suggests that PC1 and PC2 represent the effect of HS alone and HS+suramin, respectively (Fig. 2A). A clear separation between suramin treatments was notable only in M. trossulus at 35°C along PC2 (Fig. 2B). These results suggest that the effect of sirtuin inhibition by suramin (implicating SIRT1, 2 and 5) was strongest at the highest HS, and was particularly strong near the upper thermal limits of M. trossulus. When the 35°C HS treatment was removed from the M. trossulus PCA, there was a clear gradient from 13°C to 28°C to 32°C HS moving from left to right across PC1 (Fig. S1). Furthermore in M. trossulus, 32°C HS+suramin was separated along PC2 from all other treatments, indicating a unique interaction that was only observed at 32°C HS (Fig. S1). Thus, the effect of sirtuin inhibition started at 32°C and continued to increase up to 35°C in M. trossulus (Fig. 2B; Fig. S1). PCs that separated the treatments were analyzed further by identifying proteins that contributed the most to the separation based on their component loading values.
For M. galloprovincialis, the first PC separated the effect of increasing HS from left to right along the x-axis and accounted for 16.6% of the variation in protein abundance (Fig. 2A). Proteins representing actin binding (profilin), oxidative stress (peroxiredoxin, DyP-type peroxidase), stress signaling (Cdc42, major vault protein), energy metabolism (UDP-galactose 4-epimerase) and molecular chaperones (chaperonin, HSP22) showed the most positive loadings for PC1 (Table 1). Proteins contributing highly negative loadings mainly consisted of cytoskeletal proteins (α- and β-tubulin, actin) with the addition of one energy metabolism protein (ATP synthase β) (Table 1).
PC2, which accounted for 8.2% of the variation in protein abundance, separated the effect of HS+suramin from that of HS alone along the y-axis but in a manner that was dependent on the temperature of the HS (Fig. 2A). For example, 35°C HS was located below 35°C HS+suramin, while 32°C HS was located above 32°C HS+suramin (Fig. 2A). The observed effect is opposite depending on the HS treatment, detailing the complexity of the HS and sirtuin inhibition interaction occurring in M. galloprovincialis at thermal extremes. Proteins related to iron binding/oxidative stress (ferritin), cytoskeleton (α-tubulin), oxidative stress [Cu-Zn superoxide dismutase (SOD)], mRNA binding [DAZ-associated protein (DAZAP)], proteostasis [heat shock cognate (HSC)70, proteasome α-type 5] showed the most positive loadings for PC2 (Table 1). Proteins related to amino acid metabolism (adenosylhomocysteinase), the pay-off phase of glycolysis (α-enolase), Krebs cycle (citrate synthase and dihydrolipoyl dehydrogenase) and the electron transport system (cytochrome c reductase), actin filament bundling (fascin), stress signaling (major vault protein), cytoskeleton (β-tubulin) and oxidative stress (Cu-Zn SOD) contributed the most negative loadings for PC2 (Table 1).
For M. trossulus, PC1 separated the effect of 35°C HS+suramin from all other treatments along the x-axis and accounted for 23.5% of the variation in protein expression (Fig. 2B). Proteins representing proteostasis (HSP70, prohibitin, proteasome α-type 5), signaling (14-3-3, major vault protein, Rho GDP dissociation inhibitor, calcyphosine) and energy metabolism (ATP synthase) showed the most positive loadings for PC1 and contributed to the separation of 35°C HS+suramin from the remaining treatments (Table 2). Proteins contributing highly negative loadings mainly consisted of molecular chaperones (calreticulin, HSC71), cytoskeletal proteins (β-tubulin, intermediate filament, α-tubulin), amino acid metabolism (ornithine aminotransferase) and oxidative stress (DyP-type peroxidase), suggesting that the specific functions of these proteins generally decreased the most with the addition of suramin during HS (Table 2).
PC2 separated the effect of 35°C HS alone from all other treatments along the y-axis and accounted for 11.1% of the variation in protein abundance (Fig. 2B). Proteins related to signaling (Rab GDP dissociation inhibitor α, guanine nucleotide binding protein β), actin binding (fascin, F-actin capping protein α, gelsolin), energy metabolism [UDP-galactose 4-epimerase, antiquitin (or ALDH7A1, an aldehyde dehydrogenase), NADP-isocitrate dehydrogenase] and the cytoskeleton (actin, RIB43A) showed the most positive loadings for PC2, meaning that they generally decreased with 35°C HS (Table 2). Proteins related to microtubule stabilization (translationally controlled tumor protein), DNA/RNA synthesis (ribonucleoside-diphosphate reductase α), the cytoskeleton (β- and α-tubulin, actin), signaling (GTP-binding protein and 14-3-3) and amino acid metabolism (ornithine aminotransferase, OAT) showed the most negative loadings for PC2 contributing to the separation of 35°C HS from all other treatments (Table 2).
Analysis of acetylated proteins and western blot analysis of SIRT5
2DGE and western blotting with an acetyl-lysine antibody showed that 10% (of 48 spots total) and 27% (of 52 spots) of acetylated proteins were affected during HS (up to 35°C) by suramin in M. galloprovincialis and M. trossulus, respectively (two-way ANOVA, P≤0.05; Fig. 3). While both detection methods showed that ∼2.7 times as many spots were affected by suramin during HS in M. trossulus compared with M. galloprovincialis, the number of acetylated spots that changed using the acetyl-lysine antibody was much lower than the number detected with the total protein stain.
The most abundant SIRT5 band (band 2) was ∼3 times more abundant in M. trossulus than in M. galloprovincialis (species main effect; GLM, χ2=50.4, d.f.=1, P≤0.0001; Fig. 4A). In M. trossulus without suramin, the SIRT5 band abundance was significantly greater at 35°C than at either 28°C (Tukey post hoc comparison; P=0.002) or 32°C (Tukey post hoc comparison; P=0.0005) (Fig. 4B). The interaction of HS and suramin was statistically significant only in M. trossulus (GLM, χ2=7.7, d.f.=3, P=0.05), leading to higher levels of SIRT5 during inhibition at 28 and 32°C than in the non-inhibited control groups. However, pairwise comparisons showed that the SIRT5 abundance was greater at 35°C HS+suramin than at 32°C HS without suramin only, in M. trossulus (Tukey post hoc comparison; P=0.03) (Fig. 4B).
Protein chaperoning and degradation (proteostasis)
A comparison of identified proteins by functional categories showed that changes in the abundance of molecular chaperones in response to HS were more affected by sirtuin inhibition (henceforth referred to as an interaction effect) in M. trossulus than in M. galloprovincialis (Figs 5 and 6). In general, molecular chaperones function to maintain protein conformation while proteins related to the proteasome function in degradation (Ciechanover, 2005; Glickman and Ciechanover, 2002; Kim et al., 2013). In M. galloprovincialis, we identified two HSP/HSC70s, one isoform of the chaperonin T-complex protein 1 (TCP1), which catalyzes the folding of cytoskeletal proteins (Sternlicht et al., 1993), one endoplasmic reticulum (ER) chaperone (malectin) and one proteasome subunit, all of which were significant for interaction effects (Figs 5A and 6A). Three other chaperones, a chaperonin (spot 33) isoform, and small heat shock protein (sHSP)21/22 functioning specifically in cytoskeletal stability (Benndorf et al., 2001; Haslbeck et al., 2005; Sternlicht et al., 1993), showed only a temperature effect (Figs 5A and 6A). However, one small HSP (HSP22), two mitochondrial chaperones (prohibitin and HSP60), three constitutively expressed HSC70s, two ER chaperones (protein disulfide isomerase PDI), and three proteasome isoforms were identified but did not change significantly (Fig. 5A). The protein abundance changes of molecular chaperones and proteasome isoforms may indicate partial disruption of the gill tissue cytoskeleton upon HS to 35°C, which is consistent with a previous proteomic study (Tomanek and Zuzow, 2010), but given the non-significance of the changes, it seems likely that the disruption of proteostasis was only partial and did not expand to the mitochondria in M. galloprovincialis.
Of the proteins in M. trossulus showing interaction effects, 11 were molecular chaperones (HSP/HSC70s, two isoforms of TCP1, three sHSPs, HSP40, HSP90, the ER PDI and calreticulin; Figs 5B and 6B). Changes in only one HSP70 and two proteasome isoforms were non-significant (Figs 5B and 6B). While the majority of these chaperones are involved in stabilizing denaturing proteins and facilitating the proper formation of protein aggregates (Kim et al., 2013; Tyedmers et al., 2010), PDI and calreticulin play a role in the maturation of secreted ER proteins (Araki and Nagata, 2012). Several chaperones showed a decrease in one isoform at 35°C HS+suramin (cluster I) while other chaperones simultaneously showed an increase under the same conditions (cluster III) in M. trossulus but not in M. galloprovincialis (Fig. 5). Specifically, two HSC70s decreased while two HSP70s increased at 35°C HS+suramin in comparison to 35°C HS alone, and they were significant for interaction effects (Figs 5B and 6B). In comparison, only one HSC70 decreased at 28°C with sirtuin inhibition in M. galloprovincialis, a much lower temperature. In general, sirtuin inhibition (indicated through an interaction) affected HSP70 and HSC70 equally in both species (one HSC70 out of two and two HSC70 out of four HSC/P70 isoforms in M. galloprovincialis and M. trossulus, respectively) (Figs 5 and 6).
We interpret the increases in protein abundance observed in cluster III in comparison to cluster I in M. trossulus as indicating the cellular response of gills to compensate for the loss of the cellular stress response under sirtuin inhibition. This is evident in the increase in protein abundance of proteasome α-type 5 (spot 79) and the two HSP70s in cluster III, which showed an interaction effect (Figs 5B and 6B). This corresponds to the compensatory response as a consequence of sirtuin inhibition seen for SIRT5 levels in M. trossulus (Fig. 4). Conversely, the lack of a difference in protein abundance with and without suramin in M. galloprovincialis may indicate that sirtuins play a less important role in its stress response. Alternatively, the signal for inducing sirtuin activity on the proteostasis proteome may not have yet been triggered at 35°C HS in M. galloprovincialis, although this is unlikely. Interestingly, SIRT1 is known to regulate the inducible heat shock response and transcription of HSP70 by deacetylating HSF1, which prolongs its binding to the heat shock promoter element in humans (Westerheide et al., 2009). Our findings indicate that the gill proteome is highly heat sensitive and that sirtuins play a global role in regulating molecular chaperones and protein degradation, and therefore protein homeostasis during heat shock, more so in M. trossulus than in M. galloprovincialis.
Sirtuin inhibition during HS treatment also affected mRNA processing and translation (Figs 5 and 7). We identified two isoforms of the translation initiation activator DAZAP (Smith et al., 2011) and the elongation factor eIF5a in M. galloprovincialis mussels as showing interaction effects under 35°C HS+suramin (Figs 5A and 7A). Interestingly, while one isoform of DAZAP (spot 175) increased in protein abundance under HS, another isoform (spot 174) increased in protein abundance almost 10 times under 35°C HS+suramin compared with 35°C HS alone, suggesting a direct role of sirtuins in regulating translation initiation in M. galloprovincialis (Figs 5A and 7A).
In M. trossulus, we identified two isoforms of Musashi1, a translation inhibitor (Gunter and McLaughlin, 2011), and serine/arginine (SR)-rich protein, a splicing repressor during HS (Shi and Manley, 2007), both of which also showed interaction effects under 35°C HS+suramin (Figs 5B and 7B). SR-rich protein (cluster II) clustered with other chaperones found to increase in abundance under 32°C HS+suramin (Figs 5B and 7B). RNA-binding protein Musashi (spot 159; cluster I) showed decreased protein abundance under 35°C HS+suramin while an additional isoform (spot 119) showed increased protein abundance under the same treatment in cluster III (Figs 5B and 7B). The occurrence of two isoforms with complementary patterns at the same temperature may represent a direct effect of sirtuins on Musashi1 and thus a difference in how sirtuins affect the proteomic responses to acute HS differently in the two Mytilus congeners. Specifically, sirtuin inhibition affected translation initiation in M. galloprovincialis while it affected translation inhibition and splicing in M. trossulus during HS.
Several enzymes of glycolysis, the pentose phosphate pathway (PPP), Krebs cycle, the ETS, oxidative phosphorylation and fatty acid β-oxidation responded to HS and suramin treatment in species-specific patterns.
In the more heat-tolerant M. galloprovincialis, we identified metabolic proteins that showed an interaction effect but without a general trend of how HS+suramin affected protein abundance (Fig. 8A). The changes involved cytosolic and mitochondrial isoforms of malate dehydrogenase (c/mMDH), dihydrolipoyl dehydrogenase, which is part of the pyruvate dehydrogenase complex, and cytochrome c reductase (Figs 8A and 9A). In contrast, UDP-galactose 4-epimerase showed a significant increase in protein abundance under 35°C HS in comparison to lower temperatures, but a decrease under 35°C HS+suramin relative to 35°C HS alone, suggesting a role for sirtuins in regulating its abundance (Figs 8A and 9A). UDP-galactose 4-epimerase functions in the conversion of UDP-galactose to UDP-glucose in the Leloir pathway (Frey, 1996), and the observed protein abundance changes in this species suggest that a greater energy demand is fulfilled by glycolysis under HS, or a shunting towards the PPP for the production of NADPH for ROS scavenging (Chandel, 2015). As changes in the glycolytic triose phosphate isomerase were non-significant and enolase showed a significant decrease with increasing temperature, we can probably rule out an increase in the rate of glycolysis and infer that the shunting towards the production of reducing equivalents in the form of NADPH via the PPP was up-regulated instead (Figs 8A and 9A).
In M. trossulus, proteins with decreased abundance under HS+suramin (cluster I) included amino acid metabolism proteins caffeoyl-CoA O-methyltransferase, ornithine aminotransferase (OAT), and the DNA/RNA synthesis protein ribonucleoside-diphosphate reductase α (spot 194) (Figs 8B and 9B). OAT functions as an enzyme in ornithine degradation of the urea cycle and studies have shown that there is a functional ornithine cycle in bivalves (Andrews and Reid, 1972; Meijer et al., 1990). Interestingly, SIRT5 is known to deacetylate and thereby activate the first enzyme of the urea cycle (carbamoyl phosphate synthetase 1) (Nakagawa et al., 2009). Therefore, the decrease in protein abundance of OAT observed under 35°C HS+suramin is most likely linked to inhibition of deacetylase activity of SIRT5 with respect to carbamoyl phosphate synthetase 1, which thereby down-regulates urea cycling and ornithine synthesis. Proteins with a general increase in abundance under HS+suramin (cluster III) included phenylalanyl-tRNA synthetase β, the malate–aspartate shuttle protein cMDH, and the ETS protein NADH dehydrogenase (Figs 8B and 9B). Finally, SIRT5 is known to hypersuccinylate lysine residues in proteins involved in the metabolic pathway for β-oxidation (Rardin et al., 2013), and we identified two proteins representing β-oxidation [electron transfer flavoprotein β (ETF) and enoyl CoA hydratase (EHHADH)] that showed an interaction and one [D-3 hydroxybutyrate dehydrogenase (BDH)] that showed a temperature effect in M. trossulus (Figs 8B and 9B).
The identification of an isoform of ATP synthase in M. trossulus (cluster III) that had a temperature main effect with increased abundance under 35°C HS and 35°C HS+suramin may indicate an increase in the activity of the fermentation of glucose to succinate and propionate, especially in combination with a sirtuin-dependent increase in NADH dehydrogenase and cMDH (Fig. 8B, cluster III). Mytilus anaerobic energy metabolism relies on malate dismutation and cMDH to reduce malate to succinate, which is further converted to propionate (Müller et al., 2012). Furthermore, the reducing equivalents produced during the dismutation of malate can be oxidized by NADH dehydrogenase and the resulting proton gradient can be used by ATP synthase to produce ATP (Müller et al., 2012). These protein changes thus represent the alternative metabolic pathway occurring during hypoxia under aerial emersion and prolonged anaerobiosis in Mytilus (Fields et al., 2014; Müller et al., 2012). This could suggest that M. trossulus relies more on this alternative anaerobic pathway with increasing HS than M. galloprovincialis and that sirtuins are involved in the regulation of this pathway during heat shock.
NADP-isocitrate dehydrogenase (IDH), which contributes to ROS scavenging by producing NADPH to reduce glutathione and thioredoxin (Go and Jones, 2008; Tomanek, 2015), and succinate dehydrogenase (SDH) showed only a temperature effect in M. trossulus (cluster II) (Figs 8B and 9B). Both NADP-IDH and SDH may have shown a temperature effect only because the sirtuin that most likely affects their acetylation (SIRT3) is not inhibited by suramin (Guan and Xiong, 2011); however, SIRT5 is also known to affect SDH through desuccinylation (Rardin et al., 2013). Phospho-glucono-lactonase, an enzyme of the PPP, which can contribute NADPH for ROS scavenging, showed individual temperature and suramin main effects but not an interaction effect in M. trossulus (Fig. 8B). These data suggest that there is a greater demand for NADPH to reduce ROS and limit oxidative damage during acute heat shock, which is consistent with previous findings in Mytilus (Tomanek, 2015; Tomanek and Zuzow, 2010).
Suramin-dependent changes in EHHADH, an enzyme that catalyzes the first reaction of fatty acid β-oxidation, and ETF, the flavoprotein that carries reducing equivalents (FADH2) to coenzyme Q of the ETS, were only observed in M. trossulus (Figs 8B and 9B). However, while sirtuin inhibition increased the abundance of EHHADH relative to the 32°C HS control, ETF levels decreased (Figs 8B and 9B). A metabolite of β-oxidation, acetoacetyl-CoA, can be converted first to acetoacetate and then to β-hydroxybutyrate. It is known that β-hydroxybutyrate is an inhibitor of histone deacetylases which provides protection against oxidative stress (Shimazu et al., 2013), and we observed a protein abundance decrease with HS (Figs 8B and 9B). Although the two β-oxidation enzymes showed opposite effects in response to sirtuin inhibition, ETF and BDH clustered together (cluster I) (Fig. 8B). Given that all three enzymes are part of the two connected β-oxidation and ketone synthesis pathways, it may be interesting to ask how the combination of these pathways may contribute to ameliorating heat stress. We offer the following possibility. Based on studies assessing the production of ROS in isolated mitochondria with different substrates and inhibitors, the transfer of reducing equivalents (FADH2) produced during β-oxidation by ETF to coenzyme Q produces significantly fewer ROS at complex I (NADH dehydrogenase) compared with when FADH2 enters the ETS through SDH and succinate is the substrate, possibly because of reverse electron transfer (Schönfeld et al., 2010). Given the greater loss of an oxidative stress response with increasing temperatures in M. trossulus (Tomanek and Zuzow, 2010), and the suramin-dependent response of oxidative stress proteins to heat in this Mytilus congener (see below), we think a shift from succinate to fatty acids (carnitine esters) as the main substrate for energy metabolism, possibly regulated by sirtuins, could reduce the production of ROS in this species as a last-ditch effort to avoid greater oxidative damage. We postulated a similar role for SIRT5 in switching from pro-oxidant NADH-producing pathways to anti-oxidant NADPH-producing pathways previously (Tomanek and Zuzow, 2010). The occurrence of a temperature effect without a suramin interaction effect for SDH may be explained by the fact that suramin does not affect SIRT3 (but see Rardin et al., 2013). Thus, we hypothesize that, in M. trossulus, the mitochondrial SIRT5 may be partly responsible for a shift in energy metabolism (e.g. β-oxidation versus Krebs cycle), mainly to reduce ROS production from the ETS. In contrast, M. galloprovincialis seems to respond more robustly to HS, and this may be best shown by its different strategy of responding to oxidative stress.
Several oxidative stress proteins, including two Cu-Zn SOD isoforms and one ALDH showed interaction effects in M. galloprovincialis (Figs 8A and 9A). Cu-Zn SOD may reside in the intermembrane space of the mitochondria (Tomanek, 2015) and we identified two isoforms (spots 163 and 167) that showed complementary and opposite abundance changes with the addition of suramin (Figs 8A and 9A). These two isoforms differ in isoelectric point by 0.08 pI units and molecular mass by 2 kDa (Table S1A), which suggests a PTM, possibly several lysine-acylations, which have been shown to affect the activity of Mn SOD (Tao et al., 2010). Therefore, while one isoform decreased in protein abundance (spot 163) the other (spot 167) increased and these changes were related to the interaction of suramin with HS at 35°C and the resulting inhibition of SIRT1, 2 and 5 (Figs 8A and 9A). A switch in abundance between two isoforms of Cu-Zn SOD with increasing heat stress was shown for M. galloprovincialis in a prior study (Tomanek and Zuzow, 2010).
In M. galloprovincialis, peroxiredoxin 4 (PRX4), which reduces hydrogen peroxide in the cytosol and ER to prevent ROS-induced cellular damage (Hanschmann et al., 2013), showed a temperature effect with increasing abundance under 35°C HS in comparison to the control temperature (Figs 8A and 9A). Given the steady increase of PRX4 with increasing HS, it is likely to be the main inducible antioxidant during HS exposure in M. galloprovincialis. With the addition of suramin under HS, increased protein abundance of a Cu-Zn SOD isoform indicates use of an additional antioxidant in conjunction with PRX4, possibly through a PTM, illustrating the ability of M. galloprovincialis to utilize Cu-Zn SOD as a second line of defense in addition to the already present PRX4. These findings are consistent with previous work identifying antioxidants as a method to combat ROS-induced cellular damage during heat shock and interspecific differences between Mytilus species (Lockwood et al., 2010; Tomanek and Zuzow, 2010). Interestingly, changes in the mitochondrial PRX5 were not significant, possibly suggesting that mitochondrial ROS production did not increase significantly or that it was sufficiently controlled by Cu-Zn SOD in M. galloprovincialis (Fig. 8A).
Mytilus trossulus showed interaction effects for Cu-Zn SOD, DyP-type peroxidase, two ALDHs and PRX4 (Figs 8B and 9B; cluster I). Soma ferritin, which binds iron to prevent the Fenton reaction with hydrogen peroxide and the subsequent formation of hydroxyl radicals (Tomanek, 2015), showed an interaction effect in M. trossulus but only a temperature effect in M. galloprovincialis (Figs 8B and 9B). In M. trossulus, proteins in cluster I (Fig. 8B) showed a trend towards decreasing abundance when under 32°C HS+suramin and 35°C HS with and without suramin, indicating the upper thermal limit of an antioxidant response in this species, which is consistent with previous findings (Tomanek and Zuzow, 2010). We observed a greater decrease in abundance of Cu-Zn SOD under 35°C HS+suramin in comparison to abundance under 13°C and suramin, or 32°C HS+suramin, and this was observed for PRX4, DyP-type peroxidase and ALDH (Figs 8B and 9B). Additionally, glutathione S-transferase (GST; cluster I) showed a temperature effect with protein abundance decreasing under 35°C HS (Figs 8B and 9B). Our data suggest that in M. trossulus, the main inducible antioxidants under 32°C HS are GST, PRX4 and Cu-Zn SOD. The lower induction temperature of these proteins in comparison to findings for M. galloprovincialis parallels the lower thermal tolerance in M. trossulus. In addition, under 35°C HS+suramin, we did not observe increased protein abundance of complementary isoforms of Cu-Zn SOD, as we did in M. galloprovincialis.
DyP-type peroxidase, a possible scavenger of peroxynitrite in the peroxisome (Trujillo et al., 2008), showed a decrease in abundance under 35°C HS+suramin compared with 35°C HS alone in M. trossulus (Figs 8B and 9B). These data suggest a decrease in protein abundance of DyP-type peroxidase is influenced by increasing acute HS (mainly at 35°C) while the addition of suramin further enhances this response (Figs 8B and 9B).
Of the proteins related to oxidative stress that were identified, only SOD and ALDH are known to have a relationship with sirtuins. Specifically, Mn SOD, a SOD distantly related to Cu-Zn SOD, is deacetylated by SIRT3, which enhances the antioxidant activity of Mn SOD (Chen et al., 2011; Tao et al., 2010). ALDH functions in the oxidation of aldehydes and as an antioxidant against lipid peroxidation (Somero et al., 2017; Tomanek, 2015), and deacetylation of ALDH by SIRT3 inhibits ALDH activity (Xue et al., 2012). However, the sirtuin inhibitor used in our study does not affect the function of SIRT3 and it is unclear whether the protein abundance changes observed are a direct result of suramin-induced sirtuin inhibition or whether there may be a possible feedback loop inducing SIRT3 activity in response to suramin.
Thus, sirtuin inhibition interferes more with the response to heat-induced oxidative stress in M. trossulus than in M. galloprovincialis, with the latter able to induce a second line of antioxidant defense (Cu-Zn SOD), in addition to PRX4, under sirtuin inhibition. In contrast, M. trossulus experienced significant decreases in antioxidant proteins (PRX4 and Cu-Zn SOD) under HS+suramin, with no secondary response observed.
Cytoskeletal and actin-binding proteins
Cytoskeletal and actin-binding proteins also showed species-specific suramin-dependent responses to HS (Figs 10 and 11). In M. galloprovincialis, profilin, a protein that binds actin and regulates the rate of actin filament growth (Witke, 2004), and Na+/H+ exchange regulatory cofactor (NHERF1), which functions as a scaffolding protein (Lee et al., 2007), showed interaction effects (Figs 10A and 11A). Nine out of the top 10 most negative loading values for M. galloprovincialis along PC1 were cytoskeletal proteins, mainly α/β-tubulin and actin (Table 1). This suggests an overall decrease in protein abundance for cytoskeletal proteins with increasing heat shock exposure, possibly indicating cytoskeletal rearrangement, which may stabilize the cytoskeleton (Tomanek and Zuzow, 2010). However, it is unclear whether the protein abundance changes observed indicate a loss of or greater cytoskeletal stability, which could be confirmed in future studies by analyzing the gill tissue structure.
In contrast, the majority of the cytoskeletal and actin-binding proteins in M. trossulus showed interaction effects, especially at 32°C and 35°C HS (14 tubulin isoforms, eight actins, five actin-regulatory proteins, including F-actin capping protein, two profilin isoforms, gelsolin, tropomyosin and fascin, and a microtubule-associated protein known as tektin-4) (Figs 10B and 11B). In addition, proteins in cluster II showed a general decrease in protein abundance with increasing HS, especially at 35°C HS, similar to M. galloprovincialis (Figs 10B and 11B). This decrease became more pronounced with the addition of suramin (Figs 10B and 11B). In comparison to M. trossulus, we identified nine α/β-tubulins and fascin as not changing significantly in M. galloprovincialis (Fig. 10A).
SIRT1 is known to deacetylate and activate the actin-binding protein cortactin (Zhang et al., 2009), while SIRT2 can interact with α-tubulin (North et al., 2003). As SIRT1 and 2 were inhibited by suramin in our study, the observed protein abundance changes demonstrate a widespread effect of sirtuins on cytoskeletal elements and actin dynamics during heat shock in M. trossulus. However, it is also possible that these effects of sirtuins are mediated through a signaling cascade, e.g. MAPK pathways (Evans and Somero, 2010).
The abundance of a number of signaling proteins changed in a suramin-dependent way during HS, much more so in M. trossulus than in M. galloprovincialis, especially at 32 and 35°C (Figs 12 and 13). In M. galloprovincialis, we identified the small G-proteins ADP-ribosylation factor (ARF) and Cdc42 (chain A), which function in regulating intracellular vesicular trafficking and actin rearrangements, respectively (Cerione, 2004; Erickson and Cerione, 2001), and both exhibited a temperature effect with increased abundance at 35°C HS (Figs 12A and 13A). In contrast, major vault protein, which functions in regulating extracellular signal-regulated kinase (ERK) signaling pathways (Liang et al., 2010), had a significant interaction effect with a trend towards decreased abundance with increasing HS+suramin (Figs 12A and 13A).
Mytilus trossulus showed interaction effects for the two small G-proteins ARF and Rab-GDP dissociation inhibitor (GDI), with increased protein abundance under HS+suramin (Figs 12B and 13B). This suggests that sirtuins inhibit the formation and transport of intracellular vesicles during heat shock (Marks et al., 2009). Several heterotrimeric and small G-proteins (Rho-GDI), isoforms of the phosphoprotein chaperone 14-3-3, and MAPKs, all regulating cytoskeletal dynamics (Marks et al., 2009), showed interaction effects at 35°C in M. trossulus, but not in M. galloprovincialis (e.g. Cdc42 and ERK1/2) (Figs 12 and 13). Of those proteins showing interaction effects in M. trossulus, 14-3-3, major vault protein, Rho-GDI and calcyphosin showed an increase in protein abundance under 35°C HS (Figs 12B and 13B). Conversely, Rab-GDI α and guanine nucleotide binding protein β also showed interaction effects for HS+suramin in which we observed a decrease in protein abundance under 35°C HS, but this response changed with the addition of suramin when protein abundances increase (Figs 12B and 13B). It is thus possible that sirtuins inhibit cytoskeletal modifications during heat shock more so in M. trossulus than in M. galloprovincialis. More generally, these results demonstrate a greater role of sirtuins in affecting stress signaling pathways that modify vesicle transport and the cytoskeleton in M. trossulus than in M. galloprovincialis.
Analysis of acetylated proteins
While the interspecific difference in the number of acetylated proteins affected by HS and sirtuin inhibition obtained using an acetyl-lysine antibody and 2DGE was the same as that with the total protein stain, we detected overall fewer acetylated proteins being affected during HS (35°C) by suramin (Fig. 3). We believe that this suggests that PTMs other than deacetylation, e.g. desuccinylation of lysines, as they were not detected by the acetyl-lysine-specific antibody, were responsible for several of the interaction effects observed (Park et al., 2013; Rardin et al., 2013). Also, it is possible that sirtuin inhibition modified signaling cascades that affected protein abundance further downstream without directly deacylating the proteins. However, these results do further support our inferences that sirtuins play distinct roles in the stress responses of the Mytilus congeners.
Western blot analysis of SIRT5
Western blot analysis of SIRT5 found 3 times greater levels in M. trossulus than in M. galloprovincialis (Fig. 4). In addition, sirtuin inhibition by suramin led to an increase in SIRT5 in M. trossulus but not M. galloprovincialis, suggesting that a feedback mechanism to sirtuin inhibition exists in M. trossulus. These results support the possibility that other observed interspecific differences in the response of the proteome to heat shock are linked to the response of SIRT5 to heat shock. For example, as a demalonylase, the activity of SIRT5 may be linked to malonyl-CoA and its ability to inhibit carnitine acyltransferase and thereby the production of carnitine esters for β-oxidation (Derrick and Ramsay, 1989). As a desuccinylase, SIRT5 has been shown to affect several enzymes of β-oxidation (Park et al., 2013; Rardin et al., 2013). We thus hypothesize that the interaction effects of two enzymes of β-oxidation are linked to the activity of SIRT5 (Du et al., 2011; Peng et al., 2011), possibly to decrease ROS production by NADH dehydrogenase (Schönfeld et al., 2010) by simultaneously affecting SDH, which changed in response to the acute HS in M. trossulus (Fig. 4B). This conjecture is consistent with the observation that the abundance of antioxidant enzymes decreases with HS, especially in M. trossulus (Tomanek and Zuzow, 2010), and that sirtuin inhibition affects the change in abundance of several oxidative stress proteins during heat shock (Fig. 8B).
In summary, heat shock in conjunction with sirtuin inhibition through suramin affected a greater number of molecular chaperones (specifically, mitochondrial and ER chaperones) in the cold-adapted M. trossulus than in the warm-adapted M. galloprovincialis. Also, we observed increasing levels of a translation initiation protein (DAZAP) in M. galloprovincialis and decreasing levels of a translation repressor (Musashi 1) in M. trossulus under HS+suramin treatment in comparison to HS alone, and therefore observed divergent effects of sirtuin inhibition on proteins involved in translation initiation and repression between the congeners. Metabolic enzymes, specifically of β-oxidation, that could affect ROS production by the ETS changed under HS+suramin in M. trossulus only. Sirtuin inhibition affected the abundance of two Cu-Zn SOD isoforms in M. galloprovincialis, while it affected several antioxidant proteins (Cu-Zn SOD, DyP-type peroxidase, PRX4 and ferritin) capable of scavenging ROS in M. trossulus. Few cytoskeletal proteins changed with HS+suramin in M. galloprovincialis while 28 showed interaction effects and decreased protein abundance in M. trossulus, especially at 32°C. Given the higher abundance of SIRT5 in M. trossulus observed through western blotting even under non-stressful conditions, and its increase in response to suramin in M. trossulus during heat shock, SIRT5 may play a greater role in the stress response of this comparatively heat-sensitive species (Du et al., 2011; Peng et al., 2011). Additionally, recent surveys of the SIRT5-dependent changes in the succinylome in mice showed an overlap with reactions and pathways that we found to change with heat stress in Mytilus, e.g. ATP synthase, TCP1, GST, respiratory chain complex, pyruvate dehydrogenase, β-oxidation and ketogenesis (Park et al., 2013; Rardin et al., 2013). These previous studies further support our hypothesis that the interspecific differences in the proteomic response to heat shock are due in part to SIRT5, which we initially identified as responding differently to an acute heat shock between the congeners (Tomanek and Zuzow, 2010).
Our results demonstrate that the higher thermal tolerance of M. galloprovincialis may come from a more robust ER and mitochondrial proteome, a metabolism that produces fewer ROS or has a constitutively more efficient ROS scavenging system, possibly based on Cu-Zn SOD and PRX4. The decrease of key molecular chaperones and antioxidants and differences in regulating translation under heat shock and sirtuin inhibition highlight the importance of sirtuins in regulating the cellular stress response in M. trossulus, more so than in M. galloprovincialis. Perhaps most importantly, our results expand our understanding of the role of sirtuins to broadly regulating the proteomic response to heat shock, demonstrate their potential role in setting interspecific differences in thermal tolerance, and suggest that one function of SIRT5 is to shift metabolic pathways towards producing fewer ROS under environmental stress.
We thank the National Science Foundation for supporting this work through grants to M.C.V. and L.T.
Conceptualization: M.C.V., L.T.; Methodology: L.T., M.J.Z.; Formal analysis: M.C.V., L.T., M.B., M.J.Z.; Investigation: M.C.V., M.B., M.J.Z.; Resources: L.T.; Writing - original draft: M.C.V., L.T., S.B.; Writing - review & editing: M.C.V., L.T., S.B.; Visualization: M.C.V., L.T., M.J.Z.; Supervision: L.T.; Project administration: L.T.; Funding acquisition: L.T.
This work was supported by the National Science Foundation (grant IOS-0717087 to L.T.). During writing of the manuscript, M.C.V. was supported by funding from the National Science Foundation (PRFB DBI-1401357).
The authors declare no competing or financial interests.