"Polarizable force field has been successfully used in molecular modeling for years, especially in biological and protein simulations. In this research thesis, development of a new polarizable force field ―POSSIM (POlarizable Simulations with Second order Interaction Model) involving electrostatic polarization is described and parameters for several protein residues were produced. In this research thesis, the POSSIM force field was extended to the side chains of the following residues: lysine, glutamic acid, prontonated hisidine, phenylalanine and tryptophan. This work involved producing parameters for methyl ammonium, acetate ion, imidazolium cation, benzene and pyrrole molecules. The parameters fitting procedure starts from the molecular complex with dipolar electrostatic probes of a many-body system to produce polarizabilities, compute the energies, then charges and Lennard-Jones parameters are produced by fitting to gas-phase dimerization calculations, followed by the torsional parameters fitting and end up with the pure liquid simulations. In all the cases, three-body energies, dimerization energies and distances agree well to the accurate quantum mechanical results. The final parameters obtained assured the error of less than 2% in the heat of vaporization and average volume results compared with the available experimental data. Unlike the quantum mechanical calculations, the polarizable force field computations require a relatively small amount of computational resources. Moreover, compared to fixed-charges empirical force fields, polarizable force fields are much more accurate in a number of energy calculations. In the following chapters, the results obtained with this particular polarizable force field are discussed."

Accurate prediction of protein-ligand binding affinity is essential to computational drug discovery. Current approaches are limited by the accuracy of the underlying potential energy model that describes atomic interactions. A more rigorous physical model is critical for evaluating molecular interactions to chemical accuracy. The objective of this thesis research is to develop a polarizable force field with an accurate representation of electrostatic interactions, and apply this model to protein-ligand recognition and to ultimately solve practical problems in computer aided drug discovery. By calculating the hydration free energies of a series of organic small molecules, an optimal protocol is established to develop the electrostatic parameters from quantum mechanics calculations. Next, the systematical development and parameterization procedure of AMOEBA protein force field is presented. The derived force field has gone through extensive validations in both gas phase and condensed phase. The last part of the thesis involves the application of AMOEBA to study protein-ligand interactions. The binding free energies of benzamidine analogs to trypsin using molecular dynamics alchemical perturbation are calculated with encouraging accuracy. AMOEBA is also used to study the thermodynamic effect of constraining and hydrophobicity on binding energetics between phosphotyrosine(pY)-containing tripeptides and the SH2 domain of growth receptor binding protein 2 (Grb2). The underlying mechanism of an "entropic paradox" associated with ligand preorganization is explored. / text

"The accurate, fast and low cost computational tools are indispensable for studying the structure and dynamics of biological macromolecules in aqueous solution. The goal of this thesis is development and validation of continuum Fuzzy-Border (FB) solvation model to work with the Polarizable Simulations Second-order Interaction Model (POSSIM) force field for proteins developed by Professor G A Kaminski. The implicit FB model has advantages over the popularly used Poisson Boltzmann (PB) solvation model. The FB continuum model attenuates the noise and convergence issues commonly present in numerical treatments of the PB model by employing fixed position cubic grid to compute interactions. It also uses either second or first-order approximation for the solvent polarization which is similar to the second-order explicit polarization applied in POSSIM force field. The FB model was first developed and parameterized with nonpolarizable OPLS-AA force field for small molecules which are not only important in themselves but also building blocks of proteins and peptide side chains. The hydration parameters are fitted to reproduce the experimental or quantum mechanical hydration energies of the molecules with the overall average unsigned error of ca. 0.076kcal/mol. It was further validated by computing the absolute pKa values of 11 substituted phenols with the average unsigned error of 0.41pH units in comparison with the quantum mechanical error of 0.38pH units for this set of molecules. There was a good transferability of hydration parameters and the results were produced only with fitting of the specific atoms to the hydration energy and pKa targets. This clearly demonstrates the numerical and physical basis of the model is good enough and with proper fitting can reproduce the acidity constants for other systems as well. After the successful development of FB model with the fixed charges OPLS-AA force field, it was expanded to permit simulations with Polarizable Simulations Second-order Interaction Model (POSSIM) force field. The hydration parameters of the small molecules representing analogues of protein side chains were fitted to their solvation energies at 298.15K with an average error of ca.0.136kcal/mol. Second, the resulting parameters were used to reproduce the pKa values of the reference systems and the carboxylic (Asp7, Glu10, Glu19, Asp27 and Glu43) and basic residues (Lys13, Lys29, Lys34, His52 and Lys55) of the turkey ovomucoid third domain (OMTKY3) protein. The overall average unsigned error in the pKa values of the acid residues was found to be 0.37pH units and the basic residues was 0.38 pH units compared to 0.58pH units and 0.72 pH units calculated previously using polarizable force field (PFF) and Poisson Boltzmann formalism (PBF) continuum solvation model. These results are produced with fitting of specific atoms of the reference systems and carboxylic and basic residues of the OMTKY3 protein. Since FB model has produced improved pKa shifts of carboxylic residues and basic protein residues in OMTKY3 protein compared to PBF/PFF, it suggests the methodology of first-order FB continuum solvation model works well in such calculations. In this study the importance of explicit treatment of the electrostatic polarization in calculating pKa of both acid and basic protein residues is also emphasized. Moreover, the presented results demonstrate not only the consistently good degree of accuracy of protein pKa calculations with the second-degree POSSIM approximation of the polarizable calculations and the first-order approximation used in the Fuzzy-Border model for the continuum solvation energy, but also a high degree of transferability of both the POSSIM and continuum solvent Fuzzy Border parameters. Therefore, the FB model of solvation combined with the POSSIM force field can be successfully applied to study the protein and protein-ligand systems in water. "

Genetic sequences are being collected at an ever increasing rate due to rapid cost reductions; however, experimental approaches to determine the structure and function of the protein(s) each gene codes are not keeping pace. Therefore, computational methods to augment experimental structures with comparative (i.e. homology) models using physics-based methods for building residues, loops and domains are needed to thread new sequences onto homologous structures. In addition, even experimental structure determination relies on analogous first principles structure refinement and prediction algorithms to place structural elements that are not defined by the data alone. Computational methods developed to find the global free energy minimum of an amino acid sequence (i.e. the protein folding problem) are increasingly successful, but limitations in accuracy and efficiency remain. Optimization efforts have focused on subsets of systems and environments by utilizing potential energy functions ranging from fixed charged force fields (Fiser, Do, & Sali, 2000; Jacobson et al., 2004), statistical or knowledge based potentials (Das & Baker, 2008) and/or potentials incorporating experimental data (Brunger, 2007; Trabuco, Villa, Mitra, Frank, & Schulten, 2008). Although these methods are widely used, limitations include 1) a target function global minimum that does not correspond to the actual free energy minimum and/or 2) search protocols that are inefficient or not deterministic due to rough energy landscapes characterized by large energy barriers between multiple minima. Our Global Optimization Using Metadynamics and a Polarizable Force Field (GONDOLA) approach tackles the first limitation by incorporating experimental data (i.e. from X-ray crystallography, CryoEM or NMR experiments) into a hybrid target function that also includes information from a polarizable molecular mechanics force field (Lopes, Roux, & MacKerell, 2009; Ponder & Case, 2003). The second limitation is overcome by driving the sampling of conformational space by adding a time-dependent bias to the objective function, which pushes the search toward unexplored regions (Alessandro Barducci, Bonomi, & Parrinello, 2011; Zheng, Chen, & Yang, 2008). The GONDOLA approach incorporates additional efficiency constructs for search space exploration that include Monte Carlo moves and fine grained minimization. Furthermore, the dimensionality of the search is reduced by fixing atomic coordinates of known structural regions while atoms of interest explore new coordinate positions. The overall approach can be used for optimization of side-chains (i.e. set side-chain atoms active while constraining backbone atoms), residues (i.e. side-chain atoms and backbone atoms active), ligand binding pose (i.e. set atoms along binding interface active), protein loops (i.e. set atoms connecting two terminating residues active) or even entire protein domains or complexes. Here we focus on using the GONDOLA general free energy driven optimization strategy to elucidate the structural details of missing protein loops, which are often missing from experimental structures due to conformational heterogeneity and/or limitations in the resolution of the data. We first show that the correlation between experimental data and AMOEBA (i.e. a polarizable force field) structural minima is stronger than that for OPLS-AA (i.e. a fixed charge force field). This suggests that the higher order multipoles and polarization of the AMOEBA force field more accurately represented the true crystalline environment than the simpler OPLS-AA model. Thus, scoring and optimization of loops with AMOEBA is more accurate than with OPLS-AA, albeit at a slightly increased computational cost. Next, missing PDZ domain protein loops and protein loops from a loop decoy data set were optimized for 5 ns using the GONDOLA approach (i.e. under the AMOEBA polarizable force field) as well as a commonly used global optimization procedure (i.e. simulated annealing under the OPLS-AA fixed charge force field). The GONDOLA procedure was shown to provide more accurate structures in terms of both experimental metrics (i.e. lower Rfree values) and structural metrics (i.e. using the MolProbity structure validation tool). In terms of Rfree, only one out of seven simulated annealing results was better than the Gondola global optimization. Similarly, one simulated anneal loop had a better MolProbity score, but none of the simulated annealing loops were better in both categories. On average, GONDOLA achieved an Rfree value 19.48 and simulated annealing saw an average Rfree value of 19.63, and the average MolProbity scores were 1.56 for GONDOLA and 1.75 for simulated annealing. In addition to providing more accurate predictions, GONDOLA was shown to converge much faster than the simulated annealing protocol. Ten separate 5 ns optimizations of the 4 residue loop missing from one of the PDZ domains were conducted. Five were done using GONDOLA and five with the simulated annealing protocol. The fastest four converging results belonged to the GONDOLA approach. Thus, this work demonstrates that GONDOLA is well-suited to refine or predict the coordinates of missing residues and loops because it is both more accurate and converges more rapidly.

Specific recognition of ligands including metal ions by proteins is the key of many crucial biological functions and systems. Accurate prediction of the binding strength not only sheds light on the mechanism of the molecular recognition but also provides the most important prerequisite of drug discovery. Computational modeling of molecular binding has gained a great deal of attentions in the last few decades since the advancement of computer power and availability of high-resolution crystal structures. However there still exist two major challenges in the field of molecular modeling, i.e. sampling issue and accuracy of the models. In this work, I have dedicated to tackling these two problems with a noval polarizable force field which is believed to produce more accurate description of molecular interactions than classic non-polarizable force fields. We first developed the model for divalent cations Mg²⁺ and Ca²⁺, deriving the parameters from quantum mechanics. To understand the hydration thermodynamics of these ions we have performed molecular dynamics simulations using our AMOEBA force field. Both the water structures around ions and the solvation free energies were in great accordance with experiment data. We have also simulated and calculated the binding free energies of a series of benzamidine-like inhibitors to trypsin using explicit solvent approach by free energy perturbation. The calculated binding free energies are well within the accuracy of experimental measurement and the direction of change is predicted correctly in all cases. Finally, we computed the hydration free energies of a few organic molecules and automated the calculation procedure. / text

In this work, the development of parameters for a series of imidazolium-based ionic liquids molecules, now included in the AMOEBA force field, is discussed. The quality of obtained parameters is tested in a variety of calculations to reproduce structural, thermodynamic, and transport properties. First, it is proposed a novel method to parameterize in a faster, and more efficient way parameters for the AMOEBA force field that can be applied to any imidazolim-based cation. Second, AMOEBA-IL polarizable force field is applied to study the N-tert-butyloxycarbonylation of aniline reaction mechanism in water/[EMIM][BF4] solvent via QM/MM approach and compared with the reaction carried out in gas-phase and implicit solvent media. Third, AMOEBA-IL force field is applied in alchemical calculations. Free energies of solvation for selected solutes solvated in [EMIm][OTf] are calculated via BAR method implemented in TINKER considering the effect of polarization as well as the methodology to perform the sampling of the alchemical process. Finally, QM/MM calculations using AMOEBA to get more insights into the catalytic reaction mechanism of horseradish peroxidase enzyme, particularly the structures involved in the transition from Cp I to Cp II.

Yes / Materials-binding peptides offer promising routes to the production of tailored Pd nanomaterials in aqueous media, enabling the optimization of catalytic properties. However, the atomic-scale details needed to make these advances are relatively scarce and challenging to obtain. Molecular simulations can provide key insights into the structure of peptides adsorbed at the aqueous Pd interface, provided that the force-field can appropriately capture the relevant bio-interface interactions. Here, we introduce and apply a new polarizable force field, PdP-CHARMM, for the simulation of biomolecule–Pd binding under aqueous conditions. PdP-CHARMM was parametrized with density functional theory (DFT) calculations, using a process compatible with similar polarizable force-fields created for Ag and Au surfaces, ultimately enabling a direct comparison of peptide binding modes across these metal substrates. As part of our process for developing PdP-CHARMM, we provide an extensive study of the performance of ten different dispersion-inclusive DFT functionals in recovering biomolecule–Pd(111) binding. We use the functional with best all-round performance to create PdP-CHARMM.We then employ PdP-CHARMM and metadynamics simulations to estimate the adsorption free energy for a range of amino acids at the aqueous Pd(111) interface. Our findings suggest that only His and Met favor direct contact with the Pd substrate, which we attribute to a remarkably robust interfacial solvation layering. Replica-exchange with solute tempering molecular dynamics simulations of two experimentally-identified Pd-binding peptides also indicate surface contact to be chiefly mediated by His and Met residues at aqueous Pd(111). Adsorption of these two peptides was also predicted for the Au(111) interface, revealing distinct differences in both the solvation structure and modes of peptide adsorption at the Au and Pd interfaces. We propose that this sharp contrast in peptide binding is largely due to the differences in interfacial solvent structuring. / Air Force Office for Scientfi c Research (Grant #FA9550-12-1-0226)

With recent advances and cost reductions in next generation sequencing (NGS), the amount of genetic sequence data is increasing rapidly. However, before patient specific genetic information reaches its full potential to advance clinical diagnostics, the immense degree of genetic heterogeneity that contributes to human disease must be more fully understood. For example, although large numbers of genetic variations are discovered during clinical use of NGS, annotating and understanding the impact of such coding variations on protein phenotype remains a bottleneck (i.e. what is the molecular mechanism behind deafness phenotypes). Fortunately, computational methods are emerging that can be used to efficiently study protein coding variants, and thereby overcome the bottleneck brought on by rapid adoption of clinical sequencing. To study proteins via physics-based computational algorithms, high-quality 3D structural models are essential. These protein models can be obtained using a variety of numerical optimization methods that operate on physics-based potential energy functions. Accurate protein structures serve as input to downstream variation analysis algorithms. In this work, we applied a novel amino acid side-chain optimization algorithm, which operated on an advanced model of atomic interactions (i.e. the AMOEBA polarizable force field), to a set of 164 protein structural models implicated in deafness. The resulting models were evaluated with the MolProbity structure validation tool. MolProbity “scores” were originally calibrated to predict the quality of X-ray diffraction data used to generate a given protein model (i.e. a 1.0 Å or lower MolProbity score indicates a protein model from high quality data, while a score of 4.0 Å or higher reflects relatively poor data). In this work, the side-chain optimization algorithm improved mean MolProbity score from 2.65 Å (42nd percentile) to nearly atomic resolution at 1.41 Å (95th percentile). However, side-chain optimization with the AMOEBA many-body potential function is computationally expensive. Thus, a second contribution of this work is a parallelization scheme that utilizes nVidia graphical processing units (GPUs) to accelerate the side-chain repacking algorithm. With the use of one GPU, our side-chain optimization algorithm achieved a 25 times speed-up compared to using two Intel Xeon E5-2680v4 central processing units (CPUs). We expect the GPU acceleration scheme to lessen demand on computing resources dedicated to protein structure optimization efforts and thereby dramatically expand the number of protein structures available to aid in interpretation of missense variations associated with deafness.

Biological molecules are essential parts of organisms and participate in a variety of biological processes within cells. Understanding the relationship between sequence, structure, and function of biological molecules are of fundamental importance in life science and the health care industry. In this dissertation, a multi-scale approach was utilized to develop coarse-grained molecular models for protein and RNA simulations. By simplifying the atomistic representation of a biomolecular system, the coarse-grained approach enables the molecular dynamics simulations to reveal the biological processes, which occur on the time and length scales that are inaccessible to the all-atom models. For RNA, an "intermediate" coarse-grained model was proposed to provide both accuracy and efficiency for RNA 3D structure modeling and prediction. The overall potential parameters were derived based on structural statistics sampled from experimental structures. For protein, a general, transferable coarse-grain framework based on the Gay-Berne potential and electrostatic point multipole expansion was developed for polypeptide simulations. Next, an advanced atomistic model was developed to model electrostatic interaction with high resolution and incorporates electronic polarization effect that is ignored in conventional atomistic models. The last part of my thesis work involves applying all-atom molecular simulations to address important questions and problems in biophysics and structural biology. For example, the interaction between protein and miRNA, the recognition mechanism of antigen and antibody, and the structure dynamics of protein in mixed denaturants. / text

Amyloidogenic diseases, such as Alzheimer's disease (AD) and Type II Diabetes (T2D), are characterized by the accumulation of amyloid aggregates. Despite having very different amino-acid sequences, the underlying amyloidogenic proteins form similar supramolecular fibril structures that are highly stable and resistant to physical and chemical denaturation. AD is characterized by two toxic lesions: extracellular amyloid β-peptide (Aβ) plaques and intracellular neurofibrillary tangles composed of microtubule-associated protein tau. Similarly, a feature of T2D is the deposition of islet amyloid polypeptide (IAPP) aggregates in and around the pancreas. The mechanisms by which Aβ, tau, and IAPP aggregate, and cause cell death is unknown; thus, gaining greater insight into the stabilizing forces and initial unfolding events is crucial to our understanding of these amyloidogenic diseases. This work uses molecular dynamics (MD) simulations to study the secondary, tertiary, and quaternary structure of Aβ, tau, and IAPP. Specifically, this work used the Drude polarizable force field (FF), which explicitly represents electronic polarization allowing charge distributions to change in response to perturbations in local electric fields. This model allows us to describe the role charge plays on protein folding and stability and how perturbations to the charge state drive pathology. Studies were conducted to address the following questions: 1) What are the stabilizing forces of fibril and oligomeric structures? 2) How do charge-altering mutations modulate the conformational ensemble and thermodynamic properties of Aβ? 3) How do charge-altering post-translational modifications of Aβ and tau modulate changes in the conformational ensembles? These studies establish that shifts in local microenvironments play a role in fibril and oligomer stability. Furthermore, these studies found that changes in protein sequence and charge are sufficient to disrupt and change the secondary and tertiary structure of these amyloidogenic proteins. Overall, this dissertation describes how charge modulates protein unfolding and characterizes the mechanism of those changes. In the long term, this work will help in the development of therapeutics that can target these changes to prevent protein aggregation that leads to cell death. / Doctor of Philosophy / Protein aggregation is the hallmark of many chronic diseases, such as Alzheimer's disease (AD) and Type II Diabetes (T2D). The formation of two toxic aggregates: amyloid β-peptide (Aβ) plaques and neurofibrillary tangles composed of microtubule-associated protein tau are some of the key characteristics of AD. In addition, the formation of islet amyloid polypeptide (IAPP) aggregates in the pancreas is thought to play a role in the development of T2D. The pathways by which the proteins Aβ, tau, and IAPP aggregate are unknown; thus, gaining a greater insight into the properties that may cause these diseases is necessary to develop treatments. By studying these proteins at the atomistic level, we can understand how small changes to these proteins alter how they misfold in a way that promotes toxicity. Herein, we used a computational technique called molecular dynamics (MD) simulations to gain new insights into how protein structure changes. We explored the dynamics of these proteins and investigated the role that charge plays in protein folding and described how charge modulates protein folding and characterized the mechanism of those changes. This work serves as a characterization of protein folding and sets the ground for future structural studies and drug development.