Advances in Modeling the Stability of Noncovalent Complexes in

Perspective pubs.acs.org/ac

Advances in Modeling the Stability of Noncovalent Complexes in Charged Droplets with Applications in Electrospray Ionization-MS Experiments Styliani Consta,*,†,¶ Mahmoud Sharawy,† Myong In Oh,† and Anatoly Malevanets‡ †

Department of Chemistry, The University of Western Ontario, London, Ontario N6A 5B7, Canada Department of Electrical and Computer Engineering, The University of Western Ontario, London, Ontario N6A 5B9, Canada ¶ Department of Chemistry, University of Cambridge, Lensfield Road, Cambridge, CB2 1EW, United Kingdom ‡

S Supporting Information *

ABSTRACT: Electrospray ionization mass spectrometry is used extensively to measure the equilibrium constant of noncovalent complexes. In this Perspective, we attempt to present an accessible introduction to computational methodologies that can be applied to determine the stability of weak noncovalent complexes in their journey from bulk solution into the gaseous state. We demonstrate the usage of the methods on two typical examples of noncovalent complexes drawn from a broad class of nucleic acids and transient protein complexes found in aqueous droplets. We conclude that this new emerging direction in the use of simulations can lead to estimates of equilibrium constant corrections due to complex dissociation in the carrier droplet and finding of agents that may stabilize the protein interfaces.

T

reliability of the detected structures of protein assemblies is ensured by using a combination of methods that include native mass spectrometry, X-ray crystallography, NMR spectroscopy, and (cryo-)electron and computational methods.23,26 For protein dimers and nucleic acids, the reliability of ESI-MS in measuring equilibrium constants has been validated by comparisons mainly with isothermal titration calorimetry.10,12,27−29 Despite the improvements, the reliable detection of the equilibrium constant of a class of weak noncovalent protein complexes and that of nucleic acids is still vigorously debated.6,7,30 It has also been argued that analysis of every bimolecular complex requires an individual experimental protocol11 that, as expected, reduces the efficiency of the analysis. To our knowledge, the role of the droplet environment in the stability of complexes has not been investigated yet. We think that this knowledge may lead to a rational selection of the experimental conditions and reduction of the experimentation required to establish the best analytical chemistry protocol. The effect of the carrier droplet cannot be readily detected in experiments due to the small sizes of the droplets and the high rate at which the various processes take place. To understand the droplet environment, experiments have attempted to profile the

he relation between the gas phase and the solvation phase structures of protein assemblies and nucleic acids detected by native mass spectrometry has been a recurrent question for over two decades.1−5 The same question arises when the equilibrium constants of noncovalent complexes are determined by using electrospray ionization mass spectrometry (ESIMS).6−17 The noncovalent complexes may be nucleic acids, proteins, and their combination, as well as protein−ligand complexes.18 The ligands include drug molecules, peptides, oligonucleotides, carbohydrates, and lipids. The noncovalent complexes of biological molecules play a fundamental role in virtually all cellular functions. For instance, transient protein complexes are responsible for signaling and regulatory mechanisms in a cell,19 and nucleic acids store genetic information. The dissociation constant of a complex ranges from nanomolar concentration for drug−protein interactions to millimolar concentration for protein complexes participating in signaling pathways. Using double-stranded DNA with different nucleotide sequences, 6 one may design and investigate complexes with different values of the equilibrium constant. The experimental conditions that lead to stable protein assemblies and protein or nucleic acid dimers have been thoroughly investigated over 20 years.4,5,20 These efforts resulted in considerable improvement in maintaining weak protein− ligand, protein−protein interfaces, and nucleic acid pairing throughout the complex transfer from bulk solution into the gaseous state.7,9,14,15,21−25 In the field of structural biology, the © XXXX American Chemical Society

Received: May 22, 2017 Accepted: July 19, 2017 Published: July 19, 2017 A

DOI: 10.1021/acs.analchem.7b01941 Anal. Chem. XXXX, XXX, XXX−XXX

Perspective

Analytical Chemistry pH of the shrinking droplets,31−33 but the interactions of the complexes within the droplets are yet unknown. Our computational efforts aim at establishing the factors that determine the stability of the complexes in the droplet environment, using atomistic molecular dynamics (MD). Simulations have revealed several distinct features that differentiate the droplet environment from its bulk analogue. The large surface-to-volume ratio is a natural characteristic of finite-sized systems that makes the surface, in proportion, a significantly more active region in the overall droplet chemistry. One consequence of the existence of a free surface is that there are additional interfaces (vapor−liquid and vapor−macromolecule) relative to bulk solution. In addition, for a droplet with linear dimensions in the nanometer range, the shape fluctuations are large relative to the spherical shape. The second feature is that electrosprayed droplets undergo many dynamic processes, which include solvent evaporation and fragmentation by ejection of solvated ions. The third feature is that the droplets are highly charged due to the presence of the macroions, which in turn induce particular structure in the solvent and in the overall droplet shape.34,35 The combination of all these factors may result in different mechanisms that stabilize the interactions within complexes from those in the bulk solution analogues.35,36 Through this Perspective, we aim at introducing the computational methods used in the study of complexes and the key findings of the studies. We demonstrate the methods using two examples drawn from a broad class of nucleic acids and transient protein complexes. We discuss how the droplet evaporation may affect the conformations and the charge states of double-stranded oligodeoxynucleotide (dsDNA) and how one may determine computationally the dissociation constant of a protein complex in a droplet at variable pH. An additional outcome of the latter study is the mechanism of protein complex dissociation. Our purpose here is 2-fold. On the one hand, we aim at stimulating the synergy of computational and experimental research aiming at systematic development of experimental protocols that stabilize weakly bound protein interfaces in droplets. On the other hand, we demonstrate the manner in which the droplet evaporation and the late desolvation may affect the conformations and charge states of complexes.

droplet is below the Rayleigh limit when the surface forces, which tend to keep the droplet connected, overcome the electrostatic repulsions among the ions that tend to fragment the droplet. In the following discussion, the Rayleigh limit is used when we assign the charge on the protein, but it is not used in the analysis of the molecular dynamics trajectories. To find the relation between the concentration of charge and the radius of the droplet, eq 1 is recast in the following form [Z ] =

1000NAVe

=

1 1000NA

36ε0γ r0 3e 2

(2)

4

where V = 3 πr0 3 is the volume of a spherical droplet, e is the positive elementary charge, and NA is Avogadro’s number. The plot of [Z] vs r0 is shown in Figure 1. We note that Qr is the net

Figure 1. Log−log plot of the charge concentration [Z] vs droplet radius given by eq 2 (solid line) and the average charge density of a macroion vs its radius37,38 (dashed line) (see text for details). IEM in the legend stands for ion-evaporation mechanism. The snapshot at the upper right corner of the plot is from simulations of an aqueous droplet charged with Na+ ions (depicted by blue spheres). The oxygen sites of the H2O molecules are represented by red spheres and the hydrogen sites by white spheres. The sizes of the water molecules and Na+ ions have been rescaled for clarity.

■

METHODOLOGY Macroscopic Stability Model of a Charged Droplet. The fundamental model that has been used in the field of MS to understand the amount of charge that a droplet can hold is that of Lord Rayleigh developed in 1882.39 The Rayleigh model is macroscopic. It assumes a charged spherical conductor that may undergo small shape fluctuations relative to the spherical shape.35,40 The different shapes have a constant volume, but the surface area may change. The model expresses the total energy of the system as the sum of surface energy and electrostatic energy. Linear stability analysis of the total energy of the droplet40,41 leads to the Rayleigh limit for the stability of the droplet, which is given by Q r 2 = 64π 2ε0γr0 3

Qr

charge of the droplet. To obtain a sense of the increase in the ion concentration in a droplet, we may compare the Na + concentration in physiological saline solution which is approximately 150 mM to that in a shrinking droplet. We assume that an electrosprayed droplet will have an amount of NaCl as in the parent bulk solution and an excess amount of Na+ ions because of its enrichment on positive charge. In a droplet of 5000 H2O molecules (corresponding approximate radius of 3.3 nm), the excess Na+ concentration becomes approximately 265 mM and in a 1000 H2O molecule droplet (corresponding approximate radius of 1.9 nm) becomes 592 mM. The Rayleigh limit shows that a 100 times decrease in a droplet radius results in a 1000 times increase in the net charge and, therefore, the ion concentration that carries the excess charge. We note that eq 1 holds for separable charges. For instance, droplets with alkali metal ion, halogen, ammonium, and acetate ions are typically found in experiments. Eq 1 is also used for droplets that contain a macroion and simpler ions. Exper-

(1)

where Qr denotes the charge of the droplet at the Rayleigh limit, ε0 is the vacuum permittivity, γ is the surface tension, and r0 is the droplet radius. Hereafter, we will refer to the droplet being found “below the Rayleigh limit” when Qr2 < 64π2ε0γr03. When the relation is reversed, a droplet is referred to be found “above the Rayleigh limit”. Another way to express this criterion is that a B


Perspective

Analytical Chemistry

where rij is the distance between two atomic sites i and j, ε0 is the vacuum permittivity, qi is the partial charge on the atomic site i, ϵ and σ are LJ parameters that model the strength of the interatomic interactions and the diameter of the atomic sites, respectively. The first term in the summation of eq 3 is the electrostatic energy, and the second term is the LJ energy. The intramolecular interactions are augmented by potential energy functions that account for vibrations of bonds and angles as well as torsions. The fast vibrational motions may also be kept frozen by applying constraints in the Newton’s equations of motion.57 There is a variety of force fields that are generally used in simulations. For example, two force fields that we use here are the CHARMM General Force Field58,59 and the AMBER99SBILDN.60 Here, the CHARMM force field is used for modeling the ubiquitin−ubiquitin associated domain protein complex (RCSB PDB61 code 2MRO62) in aqueous droplets, where the water molecules are represented by the TIP3P (three-site transferable intermolecular potential)63 model. The AMBER99SB-ILDN force field features explicit hydrogen atomic sites and is used here to model a DNA duplex in TIP3P water. A typical snapshot of an aqueous droplet containing the DNA and ions modeled by the AMBER99SB-ILDN is shown in Figure 2.

imentally, it has been found that a droplet may release its charge below the Rayleigh limit by the ion-evaporation mechanism (IEM).37,38,42−51 However, the size of the droplet does not change significantly relative to its size at the Rayleigh limit, because the ions are released while solvated with a small number of water molecules. In Figure 1, we compare the Rayleigh predictions to that of the average charge density of a macroion obtained via the ion concentration kinetic model of Labowsky et al.37 implemented for protein containing droplets by Hogan and de la Mora.38 In generating the plot, we used eq 4 in ref 38. Reference 38 provides evidence that the charge on protein ions is determined by triethylammonium buffer ion evaporation from electrosprayed drops. Modeling and Simulation Methods. Determining the Charge State of the Macroion. Molecular simulations of droplets, in general, are based to a great extent on the usage of MD.52 Simulations of a complex in an evaporating droplet require the use of a multiscale approach of computational methods. Experimental evidence demonstrates that the charge state of a complex differs from that in the parent solution.6,53 Therefore, the simulation approach should accommodate changes in the charge state of a macromolecule throughout its journey from the bulk solution to the gaseous state. The choice of the charge state of a macromolecule depends on the droplet evolution as well as the composition of the solvent. The use of quantum mechanical methods coupled to MD, such as CarParrinello MD54 or Born-Oppenheimer MD55 that would allow changes in the molecular structure in the course of simulations, are prohibitive even for the smallest droplets. Hence, one must carefully select the molecular structures and their protonation states for each stage of the droplet evolution. On a macroscopic level, the macromolecular charge state ought to be deduced by solving the Poisson-Boltzmann equation self-consistently coupled to a system of chemical equilibrium equations that involve the reacting species.56 As one determines the apparent pH level in the droplet, one may protonate the macromolecule based on the pKa of the titratable groups.35,56 In addition, the varying pH of a droplet alters the chemical composition of the solvent. Therefore, for every chosen protonation state, a separate molecular simulation study has to be carried out. In one of the following sections, we present the equations that are used to determine the protonation state of a protein complex. Once the droplet composition and the macroion’s charge state are decided, one may perform MD on the system in order to monitor its time evolution. MD Protocol. To perform MD simulations, we first define the molecular composition of a droplet and the detailed molecular structure of every macromolecule in the system. For the given topology and composition, we use a force field that describes the interactions between the atomic sites. In atomistic simulations, the molecules in the droplets are often modeled by molecular mechanics force fields. At this level of modeling, the nucleus and the electrons of an atom are represented as a single entity by a spherical atomic site that carries a partial charge. The molecular mechanics force fields include two basic types of interatomic interactions, Lennard-Jones (LJ) and Coulomb. The potential energy of interaction (Uij) between two atomic sites i and j is expressed as ⎡⎛ ⎞12 ⎛ ⎞6 ⎤ 1 qiqj σ σ Uij = + 4ϵ⎢⎢⎜⎜ ⎟⎟ − ⎜⎜ ⎟⎟ ⎥⎥ 4πε0 rij r r ⎝ ij ⎠ ⎦ ⎣⎝ ij ⎠

Figure 2. A typical snapshot of a water droplet composed of the 11-mer A·T (adenine-thymine) [dsDNA]20−. The ribbons represent the backbone of [dsDNA]20−. The sizes of the water molecules, Na+ (purple), and Cl− (green) ions are rescaled for clarity. The oxygen sites are represented by red spheres and the hydrogen sites by white spheres.

Here, we emphasize an important point in the treatment of the electrostatic interactions. The electrostatics interactions are longranged, and for this reason, an arbitrary cutoff should not be applied in the simulations of a droplet. It is still possible to apply a cutoff scheme, but its distance has to be larger than any of the linear dimensions of the shape fluctuations of the droplet. The cutoff distance should be longer than the elongated shapes due to fragmentation of the droplet and release of ions. An application of a short cutoff distance in an strong electrolyte system may lead to artifacts. An alternate method to the direct treatment of the electrostatic interactions is the multilevel summation method implemented in the molecular simulation software NAMD.64

(3) C


Perspective

Analytical Chemistry Molecular dynamics52 generates the trajectories of the system by integrating Newton’s equation of motion for each atomic site. Common algorithms to do the time evolution of the positions of the atomic sites are the Verlet65 and velocity Verlet66 algorithms.52 The time step used in all of our simulations has been 0.5 or 1.0 fs. Molecular dynamics of atomistically modeled systems is very time-consuming, and for this reason, it is limited to the study of droplets with linear dimensions of several nanometers. Indicative sizes of atomistically modeled droplets in our studies are those with 7000 H2O molecules (≈ diameter of 8 nm)67,68 and 11 500 CH3CN molecules69 (≈ diameter of 12 nm) down to a few hundred solvent molecules.68 To our knowledge, these aqueous and CH3CN droplets are among the largest that have been simulated so far by atomistic models for typical simulation times of tens to hundreds of nanoseconds. Nanodrops with linear dimensions of several nanometers correspond to the latest stage of solvent evaporation in an ESI-MS experiment. We note that these simulated droplet sizes are small relative to the initially generated droplets by ESI where droplet dimensions of nanospray are approximately 100−150 nm. Atomistic modeling can not currently explore the chemical processes in droplets with a diameter of tens of nanometers or larger. The processes in these larger droplets might be similar to those in the bulk solvent or bulk vapor−liquid interface; hence, one may consider performing bulk simulations. Atomistic modeling of the smaller droplets is a feasible and valuable tool because conformational changes and charging of the macroions may still take place in these tiny droplets.36,67,68,70−72 On the other end of modeling, an analytical theory can be applied to any system size. Such models have been developed to capture the charge-induced instabilities and the release from droplets of linear macroions that do not undergo chemical modifications.35,56,73 In the study of the complexes of macromolecules in droplets, MD can be used to perform simulations of the evaporation of the droplets in order to find out (a) the manner in which an analyte emerges from a droplet and (b) the correction to the bulk equilibrium constant due to the complex dissociation in the droplet. Nonequilibrium Simulations. Nonequilibrium simulations intend to explore complex dynamic processes that occur in the evaporating droplets. Hereafter, the term “evaporation” will be used for the removal of solvent molecules from the droplet and the term “desolvation” to describe detachment of the few solvent molecules remaining on a macromolecule or the complex already emerged from a droplet. The exact experimental conditions of temperature and pressure under which the evaporation and desolvation of the gas phase ions take place are still not known. Therefore, one may consider comparing simulations at constant energy and constant temperature conditions. As we mentioned in the previous section, MD is often performed by using the Verlet or the velocity Verlet integration algorithm.52 In the nonequilibrium simulations, the droplet is placed in vacuo. Constant temperature simulations are usually carried out by coupling the integration algorithm to a heat bath. The temperature control in our simulations is performed either by Nosé-Hoover dynamics or Lowe-Andersen thermostat.74 When constant energy simulations are performed, the connected part of the droplet cools down because of solvent evaporation. In many situations, constant temperature and constant energy MD lead to the same outcome with respect to the charge state and the conformation of the macroion.67 However, there are cases where the outcomes may be different. Figure 3a shows a typical snapshot of an aqueous droplet comprising a partially sodiated

Figure 3. (a) Aqueous droplet containing a PEG64 (64 denotes the number of −CH2−O−CH2− monomeric units in the chain), Na+ (purple spheres), and Cl− (green spheres); (b) PEG64-5Na+ and a NaCl crystallite resulting from evaporation of (a) via constant energy simulations; (c) release of PEG64-4Na+ from evaporation of (a) via constant temperature simulations.

PEG64, chloride ions, and excess of sodium ions relative to the chloride ions. This is a snapshot taken from a long MD run starting from a larger droplet that decreases in size due to solvent evaporation and solvated ion ejections. When the droplet in Figure 3a continues to evolve by constant energy MD, it leads to the release of a sodiated PEG64 at the charge state of +4 e (Figure 3b).72 The Na+ and Cl− ions remain in the droplet, where depending on their concentration, they may form ion complexes such as [NaCl2−] (when there is an excess of Cl− ions) or [Na2Cl−] (when there is an excess of Na+ ions) or larger aggregates. When the droplet in Figure 3a evolves by constant temperature MD, it leads to the formation of a charged [Nan+1Cln] aggregate, which is attached to the PEG64 chain. The overall charge of PEG64 is +5 e. We attribute the formation of a crystallite on the PEG to faster evaporation of the H2O molecules in a constant temperature run compared to that in a constant energy. Here, we note that the constant energy simulation in vacuo better represents the experimental conditions at the latest stage of droplet desolvation under conditions of very low pressure, where the droplet is not constantly thermalized due to collisions with the gas. In performing the simulations, there are two points to be considered: (a) Depending on the software used, one might have to remove the evaporated solvent molecules when they are at a certain distance (≈70 nm) from the connected body of the droplet. The removal of the evaporated molecules (solvent molecules and ejected solvated charged ions) is required if GROMACS is used, while NAMD software handles itself the evaporated H 2O molecules. (b) In the nonequilibrium simulations, a number of carefully prepared initial conditions have to be used. The same simulation protocol that was outlined for single macromolecules in droplets is also used for complexes of macromolecules. We present MD simulations of the evaporation of aqueous droplets comprising a [dsDNA]20− with various concentrations of sodium and chloride ions. Details of the simulations are found in Sharawy and Consta.67

■

SIMULATIONS OF NUCLEIC ACIDS IN DROPLETS We used molecular modeling to study the manner in which a [dsDNA]20− emerges from an aqueous droplet that also contains Na+ and Cl− ions. We expect that the polyion [dsDNA]20− sprayed from the initial solution will carry an ion atmosphere. The ionic atmosphere carries an excess of counterions (e.g., Na+ or NH4+) and most likely a smaller concentration of co-ions (e.g., Cl− or CH3COO−). D


Perspective


Figure 4. Latest stages of aqueous droplet evaporation. The initial configuration, shown in ref 67, comprised [dsDNA]20−, 1500 H2O molecules, and Na+ and Cl− ions. The color coding is as in Figure 2. The initial nanodroplet carries an overall negative charge because of excess Cl− ions. The droplet charge state changes its value during evaporation because of emission of Cl− ions. The snapshots in each column show (i) an intermediate state and (ii) the [dsDNA]20− forming adducts with Na+ and Cl− ions at the latest stage of the evaporation, where the complex is almost bare of water molecules. (a) Fragmentation of an aqueous [dsDNA + 5Na+]15− droplet. In (ii), “A” stands for the adenine strand strand and “T” for the thymine strand (see also Figure 2). (b) Drying-out of an aqueous [dsDNA + 10Na+]10− droplet. (c) Drying-out of an aqueous [dsDNA + 20Na+ + yCl−]y− droplet where y = 9, 7. (d) Drying-out of an aqueous [dsDNA + 25Na+ + 10Cl−]5− droplet.

The form of the charge-induced instability is a key point to note in droplets that contain a polyion. An example of the instability is shown in Figure 5a,b that shows a [dsDNA·5Na]15− complex in an aqueous droplet. During the evaporation of water, the droplet passes through states that are found above the Rayleigh limit. These states are found in the instability regime, but they cannot fragment and emit ions in order to lower their charge.50,51 Here, the ionic atmosphere is strongly bound to the [dsDNA]20−; therefore, the macroion and its ionic atmosphere behave like a single macroion with fixed charges. The droplet has to change its shape to accommodate this type of instability. The new shape is characterized by conical protrusions (hereafter, we will call them “spikes”) of the solvent surrounding the macroion. We first found this form of instability along with other instabilities in 2010,34 and we call it a “star”-like shape because when the central ion is spherical, it literally looks like a “star”. The structure of the star shapes and the equations that govern their formations are distinct. Details on their characterization and their physical chemistry are found in Sharawy and Consta69 and Consta et al.35 A more subtle “spiky” form of instability also appears in charged droplets that carry less net negative charge (Figure 5c,d). In Figure 5c,d, the spikes on the surface of the droplet containing a [dsDNA·10Na]10− are less pronounced than in Figure 5 because the droplet here carries a less negative charge. Despite having a less net negative charge, the droplet is still above the Rayleigh limit. The consequences of the star-shaped instability in maintaining the binding of the complexes have been discussed in refs 36 and 67. One of the findings is that the spiky droplets shapes increase the solvent evaporation rate due to the increase of the surface area. Consequently, the star-shaped droplet may carry a noncovalently bound complex even when the macroion is surrounded by only a single solvation shell. We think that the associated complex is maintained because the solvent evaporation rate becomes higher than the complex dissociation rate.

In charged droplets containing mobile charges, the total charge of the droplet is assumed to be given by the Rayleigh expression, eq 2. Given the size of the droplet and the composition of the solvent, we can deduce the charge on the macroion and vice versa. We note that the macroion charge in the intermediate sized droplets of 1 nm < radius 1. Thus, the overall droplet has an excess of negative charge which mimics the conditions of the negative-ion mode in an ESI experiment. The initial conditions are presented in ref 67. In all the systems that we examined, if the number of Na+ ions was less than half the charge on the [dsDNA]20− then the Cl− ions were ejected, and the droplets were left only with the [dsDNA]20− and Na+ ions. Typical snapshots of the evolution of droplets that hold 5 and 10 Na+ ions are shown in Figure 4a,b, respectively. In contrast, when the Na+ ions almost compensate the charge of the [dsDNA]20−, several Cl− ions remain in the droplet. Moreover, as the system dries out, salt-aggregates condensate on the [dsDNA]20−. The Rayleigh limit determines the number of Na+ and Cl− that will remain in a droplet. For instance, the Rayleigh limit for a droplet of 1000 H2O molecules is at a charge of −9 e. This estimate is made by applying eq 1 to a droplet of 1000 H2O molecules with surface tension 0.053 N/m for the TIP3P water model at T = 300 K.75 If the droplet contains a dsDNA with a charge of −30 e, it will maintain at most 21 Na+ ions and it will expel all the Cl− ions. E


Perspective


Figure 5. Two sets (a, b) and (c, d) of typical snapshots of the evaporation of an aqueous droplet from the simulation of the systems, Figure 4a,b, respectively. The color coding of the atomic sites is the same as in Figure 2. The initial number of Cl− ions was smaller than that of the Na+ ions, but the droplets are overall negatively charged. The initial configurations that contained Cl− ions are not shown because all the initial Cl− ions evaporated from the droplet. (a, b) comprised initially 5 Na+ ions and (c, d) comprised 10 Na+ ions. Simulation with 5 Na+ ions: (a) is composed of 1240 H2O molecules and (b) 960. Simulation with 10 Na+ ions: (c) is composed of 990 H2O molecules and (d) 570.

Figure 6. (a) Snapshot Figure 4b-ii. (b) Schematic picture of the [dsDNA + 10 Na+]10− conformation and charge distribution in the gas phase. During the transition of [dsDNA + 10 Na+]10− from the droplet into the gas phase, the two strands slide in opposite direction about half their length. The two strands twist at the points indicated by two brackets such that their backbones are held together by the excess of Na+ ions. There is no net charge in the central region (highlighted by the bubble), where the two halves of the strands are held together by 10 Na+ ions. The exposed halves contribute a net charge of −10 e.

We continue to explore the consequences of the star-shaped structures in the stability of the complexes. [dsDNA + xNa+ + yCl−]−20+x−y in the Gaseous State. Another interesting finding coming from the simulations of a [dsDNA]20− in a charged droplet is the effect of the charge on the conformations of the dsDNA in the gas phase. We realized three tertiary structures of the duplex in the gas phase: denatured, double-stranded, and inverted zipped-up DNA.67 These final structures are dependent on the net charge that the droplet carries in the gas phase. In the systems that we investigated when the total charge is more negative than half of the fixed charge on the [dsDNA]20−, the duplex dissociates in the gas phase (Figure 4a). In contrast, when the total charge is less negative than half of the charge on the [dsDNA]20−, the dsDNA maintains its duplex conformation in the gas phase (Figure 4c,d). An interesting structure of the DNA duplex in the gas phase is generated when the total negative charge of the droplet is half of the negative charge on the DNA (Figure 4b). In this charge state, the hydrogen bonds between the complementary nucleotide pairs break, and the two strands slide in the opposite directions. In the gas phase, the bases of the nucleotides are inverted (about their pairing orientation), and the two strands are held together by the charge of the 10 Na+ ions. It is intriguing that the +10 e contribution of the Na+ is the charge required to compensate for the −10 e of the DNA in the pseudodouble strand region (Figure 6). The location of the Na+ ions in the center of DNA and the location of the negative charge of DNA at the termini agrees with the experimental findings of Favre et al.76,77 Electrospray experiments78−80 have suggested that the two strands remained bound in the gas phase. However, for the question of the negative charge limit to maintain the DNA duplex

in the gas phase, limited information was found in the literature. Gabelica et al.80 generated data on a 12-mer dsDNA which suggested that the duplex state is maintained for a net charge that is less negative than −5 e, ∼25% of the fully deprotonated 12-mer dsDNA charge. Because of the complexity of real systems, we cannot argue that our computational study replicates the experiment. However, our findings generally agree with those obtained in ESI-MS that the dsDNA can maintain its doublestranded state if its net charge in the gas phase is less than half its isolated charge. In the simulations, we were assisted by the fact that the charge carriers Na+ and Cl− were forming ion-adducts. If the simple charge carriers were to participate in protonation reactions, such as the commonly used ammonium acetate, then the simulations would require an approach that involves the quantum mechanical modeling of protonation reactions. For instance, in the experimental study by Barylyuk et al.6 on dsDNA−ssDNA (ss stands for single-stranded DNA) equilibrium, the ratio of the final charge state for ssDNA is approximately 2 = VdsDNA /VssDNA (where VdsDNA and VssDNA are the volumes of the ssDNA and dsDNA, respectively) and is determined solely by the final volume of the macroion. We estimated this ratio by simply applying eq 2 in the data provided in ref 6. The fact that in the positive ESI mode the DNA charge is positive indicates proton transfer to the nucleotide bases during droplet evolution. Interestingly, the final positive charge on the DNA necessitates that the dsDNA would have significantly higher charge in the intermediate droplets. This statement F


Perspective


simulations of different droplet sizes at the corresponding pH and complex protonation state. Since the simulations are at equilibrium, the reactivity of NH4+ or CH3COO− ions do not change the pH. The droplet is thermalized at the desired temperature by using the Lowe-Andersen thermostat.74 One may turn off the thermostat and continue the simulations in the microcanonical ensemble (constant energy, volume, and number of molecules). The question that arises now is how long should the MD run be? In order to determine the length of the MD run, one estimates the rate of solvent evaporation from a separate run. From the evaporation rate, one estimates the time required for the droplet to change its size. This is the cutoff time of the equilibrium MD runs in the cavity. We monitor whether there is protein complex dissociation and subsequent droplet fragmentation within the MD time period. If the dissociation rate is comparable or higher than the evaporation rate, it is deduced that protein dissociation may occur in the droplet and it may alter the value of the equilibrium constant of the complex. The protein dissociation times are fit to a decaying exponential function from which we determine the complex dissociation rate. By assuming first-order kinetics, we can estimate the number of dissociated complexes and, thus, the error in the experimentally measured equilibrium constant due to the complex dissociation in the droplet carrier. Determination of the Protonation State of a Protein in a Droplet. The assignment of the charge state of each ionizable group depending on its environment has been a challenging problem in biochemistry.81−91 Here, we use macroscopic equations to find the average protonation state of a protein in a droplet. The details of the modeling are described in Consta et al.35 The charge state of the protein, the pH, and the charge of a droplet depend on each other. The set of equations that show this dependency is described below. The average total droplet charge is given as the sum of all the charges of the chemical species. For example, if the free charge carriers are ammonium and acetate ions, the total droplet charge is given by the following expression

follows from the fact that the macroion has the highest charge when the charge on the ion is the same as the remaining charge in the droplet (see Malevanets and Consta56). Tantalizingly, this may explain why the response coefficient differs for dsDNA in the positive and negative modes as shown in Figure 5a of Barylyuk et al.6 Simulations of the Rate of Protein Dimer Dissociation. Typically, ESI protein complexes are found in acidic aqueous droplets. The smaller the droplet, the lower is the droplet pH. The droplet acidity is determined, among other factors that will be discussed later, by the concentrations of NH4+ and CH3COO− ions. Even in the absence of CH3COONH4, the droplet environment is still acidic because of electrolytic processes that occur in the ESI instrument. In order to study the stability of a protein complex, in principle, one may follow the simulation protocol described in the previous section. Thus, one may monitor the evaporation in a large number of simulation runs starting from different initial conditions of the droplet− complex system and analyze the final state of the complex. A major bottleneck in performing those simulations is that the protonation state of the protein may change during the droplet evaporation since the pH of the droplet decreases. The protonation state of a protein at fixed pH has been a major computational problem that is being studied over several decades in bulk solution.81−91 In order to consider the protonation state of a protein or its complex in a shrinking droplet where the pH decreases, we devised the following method: Equilibrium MD runs of a droplet of specific size are performed instead of nonequilibrium runs. The equilibrium is achieved by placing the droplet in a cavity as shown in Figure 7. The cavity is a spherical

Zr = Zp(pH) + NA × V × 1000 × ([NH4 +] − [CH3CO2−] − [OH−] + [H+])

(4)

where Zr × e = Qr is the droplet charge at the Rayleigh limit (eq 1), NA is Avogadro’s number, V is the droplet volume, Zp is the macromolecule (protein) charge at a given pH, and the square brackets denote the molarities of the components. Eq 4 holds on average; therefore, one considers an ensemble of droplets. Eq 1 is recast in the following form (Zre)2 = 48π ϵγV

Figure 7. Protein complex (RCSB PDB code 2MRO) with net charge of +14 e in a droplet of ≈1650 H2O molecules. The complex is colored in blue, and only the oxygen site of water is shown in red for clarity. In order to perform equilibrium simulations, the entire droplet system is enclosed in a cavity (colored in green). The entire system droplet and vapor are composed of 2000 H2O molecules. The radius of the cavity used in this study is 20 nm.

(5)

We aim to solve numerically the system of eqs 4−8 coupled with the Rayleigh criterion (5). To this end, we need to derive an expression of the charge on a protein Zp(pH). The species concentrations are not independent but they are calculated from the chemical equilibria.

potential that confines evaporating molecules within a volume. Water evaporation from the droplet in the vacant space of the cavity establishes an equilibrium vapor pressure. For the particular droplet size, the pH and the corresponding protonated state of the protein are determined on the basis of the method presented in ref 36. For completeness of the discussion, the basic equations are presented in the following section. We can perform

[OH−]a × [H+]a = K w

(6)

10−pH = [H+]a

(7)

[NH4 +]a = K NH4+ [H ]a × [NH3]

(8)

+

G


Perspective

Analytical Chemistry [CH3COOH] = K CH3COO− [CH3COO−]a × [H+]a

Figure 8a shows the change in the droplet pH with the droplet charge in two situations: the red line is for an aqueous droplet of 150 mM CH3COONH4 (no protein present), and the blue line is for a droplet of 150 mM CH3COONH4 and the 2MRO complex. The overall droplet charge is given by eq 4, and it is determined by the Rayleigh limit (eq 1). We will trace the red line first by examining two points. We will move from the smaller droplet toward the larger ones, which is opposite to the natural evaporation process, where the droplet evolves from larger to smaller. A droplet with 10 charges (879 H2O) will have two pairs of CH3COO−-NH4+. The charge +10 e is obtained by 2 NH4+ and 8 H3O+ ions. The two CH3COO− ions become CH3COOH. The abrupt change in pH appears when the droplet charge increases approximately from 20 e to 30 e. When Zr = 20 (3517 H2O molecules), there are going to be 10 CH3COO−-NH4+ pairs, which result in 10 NH4+ ions, 10 CH3COOH, and 10 H3O+ ions. The blue line in Figure 8a shows that the protein plays a buffering role. At pH = 3 (lowest pH), the protein attains its highest charge state (+10 e) as expected. In Figure 8, the final states are the same regardless of the conditions. However, in the case of low buffer concentration, the protein complex in the late droplets may take high charge leading to the complex dissociation. A given charge state of a protein comprises all possible protonation substates where the charge is distributed among the titratable groups of the protein. The relative sampling weight of each substate is determined by the local electric potential as well as the protein conformation ensemble of the substate. With some degree of scientific confidence, this can only be achieved by using quantum mechanics/molecular mechanics (QM/MM) methods where one allows the transfer of a proton between different groups. Therefore, a direct sampling of various protein charge states is rather impractical. In the presented method, we find the average charge of the protein using the macroscopic equations and then we assign the charge in the protein by first protonating the residues with the highest pKa values. If the average charge is fractional, one has to assign various charge states. Also, if the pKa values of different ionizable groups are similar, the role of different charge distributions on complex dissociation has to be examined. Moreover, detailed studies should be carried out to find whether a specific protonation state can strongly influence the dissociation of the complex.

(9)

where Kw is the equilibrium constant for the water autoprotolysis reaction at temperature T and we equate the logarithm of the molarity of the hydronium ions with the apparent pH. The subscript a denotes activities of the ionic species. We will be mixing concentrations and activities in the discussion that follows. In the numerical calculations and the plots that we present, we have used the Debye-Hückel correction for the activity coefficient for a point charge. For hydronium and hydroxide ions, we use concentrations instead of the activity. Using the derived expression for Zp(pH), one solves the system of eqs 4−9 coupled with the Rayleigh criterion (5). The solution of the system of equations provides the charge of the protein that is subsequently used to study the dissociation of the protein complex. Let us consider CH3COONH4 as a buffer compound. CH3COOH has Ka = 1.79 × 10−5 (pKa = 4.75); NH4+ has Ka = 5.6 × 10−10 (pKa = 9.25 and pKb = 4.75 of NH3). In Figure 8, we show the dependence of the protein charge and the droplet pH in the course of quasi-static evaporation. Since the protein charge is known, the groups with the highest pKa are protonated first. It is interesting to discuss the relation between the protein complex (or single protein) charge state, the pH of the droplet, and the concentration of NH4+ and CH3COO− ions.

■

STABILITY OF PROTEIN COMPLEXES We study the complex dissociation in a droplet as the only possible process that may alter the equilibrium constant measurements in nanosprayed droplets. Association of protomers may also change the equilibrium constant; however, this process is likely only in microdrops. The study of the complex dissociation has many parameters to examine with respect to the charge state of the protein and the size of the droplet. In our studies, we focus on the late stage of the droplet lifetime, where we performed the majority of the simulations in a cavity containing 2000 H2O molecules and a 2MRO at a charge state of +14 e at temperatures in the range of 370−390 K. Because of the equilibrium established between the vapor and the droplet in the cavity, the droplet was composed of ≈1350−1650 H2O molecules. The rest of the molecules were found in the vapor phase. We simulated 14 systems of the protein complex protonated at +14 e where +11 e is found in ubiquitin and +3 e in UbA. The protein complex is embedded in an aqueous droplet. The charge state of +14 e was selected because we expect complexes at the higher charge states to be more susceptible to

Figure 8. (a) Droplet pH vs droplet charge. The red line is for an aqueous droplet that contains only CH3COO− and NH4+ ions (free) at a concentration of 150 mM. The blue line corresponds to the same system but containing 2MRO in addition. (b) The red line shows the droplet pH (red labeled axis on the right), and the blue line is the protein’s charge (Zp) (blue labeled axis on the left). The concentration of CH3COONH4 is 150 mM. Details are discussed in the text. H


Perspective

Analytical Chemistry dissociation than those with the lowest charge states. The protein charge state is assigned on the basis of the methodology described in the previous section. In the sample of the 14 realizations in the parameter space of acidity and temperature that we perform the simulations, we found three possible outcomes in the runs within the simulation time: (i) the 2MRO maintains its conformation and its interface throughout the run, (ii) the two protomers reorient and form a new interface that does not lead to dissociation, and (iii) the protomers reorient and form a new interface that leads to dissociation and at a later time to the droplet fragmentation. Typical snapshots of the protein complex evolution corresponding to scenarios (ii) and (iii) are shown in Figure 9. Snapshots of the droplet fragmentation that follows the complex dissociation is shown in Figure 10.

Figure 10. Typical snapshots of the droplet fragmentation following 2MRO dissociation in an aqueous droplet. The entire system is found in a cavity as shown in Figure 7. The water molecules are shown in light red color (the hydrogen and oxygen sites are not distinguished). (a) Bound 2MRO complex at a charge state of +14 e in an aqueous droplet of approximately 1650 H2O molecules. (b) 2MRO dissociation within the droplet. (c) Droplet fragmentation.

Figure 9. 2MRO protein complex undergoing (a) no dissociation and (b) dissociation after protomer reorientation. The protein complex has charge +14 e and it is solvated in a droplet of ≈1350−1650 H2O molecules. Ubiquitin (Ub) is colored in red and the ubiquitin-associated domain (UbA) in blue. The purple arrows are vectors pointing from the center of mass (COM) of a protomer to the center of the hydrophobic patch of the same protomer. The dashed lines connect the COMs of the two protomers. The time elapsed for each snapshot from an arbitrarily chosen initial time (t = 0) is indicated. For the purpose of clear visualization, the water and ions surrounding the 2MRO are not shown.

The complex dissociation times and the droplet fragmentation times are presented in Table S1. We found that the time required for the protein complex to dissociate, when it occurs, varies between 4 and 6 ns. The droplet fission follows within 2 ns. It is important to note here that the complex dissociation may occur much earlier than 4 ns, depending on the initial state of the complex conformation in the parent bulk solution. The solvent evaporation time shown in Figure 11 for a droplet that contains 2MRO with a charge state of +14 e indicates that the almost complete desolvation of the complex occurs within 6 ns. The comparison between the two times indicates that complex dissociation is likely to occur for the examined set of parameters of temperature of the droplet and protein charge state.

Figure 11. Number of H2O molecules (NH2O) vs time (t) (light red line) of an aqueous droplet containing 2MRO. Typical snapshots of the droplet state are depicted in the course of solvent evaporation (the same color scheme is used as in Figure 9). The total charge of the system is +14 e. The droplet has initially 2000 H2O molecules, and the temperature is 370 K.

in a droplet is surrounded by counterions such as Na+ (or NH4+) ions that form an ionic atmosphere around the macroion. This is an expected analogy to the behavior of dsDNA in the bulk solution. Depending on the charge-squared-to-volume ratio, the droplet containing the dsDNA and the ion atmosphere may be found in the instability regime. The instability manifests itself by the formation of conical solvent protrusions on the surface of the droplet. Even though we have characterized this form of instability in terms of its structure and dynamics, its consequences in the stability of a complex still need to be investigated. The Na+ ions are not ionically bound on the dsDNA backbone, until the latest stage of desolvation. When there is almost no solvent, the ionic atmosphere collapses on the dsDNA affecting its conformation.

■

CONCLUSION AND OUTLOOK In the last three years, we have made significant progress toward understanding the stability of noncovalent complexes in droplets by using atomistic MD. Using nucleic acids as an example, we have demonstrated how the presence of counterions and co-ions in a droplet determines the conformation and the charge state of the complex in the gaseous state. It has been found that dsDNA I


Perspective


Department of Chemistry, University of California at Berkeley, that took place at the Gordon Conference on Gaseous Ions: Structures, Energetics and Reactions, 2017.

We developed a feasible method to study the rate of complex dissociation where we account for the change in the pH of the droplet with its size. We implemented the method in the study of 2MRO in an aqueous droplet. Ultimately, we can find the correction in the equilibrium constant because of complex dissociation in the droplet. An additional result of the study is the finding of the dissociation mechanisms of complexes, which will contribute to understanding critical interactions in biological systems. An interesting and ambitious direction to study in the future is the influence that fissions can have on biomolecule containing drops. Many ESI-MS systems utilize “Z-spray” in which most electrosprayed droplets are not sampled, and it is presumably only droplets resulting from Coulombic fission whose components eventually get measured. Because the droplets derive from the surfaces of evaporating parent droplets, experimental findings suggest that they are often enriched in analytes and solutes.47 Another ambitious goal that we consider feasible in the near future is to find out how ions affect the stability of a complex in a droplet. This emerging use of simulations can assist in a rational selection of the experimental conditions for the analysis of complexes. The current progress of the computational methodologies shows that we are in the path of determining the corrections in the equilibrium constants due to complex dissociation in the droplet. Currently, we study the stability of the assemblies of proteins.

■

■

(1) Balthasart, F.; Plavec, J.; Gabelica, V. J. Am. Soc. Mass Spectrom. 2013, 24, 1−8. (2) Hall, Z.; Robinson, C. V. J. Am. Soc. Mass Spectrom. 2012, 23, 1161−1168. (3) Ogorzalek Loo, R. R.; Lakshmanan, R.; Loo, J. A. J. Am. Soc. Mass Spectrom. 2014, 25, 1675−1693. (4) Lomeli, S. H.; Yin, S.; Ogorzalek Loo, R. R.; Loo, J. A. J. Am. Soc. Mass Spectrom. 2009, 20, 593−596. (5) Smith, R. D.; Light-Wahl, K. J.; Winger, B. E.; Loo, J. A. Org. Mass Spectrom. 1992, 27, 811−821. (6) Barylyuk, K.; Gülbakan, B.; Xie, X.; Zenobi, R. Anal. Chem. 2013, 85, 11902−11912. (7) Hilton, G. R.; Benesch, J. L. P. J. R. Soc., Interface 2012, 9, 801−816. (8) Mehmood, S.; Allison, T. M.; Robinson, C. V. Annu. Rev. Phys. Chem. 2015, 66, 453−474. (9) Loo, J. A. Mass Spectrom. Rev. 1997, 16, 1−23. (10) Erba, E. B.; Zenobi, R. Annu. Rep. Prog. Chem., Sect. C: Phys. Chem. 2011, 107, 199−228. (11) Chen, F.; Gülbakan, B.; Weidmann, S.; Fagerer, S. R.; Ibáñez, A. J.; Zenobi, R. Mass Spectrom. Rev. 2016, 35, 48−70. (12) Erba, E. B.; Barylyuk, K.; Yang, Y.; Zenobi, R. Anal. Chem. 2011, 83, 9251−9259. (13) Gülbakan, B.; Barylyuk, K.; Zenobi, R. Curr. Opin. Biotechnol. 2015, 31, 65−72. (14) Ngounou Wetie, A. G.; Sokolowska, I.; Woods, A. G.; Roy, U.; Loo, J. A.; Darie, C. C. Proteomics 2013, 13, 538−557. (15) Yamaguchi, K. J. Mass Spectrom. 2003, 38, 473−490. (16) Barylyuk, K.; Balabin, R. M.; Grünstein, D.; Kikkeri, R.; Frankevich, V.; Seeberger, P. H.; Zenobi, R. J. Am. Soc. Mass Spectrom. 2011, 22, 1167−1177. (17) Yao, Y.; Richards, M. R.; Kitova, E. N.; Klassen, J. S. J. Am. Soc. Mass Spectrom. 2016, 27, 498−506. (18) Schmidt, C.; Robinson, C. V. FEBS J. 2014, 281, 1950−1964. (19) Acuner Ozbabacan, S. E.; Engin, H. B.; Gursoy, A.; Keskin, O. Protein Eng., Des. Sel. 2011, 24, 635. (20) Light-Wahl, K. J.; Springer, D. L.; Winger, B. E.; Edmonds, C. G.; Camp, D. G.; Thrall, B. D.; Smith, R. D. J. Am. Chem. Soc. 1993, 115, 803−804. (21) Heck, A. J. Nat. Methods 2008, 5, 927. (22) Rose, R. J.; Damoc, E.; Denisov, E.; Makarov, A.; Heck, A. J. Nat. Methods 2012, 9, 1084−1086. (23) Lössl, P.; van de Waterbeemd, M.; Heck, A. J. EMBO J. 2016, 35, 2634−2657. (24) Robinson, C. V.; Chung, E. W.; Kragelund, B. B.; Knudsen, J.; Aplin, R. T.; Poulsen, F. M.; Dobson, C. M. J. Am. Chem. Soc. 1996, 118, 8646−8653. (25) Allen, S. J.; Schwartz, A. M.; Bush, M. F. Anal. Chem. 2013, 85, 12055−12061. (26) Robinson, C. V.; Sali, A.; Baumeister, W. Nature 2007, 450, 973− 982. (27) Pierce, M. M.; Raman, C.; Nall, B. T. Methods 1999, 19, 213−221. (28) Jecklin, M. C.; Schauer, S.; Dumelin, C. E.; Zenobi, R. J. Mol. Recognit. 2009, 22, 319−329. (29) Velazquez-Campoy, A.; Leavitt, S. A.; Freire, E. In Protein-Protein Interactions: Methods and Applications; Fu, H., Ed.; Humana Press: Totowa, NJ, 2004; pp 35−54. (30) Kitova, E. N.; El-Hawiet, A.; Schnier, P. D.; Klassen, J. S. J. Am. Soc. Mass Spectrom. 2012, 23, 431−441. (31) Wortmann, A.; Kistler-Momotova, A.; Zenobi, R.; Heine, M. C.; Wilhelm, O.; Pratsinis, S. E. J. Am. Soc. Mass Spectrom. 2007, 18, 385− 393. (32) Zhou, S.; Prebyl, B. S.; Cook, K. D. Anal. Chem. 2002, 74, 4885− 4888.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.analchem.7b01941. Dissociation rate data (PDF)

■

REFERENCES

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Styliani Consta: 0000-0001-8869-4155 Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS S.C. acknowledges a Marie Curie International Incoming Fellowship Grant Number 628552 of the European Commission held at the Department of Chemistry, University of Cambridge, and NSERC-Discovery grant for supporting this research. M.I.O. acknowledges financial support from the Alexander Graham Bell Canada Graduate Scholarships-Doctoral Program (CGS D), NSERC. Sci-Net, SHARCNET, and Compute Canada are acknowledged for providing the computing facilities. S.C. and M.I.O. thank Prof. D. Frenkel, Department of of Chemistry, University of Cambridge, UK, for the very insightful discussions on the association constants of complexes of macromolecules under confinement and discussions on the stability of charged systems. S.C. and M.I.O. acknowledge stimulating discussions with Prof. D. Russell, Department of of Chemistry, Texas A& M University, Prof. M. Bush, Department of Chemistry, University of Washington, Prof. C. J. Hogan Jr., Department of Mechanical Engineering, University of Michigan, Prof. E. R. Williams, J


Perspective

Analytical Chemistry (33) Zhou, S.; Edwards, A. G.; Cook, K. D.; van Berkel, G. J. Anal. Chem. 1999, 71, 769−776. (34) Consta, S. J. Phys. Chem. B 2010, 114, 5263−5268. (35) Consta, S.; Oh, M. I.; Malevanets, A. Chem. Phys. Lett. 2016, 663, 1−12. (36) Oh, M. I.; Consta, S. J. Phys. Chem. Lett. 2017, 8, 80−85. (37) Labowsky, M.; Fenn, J.; de la Mora, J. F. Anal. Chim. Acta 2000, 406, 105−118. (38) Hogan, C. J.; de la Mora, J. F. J. Am. Soc. Mass Spectrom. 2011, 22, 158−172. (39) Rayleigh, L. Philos. Mag. 1882, 14, 184−186. (40) Consta, S.; Malevanets, A. Mol. Simul. 2015, 41, 73−85. (41) Hendricks, C.; Schneider, J. Am. J. Phys. 1963, 31, 450−453. (42) Iribarne, J. V.; Thomson, B. A. J. Chem. Phys. 1976, 64, 2287− 2294. (43) Thomson, B.; Iribarne, J. J. Chem. Phys. 1979, 71, 4451−4463. (44) Higashi, H.; Tokumi, T.; Hogan, C. J.; Suda, H.; Seto, T.; Otani, Y. Phys. Chem. Chem. Phys. 2015, 17, 15746−15755. (45) Hogan, C. J., Jr; Carroll, J. A.; Rohrs, H. W.; Biswas, P.; Gross, M. L. J. Am. Chem. Soc. 2008, 130, 6926−6927. (46) Allen, S. J.; Schwartz, A. M.; Bush, M. F. Anal. Chem. 2013, 85, 12055−12061. (47) Haddrell, A. E.; Agnes, G. R. Anal. Chem. 2004, 76, 53−61. (48) Gamero-Castano, M.; de la Mora, J. F. Anal. Chim. Acta 2000, 406, 67−91. (49) Hautreux, M.; Hue, N.; de Kerdaniel, A. D. F.; Zahir, A.; Malec, V.; Laprévote, O. Int. J. Mass Spectrom. 2004, 231, 131−137. (50) Consta, S. J. Mol. Struct.: THEOCHEM 2002, 591, 131−140. (51) Consta, S.; Mainer, K. R.; Novak, W. J. Chem. Phys. 2003, 119, 10125−10132. (52) Allen, M.; Tildesley, D. Computer Simulation of Liquids; Clarendon Press: Oxford, 1987. (53) Servage, K. A.; Silveira, J. A.; Fort, K. L.; Clemmer, D. E.; Russell, D. H. J. Phys. Chem. Lett. 2015, 6, 4947−4951. (54) Hutter, J. Wiley Interdisciplinary Reviews: Computational Molecular Science 2012, 2, 604−612. (55) Niklasson, A. M.; Tymczak, C.; Challacombe, M. Phys. Rev. Lett. 2006, 97, 123001. (56) Malevanets, A.; Consta, S. J. Chem. Phys. 2013, 138, 184312. (57) Ryckaert, J.-P.; Ciccotti, G.; Berendsen, H. J. J. Comput. Phys. 1977, 23, 327−341. (58) Vanommeslaeghe, K.; Hatcher, E.; Acharya, C.; Kundu, S.; Zhong, S.; Shim, J.; Darian, E.; Guvench, O.; Lopes, P.; Vorobyov, I.; Mackerell, A. D. J. Comput. Chem. 2010, 31, 671−690. (59) Yu, W.; He, X.; Vanommeslaeghe, K.; MacKerell, A. D. J. Comput. Chem. 2012, 33, 2451−2468. (60) Lindorff-Larsen, K.; Piana, S.; Palmo, K.; Maragakis, P.; Klepeis, J. L.; Dror, R. O.; Shaw, D. E. Proteins: Struct., Funct., Genet. 2010, 78, 1950−1958. (61) Berman, H. M.; Westbrook, J.; Feng, Z.; Gilliland, G.; Bhat, T. N.; Weissig, H.; Shindyalov, I. N.; Bourne, P. E. Nucleic Acids Res. 2000, 28, 235−242. (62) Nowicka, U.; Zhang, D.; Walker, O.; Krutauz, D.; Castañeda, C. A.; Chaturvedi, A.; Chen, T. Y.; Reis, N.; Glickman, M. H.; Fushman, D. Structure 2015, 23, 542−557. (63) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926−935. (64) Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R. D.; Kalé, L.; Schulten, K. J. Comput. Chem. 2005, 26, 1781−1802. (65) Verlet, L. Phys. Rev. 1967, 159, 98. (66) Andersen, H. C. J. Comput. Phys. 1983, 52, 24−34. (67) Sharawy, M.; Consta, S. Phys. Chem. Chem. Phys. 2015, 17, 25550−25562. (68) Soltani, S.; Oh, M. I.; Consta, S. J. Chem. Phys. 2015, 142, 114307. (69) Sharawy, M.; Consta, S. J. Phys. Chem. A 2016, 120, 8871−8880. (70) Consta, S.; Chung, J. K. J. Phys. Chem. B 2011, 115, 10447−10455. (71) Chung, J. K.; Consta, S. J. Phys. Chem. B 2012, 116, 5777−5785. (72) Sharawy, M.; Consta, S. J. Chem. Phys. 2014, 141, 104321.

(73) Consta, S.; Malevanets, A. Phys. Rev. Lett. 2012, 109, 148301. (74) Koopman, E. A.; Lowe, C. P. J. Chem. Phys. 2006, 124, 204103. (75) Vega, C.; de Miguel, E. J. Chem. Phys. 2007, 126, 154707. (76) Favre, A.; Gonnet, F.; Tabet, J.-C. Int. J. Mass Spectrom. 1999, 190−191, 303−312. (77) Favre, A.; Gonnet, F.; Tabet, J.-C. Eur. Mass Spectrom. 2000, 6, 389. (78) Gale, D. C.; Smith, R. D. J. Am. Soc. Mass Spectrom. 1995, 6, 1154− 1164. (79) Hofstadler, S. A.; Griffey, R. H. Chem. Rev. 2001, 101, 377−390. (80) Gabelica, V.; de Pauw, E.; Rosu, F. J. Mass Spectrom. 1999, 34, 1328−1337. (81) Anandakrishnan, R.; Aguilar, B.; Onufriev, A. Nucleic Acids Res. 2012, 40, W537. (82) Antosiewicz, J.; McCammon, J.; Gilson, M. Biochemistry 1996, 35, 7819−7833. (83) Baptista, A.; Teixeira, V.; Soares, C. J. Chem. Phys. 2002, 117, 4184−4200. (84) Bas, D.; Rogers, D.; Jensen, J. Proteins: Struct., Funct., Genet. 2008, 73, 765−783. (85) Olsson, M.; SØndergaard, C.; Rostkowski, M.; Jensen, J. J. Chem. Theory Comput. 2011, 7, 525−537. (86) Chen, J.; Brooks, C., III; Khandogin, J. Curr. Opin. Struct. Biol. 2008, 18, 140−148. (87) Cramer, C.; Truhlar, D. Acc. Chem. Res. 2008, 41, 760−768. (88) Gilson, M. Curr. Opin. Struct. Biol. 1995, 5, 216−223. (89) Khandogin, J.; Brooks, C., III Biophys. J. 2005, 89, 141−157. (90) Schreiber, G.; Haran, G.; Zhou, H.-X. Chem. Rev. 2009, 109, 839− 860. (91) Warshel, A.; Sharma, P.; Kato, M.; Parson, W. Biochim. Biophys. Acta, Proteins Proteomics 2006, 1764, 1647−1676.

K


Advances in Modeling the Stability of Noncovalent Complexes in

Recommend Documents