ORIGINAL RESEARCH article
Thermodynamic Origin of Differential Excipient-Lysozyme Interactions
- 1Manchester Institute of Biotechnology, The University of Manchester, Manchester, United Kingdom
- 2Department of Chemistry, The University of Manchester, Manchester, United Kingdom
- 3Departments of Chemical Engineering and Analytical Science, The University of Manchester, Manchester, United Kingdom
- 4Division of Molecular and Cellular Function, School of Biological Sciences, Faculty of Biology, Medicine and Health, The University of Manchester, Manchester, United Kingdom
- 5Faculty of Medicine and Health, Sydney Medical School, The University of Sydney, Sydney, NSW, Australia
Understanding the intricate interplay of interactions between proteins, excipients, ions and water is important to achieve the effective purification and stable formulation of protein therapeutics. The free energy of lysozyme interacting with two kinds of polyanionic excipients, citrate and tripolyphosphate, together with sodium chloride and TRIS-buffer, are analysed in multiple-walker metadynamics simulations to understand why tripolyphosphate causes lysozyme to precipitate but citrate does not. The resulting multiscale decomposition of energy and entropy components for water, sodium chloride, excipients and lysozyme reveals that lysozyme is more stabilised by the interaction of tripolyphosphate with basic residues. This is accompanied by more sodium ions being released into solution from tripolyphosphate than for citrate, whilst the latter instead has more water molecules released into solution. Even though lysozyme aggregation is not directly probed in this study, these different mechanisms are suspected to drive the cross-linking between lysozyme molecules with vacant basic residues, ultimately leading to precipitation.
Protein therapeutics are increasingly being developed in the biopharmaceutical industry to combat a wide range of diseases (Dimitrov, 2012; Faber and Whitehead, 2019). Compared with traditional small drug molecules, the large and complex structures and marginal stability of biomolecules make necessary the development of sophisticated formulations which can stabilise the structure of such therapeutics to ensure their safe administration and efficacy (Wang, 2015). Without the correct formulation conditions, proteins are prone to aggregation, precipitation or phase separation (Moussa et al., 2016). The modulation of protein-protein interactions (PPIs) in formulations is commonly achieved by the addition of small molecules, termed excipients. However, a lack of understanding of how excipients operate is hampering further development because such systems comprise multiple, transient, weak interactions between proteins, excipients, other ions and water in solution (Falconer, 2019). This is an extension of the difficulties in understanding aqueous electrolytes, where behaviour even of simple salt solutions is not well explained, particularly at higher concentrations (Collins, 1997). It is therefore important to develop new strategies to identify the effects of excipients on protein stability to improve therapeutic formulations in the biopharmaceutical industry.
A particularly intriguing phenomenon relating to the mechanism of excipients operation was recently revealed for the case of the protein lysozyme interacting with two excipients, namely tripolyphosphate (TPP) and citrate (CIT). Lysozyme is a small, stable protein that is widely studied and known to remain largely folded in both experiments and simulations and so does not require extensive sampling of protein tertiary structural changes. Experimental data, as illustrated schematically in Figure 1, has shown how TPP excipients cause lysozyme to precipitate out of solution and then resolubilise at higher TPP concentration, but similarly charged CIT excipients do not cause lysozyme precipitation at any concentration (Bye and Curtis, 2019).
FIGURE 1. Lysozyme precipitation at intermediate tripolyphosphate (TPP) concentrations and re-solubilisation at higher concentration.
The reason for these trends is not clearly understood, but it is hypothesised that TPP cross-links lysozyme due to more polarised, multi-directional P-O- bonds which form stronger interactions with basic residues. Conversely, CIT may have sterically hindered hydrogen bonding with donors on the protein surface, which prevents cross-linking between multiple lysozyme molecules (Bye and Curtis, 2019). In this work we concentrate on the differences between TPP and CIT at the excipient concentration where TPP causes lysozyme precipitation but CIT does not.
Computational methods can contribute much to understanding protein-excipient binding, given the high level of detail that they provide. However, conventional methods for molecular binding are less suitable for such problems because they tend to be aimed at systems in which there are specific binding sites between a substrate and its target, owing to the long-standing influence of structure-based drug design. Rather, it is typically the case for excipients that they bind non-specifically to protein surfaces and in some cases not even directly to proteins at all (Zalar et al., 2020). This diminishes the usefulness of computational techniques such as docking studies and alchemical free-energy methods in favour of ensemble-based methods such as molecular dynamics (MD) or Monte Carlo (MC) simulations that are able to sample the wide range of relevant configurations (Schames et al., 2004; Spyrakis et al., 2015; Yu et al., 2018; Barata et al., 2016). Given the large number of possible configurations that need to be sampled, more efficient enhanced-sampling simulation methods are required. Most enhanced sampling methods for proteins have been conducted to explore protein folding, allostery or protein-ligand binding (Singh et al., 2017; Du et al., 2017; Barducci et al., 2013; Verona et al., 2017). Studies of multiple ligands interacting with a protein surface have shown how increased sampling produces binding affinities comparable to experiment (Troussicot et al., 2015) as well as how co-solutes can affect dissociation of proteins (Banerjee et al., 2019). The use of multiple simulations starting from different poses of two proteins has been shown to be effective in reproducing native protein-protein association structures and in the prediction of new bound configurations (Plattner et al., 2017). Such methods may yield the extent of excipient binding via brute-force probabilities between bound and unbound but are unable to explain the calculated stability for different excipients. Achieving this requires more detailed strategies such as determining how all of the molecules in a system contribute to the total Gibbs free energy. However, most such studies on proteins have focused on the contribution of water molecules or the protein itself (Tarek and Tobias, 2000; Gerogiokas et al., 2016; Chong and Ham, 2017).
Here we seek to understand how lysozyme is differentially stabilised by the excipients TPP and CIT to help explain its aggregation behaviour by applying a free energy method that calculates energy and entropy directly from a simulation called EE-MCC (Energy-Entropy Multiscale Cell Correlation). Energy is calculated directly from the system Hamiltonian by summing over per-atom energies. Entropy is calculated for all molecules in the system from forces and coordinates at multiple length scales using MCC, which has been applied to liquids (Higham et al., 2018; Ali et al., 2019), chemical reactions (Ali et al., 2020), and proteins (Chakravorty et al., 2020). The three length scales employed here from smallest to largest are 1) water and monatomic ions, 2) excipients and residues, and 3) the whole protein, which are classified here as united-atom, monomer, and polymer levels, respectively. We examine the Gibbs free energy for mixing a lysozyme dimer with TRIS buffer and counterions with five excipient molecules, either TPP or CIT, together with counterions, corresponding to the concentration at which TPP-induced aggregation is detected experimentally. Protein-excipient configurations are sampled using metadynamics and multiple-walker simulations to allow the exploration of more transient interactions than would be possible using a conventional molecular dynamics simulation. Although, the large phase space required to sample solute interactions remains a bottle-neck, we observe several important differences between CIT and TPP containing solutions that may help explain their effects on lysozyme aggregation.
2.1 Multiscale Cell Correlation Entropy Theory
In MCC, entropy S is calculated in a multiscale fashion in terms of cells of correlated units. The total entropy is calculated as a sum of components
In this equation, S is calculated for each kind of molecule i, at the appropriate length scales j for each molecule, in terms of translational or rotational motion l over all units at that level, and in terms of vibration or topography k for each type of motion. Vibrational entropy relates to the average size of energy wells for that unit while topographical entropy relates to the probability distribution of the energy wells. Length scales are defined at the united-atom (UA), monomer (M) and polymer (P) levels for the protein, at the UA and M levels for excipients, and the UA level for water and monatomic ions.
2.1.1 Vibrational Entropy
The entropy of internal molecular vibrations of a unit in the
where h is Planck’s constant,
The forces and torques are rotated into appropriate reference frames defined by the type of unit as shown in Supplementary Table S1. In the mean-field approximation, forces at the polymer level and torques at all levels are halved. For all but the highest length scale in each molecule, vibrational entropies associated with the six smallest frequencies in the force covariance matrix are removed to avoid double counting translation and rotation at the higher level.
2.1.2 Conformational Entropy
Translational topographical entropy at the united-atom level is otherwise known as conformational entropy. It is calculated from the probabilities
2.1.3 Orientational Entropy
The entropy arising from water molecule orientations with respect to neighbouring united atoms is a generalisation from before (Higham et al., 2018) that captures anisotropy due to the presence of solutes. Orientational entropy is calculated using
where c is the coordination-shell type of a water molecule and σ is the symmetry number of water, equal to 2. The average bias in hydrogen bonds (HBs)
Neighbouring water molecules are identified by the solute they are closest to when in the solute coordination shell and labeled as bulk otherwise. Neighbouring solutes are identified as their particular united atoms that are proximal to the water molecule. Contacts between united atoms are defined using the relative angular distance (RAD) algorithm (Higham and Henchman, 2016). HBs are defined topologically (Henchman and Irudayam, 2010; Irudayam et al., 2010; Irudayam and Henchman, 2011) based on the acceptor with the most negative
2.2 Molecular Energy
The potential (U) and kinetic (K) energies per atom i from simulation are summed to give the enthalpy H of a given unit
ignoring the negligible pressure-volume term.
2.3 Water Molecule Assessment Based on Coordination Shell Neighbours
The immediate environment of a water molecule is defined by what molecules are in its first coordination shell. This is used to assess the free energy of water molecules depending on what solutes they interact with. As there are many possible local configurations of a coordination shell, we simplify the definition of local water environments by defining five categories of water contacts with various combinations of excipients and counterion or lysozyme molecules as follows:
2.4 Change in Free Energy for Excipient Binding
The free-energy change for binding five TPP or CIT excipient molecules to the buffered protein relative to an infinitely dilute excipient with their neutralising counterions is calculated using
where X is each type of solute,
FIGURE 2. Schematic of the binding process between a set of five excipients, each with neutralising
2.5 Well-Tempered Metadynamics
Metadynamics allows for enhanced sampling compared to standard molecular dynamics (MD) simulations by forcing simulations to sample unexplored regions of free energy landscapes (Laio and Parrinello, 2002). Sampling is directed through the use of collective variables (CVs), whose choice depends on the system studied and the problem under consideration. Each CV is selected to most efficiently sample regions of interest, such as sampling dihedrals on a flexible peptide to find new and metastable conformations. To bias the system to previously unexplored regions, a biasing potential is used to add Gaussian functions to the potential, which depends on what regions of the potential energy landscape have been historically explored. This potential has the form
The CVs chosen are to efficiently fill metastable states and overcome barriers to neighbouring unexplored metastable states. To improve the convergence of metadynamics simulations, the height of the Gaussian is decreased with longer simulation time using well-tempered metadynamics (WT-MTD) (Barducci et al., 2008)
This is the ratio between the temperature of the CVs and the system temperature. WT-MTD helps with free-energy barrier crossing by simulating CVs at higher temperatures and reduces noise by gradually reducing Gaussian heights.
2.5.1 Reweighting Simulations From the Bias Potential
Because the bias potential in WT-MTD simulations is time-dependent, the weighting w for each simulation frame is calculated from the saved history-dependent bias potential applied at a given time frame t (Tiwary and Parrinello, 2015)
Statistical values (A) used to calculate the free energy of excipients and water molecules at frame t are weighted (
2.6 Selection of Collective Variables and Metadynamics Variables
CVs are selected based on how efficiently the interactions between solutes can be sampled without having to sample more expensive molecular conformations. We therefore consider the number of contacts made between molecule types in the system. Three possible contact CVs are sampled: 1) protein-protein contacts of any side-chain oxygen or nitrogen between each protein, 2) polyanion-protein contacts of any oxygen on the polyanion with any oxygen or nitrogen on residue side chains and 3) any contacts between water oxygens and any O or N atom on a protein or polyanion. We therefore efficiently sample protein-protein interactions, protein-polyanion interactions and water-solute interactions. For each CV, the rational switching function,
where a contact
The Gaussian width σ for each CV is selected based on the standard deviation observed from unbiased simulations. Contacts between solutes are set with
2.7 Multiple Walkers
Multiple walkers are multiple simulations running independently, but sharing information about already visited configurations along several CVs. Because we want to study non-specific interactions, which is characterised by many possibilities of weak or indirect interactions with the protein, there is no precise starting configuration to use. We thus consider multiple starting structures for sampling of interactions between proteins and polyanions. We use multiple walkers of 25 possible starting structures of various orientations between two proteins as shown in Figure 3. Starting poses are generated in VMD by sampling 90° orientations about the protein principal axes, giving 21 starting dimers with different protein surfaces facing each other. An additional four poses are included of the same protein surfaces facing each other but rotated about the dimer principal axes. While more starting poses and simulations would be required to achieve full protein-excipient sampling, this number of structures was chosen to resolve the different effects of the two excipients in a computationally efficient manner. A similar method to sample between multiple starting configurations is bias-exchange metadynamics, where many collective variables are sampled and exchanged between multiple simulations (Piana and Laio, 2007). Here we use multiple walkers instead of bias-exchange for a few reasons. Multiple walkers allow for simulations to start and run at different times, unlike with the bias-exchange method, where information is swapped every step between simulations. For multiple walkers, previously visited configurations over all steps and simulations are read in and biased toward unvisited configurations. Given that we use a high performance computing (HPC) cluster, our simulations can run independently without the need to wait for all 25 simulations to be simultanously running, therefore efficiently using the CPU resources. However, a disadvantage of multiple walkers is that only a few collective variables can be sampled. Guided by the previous experimental work (Bye and Curtis, 2019), we assume that sampling contacts between protein and excipients is sufficient to determine possible mechanisms for solute interactions. Lastly, because we use LAMMPS for per-atom energy, we are restricted to using multiple walkers rather than replica exchange, which is not currently supported for use between PLUMED and LAMMPS, as the latter uses temperature swaps, while PLUMED uses coordinate swaps between biased simulations.
FIGURE 3. The 25 starting poses for two lysozyme proteins for each multiple-walker simulation. First and last protein residues in the sequence are highlighted in orange and cyan, respectively. The
2.8 System Setup
Structures of systems containing two lysozyme proteins surrounded by 14,400 water molecules to represent 100 mM protein concentration are created for three different conditions of excipient: tris(hydroxymethyl)aminomethane (TRIS) buffer only, citrate (CIT) with buffer and tripolyphosphate (TPP) with buffer. The two polyanion-containing systems have five polyanions to represent a concentration of approximately 20 mM and three TRIS molecules for approximately 10 mM buffer concentration for all systems. These conditions are comparable to experiment at the concentration where TPP precipitates lysozyme, while CIT does not (Bye and Curtis, 2019). Each system is created with 25 replicates of different starting poses of lysozymes with respect to each other. CIT or TPP are randomly assigned within a 30 Å radius sphere centered within a box of approximately 80
2.9 Simulation Protocol
Each system is treated as a periodic box and initially minimized for 5,000 steepest descent steps with sander in AMBER18 (Pearlman et al., 1995). Topology and initial coordinates are converted to LAMMPS formatted files with InterMol (Shirts et al., 2017). A short 0.2 ns equilibration is conducted at NVT (number, volume, temperature) conditions to slowly heat the systems to 298 K and then equilibrated for 12 ns at constant NPT (number, pressure, temperature) conditions at 1 atm pressure with 1 fs time-steps in the LAMMPS simulation package (Plimpton, 1995). Restraints on protein Cα atoms are placed with a 500
3.1 Free Energy Change for the Formation of Each Lysozyme-Polyanion System
The total change in Gibbs energy, enthalpy and entropy for mixing each lysozyme-polyanion system from the individually separated polyanions and the protein-TRIS systems are shown in Table 1.
Also included is a decomposition over all five kinds of solute in the system, namely protein lysozyme, polyanion TPP or CIT, buffer TRIS, counterions
FIGURE 4. Change in free energy (orange line), enthalpy (grey bar) and entropy (bars in shades of blue) of (A) solutes and (B) water molecules around solutes in the TPP and CIT systems. (C) Change in free energy for CIT vs. TPP of each residue classified by type (top) and of water around each residue (bottom). The dashed line represents
When interpreting the results, we refer to
Also of note is the number of contacts between each kind of species. Table 2 gives the number of contacts for the unbound systems, namely the separated excipient with
For proteins, there is a lower energy for both polyanions due to favorable electrostatic interactions between oppositely charges molecules, in particular, involving charged patches on lysozyme as shown in Figure 5. However, protein free energy is seen to be more stable with TPP than CIT. This is due to both energy and entropy (Table 1), with contributions from all components of protein entropy (Figures 4A,B). Interestingly, both protein conformational and surrounding water orientational entropy are larger in the presence of polyanions, which suggest weaker interactions of the protein with its environment. Correspondingly, an increase in residue RMSD is observed in the CIT and TPP systems in Figure 5C. Water around proteins is destabilised in both energy and entropy in the presence of polyanions, bringing about slightly more stabilisation for CIT than TPP. Further decomposition of the change in free energies for each kind of residue according to charge, hydrophilicity and hydrophobicity and for their surrounding waters are plotted in Figure 4C. The free energy of basic residues, particularly Arg, is reduced for both polyanions but more so in the TPP system than CIT system. The free energy of acidic residues is generally more destabilised with TPP than with CIT.
FIGURE 5. A) Charged patches on lysozyme at pH 9 for two opposing orientations (top). Average percentage of residue contacts (darker regions have more contacts) of polyanions mapped onto the protein surface for the same orientations (middle and bottom). Regions with no contacts are represented as CPK structures. (B) Percentages of solute-solute interactions in TPP (left) and CIT (right) systems. C Change of per-residue RMSD compared to the protein in buffer only plotted against percentages of polyanion interactions with each residue.
3.2 Interactions Between Solutes
Given the free-energy destabilisation for TPP but stabilisation for CIT, we assess further how interactions with other solute molecules may cause these differences in free energy. In dilute solution TPP has more UA contacts with
From Figure 5A, both polyanions have similar numbers of contacts in similar regions of positively charged patches on the protein surface. Patches are generated on the protein using Poisson-Boltzmann electrostatics as described in previous work (Kalayan et al., 2020). CIT has overall weaker interactions over the protein surface, while TPP interacts with fewer regions of the protein but with higher occupancy. From Figure 5B, TPP still has a high percentage of
As to why TPP remains bound to the lysozyme surface even though it is destabilised, we observe stabilisation of residues that TPP interacts with most frequently from Figures 5C, 4C, where highest percentage occupied residues are also reduced the most in free energy. Therefore protein stabilisation may prevent TPP from being released into solution. Another possible reason for TPP remaining in the bound state may be attributed to surrounding water molecules as described in more detail next.
3.3 Specific Residue Interactions of Water Molecules
In Figure 6, water molecules are analysed based on the solutes in their first coordination shell. Total numbers of each water type are given in Table 4. By studying interactions of water with solutes, we can also infer which solutes are in close proximity to one another from the reduction in water coordination. Representation of each water type is shown in part A. In part B, distributions of water energy and entropy for each water type are shown, where water molecules interacting with just protein residues (
FIGURE 6. A) Water types (
Part C in Figure 6 shows the water environments based on their interactions with the nearest two residues, ordered by acidic, basic, uncharged polar and non-polar from left to right. More destabilised
Protein precipitation is often characterised by the propensity for ions to “salt-out” a protein. Furthermore, the ion-specificity to cause precipitation can be ranked in the Hofmeister series (Hofmeister, 1888). The mechanism by which protein precipitation occurs is generally considered to be caused by ions forming stronger interactions than proteins with surrounding solvent. It is assumed that stronger ion-solvent interactions leads to a more structured water HB network, explaining the designation of such ions as kosmotropes (Collins, 1997). It is the loss of water at the protein surface that then drives PPIs between hydrophobic regions. Suggested mechanisms for protein precipitation assume that ions do not interact preferentially with the protein surface, and instead remain fully hydrated in solution due to strong interactions with water. In the case of multivalent anions or cations, however, these ions bind onto the protein surface as shown by the change in net protein charge in experimental measurements of their zeta potentials (Matsarskaia et al., 2016; Roosen-Runge et al., 2013; Matsarskaia et al., 2018). At the concentration when the net charge of a globular protein is neutral, the propensity to precipitate is high because repulsive interactions are screened by bound ions. This phenomenon is observed for polyvalent cations (
To understand why ion specific protein precipitation occurs, the change in both the stability and contacts between solutes molecules is studied for conditions where TPP causes lysozyme to precipitate but CIT does not. Stability is assessed by the change in Gibbs free energy for mixing lysozyme with each excipient and the contributions of all molecules involved. This analysis shows several differences between TPP and CIT-containing systems. First we observe TPP destabilisation upon interaction with lysozyme, caused by the release of half of the strongly bound
The stability and contacts of water surrounding each solute reveals the stabilisation of water around TPP, but overall destabilisation around CIT. Upon binding to the protein, TPP does not lose any water-molecule contacts. Instead, water-molecule energy is stabilised but entropy is lost, therefore forming a more structured HB network around TPP and agreeing with the description of TPP as a kosmotropic ion. Thus instead of water being stabilised upon release into solution as is usually the case for PPIs, water is stabilised by forming strong interactions with protein-bound TPP when replacing released
Although water is slightly destabilised around lysozyme in the presence of TPP, it does not seem to be destabilised enough to cause precipitation via hydrophobic interactions. Instead the combination of overall solvent and protein stabilisation upon TPP binding is expected to drive binding via cross-linking with other proteins at moderate TPP concentration. In the conditions simulated here, five polyanion molecules per lysozyme pair closely represents a 20 mM concentration. There are a total of 34 basic residues on the lysozyme pair at pH 9 with a net charge of 16+, therefore giving many vacancies for TPP cross-linking. At 100 mM TPP concentration, experiment shows that lysozyme is fully resolublised into solution (Bye and Curtis, 2019) and there are approximately 29 TPP molecules per lysozyme pair. Each basic residue can be bound individually at 100 mM TPP, and therefore no cross-linking between proteins would be required. To more conclusively explain the experimental data, much longer simulations with more starting poses would be required that can explicitly account for direct lysozyme-lysozyme binding and cross-linking of lysozyme molecules by excipients at varying concentrations of each excipient. Adequately addressing such questions may require coarse-grain simulations to achieve the necessary sampling. At the same time, greater model accuracy may also be required by accounting for the possibility of variable protonation states of amino acids and excipients arising from their changing environments using constant pH simulations. This makes clear that much work remains in order to fully understand the effects of excipients on protein behaviour.
Lysozyme-excipient systems have been studied to understand how excipients differentially interact with lysozyme using a statistical-mechanics based free-energy method called Energy-Entropy Multiscale Cell Correlation that calculates the energy and entropy of each solute and hydration shell water molecules. Simulations are conducted using multiple walker metadynamics simulations followed by reweighting to give a Boltzmann distribution. We observe different contacts and thermodynamic properties when the excipient TPP interacts with lysozyme compared to the excipient CIT. A possible mechanism by which TPP precipitates lysozyme is suggested to be due to the stabilisation in free energy of both lysozyme and solvent surrounding TPP upon TPP interacting with basic residues and the release of bound
Data Availability Statement
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.
JK carried out the research; all authors designed the research; JK and RH devised the theory; JK and RH wrote the manuscript with support from RC and JW; RH supervised the project.
JK was supported by Centre of Doctoral Training in Emergent Macromolecular Therapy funded by the EPSRC under grant codes EP/L015218/1 and EP/N025105/1.
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.
The authors thank IT Services at the University of Manchester for access to the Computational Shared Facility.
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmolb.2021.689400/full#supplementary-material
Ali, H. S., Higham, J., de Visser, S. P., and Henchman, R. H. (2020). Comparison of Free-Energy Methods to Calculate the Barriers for the Nucleophilic Substitution of Alkyl Halides by Hydroxide. J. Phys. Chem. B. 124(31):6835–6842. doi:10.1021/acs.jpcb.0c02264
Barducci, A., Bonomi, M., Prakash, M. K., and Parrinello, M. (2013). Free-energy Landscape of Protein Oligomerization from Atomistic Simulations. Proc. Natl. Acad. Sci. 110, E4708–E4713. doi:10.1073/pnas.1320077110
Bonomi, M., Bussi, G., Camilloni, C., Tribello, G. A., Banáš, P., Barducci, A., et al. (2019). Promoting Transparency and Reproducibility in Enhanced Molecular Simulations. Nat. Methods 16, 670–673. doi:10.1038/s41592-019-0506-8
Dolinsky, T. J., Nielsen, J. E., McCammon, J. A., and Baker, N. A. (2004). PDB2PQR: an Automated Pipeline for the Setup of Poisson-Boltzmann Electrostatics Calculations. Nucleic Acids Res. 32, W665–W667. doi:10.1093/nar/gkh381
Du, H., Liu, Z., Jennings, R., and Qian, X. (2017). The Effects of Salt Ions on the Dynamics and Thermodynamics of Lysozyme Unfolding. Separation Sci. Technol. 52, 320–331. doi:10.1080/01496395.2016.1229336
Gerogiokas, G., Southey, M. W. Y., Mazanetz, M. P., Heifetz, A., Bodkin, M., Law, R. J., et al. (2016). Assessment of Hydration Thermodynamics at Protein Interfaces with Grid Cell Theory. J. Phys. Chem. B. 120, 10442–10452. doi:10.1021/acs.jpcb.6b07993
Hanwell, M. D., Curtis, D. E., Lonie, D. C., Vandermeersch, T., Zurek, E., and Hutchison, G. R. (2012). Avogadro: An Advanced Semantic Chemical Editor, Visualization, and Analysis Platform. J. Cheminform 4, 17. doi:10.1186/1758-2946-4-17
Higham, J., Chou, S.-Y., Gräter, F., and Henchman, R. H. (2018). Entropy of Flexible Liquids from Hierarchical Force-Torque Covariance and Coordination. Mol. Phys. 116, 1965–1976. doi:10.1080/00268976.2018.1459002
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
Joung, I. S., and Cheatham, T. E. (2008). Determination of Alkali and Halide Monovalent Ion Parameters for Use in Explicitly Solvated Biomolecular Simulations. J. Phys. Chem. B. 112, 9020–9041. doi:10.1021/jp8001614
Maier, J. A., Martinez, C., Kasavajhala, K., Wickstrom, L., Hauser, K. E., and Simmerling, C. (2015). ff14SB: Improving the Accuracy of Protein Side Chain and Backbone Parameters from ff99SB. J. Chem. Theor. Comput. 11, 3696–3713. doi:10.1021/acs.jctc.5b00255
Martínez, L., Andrade, R., Birgin, E. G., and Martínez, J. M. (2009). PACKMOL: A Package for Building Initial Configurations for Molecular Dynamics Simulations. J. Comput. Chem. 30, 2157–2164. doi:10.1002/jcc.21224
Matsarskaia, O., Braun, M. K., Roosen-Runge, F., Wolf, M., Zhang, F., Roth, R., et al. (2016). Cation-Induced Hydration Effects Cause Lower Critical Solution Temperature Behavior in Protein Solutions. J. Phys. Chem. B. 120, 7731–7736. doi:10.1021/acs.jpcb.6b04506
Matsarskaia, O., Roosen-Runge, F., Lotze, G., Möller, J., Mariani, A., Zhang, F., et al. (2018). Tuning Phase Transitions of Aqueous Protein Solutions by Multivalent Cations. Phys. Chem. Chem. Phys. 20, 27214–27225. doi:10.1039/C8CP05884A
Michaud-Agrawal, N., Denning, E. J., Woolf, T. B., and Beckstein, O. (2011). MDAnalysis: A Toolkit for the Analysis of Molecular Dynamics Simulations. J. Comput. Chem. 32, 2319–2327. doi:10.1002/jcc.21787
Moussa, E. M., Panchal, J. P., Moorthy, B. S., Blum, J. S., Joubert, M. K., Narhi, L. O., et al. (2016). Immunogenicity of Therapeutic Protein Aggregates. J. Pharm. Sci. 105, 417–430. doi:10.1016/j.xphs.2015.11.002
Pearlman, D. A., Case, D. A., Caldwell, J. W., Ross, W. S., Cheatham, T. E., DeBolt, S., et al. (1995). AMBER, a Package of Computer Programs for Applying Molecular Mechanics, normal Mode Analysis, Molecular Dynamics and Free Energy Calculations to Simulate the Structural and Energetic Properties of Molecules. Computer Phys. Commun. 91, 1–41. doi:10.1016/0010-4655(95)00041-D
Plattner, N., Doerr, S., De Fabritiis, G., and Noé, F. (2017). Complete Protein-Protein Association Kinetics in Atomic Detail Revealed by Molecular Dynamics Simulations and Markov Modelling. Nat. Chem. 9, 1005–1011. doi:10.1038/nchem.2785
Roosen-Runge, F., Heck, B. S., Zhang, F., Kohlbacher, O., and Schreiber, F. (2013). Interplay of pH and Binding of Multivalent Metal Ions: Charge Inversion and Reentrant Condensation in Protein Solutions. J. Phys. Chem. B. 117, 5777–5787. doi:10.1021/jp401874t
Ryckaert, J.-P., Ciccotti, G., and Berendsen, H. J. C. (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
Schames, J. R., Henchman, R. H., Siegel, J. S., Sotriffer, C. A., Ni, H., and McCammon, J. A. (2004). Discovery of a Novel Binding Trench in HIV Integrase. J. Med. Chem. 47, 1879–1881. doi:10.1021/jm0341913
Shirts, M. R., Klein, C., Swails, J. M., Yin, J., Gilson, M. K., Mobley, D. L., et al. (2017). Lessons Learned from Comparing Molecular Dynamics Engines on the SAMPL5 Dataset. J. Comput. Aided Mol. Des. 31, 147–161. doi:10.1007/s10822-016-9977-1
Spyrakis, F., Benedetti, P., Decherchi, S., Rocchia, W., Cavalli, A., Alcaro, S., et al. (2015). A Pipeline to Enhance Ligand Virtual Screening: Integrating Molecular Dynamics and Fingerprints for Ligand and Proteins. J. Chem. Inf. Model. 55, 2256–2274. doi:10.1021/acs.jcim.5b00169
Tarek, M., and Tobias, D. J. (2000). The Dynamics of Protein Hydration Water: A Quantitative Comparison of Molecular Dynamics Simulations and Neutron-Scattering Experiments. Biophysical J. 79, 3244–3257. doi:10.1016/S0006-3495(00)76557-X
Troussicot, L., Guillière, F., Limongelli, V., Walker, O., and Lancelin, J.-M. (2015). Funnel-Metadynamics and Solution NMR to Estimate Protein-Ligand Affinities. J. Am. Chem. Soc. 137, 1273–1281. doi:10.1021/ja511336z
Wang, J., Wang, W., Kollman, P. A., and Case, D. A. (2006). Automatic Atom Type and Bond Type Perception in Molecular Mechanical Calculations. J. Mol. Graphics Model. 25, 247–260. doi:10.1016/j.jmgm.2005.12.005
Yu, J., Ciancetta, A., Dudas, S., Duca, S., Lottermoser, J., and Jacobson, K. A. (2018). Structure-Guided Modification of Heterocyclic Antagonists of the P2Y14 Receptor. J. Med. Chem. 61, 4860–4882. doi:10.1021/acs.jmedchem.8b00168
Keywords: statistical mechanics, entropy, free energy methods, multiscale, metadynamics method, protein-protein binding, protein-excipient binding, protein hydration
Citation: Kalayan J, Curtis RA, Warwicker J and Henchman RH (2021) Thermodynamic Origin of Differential Excipient-Lysozyme Interactions. Front. Mol. Biosci. 8:689400. doi: 10.3389/fmolb.2021.689400
Received: 31 March 2021; Accepted: 25 May 2021;
Published: 11 June 2021.
Edited by:Fabio Trovato, Freie Universität Berlin, Germany
Reviewed by:Yandong Huang, Jimei University, China
Riccardo Nifosì, National Research Council (CNR), Italy
Giuseppe Pellicane, University of Messina, Italy
Copyright © 2021 Kalayan, Curtis, Warwicker and Henchman. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Jas Kalayan, firstname.lastname@example.org