Article pubs.acs.org/JPCC
Reconstructing Multiwell Potentials with Steep Gradients Using Stochastically Excited Spring Probes Ahmad Haider,†,‡ Daniel Potter,†,‡ and Todd A. Sulchek*,†,‡ †
Woodruff School of Mechanical Engineering, and ‡Parker H. Petit Institute for Bioengineering and Biosciences, Georgia Institute of Technology, Atlanta, Georgia 30332, United States S Supporting Information *
ABSTRACT: Measurements of free energy landscapes are critical for understanding the basis of many physical, chemical, and biological interactions. Statistical mechanics provide exact equations to calculate free energies, but are built on the assumption that all possible configurations of the system are sampled. The most pronounced limit to accurate free energy computations is therefore the imperfect sampling of a potential field, particularly in the case of interactions with steep gradients and short reaction coordinates. We show through simulations that increasing the stochastic fluctuations of a harmonic probe by active excitation results in increased sampling times of high gradient adhesive interactions and leads to the reconstruction of a more accurate energy landscape. We use Brownian dynamics simulations to test the impact of probe approach velocity, stiffness, and thermal energy to sample complex energy landscapes with multiple wells of various depths and slopes to understand the accuracy of energy surface reconstruction. We then show experimentally that through the application of optimal stochastic excitations, we are able to obtain accurate energy landscape reconstructions for different probe and landscape parameters due to improved sampling of previously poorly probed interactions.
■
INTRODUCTION One of the most complete quantitative descriptions of intermolecular and interfacial interactions is the free energy landscape. Free energy landscapes can represent the energetics of a variety of systems of interest, such as protein folding,1,2 ligand receptor binding,3,4 nucleic acid base flipping,5 complex conformational changes in macromolecules,6 and boundary lubrication and adhesion in tribology.7 Energy landscapes are also commonly employed in computational studies of processes including hydrophobic interactions,8 organic reactions in water,9 proton transfer,10 ionic permeation through membrane channels,11 and peptide12 and protein13 equilibria. Energy landscapes therefore provide a valuable framework to describe the behavior of complicated systems in a straightforward and statistically rigorous manner. The energy landscape, G(x), is defined up to a constant C by inverting the Boltzmann relation: G(x) = −kBT ln(P(x))) + C
positional values of the reaction coordinate into histograms to obtain an approximation of the true probability distribution.23−25 Alternatively, energy landscapes can also be calculated from repeated nonequilibrium measurements by determining the distribution of values of the mechanical work used to drive the system and transforming it into free energy.26−28 A limitation of experimental energy reconstructions of steep, adhesive landscapes is the lack of sufficient data for a highquality reconstruction at all points along the reaction coordinate. This limitation stems from a finite sampling time23 due to sharp energy landscape barriers and probe dynamics, which do not accurately track energy gradients in equilibrium experiments.26 In the equilibrium regime, the system is typically seen to be in a bistable state, driven by thermal fluctuations across the energy barrier between completely bound and fully unbound states.29 For many realistic interactions, such as those with energy surfaces with high barriers, the transitions across the barrier may be so rare as to be impossible to achieve sufficient equilibrium sampling in the time frame of an experiment. Nonequilibrium force ramp experiments are one way to increase sampling in a given time frame30,31 and utilize fluctuations theorems such as the Jarzynski equality27 to convert the nonequilibrium sampled
(1)
where kBT is the unit of thermal energy and P(x) is the probability distribution with respect to a reaction coordinate of interest, x. The reaction coordinate14,15 represents a variable that can be used to follow the progress of an interaction, such as intermolecular distance or angle. Experimental force measurements, for example, with atomic force microscopy (AFM)16−19 or with optical tweezers,20−22 provide the data for reconstruction of free energy landscapes. Free energy landscapes are reconstructed from experiments performed under equilibrium conditions by binning the © 2017 American Chemical Society
Received: December 24, 2016 Revised: February 28, 2017 Published: March 14, 2017 7248
DOI: 10.1021/acs.jpcc.6b12961 J. Phys. Chem. C 2017, 121, 7248−7258
Article
The Journal of Physical Chemistry C
Brownian Dynamics Simulations. Overdamped Brownian dynamics (Langevin) simulations are used to model the behavior of a fluctuating cantilever probe under the influence of a potential landscape. At any given instant, the probe is subjected to fluctuating forces due to random impulses from surrounding fluid molecules as well as the forces due to the external potential field, while being constrained by the restoring force of the probe itself. Hence, the position xi of the tip at time ti during the ith step of the simulation can be given by the following equation:38−40
work distributions to equilibrium free energies. However, the slow convergence of the Jarzynski equality with the increase in the number of samples32 limits its applicability to near equilibrium fluctuation regimes only. As a thermal probe traverses an adhesive energy landscape, it can encounter multiple energy barriers and wells. If the barrier curvature is greater than the probe stiffness, the probe will be quickly pulled to the bottom of its adjacent well without sufficient time to sample the transition region of the landscape from barrier to well bottom. Once it reaches the bottom of the energy well, the limited thermal energy fluctuations of the probe are not sufficient to hop over the adjacent energy barrier and escape the well. Hence, the probe spends most of its time at the bottom of the well in the lowest energy state, resulting in a limited sampling of the entire energy landscape. Hoh et al.33 used white noise perturbations to the thermally driven motion of an AFM cantilever to probe adhesive tip− sample interactions. External white noise excitation increases the cantilever’s magnitude of fluctuations, and thereby its reaction coordinate sampling range at any given point on the coordinate. By combining this approach of enhanced noise to increase the sampling range with a method to consistently transform fluctuations to energies for the entire reaction coordinate using the weighted histogram analysis method (WHAM),34−37 the entire energy landscape can be reconstructed.31 We have previously demonstrated in experimental studies that enhanced stochastic fluctuations can fully determine the depth and width of adhesive energy landscapes with steep gradient interactions in a silicon nitride−water−mica system. In this Article, we use damped Brownian dynamics simulations to study the influence of enhanced stochastic sampling to reconstruct multiwell free energy landscapes by performing a parametric sweep of the energy barrier height, barrier slope, probe stiffness, and probe velocity. Using the simulations, we are able to determine both the parametric trends and the optimal excitations required to best sample an interaction and minimize energy reconstruction errors. Finally, we demonstrate the applicability of our simulation results with AFM experiments conducted on a silicon nitride−water−mica system at different probe velocities and validate that optimal probe velocities for best interaction sampling indeed follow the same trend as those predicted by the simulations. This technique may be used to analyze rough energy landscapes typical in many biophysically interesting systems without incorrect energy reconstruction estimations due to finite sampling errors.
xi = xi − 1 −
(3)
Figure 1. (A) A symmetric bimodal Lorentzian energy profile with overlapping wells. (B) Force curve generated by an AFM probe tip while traversing the energy landscape. The base of the probe is translated in the positive z-direction. Inset on top shows a section of the force curve when the cantilever enters the energy well, and inset on bottom shows another section of the force curve when the cantilever jumps between the wells typically observed in an adhesion/deadhesion event.
METHODS Energy Landscape Model. We model a complex energy landscape E(x) consisting of two overlapping wells with a onedimensional Lorentzian profile, of the form ⎤ ⎥ ⎥ ⎥⎦
⎡ d E (x ) ⎤ 2DΔt wi − ⎢ ⎥ ηΔt ⎣ dx ⎦i − 1
where η = 5 × 10−5 pN·s/nm is the friction coefficient, Δt = 5 ns is the time step of the simulation, k is the cantilever’s spring constant, wi is a Gaussian random number with a mean value of 0 and a variance of 1, and D is the diffusion coefficient related to the friction coefficient and thermal energy of the cantilever by Einstein’s relation ηD = kBT. 400 time steps are used between each recorded data point, which is obtained by dividing the data collection rate of 50 kHz (typical for experimental AFM systems) by the simulation time step of 5 ns. The cantilever spring, which superimposes a harmonic biasing potential to the underlying energy landscape, sweeps the reaction coordinate from −0.5 to 0.5 μm for the duration of the simulation at a constant velocity. Figure 1 shows a representative sample of the force curve generated by the simulation as the tip sweeps over the energy
■
⎡ ⎢ A1 A2 + E (x ) = − ⎢ (x − μ1) (x − μ2 ) 1+ ⎢⎣ 1 + σ 2 σ22 1
k xi − 1Δt + η
landscape. As the tip traverses the energy landscape (from bottom to top in Figure 1A), it experiences negligible interaction forces until −80 nm at which point the tip enters into the first well. Hence, on account of negligible tip forces, the tip and base of the cantilever on average move synchronously from −250 to −80 nm as evidenced by the 45° line in Figure 1B. From −80 to −30 nm, the tip experiences a steep potential gradient, and hence the tip’s downward movement surpasses the Z-movement of the base, as seen as a sharp transition with force curve gradient greater than 45° in the upper inset of Figure 1B. Thereafter, the tip moves to the
(2)
where A1 and A2 correspond to well depths, μ1 and μ2 are the well location parameters, and σ1 and σ2 are the well scaling parameters. We simplify the calculations by making the landscape symmetric with the zero point of the reaction coordinate centrally placed between the two wells. A range of values of μ and σ are employed to simulate different landscapes. 7249
DOI: 10.1021/acs.jpcc.6b12961 J. Phys. Chem. C 2017, 121, 7248−7258
Article
The Journal of Physical Chemistry C
Figure 2. Schematic of an AFM tip approaching a symmetric bimodal energy landscape separated by a small barrier (A) and a large barrier (B,C) in the direction indicated by the arrow. The probability density function (PDF) histograms of the detected positions of the tip are shown for the no excitation case (blue), ideal excitation case (red), and over excitation case (black). Green dotted line depicts the minimum number of detected tip positions needed to sample the reaction coordinate (sampling threshold).
(WHAM)43,44 were employed to iteratively solve for the interfacial potential distributions, Pinteraction(x). Finally, the inverse Boltzmann relation was used to convert the interfacial potential probability distributions into the interfacial energy landscape, Ginteraction(x) = −kBT ln(Pinteraction(x)) + C. The reconstructions were carried out with 5000 point size windows and 0.05 nm bin widths using custom routines written in MATLAB (MathWorks Inc., Natick, MA). Because of imperfect sampling, two types of errors can be encountered in a reconstruction. The first type of error results from the probe not sufficiently equilibrating when sampling the reaction coordinate, leading to errors when the Boltzmann relation is used to transform these nonequilibrium fluctuations to free energies. In the second type of error, the reaction coordinate of the given energy landscape is not fully traversed such that there are no reconstructed energy values for certain reaction coordinates, leading to reconstruction errors when extrapolating the free energies to these undersampled regions.23 Techniques such as series expansion or bootstrapping45,46 may not be effective to reconstruct free energies in these undersampled regions because too much data are missing to enable a bootstrap to fill in the missing information.47 The first type of error is handled by using a root-mean-square error (RMSE) measure to quantify the energy difference between the reconstructed energy and the underlying energy surface. To account for the second type of error, we adjust the calculated RMSE by dividing it by a normalized range of reaction coordinate x that is being traversed by the reconstruction as shown in eq 4.
bottom of the energy well and attains a minimum energy state whereupon it fluctuates within a small range of the first well, exhibiting a shallow force curve gradient 1 N/m may be used for these studies, but choosing such stiff cantilevers, as compared to the energy landscape curvatures, will possibly result in low deflection sensitivity of the tip such that the tip fluctuations will fall below the instrument noise level, thereby obscuring the true tip sample interaction information, and producing inaccurate reconstructions. We believe stochastic excitation may provide a way to sufficiently sample the energy well of biological interactions by allowing the probe to rapidly escape the steep gradient well. We can identify two pitfalls in reconstructing such biomolecular interactions. First, biological molecules themselves are soft and deformable, and so the interaction force may not directly transfer to the cantilever probe through the soft deformable body of the molecule. We expect that using very 7256
DOI: 10.1021/acs.jpcc.6b12961 J. Phys. Chem. C 2017, 121, 7248−7258
Article
The Journal of Physical Chemistry C
■
(11) Allen, T. W.; Andersen, O.; Roux, B. Energetics of Ion Conduction through the Gramicidin Channel. Proc. Natl. Acad. Sci. U. S. A. 2004, 86, 117−122. (12) Karplus, M.; Pettitt, B. M. The Potential of Mean Force Surface for the Alanine Dipeptide in Aqueous Solution. Chem. Phys. Lett. 1985, 121, 194−201. (13) Guo, Z.; Boczko, E. M.; Brooks, C. L. Exploring the Folding Free-Energy Surface of a Three-Helix Bundle Protein. Proc. Natl. Acad. Sci. U. S. A. 1997, 94, 10161−10166. (14) Fukui, K. Formulation of the Reaction Coordinate. J. Phys. Chem. 1970, 74, 4161−4163. (15) Frauenfelder, H.; Sligar, S. G.; Wolynes, P. G. The Energy Landscapes and Motions of Proteins. Science 1991, 254, 1598−603. (16) Rief, M.; Gautel, M.; Oesterhelt, F.; Fernandez, J. M.; Gaub, H. E. Reversible Unfolding of Individual Titin Immunoglobulin Domains by Afm. Science 1997, 276, 1109−12. (17) Bustamante, C.; Macosko, J. C.; Wuite, G. J. Grabbing the Cat by the Tail: Manipulating Molecules One by One. Nat. Rev. Mol. Cell Biol. 2000, 1, 130−6. (18) Fisher, T. E.; Marszalek, P. E.; Fernandez, J. M. Stretching Single Molecules into Novel Conformations Using the Atomic Force Microscope. Nat. Struct. Biol. 2000, 7, 719−24. (19) Engel, A.; Muller, D. J. Observing Single Biomolecules at Work with the Atomic Force Microscope. Nat. Struct. Biol. 2000, 7, 715−8. (20) Ashkin, A. Optical Trapping and Manipulation of Neutral Particles Using Lasers; World Scientific: Hackensack, NJ, 2006; Vol. xxiv, 915 pp. (21) Clausen-Schaumann, H.; Seitz, M.; Krautbauer, R.; Gaub, H. E. Force Spectroscopy with Single Bio-Molecules. Curr. Opin. Chem. Biol. 2000, 4, 524−30. (22) Neuman, K. C.; Block, S. M. Optical Trapping. Rev. Sci. Instrum. 2004, 75, 2787−809. (23) Wood, R. H.; MuhlBauer, W. C. F.; Thompson, P. Systematic Errors in Free Energy Perturbation Calculations Due to a Finite Sample of Configuration Space. J. Phys. Chem. 1991, 95, 6670−6675. (24) Pohorille, A.; Jarzynski, C.; Chipot, C. Good Practices in FreeEnergy Calculations. J. Phys. Chem. B 2010, 114, 10235−53. (25) Liu, P.; Dehez, F.; Cai, W.; Chipot, C. A Toolkit for the Analysis of Free-Energy Perturbation Calculations. J. Chem. Theory Comput. 2012, 8, 2606−16. (26) Harris, N. C.; Song, Y.; Kiang, C. H. Experimental Free Energy Surface Reconstruction from Single-Molecule Force Spectroscopy Using Jarzynski’s Equality. Phys. Rev. Lett. 2007, 99, 068101. (27) Jarzynski, C. Nonequilibrium Inequality for Free Energy Differences. Phys. Rev. Lett. 1997, 78, 2690−2693. (28) Preiner, J.; Janovjak, H.; Rankl, C.; Knaus, H.; Cisneros, D. A.; Kedrov, A.; Kienberger, F.; Muller, D. J.; Hinterdorfer, P. Free Energy of Membrane Protein Unfolding Derived from Single-Molecule Force Measurements. Biophys. J. 2007, 93, 930−7. (29) Mayor, U.; Johnson, C. M.; Daggett, V.; Fersht, A. R. Protein Folding and Unfolding in Microseconds to Nanoseconds by Experiment and Simulation. Proc. Natl. Acad. Sci. U. S. A. 2000, 97, 13518−22. (30) Sulchek, T.; Friddle, R. W.; Noy, A. Strength of Multiple Parallel Biological Bonds. Biophys. J. 2006, 90, 4686−91. (31) Haider, A.; Potter, D.; Sulchek, T. Enhanced Stochastic Fluctuations to Measure Steep Adhesive Energy Landscapes. Proc. Natl. Acad. Sci. U. S. A. 2016, 113, 1−6. (32) Gore, J.; Ritort, F.; Bustamante, C. Bias and Error in Estimates of Equilibrium Free-Energy Differences from Nonequilibrium Measurements. Proc. Natl. Acad. Sci. U. S. A. 2003, 100, 12564−12569. (33) Koralek, D. O.; Heinz, W. F.; Antonik, M. D.; Baik, A.; Hoh, J. H. Probing Deep Interaction Potentials with White-Noise-Driven Atomic Force Microscope Cantilevers. Appl. Phys. Lett. 2000, 76, 2952−2954. (34) Zhu, F.; Hummer, G. Convergence and Error Estimation in Free Energy Calculations Using the Weighted Histogram Analysis Method. J. Comput. Chem. 2012, 33, 453−65.
ASSOCIATED CONTENT
S Supporting Information *
The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcc.6b12961. K-Means clustering methodology used to partition the approach speed into distinct clusters (Figure S1); and reconstructed free energies of silicon nitride−water− mica interactions at varying voltages as compared to those predicted by a DFS model to identify the optimal excitation voltage (Figure S2) (PDF)
■
AUTHOR INFORMATION
Corresponding Author
*Tel.: 4043851887. Fax: 4048948496. E-mail: todd.sulchek@ me.gatech.edu. ORCID
Ahmad Haider: 0000-0002-4996-747X Daniel Potter: 0000-0002-4913-3934 Todd A. Sulchek: 0000-0003-4196-6293 Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS
■
REFERENCES
We thank Raymond Friddle (Sandia National Laboratory), Paul Ashby (Lawrence Berkeley National Laboratory), and Hermann Gruber (University of Linz) for helpful discussion and comments; and Aaron Enten for expertise and assistance with bandpass filter design. This work was supported by the National Science Foundation (CBET-CAREER-1055437).
(1) Manuel, A. P.; Lambert, J.; Woodside, M. T. Reconstructing Folding Energy Landscapes from Splitting Probability Analysis of Single-Molecule Trajectories. Proc. Natl. Acad. Sci. U. S. A. 2015, 112, 7183−8. (2) Shea, J. E.; Brooks, C. L., 3rd From Folding Theories to Folding Proteins: A Review and Assessment of Simulation Studies of Protein Folding and Unfolding. Annu. Rev. Phys. Chem. 2001, 52, 499−535. (3) Woo, H. J.; Roux, B. Calculation of Absolute Protein-Ligand Binding Free Energy from Computer Simulations. Proc. Natl. Acad. Sci. U. S. A. 2005, 102, 6825−30. (4) Steinbrecher, T.; Labahn, A. Towards Accurate Free Energy Calculations in Ligand Protein-Binding Studies. Curr. Med. Chem. 2010, 17, 767−85. (5) Huang, N.; Banavali, N. K.; MacKerell, A. D., Jr. ProteinFacilitated Base Flipping in DNA by Cytosine-5-Methyltransferase. Proc. Natl. Acad. Sci. U. S. A. 2003, 100, 68−73. (6) Wereszczynski, J.; Andricioaei, I. On Structural Transitions, Thermodynamic Equilibrium, and the Phase Diagram of DNA and Rna Duplexes under Torque and Tension. Proc. Natl. Acad. Sci. U. S. A. 2006, 103, 16200−5. (7) Bhushan, B. Handbook of Micro/Nano Tribology, 1st ed.; CRC Press: Boca Raton, FL, 1995. (8) Perunov, N.; England, J. L. Quantitative Theory of Hydrophobic Effect as a Driving Force of Protein Structure. Protein Sci. 2014, 23, 387−99. (9) Gao, J. L. A Priori Computation of a Solvent-Enhanced Sn2 Reaction Profile in Water - the Menshutkin Reaction. J. Am. Chem. Soc. 1991, 113, 7796−7797. (10) Hinsen, K.; Roux, B. Potential of Mean Force and Reaction Rates for Proton Transfer in Acetylacetone. J. Chem. Phys. 1997, 105, 3567−3577. 7257
DOI: 10.1021/acs.jpcc.6b12961 J. Phys. Chem. C 2017, 121, 7248−7258
Article
The Journal of Physical Chemistry C (35) Roux, B. The Calculation of the Potential of Mean Force Using Computer Simulations. Comput. Phys. Commun. 1995, 91, 275−282. (36) Souaille, M.; Roux, B. Extension to the Weighted Histogram Analysis Method: Combining Umbrella Sampling with Free Energy Calculations. Comput. Phys. Commun. 2001, 135, 40−57. (37) Kumar, S.; Rosenberg, J. M.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A. The Weighted Histogram Analysis Method for FreeEnergy Calculations on Biomolecules. I. The Method. J. Comput. Chem. 1992, 13, 1011−1021. (38) Volpe, G.; Volpe, G. Simulation of a Brownian Particle in an Optical Trap. Am. J. Phys. 2013, 81, 224−230. (39) de Grooth, B. G. A Simple Model for Brownian Motion Leading to the Langevin Equation. Am. J. Phys. 1999, 67, 1248−1252. (40) La Porta, A.; Denesyuk, N. A.; de Messieres, M. Optimal Reconstruction of the Folding Landscape Using Differential Energy Surface Analysis. Phys. Rev. E 2013, 87, 87. (41) Torrie, G. M.; Valleau, J. P. Non-Physical Sampling Distributions in Monte-Carlo Free-Energy Estimation - Umbrella Sampling. J. Comput. Phys. 1977, 23, 187−199. (42) Torrie, G. M.; Valleau, J. P. Monte-Carlo Study of a PhaseSeparating Liquid-Mixture by Umbrella Sampling. J. Chem. Phys. 1977, 66, 1402−1408. (43) Kumar, S.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A.; Rosenberg, J. M. The Weighted Histogram Analysis Method for FreeEnergy Calculations on Biomolecules 0.1. The Method. J. Comput. Chem. 1992, 13, 1011−1021. (44) Roux, B. The Calculation of the Potential of Mean Force Using Computer-Simulations. Comput. Phys. Commun. 1995, 91, 275−282. (45) Kofke, D. A.; Cummings, P. T. Quantitative Comparison and Optimization of Methods for Evaluating the Chemical Potential by Molecular Simulation. Mol. Phys. 1997, 92, 973−996. (46) Kofke, D. A.; Cummings, P. T. Precision and Accuracy of Staged Free-Energy Perturbation Methods for Computing the Chemical Potential by Molecular Simulation. Fluid Phase Equilib. 1998, 150, 41− 49. (47) Lu, N.; Kofke, D. A.; Woolf, T. B. Improving the Efficiency and Reliability of Free Energy Perturbation Calculations Using Overlap Sampling Methods. J. Comput. Chem. 2004, 25, 28−39. (48) Butt, H.-J.; Jaschke, M. Calculation of Thermal Noise in Atomic Force Microscopy. Nanotechnology 1995, 6, 1−7. (49) Stark, R. W.; Drobek, T.; Heckl, W. M. Thermomechanical Noise of a Free V-Shaped Cantilever for Atomic-Force Microscopy. Ultramicroscopy 2001, 86, 207−15. (50) Ashby, P. D.; Lieber, C. M. Brownian Force Profile Reconstruction of Interfacial 1-Nonanol Solvent Structure. J. Am. Chem. Soc. 2004, 126, 16973−80. (51) Steinley, D. K-Means Clustering: A Half-Century Synthesis. Br J. Math Stat Psychol 2006, 59, 1−34. (52) Hastie, T.; Tibshirani, R.; Friedman, J. The Elements of Statistical Learning; Springer: New York, 2001. (53) Friddle, R. W.; Noy, A.; De Yoreo, J. J. Interpreting the Widespread Nonlinear Force Spectra of Intermolecular Bonds. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 13573−8.
7258
DOI: 10.1021/acs.jpcc.6b12961 J. Phys. Chem. C 2017, 121, 7248−7258
