The Poisson-Boltzmann equation for biomolecular electrostatics: a tool for structural biology

Interfaces and hydrophobic interactions in receptor-ligand systems: A level-set variational implicit solvent approach

A model nanometer-sized hydrophobic receptor-ligand system in aqueous solution is studied by the recently developed level-set variational implicit solvent model (VISM). This approach is compared to all-atom computer simulations. The simulations reveal complex hydration effects within the (concave) receptor pocket, sensitive to the distance of the (convex) approaching ligand. The ligand induces and controls an intermittent switching between dry and wet states of the hosting pocket, which determines the range and magnitude of the pocket-ligand attraction. In the level-set VISM, a geometric free-energy functional of all possible solute-solvent interfaces coupled to the local dispersion potential is minimized numerically. This approach captures the distinct metastable states that correspond to topologically different solute-solvent interfaces, and thereby reproduces the bimodal hydration behavior observed in the all-atom simulation. Geometrical singularities formed during the interface relaxation are found to contribute significantly to the energy barrier between different metastable states. While the hydration phenomena can thus be explained by capillary effects, the explicit inclusion of dispersion and curvature corrections seems to be essential for a quantitative description of hydrophobically confined systems on nanoscales. This study may shed more light onto the tight connection between geometric and energetic aspects of biomolecular hydration and may represent a valuable step toward the proper interpretation of experimental receptor-ligand binding rates.

Predicting pKa

One of the most important physicochemical properties of small molecules and macromolecules are the dissociation constants for any weakly acidic or basic groups, generally expressed as the pK(a) of each group. This is a major factor in the pharmacokinetics of drugs and in the interactions of proteins with other molecules. For both the protein and small molecule cases, we survey the sources of experimental pK(a) values and then focus on current methods for predicting them. Of particular concern is an analysis of the scope, statistical validity, and predictive power of methods as well as their accuracy.

BioDrugScreen: a computational drug design resource for ranking molecules docked to the human proteome

BioDrugScreen is a resource for ranking molecules docked against a large number of targets in the human proteome. Nearly 1600 molecules from the freely available NCI diversity set were docked onto 1926 cavities identified on 1589 human targets resulting in >3 million receptor-ligand complexes requiring >200,000 cpu-hours on the TeraGrid. The targets in BioDrugScreen originated from Human Cancer Protein Interaction Network, which we have updated, as well as the Human Druggable Proteome, which we have created for the purpose of this effort. This makes the BioDrugScreen resource highly valuable in drug discovery. The receptor-ligand complexes within the database can be ranked using standard and well-established scoring functions like AutoDock, DockScore, ChemScore, X-Score, GoldScore, DFIRE and PMF. In addition, we have scored the complexes with more intensive GBSA and PBSA approaches requiring an additional 120,000 cpu-hours on the TeraGrid. We constructed a simple interface to enable users to view top-ranking molecules and access purchasing and other information for further experimental exploration.

Electric dipole reorientation in the interaction of botulinum neurotoxins with neuronal membranes

Botulinum neurotoxins are highly potent toxins capable of rapid and specific interaction with the presynaptic membrane. We have hypothesised that: (1) these neurotoxins possess an electric dipole with the positive pole on receptor binding domain Hc-C and that (2) on approaching the negatively charged presynaptic membrane, they reorient themselves and hit the membrane surface with Hc-C; this electrostatic effect would contribute efficient binding. Electrostatic calculations confirm these hypotheses and strongly indicate that electrostatics effects can play an important role in the unique presynaptic membrane binding properties of these neurotoxins and generally on the interaction of other plasma membrane protein ligands.

Positively charged residues in DNA-binding domains of structural proteins follow sequence-specific positions of DNA phosphate groups

We study electrostatic charge complementarity along interfaces of DNA-protein complexes. We use the Protein Data Bank atomic coordinates of DNA-protein complexes for some DNA-binding proteins to study the distribution of positively charged protein residues in the close contact with DNA. We show that large structural proteins reveal a peculiar nonuniform distribution of Arg, Lys, and His amino acids in the frame of negatively charged DNA phosphate strands. We study the nucleosome core particles, DNA complexes with prokaryotic DNA-bending histone analogues, but also the basic binding motifs of small DNA-binding proteins. For large DNA-protein complexes, where extensive DNA wrapping around protein cores occurs, we show that positive amino acids on the proteins track sequence-specific positions of individual DNA phosphates. This specificity of electrostatic interactions can contribute to DNA recognition by DNA-binding proteins, which is governed for many DNA-protein complexes primarily by the hydrogen bond formation between protein residues and DNA bases.

Recent advances in computer-aided drug design

Modern drug discovery is characterized by the production of vast quantities of compounds and the need to examine these huge libraries in short periods of time. The need to store, manage and analyze these rapidly increasing resources has given rise to the field known as computer-aided drug design (CADD). CADD represents computational methods and resources that are used to facilitate the design and discovery of new therapeutic solutions. Digital repositories, containing detailed information on drugs and other useful compounds, are goldmines for the study of chemical reactions capabilities. Design libraries, with the potential to generate molecular variants in their entirety, allow the selection and sampling of chemical compounds with diverse characteristics. Fold recognition, for studying sequence-structure homology between protein sequences and structures, are helpful for inferring binding sites and molecular functions. Virtual screening, the in silico analog of high-throughput screening, offers great promise for systematic evaluation of huge chemical libraries to identify potential lead candidates that can be synthesized and tested. In this article, we present an overview of the most important data sources and computational methods for the discovery of new molecular entities. The workflow of the entire virtual screening campaign is discussed, from data collection through to post-screening analysis.

Protein-ion binding process on finite macromolecular concentration. A Poisson-Boltzmann and Monte Carlo study

Electrostatic interactions are one of the key driving forces for protein-ligands complexation. Different levels for the theoretical modeling of such processes are available on the literature. Most of the studies on the Molecular Biology field are performed within numerical solutions of the Poisson-Boltzmann Equation and the dielectric continuum models framework. In such dielectric continuum models, there are two pivotal questions: (a) how the protein dielectric medium should be modeled, and (b) what protocol should be used when solving this effective Hamiltonian. By means of Monte Carlo (MC) and Poisson-Boltzmann (PB) calculations, we define the applicability of the PB approach with linear and nonlinear responses for macromolecular electrostatic interactions in electrolyte solution, revealing some physical mechanisms and limitations behind it especially due the raise of both macromolecular charge and concentration out of the strong coupling regime. A discrepancy between PB and MC for binding constant shifts is shown and explained in terms of the manner PB approximates the excess chemical potentials of the ligand, and not as a consequence of the nonlinear thermal treatment and/or explicit ion-ion interactions as it could be argued. Our findings also show that the nonlinear PB predictions with a low dielectric response well reproduce the pK shifts calculations carried out with an uniform dielectric model. This confirms and completes previous results obtained by both MC and linear PB calculations.

Ions and inhibitors in the binding site of HIV protease: comparison of Monte Carlo simulations and the linearized Poisson-Boltzmann theory

Proteins can be influenced strongly by the electrolyte in which they are dissolved, and we wish to model, understand, and ultimately control such ionic effects. Relatively detailed Monte Carlo (MC) ion simulations are needed to capture biologically important properties of ion channels, but a simpler treatment of ions, the linearized Poisson-Boltzmann (LPB) theory, is often used to model processes such as binding and folding, even in settings where the LPB theory is expected to be inaccurate. This study uses MC simulations to assess the reliability of the LPB theory for such a system, the constrained, anionic active site of HIV protease. We study the distributions of ions in and around the active site, as well as the energetics of displacing ions when a protease inhibitor is inserted into the active site. The LPB theory substantially underestimates the density of counterions in the active site when divalent cations are present. It also underestimates the energy cost of displacing these counterions, but the error is not consequential because the energy cost is less than kBT, according to the MC calculations. Thus, the LPB approach will often be suitable for studying energetics, but the more detailed MC approach is critical when ionic distributions and fluxes are at issue.

Molecular modeling and simulation studies of ion channel structures, dynamics and mechanisms

Ion channels are integral membrane proteins that enable selected ions to flow passively across membranes. Channel proteins have been the focus of computational approaches to relate their three-dimensional (3D) structure to their physiological function. We describe a number of computational tools to model ion channels. Homology modeling may be used to construct structural models of channels based on available X-ray structures. Electrostatics calculations enable an approximate evaluation of the energy profile of an ion passing through a channel. Molecular dynamics simulations and free-energy calculations provide information on the thermodynamics and kinetics of channel function.

Properties of the Nucleic-acid Bases in Free and Watson-Crick Hydrogen-bonded States: Computational Insights into the Sequence-dependent Features of Double-helical DNA

The nucleic-acid bases carry structural and energetic signatures that contribute to the unique features of genetic sequences. Here we review the connection between the chemical structure of the constituent nucleotides and the polymeric properties of DNA. The sequence-dependent accumulation of charge on the major- and minor-groove edges of the Watson-Crick base pairs, obtained from ab initio calculations, presents unique motifs for direct sequence recognition. The optimization of base interactions generates a propellering of base-pair planes of the same handedness as that found in high-resolution double-helical structures. The optimized base pairs also deform along conformational pathways, i.e., normal modes, of the same type induced by the binding of proteins. Empirical energy computations that incorporate the properties of the base pairs account satisfactorily for general features of the next level of double-helical structure, but miss key sequence-dependent differences in dimeric structure and deformability. The latter discrepancies appear to reflect factors other than intrinsic base-pair structure.

Defining the role of salt bridges in protein stability

Although the energetic balance of forces stabilizing proteins has been established qualitatively over the last decades, quantification of the energetic contribution of particular interactions still poses serious problems. The reasons are the strong cooperativity and the interdependence ofnoncovalent interactions. Salt bridges are a typical example. One expects that ionizable side chains frequently form ion pairs in innumerable crystal structures. Since electrostatic attraction between opposite charges is strong per se, salt bridges can intuitively be regarded as an important factor stabilizing the native structure. Is that really so? In this chapter we critically reassess the available methods to delineate the role ofelectrostatic interactions and salt bridges to protein stability, and discuss the progress and the obstacles in this endeavor. The basic problem is that formation of salt bridges depends on the ionization properties of the participating groups, which is significantly influenced by the protein environment. Furthermore, salt bridges experience thermal fluctuations, continuously break and re-form, and their lifespan in solution is governed by the flexibility of the protein. Finally, electrostatic interactions are long-range and might be significant in the unfolded state, thus seriously influencing the energetic profile. Elimination of salt bridges by protonation/deprotonation at extreme pH or by mutation provides only rough energetic estimates, since there is no way to account for the nonadditive response of the protein moiety. From what we know so far, the strength of electrostatic interactions is strongly context-dependent, yet it is unlikely that salt bridges are dominant factors governing protein stability. Nevertheless, proteins from thermophiles and hyperthermophiles exhibit more, and frequently networked, salt bridges than proteins from the mesophilic counterparts. Increasing the thermal (not the thermodynamic) stability of proteins by optimization of charge-charge interactions is a good example for an evolutionary solution utilizing physical factors.

Electrostatic similarity of proteins: application of three dimensional spherical harmonic decomposition

We present a method for describing and comparing global electrostatic properties of biomolecules based on the spherical harmonic decomposition of electrostatic potential data. Unlike other approaches our method does not require any prior three dimensional structural alignment. The electrostatic potential, given as a volumetric data set from a numerical solution of the Poisson or Poisson-Boltzmann equation, is represented with descriptors that are rotation invariant. The method can be applied to large and structurally diverse sets of biomolecules enabling to cluster them according to their electrostatic features.

Accurate Molecular Polarizabilities Based on Continuum Electrostatics

A novel approach for representing the intramolecular polarizability as a continuum dielectric is introduced to account for molecular electronic polarization. It is shown, using a finite-difference solution to the Poisson equation, that the Electronic Polarization from Internal Continuum (EPIC) model yields accurate gas-phase molecular polarizability tensors for a test set of 98 challenging molecules composed of heteroaromatics, alkanes and diatomics. The electronic polarization originates from a high intramolecular dielectric that produces polarizabilities consistent with B3LYP/aug-cc-pVTZ and experimental values when surrounded by vacuum dielectric. In contrast to other approaches to model electronic polarization, this simple model avoids the polarizability catastrophe and accurately calculates molecular anisotropy with the use of very few fitted parameters and without resorting to auxiliary sites or anisotropic atomic centers. On average, the unsigned error in the average polarizability and anisotropy compared to B3LYP are 2% and 5%, respectively. The correlation between the polarizability components from B3LYP and this approach lead to a R2 of 0.990 and a slope of 0.999. Even the F2 anisotropy, shown to be a difficult case for existing polarizability models, can be reproduced within 2% error. In addition to providing new parameters for a rapid method directly applicable to the calculation of polarizabilities, this work extends the widely used Poisson equation to areas where accurate molecular polarizabilities matter.

Accurate Calculations of Binding, Folding, and Transfer Free Energies by a Scaled Generalized Born Method

The Poisson-Boltzmann (PB) equation is widely used for modeling solvation effects. The computational cost of PB has restricted its applications largely to single-conformation calculations. The generalized Born (GB) model provides an approximation at substantially reduced cost. Currently the best GB methods reproduce PB results for electrostatic solvation energies with errors at ~5 kcal/mol. When two proteins form a complex, the net electrostatic contributions to the binding free energy are typically of the order of 5 to 10 kcal/mol. Similarly, the net contributions of individual residues to protein folding free energy are < 5 kcal/mol. Clearly in these applications the accuracy of current GB methods is insufficient. Here we present a simple scaling scheme that allows our GB method, GBr⁶, to reproduce PB results for binding, folding, and transfer free energies with high accuracy. From an ensemble of conformations sampled from molecular dynamics simulations, five were judiciously selected for PB calculations. These PB results were used for scaling GBr⁶. Tests on the binding free energies of the barnase-barstar, GTPase-WASp, and U1A-U1hpII complexes and on the folding free energy of FKBP show that the effects of point mutations calculated by scaled GBr⁶ are accurate to within 0.3 kcal/mol of PB results. Similar accuracy was also achieved for the free energies of transfer for ribonuclease Sa and insulin from the crystalline phase to the solution phase at various pH's. This method makes it possible to thoroughly sample the transient-complex ensemble in predicting protein binding rate constants and to incorporate conformational sampling in electrostatic modeling (such as done in the MM-GBSA approach) without loss of accuracy.

Assessing the performance of implicit solvation models at a nucleic acid surface

Implicit solvation models are popular alternatives to explicit solvent methods due to their ability to "pre-average" solvent behavior and thus reduce the need for computationally-expensive sampling. Previously, we have demonstrated that Poisson-Boltzmann models for polar solvation and integral-based models for nonpolar solvation can reproduce explicit solvation forces in a low-charge density protein system. In the present work, we examine the ability of these continuum models to describe solvation forces at the surface of a RNA hairpin. While these models do not completely describe all of the details of solvent behavior at this highly-charged biomolecular interface, they do provide a reasonable description of average solvation forces and therefore show significant promise for developing more robust implicit descriptions of solvent around nucleic acid systems for use in biomolecular simulation and modeling. Additionally, we observe fairly good transferability in the nonpolar model parameters optimized for protein systems, suggesting its robustness for modeling general nonpolar solvation phenomena in biomolecular systems.

Continuum molecular electrostatics, salt effects, and counterion binding--a review of the Poisson-Boltzmann theory and its modifications

This work is a review of the Poisson-Boltzmann (PB) continuum electrostatics theory and its modifications, with a focus on salt effects and counterion binding. The PB model is one of the mesoscopic theories that describes the electrostatic potential and equilibrium distribution of mobile ions around molecules in solution. It serves as a tool to characterize electrostatic properties of molecules, counterion association, electrostatic contributions to solvation, and molecular binding free energies. We focus on general formulations which can be applied to large molecules of arbitrary shape in all-atomic representation, including highly charged biomolecules such as nucleic acids. These molecules present a challenge for theoretical description, because the conventional PB model may become insufficient in those cases. We discuss the conventional PB equation, the corresponding functionals of the electrostatic free energy, including a connection to DFT, simple empirical extensions to this model accounting for finite size of ions, the modified PB theory including ionic correlations and fluctuations, the cell model, and supplementary methods allowing to incorporate site-bound ions in the PB calculations.

A paramagnetic molecular voltmeter

We have developed a general electron paramagnetic resonance (EPR) method to measure electrostatic potential at spin labels on proteins to millivolt accuracy. Electrostatic potential is fundamental to energy-transducing proteins like myosin, because molecular energy storage and retrieval is primarily electrostatic. Quantitative analysis of protein electrostatics demands a site-specific spectroscopic method sensitive to millivolt changes. Previous electrostatic potential studies on macromolecules fell short in sensitivity, accuracy and/or specificity. Our approach uses fast-relaxing charged and neutral paramagnetic relaxation agents (PRAs) to increase nitroxide spin label relaxation rate solely through collisional spin exchange. These PRAs were calibrated in experiments on small nitroxides of known structure and charge to account for differences in their relaxation efficiency. Nitroxide longitudinal (R(1)) and transverse (R(2)) relaxation rates were separated by applying lineshape analysis to progressive saturation spectra. The ratio of measured R(1) increases for each pair of charged and neutral PRAs measures the shift in local PRA concentration due to electrostatic potential. Voltage at the spin label is then calculated using the Boltzmann equation. Measured voltages for two small charged nitroxides agree with Debye-Hückel calculations. Voltage for spin-labeled myosin fragment S1 also agrees with calculation based on the pK shift of the reacted cysteine.

Treatment of charge singularities in implicit solvent models

This paper presents a novel method for solving the Poisson-Boltzmann (PB) equation based on a rigorous treatment of geometric singularities of the dielectric interface and a Green's function formulation of charge singularities. Geometric singularities, such as cusps and self-intersecting surfaces, in the dielectric interfaces are bottleneck in developing highly accurate PB solvers. Based on an advanced mathematical technique, the matched interface and boundary (MIB) method, we have recently developed a PB solver by rigorously enforcing the flux continuity conditions at the solvent-molecule interface where geometric singularities may occur. The resulting PB solver, denoted as MIBPB-II, is able to deliver second order accuracy for the molecular surfaces of proteins. However, when the mesh size approaches half of the van der Waals radius, the MIBPB-II cannot maintain its accuracy because the grid points that carry the interface information overlap with those that carry distributed singular charges. In the present Green's function formalism, the charge singularities are transformed into interface flux jump conditions, which are treated on an equal footing as the geometric singularities in our MIB framework. The resulting method, denoted as MIBPB-III, is able to provide highly accurate electrostatic potentials at a mesh as coarse as 1.2 A for proteins. Consequently, at a given level of accuracy, the MIBPB-III is about three times faster than the APBS, a recent multigrid PB solver. The MIBPB-III has been extensively validated by using analytically solvable problems, molecular surfaces of polyatomic systems, and 24 proteins. It provides reliable benchmark numerical solutions for the PB equation.

A method to determine dielectric constants in nonhomogeneous systems: application to biological membranes

Continuum electrostatic models have had quantitative success in describing electrostatic-mediated phenomena on atomistic scales; however, there continues to be significant disagreement about how to assign dielectric constants in mixed, nonhomogeneous systems. We introduce a method for determining a position-dependent dielectric profile from molecular dynamics simulations. In this method, the free energy of introducing a test charge is computed two ways: from a free energy perturbation calculation and from a numerical solution to Poisson's Equation. The dielectric profile of the system is then determined by minimizing the discrepancy between these two calculations simultaneously for multiple positions of the test charge. We apply this method to determine the dielectric profile of a lipid bilayer surrounded by water. We find good agreement with dielectric models for lipid bilayers obtained by other approaches. The free energy of transferring an ion from bulk water to the lipid bilayer computed from the atomistic simulations indicates that large errors are introduced when the bilayer is represented as a single slab of low dielectric embedded in the higher-dielectric solvent. Significant improvement results from introducing an additional layer of intermediate dielectric ( approximately 3) on each side of the low dielectric core extending from approximately 12 A to 18 A. A small dip in transfer free energy just outside the lipid headgroups indicates the presence of a very high dielectric. These results have implications for the design of implicit membrane models and our understanding of protein-membrane interactions.

Atomically detailed simulations of concentrated protein solutions: the effects of salt, pH, point mutations, and protein concentration in simulations of 1000-molecule systems

An ability to accurately simulate the dynamic behavior of concentrated macromolecular solutions would be of considerable utility in studies of a wide range of biological systems. With this goal in mind, a Brownian dynamics (BD) simulation method is reported here that allows systems to be modeled that comprise in excess of 1000 protein molecules, all of which are treated in atomic detail. Intermolecular forces are described in the method using an energy function that incorporates electrostatic and hydrophobic interactions and that is calibrated to reproduce experimental thermodynamic information with a single adjustable parameter. Using the method, BD simulations have been performed over a wide range of pH and ionic strengths for three proteins: hen egg white lysozyme (HEWL), chymotrypsinogen, and T4 lysozyme. The simulations reproduce experimental trends in second virial coefficients (B(22)) and translational diffusion coefficients, correctly capture changes in B(22) values due to single amino acid substitutions, and reveal a new explanation for the difficulties reported previously in the literature in reproducing B(22) values for protein solutions of very low ionic strength. In addition, a strong correlation is found between a residue's probability of being involved in a protein-protein contact in the simulations and its probability of being involved in an experimental crystal contact. Finally, exploratory simulations of HEWL indicate that the simulation model also gives a promising description of behavior at very high protein concentrations (approximately 250 g/L), suggesting that it may provide a suitable computational framework for modeling the complex behavior exhibited by macromolecules in cellular conditions.

Monte Carlo and modified Tanford-Kirkwood results for macromolecular electrostatics calculations

The understanding of electrostatic interactions is an essential aspect of the complex correlation between structure and function of biological macromolecules. It is also important in protein engineering and design. Theoretical studies of such interactions are predominantly done within the framework of Debye-Hückel theory. A classical example is the Tanford-Kirkwood (TK) model. Besides other limitations, this model assumes an infinitesimally small macromolecule concentration. By comparison to Monte Carlo (MC) simulations, it is shown that TK predictions for the shifts in ion binding constants upon addition of salt become less reliable even at moderately macromolecular concentrations. A simple modification based on colloidal literature is suggested to the TK scheme. The modified TK models suggested here satisfactorily predict MC and experimental shifts in the calcium binding constant as a function of protein concentration for the calbindin D(9k) mutant and calmodulin.

Balancing an accurate representation of the molecular surface in generalized born formalisms with integrator stability in molecular dynamics simulations

Different integrator time steps in NVT and NVE simulations of protein and nucleic acid systems are tested with the GBMV (Generalized Born using Molecular Volume) and GBSW (Generalized Born with simple SWitching) methods. The simulation stability and energy conservation is investigated in relation to the agreement with the Poisson theory. It is found that very close agreement between generalized Born methods and the Poisson theory based on the commonly used sharp molecular surface definition results in energy drift and simulation artifacts in molecular dynamics simulation protocols with standard 2-fs time steps. New parameters are proposed for the GBMV method, which maintains very good agreement with the Poisson theory while providing energy conservation and stable simulations at time steps of 1 to 1.5 fs.

A new outer boundary formulation and energy corrections for the nonlinear Poisson-Boltzmann equation

The nonlinear Poisson-Boltzmann equation (PBE) has been successfully used for the prediction of numerous electrostatic properties of highly charged biopolyelectrolytes immersed in aqueous salt solutions. While numerous numerical solvers for the 3D PBE have been developed, the formulation of the outer boundary treatments used in these methods has only been loosely addressed, especially in the nonlinear case. The de facto standard in current nonlinear PBE implementations is to either set the potential at the outer boundaries to zero or estimate it using the (linear) Debye-Hückel (DH) approximation. However, an assessment of how these outer boundary treatments affect the overall solution accuracy does not appear to have been previously made. As will be demonstrated here, both approximations can, under certain conditions, produce completely erroneous estimates of the potential and energy salt dependencies. A related concern for calculations carried out on grids of finite extent (e.g., all current finite difference and finite element implementations) is the contribution to the energy and salt dependence from the exterior region outside the computational grid. This too is shown to be significant, especially at low salt concentration where essentially all of the contributions to the excess osmotic pressure and ion stress energies originate from this exterior region. In this paper the authors introduce a new outer boundary treatment that is valid for both the linear and nonlinear PBE. The authors also formulate energy corrections to account for contributions from outside the computational domain. Finally, the authors also consider the effects of general ion exclusion layers upon biomolecular electrostatics. It is shown that while these layers tend to increase the surface electrostatic potential, under physiological salt conditions and high net charges their effect on the excess osmotic pressure term, which is a measure of the salt dependence of the total electrostatic free energy, is weak. To facilitate presentation and allow very fine resolutions and/or large computational domains to be considered, attention is restricted to the 1D spherically symmetric nonlinear PBE. Though geometrically limited, the modeling principles nevertheless extend to general PBE solvers as discussed in the Appendix. The 1D model can also be used to benchmark and validate the salt effect prediction capabilities of existing PBE solvers.

Electrostatic properties and macroscopic electrodiffusion in OmpF porin and mutants

The bacterial porin OmpF found in the outer membrane of E. coli is a wide channel, characterized by its poor selectivity and almost no ion specificity. It has an asymmetric structure, with relatively large entrances and a narrow constriction. By applying continuum electrostatic methods we determine the ionization states of titratable amino acid residues in the protein and calculate self-consistently the electric potential 3-D distribution within the channel. The average electrostatic properties are then represented by an effective fixed charge distribution along the pore which is the input for a macroscopic electrodiffusion model. The theoretical predictions agree with measurements performed under different salt gradients and pH. The sensitivity of reversal potential and conductance to the direction of the salt gradient and the solution pH is captured by the model. The theory is also able to explain the influence of the lipid membrane charge. The same methodology is satisfactorily applied to some OmpF mutants involving slight structural changes but a large number of net charges. The correlation found between atomic structure and ionic selectivity shows that the transport characteristics of wide channels like OmpF and its mutants are mainly regulated by the collective action of a large number of residues, rather than by the specific interactions of residues at particular locations.

Hydrophobic surface burial is the major stability determinant of a flat, single-layer beta-sheet

Formation of a flat beta-sheet is a fundamental event in beta-sheet-mediated protein self-assembly. To investigate the contributions of various factors to the stability of flat beta-sheets, we performed extensive alanine-scanning mutagenesis experiments on the single-layer beta-sheet segment of Borrelia outer surface protein A (OspA). This beta-sheet segment consists of beta-strands with highly regular geometries that can serve as a building block for self-assembly. Our Ala-scanning approach is distinct from the conventional host-guest method, in that it introduces only conservative, truncation mutations that should minimize structural perturbation. Our results showed very weak correlation with experimental beta-sheet propensity scales, statistical beta-sheet propensity scales, or cross-strand pairwise correlations. In contrast, our data showed strong positive correlation with the change in buried non-polar surface area. Polar interactions including prominent Glu-Lys cross-strand pairs contribute marginally to the beta-sheet stability. These results were corroborated by results from additional non-Ala mutations. Taken together, these results demonstrate the dominant contribution of non-polar surface burial to flat beta-sheet stability even at solvent-exposed positions. The OspA single-layer beta-sheet achieves efficient hydrophobic surface burial without forming a hydrophobic core by a strategic placement of a variety of side-chains. These findings further suggest the importance of hydrophobic interactions within a beta-sheet layer in peptide self-assembly.

Catabolite Activator Protein in Aqueous Solution: A Molecular Simulation Study

The homodimeric catabolite activator protein (CAP) is a bacterial DNA binding transcription regulator whose activity is controlled by the binding of the intracellular mediator cyclic adenosine monophosphate (cAMP). Each CAP subunit consists of a cyclic nucleotide and a DNA binding domain. Here, we investigate the structural features of the ligand-bound CAP in aqueous solution by molecular dynamics simulations based on the available X-ray structures (Passner et al. J. Mol. Biol. 2000, 304, 847-859 and Chen et al. J. Mol. Biol. 2001, 314, 63-74). Our calculations suggest that the homodimer in solution assumes a symmetric arrangement in which both DNA binding domains are separated from the respective cyclic nucleotide binding domains by a cleft. This contrasts with the X-ray structure, which exhibits instead an asymmetric conformation. On the basis of electrostatics calculations, we propose that the symmetric structure in solution may be an important feature for DNA molecular recognition.

Computation of electrostatic forces between solvated molecules determined by the Poisson-Boltzmann equation using a boundary element method

A rigorous approach is proposed to calculate the electrostatic forces among an arbitrary number of solvated molecules in ionic solution determined by the linearized Poisson-Boltzmann equation. The variational principle is used and implemented in the frame of a boundary element method (BEM). This approach does not require the calculation of the Maxwell stress tensor on the molecular surface, therefore it totally avoids the hypersingularity problem in the direct BEM whenever one needs to calculate the gradient of the surface potential or the stress tensor. This method provides an accurate and efficient way to calculate the full intermolecular electrostatic interaction energy and force, which could potentially be used in Brownian dynamics simulation of biomolecular association. The method has been tested on some simple cases to demonstrate its reliability and efficiency, and parts of the results are compared with analytical results and with those obtained by some known methods such as adaptive Poisson-Boltzmann solver.

Molecular simulation of the binding of nerve growth factor peptide mimics to the receptor tyrosine kinase A

Nerve growth factor (NGF) mimics play an important role for therapies that target the receptor tyrosine kinase A (trkA). The N-terminal fragment of the NGF (N-term@NGF) was previously demonstrated to be an important determinant for affinity and specificity in the binding to trkA. Here we use a variety of computational tools (contact surface analysis and free energy predictions) to identify residues playing a key role for the binding to the receptor. Molecular dynamics simulations are then used to investigate the stability of complexes between trkA and peptides mimicking N-term@NGF. Steered molecular dynamics calculations are finally performed to investigate the process of detaching the peptide from the receptor. Three disruptive events are observed, the first involving the breaking of all intermolecular interactions except two salt bridges, which break subsequently.

Collagen Plays an Active Role in the Aggregation of β2-Microglobulin under Physiopathological Conditions of Dialysis-related Amyloidosis

Dialysis-related amyloidosis is characterized by the deposition of insoluble fibrils of beta(2)-microglobulin (beta(2)-m) in the musculoskeletal system. Atomic force microscopy inspection of ex vivo amyloid material reveals the presence of bundles of fibrils often associated to collagen fibrils. Aggregation experiments were undertaken in vitro with the aim of reproducing the physiopathological fibrillation process. To this purpose, atomic force microscopy, fluorescence techniques, and NMR were employed. We found that in temperature and pH conditions similar to those occurring in periarticular tissues in the presence of flogistic processes, beta(2)-m fibrillogenesis takes place in the presence of fibrillar collagen, whereas no fibrils are obtained without collagen. Moreover, the morphology of beta(2)-m fibrils obtained in vitro in the presence of collagen is extremely similar to that observed in the ex vivo sample. This result indicates that collagen plays a crucial role in beta(2)-m amyloid deposition under physiopathological conditions and suggests an explanation for the strict specificity of dialysis-related amyloidosis for the tissues of the skeletal system. We hypothesize that positively charged regions along the collagen fiber could play a direct role in beta(2)-m fibrillogenesis. This hypothesis is sustained by aggregation experiments performed by replacing collagen with a poly-L-lysine-coated mica surface. As shown by NMR measurements, no similar process occurs when poly-L-lysine is dissolved in solution with beta(2)-m. Overall, the findings are consistent with the estimates resulting from a simplified collagen model whereby electrostatic effects can lead to high local concentrations of oppositely charged species, such as beta(2)-m, that decay on moving away from the fiber surface.

Direct Observation of Salt Effects on Molecular Interactions through Explicit-Solvent Molecular Dynamics Simulations: Differential Effects on Electrostatic and Hydrophobic Interactions and Comparisons to Poisson−Boltzmann Theory

Proteins and other biomolecules function in cellular environments that contain significant concentrations of dissolved salts and even simple salts such as NaCl can significantly affect both the kinetics and thermodynamics of macromolecular interactions. As one approach to directly observing the effects of salt on molecular associations, explicit-solvent molecular dynamics (MD) simulations have been used here to model the association of pairs of the amino acid analogues acetate and methylammonium in aqueous NaCl solutions of concentrations 0, 0.1, 0.3, 0.5, 1, and 2 M. By performing simulations of 500 ns duration for each salt concentration properly converged estimates of the free energy of interaction of the two molecules have been obtained for all intermolecular separation distances and geometries. The resulting free energy surfaces are shown to give significant new insights into the way salt modulates interactions between molecules containing both charged and hydrophobic groups and are shown to provide valuable new benchmarks for testing the description of salt effects provided by the simpler but faster Poisson-Boltzmann method. In addition, the complex many-dimensional free energy surfaces are shown to be decomposable into a number of one-dimensional effective energy functions. This decomposition (a) allows an unambiguous view of the qualitative differences between the salt dependence of electrostatic and hydrophobic interactions, (b) gives a clear rationalization for why salt exerts different effects on protein-protein association and dissociation rates, and (c) produces simplified energy functions that can be readily used in much faster Brownian dynamics simulations.

Assessing the acid–base and conformational properties of histidine residues in human prion protein (125–228) by means of pK a calculations and molecular dynamics simulations

A thorough study of the acid-base behavior of the four histidines and the other titratable residues of the structured domain of human prion protein (125-228) is presented. By using multi-tautomer electrostatic calculations, average titration curves have been built for all titratable residues, using the whole bundles of NMR structures determined at pH 4.5 and 7.0. According to our results, (1) only histidine residues are likely to be involved in the first steps of the pH-driven conformational transition of prion protein; (2) the pK(a)'s of His140 and His177 are approximately 7.0, whereas those of His155 and His187 are < 5.5. 10-ns long molecular dynamics simulations have been performed on five different models, corresponding to the most significant combinations of histidine protonation states. A critical comparison between the available NMR structures and our computational results (1) confirms that His155 and His187 are the residues whose protonation is involved in the conformational rearrangement of huPrP in mildly acidic condition, and (2) shows how their protonation leads to the destructuration of the C-terminal part of HB and to the loss of the last turn of HA that represent the crucial microscopic steps of the rearrangement.

MM/PBSA analysis of molecular dynamics simulations of bovine beta-lactoglobulin: free energy gradients in conformational transitions?

The pH-driven opening and closure of beta-lactoglobulin EF loop, acting as a lid and closing the internal cavity of the protein, has been studied by molecular dynamics (MD) simulations and free energy calculations based on molecular mechanics/Poisson-Boltzmann (PB) solvent-accessible surface area (MM/PBSA) methodology. The forms above and below the transition pH differ presumably only in the protonation state of residue Glu89. MM/PBSA calculations are able to reproduce qualitatively the thermodynamics of the transition. The analysis of MD simulations using a combination of MM/PBSA methodology and the colony energy approach is able to highlight the driving forces implied in the transition. The analysis suggests that global rearrangements take place before the equilibrium local conformation is reached. This conclusion may bear general relevance to conformational transitions in all lipocalins and proteins in general.

A decoy set for the thermostable subdomain from chicken villin headpiece, comparison of different free energy estimators

BACKGROUND

Estimators of free energies are routinely used to judge the quality of protein structural models. As these estimators still present inaccuracies, they are frequently evaluated by discriminating native or native-like conformations from large ensembles of so-called decoy structures.

RESULTS

A decoy set is obtained from snapshots taken from 5 long (100 ns) molecular dynamics (MD) simulations of the thermostable subdomain from chicken villin headpiece. An evaluation of the energy of the decoys is given using: i) a residue based contact potential supplemented by a term for the quality of dihedral angles; ii) a recently introduced combination of four statistical scoring functions for model quality estimation (FRST); iii) molecular mechanics with solvation energy estimated either according to the generalized Born surface area (GBSA) or iv) the Poisson-Boltzmann surface area (PBSA) method.

CONCLUSION

The decoy set presented here has the following features which make it attractive for testing energy scoring functions:1) it covers a broad range of RMSD values (from less than 2.0 A to more than 12 A);2) it has been obtained from molecular dynamics trajectories, starting from different non-native-like conformations which have diverse behaviour, with secondary structure elements correctly or incorrectly formed, and in one case folding to a native-like structure. This allows not only for scoring of static structures, but also for studying, using free energy estimators, the kinetics of folding;3) all structures have been obtained from accurate MD simulations in explicit solvent and after molecular mechanics (MM) energy minimization using an implicit solvent method. The quality of the covalent structure therefore does not suffer from steric or covalent problems. The statistical and physical effective energy functions tested on the set behave differently when native simulation snapshots are included or not in the set and when averaging over the trajectory is performed.

Thermodynamic and electrostatic properties of ternary Oxytricha nova TEBP-DNA complex

Telomeres constitute the nucleoprotein ends of eukaryotic chromosomes which are essential for their proper function. Telomere end binding protein (TEBP) from Oxytricha nova was among the first telomeric proteins, which were well characterized biologically. TEBP consists of two protein subunits (alpha, beta) and forms a ternary complex with single stranded telomeric DNA containing tandem repeats TTTTGGGG. This work presents the characterization of the thermodynamic and electrostatic properties of this complex by computational chemistry methods (continuum Poisson-Boltzmann and solvent accessible surface calculations). Our calculations give a new insight into molecular properties of studied system. Based on the thermodynamic analysis we provide a rationale for the experimental observation that alpha and ssDNA forms a binary complex and the beta subunit joins alpha:ssDNA complex only after the latter is formed. Calculations of distribution of the molecular electrostatic potential for protein subunits alone and for all possible binary complexes revealed the important role of the "guiding funnel" potential generated by alpha:ssDNA complex. This potential may help the beta subunit to dock to the already formed alpha:DNA intermediate in highly steric and electrostatic favorable manner. Our pK(a) calculations of TEBP are able to explain the experimental mobility shifts of the complex in electrophoretic non-denaturating gels.

Application of MM/PBSA colony free energy to loop decoy discrimination: toward correlation between energy and root mean square deviation

Accurate free energy estimation is needed in many predictive tasks. The molecular mechanics/Poisson-Boltzmann solvent accessible surface area (MM/PBSA) approach has proven to be accurate. However, the correlation between the estimated free energy and the distance (e.g., root mean square deviation [RMSD]) from the most stable conformation is hindered by the strong free energy dependence on minor conformational variations. In this paper, a protocol for MM/PBSA free energy estimation is designed and tested on several loop decoy sets. We show that further integration of MM/PBSA free energy estimator with the colony energy approach makes the correlation between the free energy and RMSD from the native structure apparent, for the test sets on which it could be applied. Our results suggest that (1) the MM/PBSA free energy estimator is able to detect native-like structures for most decoy sets, and (2) application of the colony energy approach greatly hampers the MM/energy strong dependence on minor conformational changes.

Comparative molecular dynamics simulations of HIV-1 integrase and the T66I/M154I mutant: Binding modes and drug resistance to a diketo acid inhibitor

HIV-1 IN is an essential enzyme for viral replication and an interesting target for the design of new pharmaceuticals for use in multidrug therapy of AIDS. L-731,988 is one of the most active molecules of the class of beta-diketo acids. Individual and combined mutations of HIV-1 IN at residues T66, S153, and M154 confer important degrees of resistance to one or more inhibitors belonging to this class. In an effort to understand the molecular mechanism of the resistance of T66I/M154I IN to the inhibitor L-731,988 and its specific binding modes, we have carried out docking studies, explicit solvent MD simulations, and binding free energy calculations. The inhibitor was docked against different protein conformations chosen from prior MD trajectories, resulting in 2 major orientations within the active site. MD simulations have been carried out for the T66I/M154I DM IN, DM IN in complex with L-731,988 in 2 different orientations, and 1QS4 IN in complex with L-731,988. The results of these simulations show a similar dynamical behavior between T66I/M154I IN alone and in complex with L-731,988, while significant differences are observed in the mobility of the IN catalytic loop (residues 138-149). Water molecules bridging the inhibitor to residues from the active site have been identified, and residue Gln62 has been found to play an important role in the interactions between the inhibitor and the protein. This work provides information about the binding modes of L-731,988, as well as insight into the mechanism of inhibitor-resistance in HIV-1 integrase.

The salt-dependence of a protein-ligand interaction: ion-protein binding energetics

Using the binding of a nucleotide inhibitor (guanosine-3'-monophosphate) to a ribonuclease (ribonuclease Sa) as a model system, we show that the salt-dependence of the interaction arises due to specific ion binding at the site of nucleotide binding. The presence of specific ion-protein binding is concluded from a combination of differential scanning calorimetry and NMR data. Isothermal titration calorimetry data are then fit to determine the energetic profile (enthalpy, entropy, and heat capacity) for both the ion-protein and nucleotide-protein interactions. The results provide insight into the energetics of charge-charge interactions, and have implications for the interpretation of an observed salt-dependence. Further, the presence of specific ion-binding leads to a system behavior as a function of temperature that is drastically different from that predicted from Poisson-Boltzmann calculations.

Probing electric fields in protein cavities by using the vibrational stark effect of carbon monoxide

To determine the magnitude and direction of the internal electric field in the Xe4 cavity of myoglobin mutant L29W-S108L, we have studied the vibrational Stark effect of carbon monoxide (CO) using infrared spectroscopy at cryogenic temperatures. CO was photodissociated from the heme iron and deposited selectively in Xe4. Its infrared spectrum exhibits Stark splitting into two bands associated with CO in opposite orientations. Two different photoproduct states can be distinguished, C' and C'', with markedly different properties. For C', characteristic temperature-dependent changes of the area, shift, and width were analyzed, based on a dynamic model in which the CO performs fast librations within a double-well model potential. For the barrier between the wells, a height of approximately 1.8 kJ/mol was obtained, in which the CO performs oscillations at an angular frequency of approximately 25 cm(-1). The magnitude of the electric field in the C' conformation was determined as 11.1 MV/cm; it is tilted by an angle of 29 degrees to the symmetry axis of the potential. Above 140 K, a protein relaxation leads to a significantly altered photoproduct, C'', with a smaller Stark splitting and a more confining potential (barrier >4 kJ/mol) governing the CO librations.

The electrostatic nature of C3d-complement receptor 2 association

The association of complement component C3d with B or T cell complement receptor 2 (CR2 or CD21) is a link between innate and adaptive immunity. It has been recognized in experimental studies that the C3d-CR2 association is pH- and ionic strength-dependent. This led us to perform electrostatic calculations to obtain a theoretical understanding of the mechanism of C3d-CR2 association. We used the crystallographic structures of human free C3d, free CR2 (short consensus repeat (SCR)1-2), and the C3d-CR2(SCR1-2) complex, and continuum solvent representation, to obtain a detailed atomic-level picture of the components of the two molecules that contribute to association. Based on the calculation of electrostatic potentials for the free and bound species and apparent pK(a) values for each ionizable residue, we show that C3d-CR2(SCR1-2) recognition is electrostatic in nature and involves not only the association interface, but also the whole molecules. Our results are in qualitative agreement with experimental data that measured the ionic strength and pH dependence of C3d-CR2 association. Also, our results for the native molecules and a number of theoretical mutants of C3d explain experimental mutagenesis studies of amino acid replacements away from the association interface that modulate binding of iC3b with full-length CR2. Finally, we discuss the packing of the two SCR domains. Overall, our data provide global and site-specific explanations of the physical causes that underlie the ionic strength dependence of C3d-CR2 association in a unified model that accounts for all experimental data, some of which were previously thought to be contradictory.

Electrostatics of Ion Stabilization in a ClC Chloride Channel Homologue from Escherichia coli

The structural determinants of electrostatics of ion stabilization within EcClC, a ClC-type chloride channel homologue from Escherichia coli, are studied using a continuum dielectric approximation. Specifically, the ion occupancy is investigated in the wild-type protein and a mutant thereof, and the contribution to the electrostatic binding free energy of local and non-local interactions is characterized at the single-residue level. This analysis shows that, in spite of the desolvation cost and the strong ion-ion repulsion, all previously reported binding sites can be occupied simultaneously. The stabilizing effect of the protein arises from hydrogen bonding as well as from longer-range favorable interactions, such as with the strictly conserved Lys131 side-chain. The latter is involved in the stabilization of the conserved GSGIP motif that delimits two of the binding sites. Interestingly, an additional low-affinity binding site, mediated by a structurally analogous motif including the side-chain of Arg340, can be identified on the extracellular side of the permeation pathway. Finally, it is shown that, in contrast to K-channels, and in analogy to the SBP/PBP sulfate/phosphate-binding proteins, the contribution of helix macrodipoles to chloride binding in EcClC is only marginal.

Constructing irregular surfaces to enclose macromolecular complexes for mesoscale modeling using the discrete surface charge optimization (DISCO) algorithm

Salt-mediated electrostatics interactions play an essential role in biomolecular structures and dynamics. Because macromolecular systems modeled at atomic resolution contain thousands of solute atoms, the electrostatic computations constitute an expensive part of the force and energy calculations. Implicit solvent models are one way to simplify the model and associated calculations, but they are generally used in combination with standard atomic models for the solute. To approximate electrostatics interactions in models on the polymer level (e.g., supercoiled DNA) that are simulated over long times (e.g., milliseconds) using Brownian dynamics, Beard and Schlick have developed the DiSCO (Discrete Surface Charge Optimization) algorithm. DiSCO represents a macromolecular complex by a few hundred discrete charges on a surface enclosing the system modeled by the Debye-Hückel (screened Coulombic) approximation to the Poisson-Boltzmann equation, and treats the salt solution as continuum solvation. DiSCO can represent the nucleosome core particle (>12,000 atoms), for example, by 353 discrete surface charges distributed on the surfaces of a large disk for the nucleosome core particle and a slender cylinder for the histone tail; the charges are optimized with respect to the Poisson-Boltzmann solution for the electric field, yielding a approximately 5.5% residual. Because regular surfaces enclosing macromolecules are not sufficiently general and may be suboptimal for certain systems, we develop a general method to construct irregular models tailored to the geometry of macromolecules. We also compare charge optimization based on both the electric field and electrostatic potential refinement. Results indicate that irregular surfaces can lead to a more accurate approximation (lower residuals), and the refinement in terms of the electric field is more robust. We also show that surface smoothing for irregular models is important, that the charge optimization (by the TNPACK minimizer) is efficient and does not depend on the initial assigned values, and that the residual is acceptable when the distance to the model surface is close to, or larger than, the Debye length. We illustrate applications of DiSCO's model-building procedure to chromatin folding and supercoiled DNA bound to Hin and Fis proteins. DiSCO is generally applicable to other interesting macromolecular systems for which mesoscale models are appropriate, to yield a resolution between the all-atom representative and the polymer level.

Performance comparison of generalized born and Poisson methods in the calculation of electrostatic solvation energies for protein structures

This study compares generalized Born (GB) and Poisson (PB) methods for calculating electrostatic solvation energies of proteins. A large set of GB and PB implementations from our own laboratories as well as others is applied to a series of protein structure test sets for evaluating the performance of these methods. The test sets cover a significant range of native protein structures of varying size, fold topology, and amino acid composition as well as nonnative extended and misfolded structures that may be found during structure prediction and folding/unfolding studies. We find that the methods tested here span a wide range from highly accurate and computationally demanding PB-based methods to somewhat less accurate but more affordable GB-based approaches and a few fast, approximate PB solvers. Compared with PB solvation energies, the latest, most accurate GB implementations were found to achieve errors of 1% for relative solvation energies between different proteins and 0.4% between different conformations of the same protein. This compares to accurate PB solvers that produce results with deviations of less than 0.25% between each other for both native and nonnative structures. The performance of the best GB methods is discussed in more detail for the application for force field-based minimizations or molecular dynamics simulations.

Observation of a Short, Strong Hydrogen Bond in the Active Site of Hydroxynitrile Lyase from Hevea brasiliensis Explains a Large pK Shift of the Catalytic Base Induced by the Reaction Intermediate

The hydroxynitrile lyase from Hevea brasiliensis (HbHNL) uses a catalytic triad consisting of Ser(80)-His(235)-Asp(207) to enhance the basicity of Ser(80)-O gamma for abstracting a proton from the OH group of the substrate cyanohydrin. Following the observation of a relatively short distance between a carboxyl oxygen of Asp(207) and the N delta(1)(His(235)) in a 1.1 A crystal structure of HbHNL, we here show by (1)H and (15)N-NMR spectroscopy that a short, strong hydrogen bond (SSHB) is formed between the two residues upon binding of the competitive inhibitor thiocyanate to HbHNL: the proton resonance of H-N delta 1(His(235)) moves from 15.41 ppm in the free enzyme to 19.35 ppm in the complex, the largest downfield shift observed so far upon inhibitor binding. Simultaneously, the D/H fractionation factor decreases from 0.98 to 0.35. In the observable pH range, i.e. between pH 4 and 10, no significant changes in chemical shifts (and therefore hydrogen bond strength) were observed for free HbHNL. For the complex with thiocyanate, the 19.35 ppm signal returned to 15.41 ppm at approximately pH 8, which indicates a pK(a) near this value for the H-N epsilon(2)(His(235)). These NMR results were analyzed on the basis of finite difference Poisson-Boltzmann calculations, which yielded the relative free energies of four protonation states of the His(235)-Asp(207) pair in solution as well as in the protein environment with and without bound inhibitor. The calculations explain all the NMR features, i.e. they suggest why a short, strong hydrogen bond is formed upon inhibitor binding and why this short, strong hydrogen bond reverts back to a normal one at approximately pH 8. Importantly, the computations also yield a shift of the free energy of the anionic state relative to the zwitterionic reference state by about 10.6 kcal/mol, equivalent to a shift in the apparent pK(a) of His(235) from 2.5 to 10. This huge inhibitor-induced increase in basicity is a prerequisite for His(235) to act as general base in the HbHNL-catalyzed cyanohydrin reaction.