Molecular Basis of S100A1 Activation and Target Regulation Within Physiological Cytosolic Ca2+ Levels

The S100A1 protein regulates cardiomyocyte function through its binding of calcium (Ca²⁺) and target proteins, including titin, SERCA, and RyR. S100A1 presents two Ca²⁺ binding domains, a high-affinity canonical EF-hand (cEF) and a low-affinity pseudo EF-hand (pEF), that control S100A1 activation. For wild-type S100A1, both EF hands must be bound by Ca²⁺ to form the open state necessary for target peptide binding, which requires unphysiological high sub-millimolar Ca²⁺ levels. However, there is evidence that post-translational modifications at Cys85 may facilitate the formation of the open state at sub-saturating Ca²⁺ concentrations. Hence, post-translational modifications of S100A1 could potentially increase the Ca²⁺-sensitivity of binding protein targets, and thereby modulate corresponding signaling pathways. In this study, we examine the mechanism of S100A1 open-closed gating via molecular dynamics simulations to determine the extent to which Cys85 functionalization, namely via redox reactions, controls the relative population of open states at sub-saturating Ca²⁺ and capacity to bind peptides. We further characterize the protein's ability to bind a representative peptide target, TRKT12 and relate this propensity to published competition assay data. Our simulation results indicate that functionalization of Cys85 may stabilize the S100A1 open state at physiological, micromolar Ca²⁺ levels. Our conclusions support growing evidence that S100A1 serves as a signaling hub linking Ca²⁺ and redox signaling pathways.

1. Introduction

S100A1 is a Ca²⁺ binding protein that is implicated in cardiac and neurological functions (Wright et al., 2009b). S100A1 regulates several targets including ryanodine receptor (RyR), sarcoplasmic/endoplasmic reticulum calcium ATPase (SERCA), phosphoglucomutase, tubulin, and tumor protein p53 (Landar et al., 1996; Santamaria-Kisiel et al., 2006; Wright et al., 2008; Duarte-Costa et al., 2014) in a Ca²⁺-dependent manner. S100A1 expressed in cardiac tissue is believed to manage contractile behavior either through modulating cytosolic Ca²⁺(Kraus et al., 2009), which triggers the initiation of contraction, or through modulating properties of the contractile fibers of the myofilament. For the latter, there is evidence that Ca²⁺-loaded S100A1 disrupts interactions between actins of the thin filament and titin (Gutierrez-Cruz et al., 2001; Granzier et al., 2010; Yamasaki et al., 2011). Specifically, there are several reports that the PEVK-rich regions of titin bind actin and thereby behave as a viscous brake during extension and contraction. Curiously, in vitro assays indicate S100A1 binds targets at Ca²⁺ levels considerably above physiological Ca²⁺ concentrations, which casts doubts on the ability of wild-type (WT) S100A1 to contribute to titin's management of contractile properties.

S100A1 belongs to the EF-hand calcium binding protein family in which the Ca²⁺ ions are bound to the helix-loop-helix motif (EF-hand). Two EF-hands exist in the S100A1 monomer, the canonical hand (cEF) and the pseudo hand (pEF). The cEF hand has a higher Ca²⁺ affinity than the pEF hand, with the dissociation constant for the former is ~27–50 μM and latter is ~250–16,700 μM, respectively (Goch et al., 2005; Wright et al., 2005). Within cells, the dissociation constant for the pEF hand is magnitudes larger than cellular Ca²⁺ concentration, implying that this pEF hand does not significantly bind Ca²⁺. However, it has been shown that binding Ca²⁺ at both the pEF and cEF hands is the prerequisite to S100A1's interactions with Ca²⁺-dependent targets (Nowakowski et al., 2013). Thus, it is of great importance to understand how S100A1 is activated under non-(half)-saturating cellular Ca²⁺ conditions.

S100A1 is generally found as a homo-dimer when subject to conditions amenable to structure determination via x-ray crystallography (Melville et al., 2017) or nuclear magnetic resonance (NMR) spectroscopy (Wright et al., 2005; Nowakowski et al., 2011, 2013); its dimerization happens at picomolar monomer concentrations (Kraus et al., 2009). Similar to most members of the S100 class of Ca²⁺ binding proteins, S100A1 activation proceeds through binding of two Ca²⁺, one each at the low-affinity pEF and high affinity cEF hands. In its fully-saturated, Ca²⁺-bound (holo) state, S100A1 presents a hydrophobic patch between helices three and four (H3 and H4) that enable binding to regulatory domains of protein targets (Wright et al., 2005; Nowakowski et al., 2011), which may be accompanied by significant increases in solvent accessible surface area relative to the apo state (Chaturvedi et al., 2020). In the absence of Ca²⁺(apo state), the hydrophobic patch is concealed by closing the H3/H4 hinge region (Figure 1) (Nowakowski et al., 2011). Our recent studies suggested that the half-saturated state of the S100A1 (Scott and Kekenes-Huskey, 2016), which we characterize as the conformation with a single Ca²⁺ bound at the cEF hand, is insufficient to maintain an exposed regulatory binding region. Since the pEF Ca²⁺ affinity is reported to be in the sub-millimolar range (Goch et al., 2005), these findings suggested that S100A1 may be incapable of recognizing Ca²⁺-dependent targets (Scott and Kekenes-Huskey, 2016). Hence, native S100A1 is unlikely to dynamically regulate protein targets, such as passive tension, within physiological Ca²⁺ ranges.

FIGURE 1

Figure 1. Schematic of S100A1 interacting with Ca²⁺-dependent targets. Upon Ca²⁺ binding, the hydrophobic patch between H3 and H4 is exposed to target peptide. The Cys85 at H4 is shown in sticks, Ca²⁺ ions are shown as yellow spheres and target peptide is colored gray.

Interestingly, functionalization of Cys85 at the C-terminus of H4 has been shown to yield 10- and 10,000-fold increases in Ca²⁺ affinity at the cEF and pEF hand, respectively (Goch et al., 2005). The increasing Ca²⁺ affinity was reported to be caused by the favorable cooperativity between the binding events in the two EF-hands (Goch et al., 2005). The enhanced Ca²⁺ affinity could therefore confer the ability to activate S100A1 at micromolar Ca²⁺. Indeed, an assay by Goch et al. (2005) indicated that the modified S100A1 protein presented its peptide-binding hydrophobic patch at physiological Ca²⁺ levels in contrast to the native protein.

Since S100A1 Cys85 is a known target for in vivo nitrosylation (Živković et al., 2012) and glutathionylation (Goch et al., 2005), we explored via molecular dynamics simulation whether mutation of the cysteine to bulkier side groups could promote S100A1 domain opening in its half-saturated state. Additionally, although atomistic-resolution structures of S100A1 complexed with RyR regulatory peptide and a 12-residue peptide “TRTK12” have been determined (Wright et al., 2008), neither the motif nor the molecular structure of an S100A1/PEVK complex have been established. Under the assumption that S100A1 binds PEVK in a conformation similar to other established Ca²⁺-dependent targets, namely TRTK12, and those exhibited among other S100 family proteins, we investigate conditions necessary for the regulatory protein's activation under physiological Ca²⁺ concentrations. We specifically examined S100A1's capacity to bind regulatory peptides in WT and C85 mutants. Furthermore, we demonstrate how changes in the Ca²⁺-sensitivity of S100A1 activation could control competition between PEVK/S100A1 and PEVK/actin binding using a multi-state macroscopic model.

2. Materials and Methods

2.1. Molecular Dynamics Simulation

The starting structures are based on the NMR apo and holo structures of S100A1 [PDBs 2L0P (Nowakowski et al., 2011) and 2LP3 (Nowakowski et al., 2013), respectively]. The mutations of C85 to E85 and R85 were performed by Charmm-Gui utility (Jo et al., 2008). For the half-saturated state, the Ca²⁺ at pEF site was deleted. The system was solvated in TIP3P (Jorgensen et al., 1983) water box with 20 Å margin. 0.15 M KCl ions were added into the system to maintain a physiological ionic strength. The system was parameterized by the AMBER ff12SB force field (Case et al., 2012) with Ca²⁺ parameters adapted from the Li-Merz work (Li et al., 2003). Each system was first subjected to an energy minimization process followed by the heating stage. During the heating stage, the weak-coupling algorithm was used. After reaching equilibrium state with temperature being 300 K, a 100 ns production MD was performed by using the PMEMD.CUDA module of the AMBER 14 package (Case et al., 2014). Clustering analysis were performed on this 100 ns production trajectory via CPPTRAJ (Roe and Cheatham, 2013) using a hierarchical agglomerative (bottom-up) approach. The representative structures of the 3 least populated clusters were used as starting structures for next cycle of production runs (each run was about 400 to 700 ns long).

The accumulative sampling time for each case is around 2 μs for the apo/half-saturated state and about 0.6 μs for the fully-saturated state. The time step was 2 fs and cutoff for non-bonded interaction was set to 10 Å. During the whole MD, the SHAKE algorithm (Ryckaert et al., 1977) was used to restraint the length of bonds involving hydrogen atom. All simulations are summarized in Table S1. The CPPTRAJ (Roe and Cheatham, 2013) program from Amber was used to calculate the root mean squared fluctuations (RMSF), α-helical probability, contact map and H3/H4 inter-helix angle values. The RMSF were calculated on backbone atoms (Cα, C, N, and O atoms). The α-helical probability for residues in the H4 C-terminus (residues 85-93) were calculated using the Define Secondary Structure of Proteins (DSSP) algorithm (Kabsch and Sander, 1983). Contact map data was calculated with distance cutoff as 7 Å and only residue pairs which are at least 6 residues apart (i and i+6) in sequence were considered.

2.2. Potential of Mean Force of TRTK12 Peptide Unbinding From S100A1

TRTK12 is a 12 residue peptide that was reported to bind to S100A1 in a Ca²⁺-dependent manner (Ivanenkov et al., 1995). To simulate the dissociation PMFs of TRTK12 peptide from S100A1, we constructed the TRTK12-S100A1 complex structures using our simulated half- and fully-saturated S100A1 structures. Specifically, the MD-sampled most probable S100A1 structures (see section S1.1 in Supplementary Material for details of determining these structures) were superimposed with NMR structure of TRTK12-S100A1 complex (PDB 2KBM). We minimized the energy after superpositioning the experimentally-determined TRTK peptide structure into the hydrophobic pocket of MD-simulated S100A1 to eliminate potentially overlapping atoms. After the minimization, we equilibrated the system and further run a 60 ns long MD simulation of the TRTK12-S100A1 complex in the WT fully-saturated state to assess the binding of the TRTK12 peptide at the hydrophobic cleft (Figure S3). The reaction coordinate (RC) was defined as the distance between the center of mass (COM) of peptide Cα atoms to the COM of Cα atoms in the H3 C-terminus (residues K30-T39) and the H4 N-terminus (residues E73-A84). We note that we used NAMD in order to make use of its support for collective reaction coordinate variables (the protein COMs). This required parameterizing the system using the CHARMM36 force field (MacKerell et al., 1998, 2004). As such, all structures subjected to PMF calculations were compared using the CHARMM36 force field, whereas all other simulations used the AMBER ff14SB parameterization to ensure consistent comparisons. After obtaining the TRTK12-S100A1 complex structures, the system was then solvated in a TIP3P water box with 14 Å margin. 0.15 M KCl ions were added into the system. The CHARMM36 (MacKerell et al., 1998, 2004) force field was used. The sampling is preformed with RC ranging from 9.5 to 33.5 Å resulting total 49 simulation windows with window width as 0.5 Å. For each window, after minimization and equilibrium, 8 ns production MD was performed in the NPT ensemble at 300 K temperature. An harmonic potential was applied at the center of each window with a force constant as 18 kcal/mol/² during the simulations. Two loose angle constraints were introduced to prevent the peptide from sliding along the S100A1 surface. All PMF calculations were performed via NAMD2.11 (Phillips et al., 2005). The PMFs along reaction coordinate were constructed using the WHAM program (Grossfield). The PMF error was estimated using the Monte Carlo Bootstrap Error Analysis function in the WHAM program (with 30 MC trials). We also performed the PMF error analysis based on block averaging. Specifically, the variance of the RC in each simulation window from umbrella sampling was estimated via block averaging analysis with block size as 200 (Figures S11A,B). The RC variance was then used to calculate PMF error (Figure S11C) based on Equation (1) in Zhu and Hummer (2012). Since the magnitudes of PMF errors from MC trails and block averaging are comparable, we show in the main text the PMF with the MC error analysis.

2.3. Molecular Mechanics-Generalized Born and Surface Area Continuum Solvation (MM-GBSA) Calculation

The per residue energy contribution to the interaction energy between monomers of S100A1 was calculated via MMGBSA.

$\begin{array}{ll}\Delta G=\langle G_{dimer}\rangle -\langle G_{monomer}\rangle & (1)\end{array}$

where 〈G_dime〉 and 〈G_monomer〉 are ensemble-averaged MMGBSA calculated free energies of S100A1 dimer and monomer, respectively. The calculations were performed on a subset of MD trajectory extracted at a 2 ns frequency. The generated sub-trajectories were used as input of MMPBSA.py in Amber16 to calculate the free energies of each part. The salt concentration was set as 0.15 M with the generalized Born model option setting as igb = 5. No quasi-harmonic entropy approximation was made during the calculation. The total ΔG was decomposited to achieve per residue contribution by setting dcomp = 2 during the calculation. We excluded the internal energies (bonded terms) from the final results as these energies are nearly identical in the dimer and monomer. Thus, our final per residue energy contribution contains three energy terms: electrostatic interactions (E_EEL), van der Waals interactions (E_vdW), and solvation energy (E_solv).

3. Results and Discussion

3.1. S100A1 Structure and Dynamics

The relative weak Ca²⁺ affinity in the S100A1 pEF hand (K_d~250–16,700 μM; Goch et al., 2005; Wright et al., 2005) suggests that only the protein's cEF site is likely occupied by Ca²⁺ at physiological intracellular Ca²⁺ levels (~0.1 μM; Clapham, 2007). Previous computational studies of S100A1 have utilized the WT structure with cysteine at the 85th position (C85) (Scott and Kekenes-Huskey, 2016). That study suggested that half-saturated S100A1 tends to assume a semi-closed state that would preclude target binding. Since post-translational modification (PTM)s at the S100A1 site C85 on helix 4 (H4) have been shown to increase its Ca²⁺ affinity at both cEF and pEF hands, with the latter exhibiting a four orders of magnitude increase (Goch et al., 2005), we sought to determine potential mechanisms for those PTM-induced changes.

In this study, we introduced glutamic acid (E) or arginine (R) at site C85 to approximate the effects of post-translational modifications including glutathionylation that modulate Ca²⁺ affinity (Goch et al., 2005). We performed extensive all-atomistic molecular dynamics (MD) simulations on the WT and C85E/R S100A1 variants in the apo, half-saturated (only cEF hand has Ca²⁺ bound) and fully-saturated (both cEF and pEF have Ca²⁺ bound) state. In present work, we use the terms “fully-saturated” and “holo” interchangeably. We also performed two extra simulations of S100A1 in the apo state with W90 mutated to alanine (W90A) and in the fully-saturated state with target peptide bound to probe the interaction between H4 and the remaining helices. All simulations as well as the starting structures and the accumulated simulation lengths are summarized in Table S1.

To determine whether using glutamic acid (E) or arginine (R) sufficiently mimicked known PTMs at C85 of S100A1, we compared the MD-sampled C85E/R apo/holo structures with those PTM S100A1 structures deposited in the Protein Data Bank. Namely, structures have been reported for the WT (PDB 2L0P-apo, 2LP3-holo) and variants with PTMs at Cys85 that increase Ca²⁺ affinity. The variants include C85-mercaptoethanol (PDB 2JPT; Zhukov et al., 2008), C85-S-nitrosylation (PDB 2LLT; Živković et al., 2012), and C85M (PDB 2LLS) in the apo state and C85-cysteine (PDB 2LP2; Nowakowski et al., 2013) in the holo state. As shown in Figure 2A, in the apo state, the MD-sampled C85E/R structures both have moderate structural deviations from β-mercaptoethanol modified S100A1, as the RMSDs are around 4.0 Å. Interestingly, most of the structural difference in the C85E/R stem from the displacement of H3 helix away from H4 helix relative to the PTM structure. This suggests that C85E/R variants tend to sample a slightly more open conformation. Additionally, the C85E/R variants present similar inter-subunit contacts as PTM S100A1, as the H1/H4 helices are highly overlapped with the PTM S100A1. Next we compared the MD-sampled C85E/R variants against homocystine-modified S100A1 in the holo state (Figure 2B). Similar to the apo state comparison, although the C85E/R have moderate structural deviations (RMSD values are 4.1 and 4.8 Å, respectively), the differences are primarily due to the H2-H3 linker and H3 helix. Specifically, in the C85E/R variants, the H2-H3 linkers are less folded and are slightly more displaced from the H4 helix. Although the two C85 variants have structural differences, the Ca²⁺ coordination patterns at the EF hands were identical to that of the homocystine-modified S100A1 structure (Figure S1). As both the C85E/R variants and the PTMs have similar inter-subunit contacts in the apo and holo state, the overall structural stability of S100A1 in the C85E/R and PTMs are comparable. These data therefore suggest that the site directed variants we considered are reasonable approximations of the chemically modified proteins reported in the Protein Data Bank. We reported RMSF values before and after Ca²⁺-binding in section S3.1 (Supplementary Material), which are consistent with previous studies.

FIGURE 2

Figure 2. (A) Superposition of MD-sampled most probable C85E/R apo state structures onto β-mercaptoethanol modified C85 apo S100A1 structure (PDB 2JPT). (B) Superposition of MD-sampled holo state structure onto homocysteine modified C85 holo S100A1 structure (PDB 2LP2). The RMSD values were calculated based on Cα atoms. The chemical structures of new residues at C85 were also shown.

3.2. Hydrophobic Pocket Opening Indicated by H3/H4 Inter-helix Angle

The fundamental physiological role of S100A1 is to bind downstream targets after chelating Ca²⁺ ions. Ca²⁺-saturated S100A1 presents a hydrophobic patch between H3 and H4 that engages in target binding, similar to other Ca²⁺-binding proteins like calmodulin and troponin C. We investigated the ability of C85 variants to maintain an open conformation of the hydrophobic patch that binds regulatory targets. In Scott and Kekenes-Huskey (2016) we utilized principal components analysis (PCA) to characterize the predominant conformational motions that distinguish the S100A1 apo from holo states. The largest mode was referred to as Principal Component 1 in that study and correlated with the opening and closing of the peptide-binding pocket formed between H3 and H4. We therefore report here the angle between H3 and H4 as an indicator of pocket opening in Figure 3 for the WT and C85 mutants.

FIGURE 3

Figure 3. (A) H3/H4 inter-helix angle of S100A1 in closed state (PDB 2L0P) and open state (PDB 2LP3). The angle is defined between Cα atoms of C85/A53 and the COM of residues L61 to D62 and Q72 to E73 (shown as a black dotted sphere). (B) Average H3/H4 angle derived from MD simulations. Chain A and chain B are represented by deep and light shaded bars, respectively. The error bars are standard deviations. (C) Superimpose of MD-sampled half-saturated and fully-saturated S100A1 mutants with the NMR structure of TRTK12 bound S100A1 (PDB 2KBM). Ca²⁺ ions are represented as yellow spheres.

The average angle values are shown in Figure 3B. For the WT structures, the apo state is generally closed and half-/fully-saturated states stay open, as the average angle values are close to values measured from the apo (2.73 × 10¹°) and holo (5.43 × 10¹°) NMR structures, respectively. We also note that peptide binding at fully-saturated state has a negligible effect on pocket opening, as indicated by the comparable angle values in the WT peptide-bound case and in the WT holo case. C85R is slightly more open than WT in the apo state, as both chains have ~ 7° larger angles than WT, while C85E has comparable angle values as WT. C85R when half-saturated is asymmetric with one chain more open and the other is more closed. C85E has comparable angle values as WT in the half-saturated state. All cases maintain the open pocket in the fully saturated state. As shown in Figure 3C, we selected the most-probable simulated structures of C85 variants (see section S1.1 in Supplementary Material for details of determining these structures) and compared them with the NMR structure of S100A1 (PDB 2KBM) in which a 12-residue TRTK12 peptide was bound at the hydrophobic patch. In both the half- and fully-saturated states, the two C85 variants are able to accommodate the TRTK12 peptide as they have no structural clashes with the TKTR12 peptide. In other words, the variants can bind targets, despite only having one bound Ca²⁺. While the two C85 variants have comparable degrees of openness with that of the WT structure, the C85R in the half-saturated state has a more opened hydrophobic patch relative to WT. We attribute this in part to the onset of H3/H4 closing reported for the WT half-saturated case in Scott and Kekenes-Huskey (2016). Overall, our simulations indicate that the two C85 variants were as good as, if not better, than WT at assuming a peptide-compatible configuration, which thereby could facilitate peptide binding. We quantify this facilitation in section 3.4 via potential of mean force calculations.

3.3. H4 Terminal Helicity

The C-terminal region of H4 (residues C85 to S93) plays a vital role in the Ca²⁺-dependent activation of S100A1 targets (Landar, 1998), as experiments show that either the deletion of this terminal region or mutations of three aromatic amino acids in this region to alanine diminish Ca²⁺-dependent activation of S100A1 targets. Further, the highly divergent C-terminal region of H4 among S100 family proteins has been suggested to account for the selectivity of targets binding (Santamaria-Kisiel et al., 2006). For example, the three-fold TRTK12 affinity difference between S100A1 and S100B was attributed in part to different residues in the H4 C-terminus with regard to S100A1 (Wright et al., 2009a). Specifically, TRTK12 in S100A1 assumes a different orientation than that of S100B due to the different residues in the H4 C terminus, resulting in less optimized hydrophobic interactions between S100A1 and TRTK12 peptide.

Experimental data indicate that in the apo state, residues from N87 to W90 in H4 are in the equilibrium between helix and random coil configurations (Nowakowski et al., 2011, 2013). Upon Ca²⁺-binding H4 adopts a complete helix. This H4 helix extension in the holo state is thought to predispose the hydrophobic residues (i.e., C85 and F88) to interact with hydrophobic residues of the target peptide (Wright et al., 2009a), however, to our knowledge, its impact on Ca²⁺ binding has not been investigated. Since H4 in S100A1 undergoes appreciable rearrangements upon Ca²⁺ binding, we speculated that such interactions might counter the free energy gain upon binding Ca²⁺, which would reduce its apparent affinity relative to a system lacking H4 self-interactions. To investigate if our C85 mutations similarly affect the α-helicity in the H4 C-terminus, we measured the α-helix probability of residues C85 to S93 in H4 (Figure 4). In the apo state, WT S100A1 residues from F88 to N92 have considerably smaller α-helix probability than that of the holo state. The reduced H4 helical content in the apo state may be caused by two contacts that hinder α-helix formation: (1) contact between the H4 C-terminus with pEF loop from the other subunit (Figure S8) and (2) contact between the H4 C-terminus with H2-H3 linker (Figure S5). These two contacts are attenuated in the holo state. We therefore mutated W90 to A90 in the apo state to disrupt the first contact and thereby permit the C-terminus of H4 to adopt a folded α-helix. However, results show that W90A has comparable α-helix probabilities to the WT, implying that the first contact does not affect the helix formation in the H4 C-terminus. Thus it is likely that the contacts between H4 and H2-H3 linker in the apo state hinder the helix formation in the H4 C-terminus. Indeed, we show in Figure S6 that in WT apo S100A1, F44 from the H2-H3 linker region maintains contacts with L81 and F88 from the H4 C-terminal region, as F44 was sandwiched by the two hydrophobic residues.

FIGURE 4

Figure 4. (A–C) α-helix probability of H4 residues C85 to S93. Chain A and chain B are represented by dark and light shaded bars, respectively. (D) Difference of per residue energy contribution to the MMGBSA-calculated interaction energy between the two monomers of S100A1. The error bars represent standard deviations.

We had expected that our C85 variants would disrupt native H4 interactions and facilitate α-helix formation. For instance, in the β-mercaptoethanol-modified C85 apo-S100A1, NOE data show that residues from C85 to F89 folded into α-helix, possibly due to the hydrophobic interactions between the newly introduced β-mercaptoethanol and aromatic residues F88/F89. However, we found that both C85E and C85R variants have reduced α-helix probability in the apo state with C85E presenting a larger degree of reduction. We speculated that charged side chains of C85E/R prevent the favorable hydrophobic interactions within the C-terminus that is needed for α-helix formation. We anticipated that the large desolvation energy penalties of these solvent-exposed charged residues [as indicated by solvent-accessible surface area (SASA) in Figure S4] hinder the formation of α-helix.

In the half-saturated WT and C85E configurations, the H4 C-terminus is unfolded to similar degrees relative to the holo state, although the latter variant features one partially-folded helix. This unfolding we believe stems from significant contacts between H4 and the H2-H3 linker in the half-saturated state (Figure S5). Meanwhile, H4 for the C85R half-saturated state exhibits folded content comparable to the holo state, which we attribute to reduced H4 and H2-H3 linker interactions. In the fully-saturated state, all cases maintain high helical content in the H4 C-terminal region. Further, binding of the target peptide has a negligible effect on helical content; this suggests that the free energy gain of helix formation likely occurs during ion binding and not thereafter. Lastly, it is interesting to note (see Figure 4D) that the apo state C85R mutant presents energetically-unfavorable interactions between its monomers that disrupt H4/linker interactions and could thereby facilitate H4 formation.

3.4. Thermodynamics of TRTK12 Peptide Binding to S100A1

To determine if the C85 variants thermodynamically facilitate target binding to S001A1 at physiological Ca²⁺ concentrations, under which only the cEF has Ca²⁺ bound, we performed potential of mean force (PMF) calculations to characterize the free energy profile of TRTK12 peptide dissociation from half-saturated S100A1. Specifically, starting from the MD-simulated most probable half-saturated structures that are compatible for target-peptide binding, a TRTK12 peptide was placed at the hydrophobic patch. After energy minimization, the TRTK12 peptide was pulled away along the reaction coordinate (RC) defined as the distance between the COM of peptide Cα atoms to the COM of Cα atoms in the H3 C-terminus (residues K30-T39) and H4 N-terminus (residues E73-A84). The PMFs along the dissociation process are shown in Figure 5. To validate the accuracy of PMF calculations, we first compared the PMF of TRTK12 peptide dissociation from WT fully-saturated S100A1 (dashed line in Figure 5A). The experimental and calculated binding free energies were ΔG_expt. ≈ −6.5 (Wright et al., 2009a) and ΔG_calc. = −9.5 kcal/mol, respectively, which are in reasonable qualitative agreement.

FIGURE 5

Figure 5. PMF of TRTK12 dissociation from half-saturated S100A1. The PMF of fully-saturated WT case was also calculated to compare with experimentally measured TRTK12 binding free energy adapted from Wright et al. (2009a). The reaction coordinate was defined as the distance between the center of mass (COM) of peptide Cα atoms to the COM of Cα atoms in H3 C-terminus (residues K30-T39) and H4 N-terminus (residues E73-A84). The error bars were obtained via the Monte Carlo Bootstrap Error Analysis function in the WHAM program (with 30 MC trials).

For all half-saturated C85 variants, we found that TRTK12 peptide exhibits a minimum in the free energy profile at RC≈14 Å, similar to the WT. However, the binding free energies are considerably more favorable. In the half-saturated state, the two C85 variants have binding free energy values of −12.3 and −13.5 vs. that of −7.6 kcal/mol for the WT type. This thermodynamic advantage in the C85 variants is likely due to the stronger hydrophobic interactions between S100A1 and TRTK12 than that of the WT. This is evidenced by the hydrogen bonding and contact map analyses in Figure S12. In the bound state, the numbers of hydrogen bonds between TRTK12 peptide and S100A1 are comparable for the WT and C85 variants. However, the contacts between TRTK12 and the H2-H3 linker in C85 variants outnumbered those of WT. These data suggest that introducing glutamic acid or arginine at C85 increases peptide binding affinity. We additionally investigated the gating kinetics of the H3/H4 patch and found no significant difference between WT and the C85 variants (see section S3.2 in the Supplementary Materials).

3.5. Relating S100A1's C85 Modifications to Physiological Function: Combining S100A1-Mediated Actin Passive Model With PMF Calculations

To exemplify the potential impact of improved peptide binding on S100A1's physiological function, we relate these functions to its capacity to bind the PEVK repeats in titin in the N2B isoform (Granzier et al., 2010) as a model system for S100A1 target regulation. The elastic PEVK domain of titin contains 70% of proline, glutamate, valine, and lysine residues and exists in three conformational states: polyproline II (PPII) helix, b-turn, and random coil (Labeit and Kolmerer, 1995; Ma and Wang, 2003). The PEVK domain consists of a repeating motif of 28 residues on average with no long-range cooperativeness between motifs (Gutierrez-Cruz et al., 2001). The extension of PEVK is believed to contribute to titin's elasticity. A competition assay demonstrated that isolated PEVK fragments washed into skinned myocyte preparations reduced passive tension over physiological sarcomere lengths (Yamasaki et al., 2011). Further, S100A1 was shown to reduce F-actin-bound I27-PEVK-I84 in a dose-dependent manner, with higher rates of reduction under conditions of elevated (0.1 mM) Ca²⁺ (Yamasaki et al., 2011).

To relate changes in peptide binding due to mutations to the potential regulation of target, we proposed a competitive S100A1/actin-binding scheme in Figure 6A. This includes an equilibrium between actin interacting with titin's PEVK segment, while S100A1 also can interact with titin's PEVK segment. The former actin-titin interaction is proposed to delay the filament sliding, while the presence of S100A1 will disrupt the actin-titin interaction and modulate muscle contraction. This model is mechanistically consistent with trends reported for S100A1-dependent reductions of actin/PEVK (AP) binding at elevated Ca²⁺(Yamasaki et al., 2011). Namely, AP binding was shown to reduce as S100A1 increased, with greater efficacy demonstrated at 0.1 mM Ca²⁺. The governing equations are listed in section S1.3 in Supplementary Material.

FIGURE 6

Figure 6. (A) Two possible interaction models between S100A1 and actin-PEVK (AP) complex. In first model, half-saturated S100A1 does not bind to AP while in second model all apo S100A1, half-saturated S100A1 (Ca.S100A1) and fully-saturated S100A1 (2Ca.S100A1) bind to AP. (B) Fitting to experimental data from Yamasaki et al. (2011). (C) Comparison of predicted actin-PEVK binding as a function of wild-type (black) and glutathionylated (red) S100A1 at pCa = 7 (solid) and pCa = 4,5 (dashed). Assumed C85 glutathionylation leads to ~48 fold enhancement for S100A1's affinity to AP complex based on the PMF calculations.

We first fitted the model to experimental data from Yamasaki et al. (2011) to obtain the dissociation constants of S100A1 to actin/PEVK (AP) complex with and without Ca²⁺. As shown in Figure 6A, in the absence and presence of 0.1 mM Ca²⁺, the fitted dissociation constants of apo-S100A1, half-saturated S100A1 and fully-saturated S100A1 to AP are K_d1 = 0.52, K_d2 = 0.13 and $K_{SBd{2}^{\prime }}=0.03$ μM, respectively. Surprisingly, the corresponding binding free energies of K_d2 and $K_{SBd{2}^{\prime }}$ are −9.40 and −10.27 kcal/mol, respectively. These two values are close to calculated binding free energy −9.50 kcal/mol of TRTK12 peptide to WT fully-saturated S100A1 (Figure 5). This agreement indicates that it is reasonable to use the TRTK12 peptide as a mimic of the PEKV fragment to study the binding affinity of PEVK to WT/mutant S100A1 systems. We anticipate that the similar binding arises from PEVK generally showing highly amphoteric charge distributions and modestly positive net charges (Forbes et al., 2005), similar to TRTK12 (Wafer et al., 2013).

In general, we have found that end-point methods for computing free energy differences as well as more rigorous approaches including potentials of mean force seem to fare well in rank ordering cases by energy. However, for a variety of reasons not limited to force field inaccuracies, difficulties in estimating entropic contributions, and finite sampling of protein configurations during limited molecular dynamics simulation times, we find that the simulation approaches for the systems we have considered are unable to accurately predict the absolute energy differences between cases that have been experimentally characterized. Hence to map our PMF results (calculated K_Ds) to the experimentally-measured K_Ds, we introduce a scaling parameter λ. The value of λ was calculated based on the alignment between experimental and calculated dissociation constants (via PMF) for TRTK12 peptide binding to WT fully-saturated S100A1:

$\begin{array}{ll}{K}_{D,expt.}=\lambda K_{D,calc.}& (2)\end{array}$

$\begin{array}{ll}\lambda =e^{(\Delta G_{expt.}-\Delta G_{calc.})/RT}=157.436& (3)\end{array}$

where we used ΔG_expt. = −6.5 and ΔG_calc. = −9.5 kcal/mol, respectively. We then use the average of C85E and C85R as the binding free energy of PEVK to S100A1 after C85-glutathionylation (ΔG_glu. = −12.9 kcal/mol). Compared to the WT, the dissociation constant of S100A1 with C85-glutathionylation would thus be reduced as:

$\begin{array}{ll}{K}_{d,glu.}=\lambda K_d{e}^{(\Delta G_{glu.}-\Delta G_{WT})/RT}=0.021K_d& (4)\end{array}$

where ΔG_WT = −7.6 kcal/mol and λ is defined in Equation (2) and K_d refers to K_d2 and $K_{SBd{2}^{\prime }}$ in Figure 6B. When used in Figure 6A, the C85-glutathionylated S100A protein more rapidly reduces the PEVK-actin interaction at various increasing, physiological Ca²⁺ concentrations, as shown in Figure 6C.

4. Conclusions

Previous studies established that post-translational modifications of C85 at H4 of S100A1 increase S100A1's Ca²⁺ sensitivity of activation (Zhukov et al., 2008; Živković et al., 2012; Nowakowski et al., 2013). By using computational methods in this study, including molecular dynamics simulations and thermodynamic models of binding, we determined potential mechanisms governing how C85 modified S100A1 can bind Ca²⁺-dependent targets at sub-saturating Ca²⁺. Specifically, we used two variants (C85E/R) that have bulky side chains as steric surrogates of post-translational modifications at the C85 position in S100A1. Our data show that C85E/R variants have similar structural effects as post-translational modifications on displacing the C-terminus of H3 from H4 in the apo state of S100A1. We further focused on structures bound with one equivalent of Ca²⁺(half-saturated) in its pEF domain, which are likely to predominate (Scott and Kekenes-Huskey, 2016) at physiological Ca²⁺ concentrations (100 to 1000 nM) found within cells (Berridge et al., 2000). We found for the C85E/R variants relative to the WT that (1) the mutations disrupt the half-saturated structures by increasing the solvent exposure of its target binding domain (the hydrophobic patch) found between H3 and H4, (2) yield stronger TRTK12 binding in the half-saturated variants vs. WT and (3) for the half-saturated configuration, the H4 C-terminus in the two variants have greater alpha helical character than the WT and are consistent with levels exhibited in the target-bound configuration. Ultimately, these findings are suggestive of cysteine-targeted post-translational modifications priming S100A1 for target regulation within physiological ranges of Ca²⁺. Importantly, the computational studies further support the notion that S100A1 toggles its Ca²⁺-dependent regulation of downfield targets in response to Cys modification, as is common in redox pathways such as those using glutathionylation (Zhukova et al., 2004).

A variety of studies implicate S100Al in regulating proteins that mediate Ca²⁺ signaling or alter their mechanical properties in response to Ca²⁺. S100A1's regulatory roles are most apparent in the heart in which the protein is predominantly expressed (Kato and Kimura, 1985). In cardiac tissue, S100A1 has inotropic effects on Ca²⁺ handling, that is, it helps increase the generation of contractile force (Kraus et al., 2009). This is accomplished through priming sarcoplasmic reticulum (SR) Ca²⁺ concentration and release (Kettlewell et al., 2005), through interactions with targets including SERCA, RyR, the L-type Ca²⁺ channel, and the sodium calcium exchanger (Rohde et al., 2010; Völkers et al., 2010; Duarte-Costa et al., 2014). While S100A1 appears to dualy regulate RyR at both diastolic (~100 nM; Berridge et al., 2000) and systolic (>1 μM; Berridge et al., 2000; Yamaguchi et al., 2011), the WT S100A1 likely acts on its Ca²⁺-dependent targets only at saturating Ca²⁺ conditions (Nowakowski et al., 2013), under which both EF hands of the protein are bound with Ca²⁺. While it has been speculated that S100A1 could modulate target proteins at systolic Ca²⁺ levels such as in the case of S100A1/titin interactions facilitating myocyte contraction (Granzier et al., 2010; Yamasaki et al., 2011) we would expect that the low binding affinity of its pseudo EF hand (K_D ~250–16,700 μM) and reduced ability to maintain an open, target-peptide compatible binding site (Scott and Kekenes-Huskey, 2016) would not be sufficient for significant regulation of the intended targets. This raises the question of how S100A1 modulates its targets in vivo where cytosolic Ca²⁺ concentrations are generally far below the protein's K_Ds of Ca²⁺.

Post-translational modifications of S100A1 likely explain this enigma. Intriguingly, Cys85 is a redox sensitive residue presenting a variety of oxidizing functional groups (Nowakowski et al., 2013). Previous studies have indicated that S100A1 species with mercaptoethanol and glutathione generally increased the apparent Ca²⁺ affinity relative to the WT by up to four orders of magnitude (Goch et al., 2005). The most apparent rationale for the enhanced Ca²⁺ affinity from our simulations was that the mutations we considered disrupt the apo state H3/H4 folding, namely by compromising H4 interactions within a given monomer and its opposing monomer of the dimeric state. Since those interactions are dramatically reduced in the Ca²⁺ bound state relative to the apo, we speculate that their weakening in the apo state reduces the thermodynamic penalty they could impose upon Ca²⁺ binding observed for the WT. Nonetheless, the half-saturated S100A1 state appeared to demonstrate thermodynamically favorable, albeit weaker, TRKT12 peptide binding, which suggests that the protein may have a modest ability to bind targets in its wild-type form. We propose therefore that post-translational modification of the H3/H4 interface may constitute a general mechanism for controlling Ca²⁺-dependent activation of protein/protein interactions in S100 families, given the prevalence of H3/H4 binding patches featured in protein-protein interactions (Zimmer et al., 2003).

If S100A1 demonstrates increased activity following post-translational modifications, when would such modifications be expected in vivo? It is reasonable to assume that S100A1 oxidation would be most significant during conditions of enhanced reactive oxygen species (ROS) signaling. ROS are particularly prevalent during metabolic stress, ischemic-reperfusion, and physiological reactive oxygen species-based signaling (Jeong et al., 2012). In fact, glutathionylation of protein targets including S100A1 is a vital component of cardiovascular ROS signaling (Pastore and Piemonte, 2013) Physiological conditions including exercise for instance demonstrate significant cardiac RyR gluthionylation that enhance its activity to compensate for increased demand (Sánchez et al., 2008). Analogous modifications of the SERCA Ca²⁺ pump promote enhanced Ca²⁺ uptake and smooth muscle relaxation (Pastore and Piemonte, 2013). It is intriguing that reducing agents mitigate these effects (Volkers, unpublished from Völkers et al., 2010), which is expected if glutathionylations are prevalent. Since both RyR and SERCA are S100A1-dependent targets and S100A1 itself is subject to glutathionylation, this suggests redox regulation of inotropy may be controlled both directly and indirectly by glutathione modifications. Hence, in physiological systems, ROS signaling, especially as mediated by glutathionylation, might prime inotropic effects relative to basal or reduced conditions (Nikolaienko et al., 2018), through augmenting S100A1 stimulation of its targets.

4.1. Limitations

Our study includes several limitations of note. One, in this study we have used glutamic acid or arginine substitutions to probe how introducing larger, polar amino acids into the redox sensitive C85 site impact S100A1 function and peptide binding. Our choice for these variants was based on observations in several S100A1 structures with either redox modifications or amino acid substitution (C85M) that have been deposited in the Protein Data Bank. These structures exhibited more open-like character than the WT, which suggests that the opening behavior may be more sensitive to the size of the introduced side group than its specific chemical properties. Nonetheless, simulations that include the specific functional group in question, such as the glutathione group investigated in Goch et al. (2005), are likely to provide more fine detail into the mechanism of its effect on S100A1. We found that substitution of glutamic acid or arginine at the C85 position as a steric surrogate of post-translational modifications yielded half-saturated S100A1 structures that more closely overlapped with the peptide-bound (TKT12) protein.

Second, we also note that to map the PMF-predicted relative energy differences between S100A1 states to their experimentally-measured absolute differences we introduced a scaling parameter λ. In our experience, we have found computational free energy methods perform reasonably well at rank-ordering system configurations according to their experimentally-reported values, but absolute differences in energy have been less successful. We anticipate the further improvements in force field parameterizations and sampling techniques could potentially better align the relative energy predictions with absolute differences and thus obviate the scaling term we used here. Additionally, although the S100A1 open state appears to be necessary to bind target peptides as part of its regulatory function, it may be of interest to examine whether the half-saturated variants we considered significantly sample states resembling the apo (closed) configuration, as we demonstrated for the WT in Scott and Kekenes-Huskey (2016). This could be accomplished using the biased sampling technique described in that study as well as other enhanced sampling techniques including accelerated-MD (Hamelberg et al., 2004).

Lastly, in order to quantitatively link the effects of potential post-translation modifications of S100A1 to a physiological process, we examined its putative binding to titin PEVK fragments discussed in Yamasaki et al. (2011). In that study, a competition assay was used as a proxy for measuring the protein's impact on passive tension in muscle fibrils (Yamasaki et al., 2011). Passive tension is described as the force when muscle cells are stretched beyond their resting length, independent of Ca²⁺. Actin/titin interactions have been suggested as an important mechanism for controlling myofilament passive tension (Granato et al., 2010). Competitive binding assays conducted by Yamasaki et al. (2011) demonstrated that S100A1 interferes with actin/titin interactions by competitively binding the titin PEVK domain. However, it is important to note that myriad factors contribute to passive tension, including tubulin and collagen, nebulin/PEVK interactions (Gutierrez-Cruz et al., 2001) or even be recapitulated without changes in titin stiffness by modulating bound myosin/actin populations (Campbell, 2009). Further since titin/actin-dependent effects are more evident for the N2B isoform, while the N2BA tends to predominate in humans (Granzier et al., 2010), the significance of S100A1 in modulating titin varies across species.

Data Availability Statement

All datasets generated for this study are included in the article/Supplementary Material.

Author Contributions

BS and PK-H designed the simulations and wrote the manuscript. BS performed the simulations and analyzed the data. All authors contributed to the article and approved the submitted version.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

Research reported in this publication, release was supported by the Maximizing Investigators' Research Award (MIRA) (R35) from the National Institute of General Medical Sciences (NIGMS) of the National Institutes of Health (NIH) under grant number R35GM124977. This work used the Extreme Science and Engineering Discovery Environment (XSEDE) (Towns et al., 2014), which is supported by the National Science Foundation under grant ACI-1548562.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmolb.2020.00077/full#supplementary-material

References

Berridge, M. J., Lipp, P., and Bootman, M. D. (2000). The versatility and universality of calcium signalling. Nat. Rev. Mol. Cell Biol. 1, 11–21. doi: 10.1038/35036035

PubMed Abstract | CrossRef Full Text | Google Scholar

Campbell, K. S. (2009). Interactions between connected half-sarcomeres produce emergent mechanical behavior in a mathematical model of muscle. PLoS Comput. Biol. 5:e1000560. doi: 10.1371/journal.pcbi.1000560

PubMed Abstract | CrossRef Full Text | Google Scholar

Case, D., Darden, T., Cheatham, T. E., Simmerling, C., Wang, J., Duke, R. E., et al. (2012). Amber 12. San Francisco, CA: University of California.

Google Scholar

Case, D. A., Babin, V., Berryman, J. T., Betz, R. M., Cai, Q., Cerutti, D. S., et al. (2014). Amber 14. San Francisco, CA: University of California.

Google Scholar

Chaturvedi, N., Ahmad, K., Yadav, B. S., Lee, E. J., Sonkar, S. C., Marina, N., et al. (2020). Understanding calcium-dependent conformational changes in S100A1 protein: a combination of molecular dynamics and gene expression study in skeletal muscle. Cells 9:E181. doi: 10.3390/cells9010181

PubMed Abstract | CrossRef Full Text | Google Scholar

Clapham, D. E. (2007). Calcium signaling. Cell 131, 1047–1058. doi: 10.1016/j.cell.2007.11.028

PubMed Abstract | CrossRef Full Text | Google Scholar

Duarte-Costa, S., Castro-Ferreira, R., and Neves, J. S. (2014). S100A1: a major player in cardiovascular performance. Physiol. Res. 63, 669–681.

PubMed Abstract | Google Scholar

Forbes, J. G., Jin, A. J., Ma, K., Gutierrez-Cruz, G., Tsai, W. L., and Wang, K. (2005). Titin PEVK segment: charge-driven elasticity of the open and flexible polyampholyte. J. Muscle Res. Cell. Motil. 26, 291–301. doi: 10.1007/s10974-005-9035-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Goch, G., Vdovenko, S., Kozłowska, H., and Bierzyński, A. (2005). Affinity of S100A1 protein for calcium increases dramatically upon glutathionylation. FEBS J. 272, 2557–2565. doi: 10.1111/j.1742-4658.2005.04680.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Granato, M. A., Jorge, M., Vlugt, T. J. H., and Rodrigues, A. E. (2010). Diffusion of propane, propylene and isobutane in 13X zeolite by molecular dynamics. Chem. Eng. Sci. 65, 2656–2663. doi: 10.1016/j.ces.2009.12.044

CrossRef Full Text | Google Scholar

Granzier, H., Fukushima, H., and Chung, C. S. (2010). Titin-isoform dependence of titin-actin interaction and its regulation by S100A1/Ca2+ in skinned myocardium. J. Biomed. Biotechnol., 2010: 727239. doi: 10.1155/2010/727239

PubMed Abstract | CrossRef Full Text | Google Scholar

Grossfield, A. WHAM: The Weighted Histogram Analysis Method, Version 2.0.10. Available online at: http://membrane.urmc.rochester.edu/wordpress/?page_id=126

Google Scholar

Gutierrez-Cruz, G., Van Heerden, A. H., and Wang, K. (2001). Modular motif, structural folds and affinity profiles of the PEVK segment of human fetal skeletal muscle titin. J. Biol. Chem. 276, 7442–7449. doi: 10.1074/jbc.M008851200

PubMed Abstract | CrossRef Full Text | Google Scholar

Hamelberg, D., Mongan, J., and McCammon, J. A. (2004). Accelerated molecular dynamics: a promising and efficient simulation method for biomolecules. J. Chem. Phys. 120, 11919–11929. doi: 10.1063/1.1755656

PubMed Abstract | CrossRef Full Text | Google Scholar

Ivanenkov, V. V., Jamieson, G. A., Gruenstein, E., and Dimlich, R. V. (1995). Characterization of S-100b binding epitopes. Identification of a novel target, the actin capping protein, CapZ. J. Biol. Chem. 270, 14651–14658. doi: 10.1074/jbc.270.24.14651

PubMed Abstract | CrossRef Full Text | Google Scholar

Jeong, E. M., Liu, M., Sturdy, M., Gao, G., Varghese, S. T., Sovari, A. A., et al. (2012). Metabolic stress, reactive oxygen species, and arrhythmia. J. Mol. Cell. Cardiol. 52, 454–463. doi: 10.1016/j.yjmcc.2011.09.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Jo, S., Kim, T., Iyer, V. G., and Im, W. (2008). Charmm-gui: a web-based graphical user interface for charmm. J. Comput. Chem. 29, 1859–1865. doi: 10.1002/jcc.20945

PubMed Abstract | CrossRef Full Text | Google Scholar

Jorgensen, W. L., Chandrasekhar, J., Madura, J. D., Impey, R. W., and Klein, M. L. (1983). Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 79, 926–935. doi: 10.1063/1.445869

CrossRef Full Text | Google Scholar

Kabsch, W. W., and Sander, C. C. (1983). Dictionary of protein secondary structure: pattern recognition of hydrogen-bonded and geometrical features. Biopolymers 22, 2577–2637.

PubMed Abstract | Google Scholar

Kato, K., and Kimura, S. (1985). S100a o (αα) protein is mainly located in the heart and striated muscles. BBA Gen. Subj. 842, 146–150. doi: 10.1016/0304-4165(85)90196-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Kettlewell, S., Most, P., Currie, S., Koch, W. J., and Smith, G. L. (2005). S100A1 increases the gain of excitation-contraction coupling in isolated rabbit ventricular cardiomyocytes. J. Mol. Cell. Cardiol. 39, 900–910. doi: 10.1016/j.yjmcc.2005.06.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Kraus, C., Rohde, D., Weidenhammer, C., Qiu, G., Pleger, S. T., Voelkers, M., et al. (2009). S100A1 in cardiovascular health and disease: closing the gap between basic science and clinical therapy. J. Mol. Cell. Cardiol. 47, 445–455. doi: 10.1016/j.yjmcc.2009.06.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Labeit, S., and Kolmerer, B. (1995). Titins: giant proteins in charge of muscle ultrastructure and elasticity. Science 270, 293–296. doi: 10.1126/science.270.5234.293

PubMed Abstract | CrossRef Full Text | Google Scholar

Landar, A. (1998). S100A1 utilizes different mechanisms for interacting with calcium- dependent and calcium-independent target proteins. Biochemistry 37, 17429–17438. doi: 10.1021/bi9817921

PubMed Abstract | CrossRef Full Text | Google Scholar

Landar, A., Caddell, G., Chessher, J., and Zimmer, D. B. (1996). Identification of an S100A1/S100B target protein: phosphoglucomutase. Cell Calcium 20, 279–285. doi: 10.1016/S0143-4160(96)90033-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, M. X., Wang, X., Lindhout, D. A., Buscemi, N., Van Eyk, J. E., and Sykes, B. D. (2003). Phosphorylation and mutation of human cardiac troponin I deferentially destabilize the interaction of the functional regions of troponin I with troponin C. Biochemistry 42, 14460–14468. doi: 10.1021/bi035408y

PubMed Abstract | CrossRef Full Text | Google Scholar

Ma, K., and Wang, K. (2003). Malleable conformation of the elastic PEVK segment of titin: non-co-operative interconversion of polyproline II helix, beta-turn and unordered structures. Biochem. J. 374(Pt 3), 687–695. doi: 10.1042/bj20030702

PubMed Abstract | CrossRef Full Text | Google Scholar

MacKerell, A. D., Bashford, D., Bellott, M., Dunbrack, R. L., Evanseck, J. D., Field, M. J., et al. (1998). All-atom empirical potential for molecular modeling and dynamics studies of proteins. J. Phys. Chem. B 102, 3586–3616. doi: 10.1021/jp973084f

PubMed Abstract | CrossRef Full Text | Google Scholar

MacKerell, A. D., Feig, M., and Brooks, C. L. (2004). Improved treatment of the protein backbone in empirical force fields. J. Am. Chem. Soc. 126, 698–699. doi: 10.1021/ja036959e

PubMed Abstract | CrossRef Full Text | Google Scholar

Melville, Z., Aligholizadeh, E., McKnight, L. E., Weber, D. J., Pozharski, E., and Weber, D. J. (2017). X-ray crystal structure of human calcium-bound S100A1. Acta Crystallogr. Sect. F Struct. Biol. Commun. 73, 215–221. doi: 10.1107/S2053230X17003983

PubMed Abstract | CrossRef Full Text | Google Scholar

Nikolaienko, R., Bovo, E., and Zima, A. V. (2018). Redox dependent modifications of ryanodine receptor: basic mechanisms and implications in heart diseases. Front. Physiol. 9:1775. doi: 10.3389/fphys.2018.01775

PubMed Abstract | CrossRef Full Text | Google Scholar

Nowakowski, M., Jaremko, Ł., Jaremko, M., Zhukov, I., Belczyk, A., Bierzyński, A., et al. (2011). Solution NMR structure and dynamics of human apo-S100A1 protein. J. Struct. Biol. 174, 391–399. doi: 10.1016/j.jsb.2011.01.011

PubMed Abstract | CrossRef Full Text | Google Scholar

Nowakowski, M., Ruszczyńska-Bartnik, K., Budzińska, M., Jaremko, Ł., Jaremko, M., Zdanowski, K., et al. (2013). Impact of calcium binding and thionylation of S100A1 protein on its nuclear magnetic resonance-derived structure and backbone dynamics. Biochemistry 52, 1149–1159. doi: 10.1021/bi3015407

PubMed Abstract | CrossRef Full Text | Google Scholar

Pastore, A., and Piemonte, F. (2013). Protein glutathionylation in cardiovascular diseases. Int. J. Mol. Sci. 14, 20845–20876. doi: 10.3390/ijms141020845

PubMed Abstract | CrossRef Full Text | Google Scholar

Phillips, J. C., Braun, R., Wang, W., Gumbart, J., Tajkhorshid, E., Villa, E., et al. (2005). Scalable molecular dynamics with namd. J. Comput. Chem. 26, 1781–1802. doi: 10.1002/jcc.20289

PubMed Abstract | CrossRef Full Text | Google Scholar

Roe, D. R., and Cheatham, T. E. (2013). PTRAJ and CPPTRAJ: software for processing and analysis of molecular dynamics trajectory data. J. Chem. Theor. Comput. 9, 3084–3095. doi: 10.1021/ct400341p

PubMed Abstract | CrossRef Full Text | Google Scholar

Rohde, D., Ritterhoff, J., Voelkers, M., Katus, H. A., Parker, T. G., and Most, P. (2010). S100A1: a multifaceted therapeutic target in cardiovascular disease. J. Cardiovasc. Transl. Res. 3, 525–537. doi: 10.1007/s12265-010-9211-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Ryckaert, J.-P., Ciccotti, G., and Berendsen, H. J. (1977). Numerical integration of the cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes. J. Comput. Phys. 23, 327–341. doi: 10.1016/0021-9991(77)90098-5

CrossRef Full Text | Google Scholar

Sánchez, G., Escobar, M., Pedrozo, Z., Macho, P., Domenech, R., Härtel, S., et al. (2008). Exercise and tachycardia increase NADPH oxidase and ryanodine receptor-2 activity: possible role in cardioprotection. Cardiovasc. Res. 77, 380–386. doi: 10.1093/cvr/cvm011

PubMed Abstract | CrossRef Full Text | Google Scholar

Santamaria-Kisiel, L., Rintala-Dempsey, A. C., and Shaw, G. S. (2006). Calcium-dependent and -independent interactions of the S100 protein family. Biochem. J. 396, 201–214. doi: 10.1042/BJ20060195

PubMed Abstract | CrossRef Full Text | Google Scholar

Scott, C. E., and Kekenes-Huskey, P. M. (2016). Molecular basis of S100A1 activation at saturating and subsaturating calcium concentrations. Biophys. J. 110, 1052–1063. doi: 10.1016/j.bpj.2015.12.040

CrossRef Full Text | Google Scholar

Towns, J., Cockerill, T., Dahan, M., Foster, I., Gaither, K., Grimshaw, A., et al. (2014). Xsede: accelerating scientific discovery. Comput. Sci. Eng. 16, 62–74. doi: 10.1109/MCSE.2014.80

CrossRef Full Text | Google Scholar

Völkers, M., Rohde, D., Goodman, C., and Most, P. (2010). S100A1: a regulator of striated muscle sarcoplasmic reticulum Ca²⁺ handling, sarcomeric, and mitochondrial function. J. Biomed. Biotechnol. 2010:178614. doi: 10.1155/2010/178614

PubMed Abstract | CrossRef Full Text | Google Scholar

Wafer, L. N., Tzul, F. O., Pandharipande, P. P., and Makhatadze, G. I. (2013). Novel interactions of the TRTK12 peptide with S100 protein family members: Specificity and thermodynamic characterization. Biochemistry 52, 5844–5856. doi: 10.1021/bi400788s

PubMed Abstract | CrossRef Full Text | Google Scholar

Wright, N. T., Cannon, B. R., Wilder, P. T., Morgan, M. T., Varney, K. M., Zimmer, D. B., et al. (2009a). Solution structure of S100A1 bound to the CapZ peptide (TRTK12). J. Mol. Biol. 386, 1265–1277. doi: 10.1016/j.jmb.2009.01.022

PubMed Abstract | CrossRef Full Text | Google Scholar

Wright, N. T., Cannon, B. R., Zimmer, D. B., and Weber, D. J. (2009b). S100A1: structure, function, and therapeutic potential. Curr. Chem. Biol. 3, 138–145. doi: 10.2174/2212796810903020138

PubMed Abstract | CrossRef Full Text | Google Scholar

Wright, N. T., Prosser, B. L., Varney, K. M., Zimmer, D. B., Schneider, M. F., and Weber, D. J. (2008). S100A1 and calmodulin compete for the same binding site on ryanodine receptor. J. Biol. Chem. 283, 26676–26683. doi: 10.1074/jbc.M804432200

PubMed Abstract | CrossRef Full Text | Google Scholar

Wright, N. T., Varney, K. M., Ellis, K. C., Markowitz, J., Gitti, R. K., Zimmer, D. B., et al. (2005). The three-dimensional solution structure of Ca2+-bound S100A1 as determined by NMR spectroscopy. J. Mol. Biol. 353, 410–426. doi: 10.1016/j.jmb.2005.08.027

PubMed Abstract | CrossRef Full Text | Google Scholar

Yamaguchi, N., Prosser, B. L., Ghassemi, F., Xu, L., Pasek, D. A., Eu, J. P., et al. (2011). Modulation of sarcoplasmic reticulum Ca2+ release in skeletal muscle expressing ryanodine receptor impaired in regulation by calmodulin and S100A1. Am. J. Physiol. Cell Physiol. 300, C998–C1012. doi: 10.1152/ajpcell.00370.2010

PubMed Abstract | CrossRef Full Text | Google Scholar

Yamasaki, R., Berri, M., Wu, Y., Trombitá, K., Mcnabb, M., Kellermayer, M. S. Z., et al. (2011) Titin–actin interaction in mouse myocardium: passive tension modulation its regulation by calcium/S100A1. Biophys. J. 81, 2297–2313. doi: 10.1016/S0006-3495(01)75876-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhu, F., and Hummer, G. (2012). Convergence and error estimation in free energy calculations using the weighted histogram analysis method. J. Comput. Chem. 33, 453–465. doi: 10.1002/jcc.21989

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhukov, I., Ejchart, A., Bierzynski, A., and Bierzyński, A. (2008). Structural and motional changes induced in apo-S100A1 protein by the disulfide formation between its Cys 85 residue and β-mercaptoethanol. Biochemistry 47, 640–650. doi: 10.1021/bi701762v

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhukova, L., Zhukov, I., Bal, W., and Wyslouch-Cieszynska, A. (2004). Redox modifications of the C-terminal cysteine residue cause structural changes in S100A1 and S100B proteins. Biochim. Biophys. Acta Mol. Cell Res. 1742, 191–201. doi: 10.1016/j.bbamcr.2004.10.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Zimmer, D. B., Wright Sadosky, P., and Weber, D. J. (2003). Molecular mechanisms of S100-target protein interactions. Microsc. Res. Tech. 60, 552–559.

PubMed Abstract | Google Scholar

Živković, M. L., Zarȩba-Kozioł, M., Zhukova, L., Poznański, J., Zhukov, I., and Wysłouch-Cieszyńska, A. (2012). Post-translational S-nitrosylation is an endogenous factor fine tuning the properties of human S100A1 protein. J. Biol. Chem. 287, 40457–40470. doi: 10.1074/jbc.M112.41839228

PubMed Abstract | CrossRef Full Text | Google Scholar