Journal of Theoretical and Computational Chemistry

Fusion of HIV-1 viral and host cellular membranes is an important step for HIV infection. The HIV-1 envelope glycoprotein mediating the membrane fusion consists of subunits gp120 and gp41 whereas gp120 recognizes the cell-surface receptors and gp41 promotes viral-cell membrane fusion. The trimeric helical complex composed of heterodimer of N-terminal and C-terminal extraviral segments has been used for the gp41 function study, and the trimeric N-terminal teptad repeat (NHR) is considered as an antiviral drug target for developing HIV-1 membrane fusion inhibitors. By using computational solvent probe mapping, we have explored druggable sites on the trimeric NHR peptides, and identified residues K574 and R579 as the hot spots for inhibitor designing. We further demonstrated that although NB-2 and NB-64 are all N-substituted Pyrrole derivatives and have very similar chemical structures, it is possible that diverse inhibitory mechanisms targeting different negative electrostatic residues (K574 and R579) exist. Results from fragment-based virtual screening identified series of potential lead compounds which could be used for further design of fusion inhibitors.

The deamination process of isoxanthopterin catalyzed by isoxanthopterin deaminase was determined using the combined QM(PM3)/MM molecular dynamics simulations. In this paper, the updated PM3 parameters were employed for zinc ions and the initial model was built up based on the crystal structure. Proton transfer and following steps have been investigated in two paths: Asp336 and His285 serve as the proton shuttle, respectively. Our simulations showed that His285 is more effective than Aap336 in proton transfer for deamination of isoxanthopterin. As hydrogen bonds between the substrate and surrounding residues play a key role in nucleophilic attack, we suggested mutating Thr195 to glutamic acid, which could enhance the hydrogen bonds and help isoxanthopterin get close to the active site. The simulations which change the substrate to pterin 6-carboxylate also performed for comparison. Our results provide reference for understanding of the mechanism of deaminase and for enhancing the deamination rate of isoxanthopterin deaminase.

Inhibition of Aurora A kinase interaction is considered to be a promising approach for the discovery of new molecularly targeted cancer therapeutics. In this study, the binding mechanisms of two different inhibitors with a contrasting binding affinity to Aurora A were investigated by long time scale GPU molecular dynamics (MD) simulations coupled with molecular mechanics-Poisson–Boltzmann/generalized Born surface area (MM-PB/GBSA) method. The results showed that the predicted binding free energies of these two complexes were consistent with the experimental data. Through analyzing the individual energy components of binding free energy, we found that the van der Waals contribution was the main force to drive the inhibitor–protein binding and the electrostatic contribution was also a crucial factor for the inhibitor–Aurora A binding. The structural analysis demonstrated that the inhibitor HPM could produce more hydrophobic interaction contacts with Aurora A than that of 2JZ, and the loss of key hydrogen bonds between the inhibitor and residue Arg137 in the hinge region of Aurora A was another important reason for the weaker binding affinity of 2JZ to Aurora A. This study sheds more light on the development of the efficient inhibitors targeting the Aurora A.

Prions are associated with neurodegenerative diseases induced by transmissible spongiform encephalopathies. The infectious scrapie form is referred to as PrP^Sc, which has conformational change from normal prion with predominant α-helical conformation to the abnormal PrP^Sc that is rich in β-sheet content. Neurodegenerative diseases have been found from both human and bovine sources, but there are no reports about infected by transmissible spongiform encephalopathies from rabbit, canine and horse sources. Here we used coarse-grained Gō model to compare the difference among human, bovine, rabbit, canine, and horse normal (cellular) prion proteins. The denatured state of normal prion has relation with the conversion from normal to abnormal prion protein, so we used all-atom Gō model to investigate the folding pathway and energy landscape for human prion protein. Through using coarse-grained Gō model, the cooperativity of the five prion proteins was characterized in terms of calorimetric criterion, sigmoidal transition, and free-energy profile. The rabbit and horse prion proteins have higher folding free-energy barrier and cooperativity, and canine prion protein has slightly higher folding free-energy barrier comparing with human and bovine prion proteins. The results from all-atom Gō model confirmed the validity of C_α-Gō model. The correlations of our results with previous experimental and theoretical researches were discussed.

RNA polymerase is the enzyme that synthesizes RNA during the transcription process. To understand its mechanism, structural studies have provided us pictures of the series of steps necessary to add a new nucleotide to the nascent RNA chain, the steps altogether known as the nucleotide addition cycle (NAC). However, these static snapshots do not provide dynamic information of these processes involved in NAC, such as the conformational changes of the protein and the atomistic details of the catalysis. Computational studies have made efforts to fill these knowledge gaps. In this review, we provide examples of different computational approaches that have improved our understanding of the transcription elongation process for RNA polymerase, such as normal mode analysis, molecular dynamic (MD) simulations, Markov state models (MSMs). We also point out some unsolved questions that could be addressed using computational tools in the future.

This work extends our earlier two-domain formulation of a differential geometry based multiscale paradigm into a multidomain theory, which endows us the ability to simultaneously accommodate multiphysical descriptions of aqueous chemical, physical and biological systems, such as fuel cells, solar cells, nanofluidics, ion channels, viruses, RNA polymerases, molecular motors, and large macromolecular complexes. The essential idea is to make use of the differential geometry theory of surfaces as a natural means to geometrically separate the macroscopic domain of solvent from the microscopic domain of solute, and dynamically couple continuum and discrete descriptions. Our main strategy is to construct energy functionals to put on an equal footing of multiphysics, including polar (i.e. electrostatic) solvation, non-polar solvation, chemical potential, quantum mechanics, fluid mechanics, molecular mechanics, coarse grained dynamics, and elastic dynamics. The variational principle is applied to the energy functionals to derive desirable governing equations, such as multidomain Laplace–Beltrami (LB) equations for macromolecular morphologies, multidomain Poisson–Boltzmann (PB) equation or Poisson equation for electrostatic potential, generalized Nernst–Planck (NP) equations for the dynamics of charged solvent species, generalized Navier–Stokes (NS) equation for fluid dynamics, generalized Newton's equations for molecular dynamics (MD) or coarse-grained dynamics and equation of motion for elastic dynamics. Unlike the classical PB equation, our PB equation is an integral-differential equation due to solvent–solute interactions. To illustrate the proposed formalism, we have explicitly constructed three models, a multidomain solvation model, a multidomain charge transport model and a multidomain chemo-electro-fluid-MD-elastic model. Each solute domain is equipped with distinct surface tension, pressure, dielectric function, and charge density distribution. In addition to long-range Coulombic interactions, various non-electrostatic solvent–solute interactions are considered in the present modeling. We demonstrate the consistency between the non-equilibrium charge transport model and the equilibrium solvation model by showing the systematical reduction of the former to the latter at equilibrium. This paper also offers a brief review of the field.

We theoretically study the second-order nonlinear optical properties of six fluorescent proteins (FPs), such as green fluorescent protein (GFP), BFP, enhanced BFP (eBFP), CFP, YFP, and DsRed. To begin with, the geometries of all these FP chromophores are optimized at B3LYP/6-311++G** level in a water medium and the polarized continuum model (PCM in water) method is adopted. Using a time-dependent density functional theory (TDDFT) method, the electronic structures and excited-state properties of chromophores are determined. Here we employ TDDFT combining with the sum-over-states (SOS) method to calculate the first-order hyperpolarizability for second-harmonic generation (SHG) optical process. Moreover, we discuss the origin of the nonlinear optical response and determine what caused the variation of first-order hyperpolarizability. Our calculations show that the charge transfers of π → π* in the central conjugated structure and p → π* charge transfers from the side chain R1 to conjugated structure of chromophores markedly affect the first-order hyperpolarizability.

A multi-domain protein is able to exist as equilibrium of different conformations in solution, which may be critical to its biological function. Besides experimental techniques, computational methods like molecular dynamics (MD) simulations are suitable to study inter-domain motions of the protein and sample different conformational states. A MD simulation usually generates a trajectory containing large amount of protein structures, and a post-processing cluster analysis would be necessary to group similar structures into clusters and identify these typical conformations of the multi-domain protein. In this paper, the widely used k-means clustering algorithm is implemented in the protein essential dynamics (ED) subspace defined by principal component analysis on the MD trajectory. Cluster analysis of the formin binding protein 21 (FBP21) tandem WW domains demonstrate that the k-means clustering results by measuring distances between structures in the ED subspace are superior to those by using other metrics like pairwise inter-domain residue distances.

Experiments revealed that cytochrome P450 2C8 enzyme (CYP2C8) has two distinct substrate binding sites to the physiologically important molecules, retinoic acids, and the main difference between these two binding sites is whether there is a salt bridge interaction between the anionic carboxylate tail of retinoic acids and the surrounding protein environment. However, the influence of such salt bridge interaction toward catalysis is still elusive. In the present paper, density functional theory (DFT) calculations were employed to research the reaction mechanism of all-trans-retinoic acid (atRA) 4-hydroxylation mediated by CYP2C8. Our DFT calculations revealed that such salt bridge interaction has obvious effects on the reaction mechanism of atRA 4-hydroxylation. In the binding site containing a salt bridge interaction between the anionic carboxylate tail of atRA and the cationic guanidine group of Arg241, C–H bond activation proceeds via a normal hydrogen atom transfer (HAT) mechanism; in the other site without this salt bridge interaction, however, C–H bond activation is achieved via a stepwise electron transfer and hydrogen atom transfer, thus, a novel ET/HAT mechanism. These findings enrich the mechanism patterns of C–H bond activation catalyzed by metalloenzymes and their biomimetics. Meanwhile, the self-interaction error (SIE) problem encountered during our calculations in vacuum was affected and removed by the inclusion of an external electric field in the calculations.

A central theme in systems biology is to reveal the intricate relationship between structure and dynamics of many complex biological networks. Using Boolean models that describe yeast cell cycle process, we developed a unique logic-based theoretical framework to quantitatively determine the structure-dynamics mapping, also known as genotype–phenotype mapping. Moreover, under the dominant inhibition condition, we used a superposition property to show rigorously that the neutral network — the network of all possible structures that encode the same dynamics and are connected via single interaction mutations — forms one giant connected and conductive component. This may help shed light on the evolution landscape of biological networks based on the distance and speed a network can evolve on this neutral network.

In this paper the process of α-tocopherol (TCP) passing across DPPC membrane was simulated using both the potential of mean force (PMF) and the accelerated molecular dynamics (aMD) methods, respectively. Energy properties, hydrogen bonds and orientation have been compared between these two methods and several conclusions have been obtained. The results indicate that TCP tends to stay at z = 1.2 nm of lipid bilayer. The binding free energy profiles of these two methods are alike. All these show that aMD could obtain comparable results as PMF method, while needs less computation time and resources. Therefore, aMD method could be used as an alternative method for prediction of transport properties of drug-lipid system.

Polymorphism is widely observed in amyloid fibrils associated with many neurodegenerative diseases. Recent experimental study reported that fibrils formed by the segment NNQNTF of elk prion protein exhibited facial polymorphism with the two β-sheets either in back-to-back (BB) or in face-to-face (FF) packing arrangement. In the BB packing, the side chains of N2, N4 and F6 are interdigitated to form steric zipper, while in the FF packing, the side chains of N1, Q3 and T5 form the interdigitated interface. In this study, we investigate the water-mediated assembly of two preformed β-sheets and the physical interactions that stabilize the two different fibrils using all-atom molecular dynamics (MD) simulations. Multiple MD simulations have been carried out by starting from FF or BB packing of two β-sheets according to the facial polymorphism revealed by X-ray microcrystallography. For both packing patterns, we observe that the assembly of β-sheets is mediated by water molecules in the interface between β-sheets, leading to a long-lived protofibrils with wet interface prior to the formation of dry amyloid fibrils. Detailed structural analysis shows that besides the side chain steric zipper interactions, intra-sheet hydrogen bonding and aromatic stacking interactions play an important role on the stabilization of the protofibril with BB packing, while the intra-sheet and inter-sheet hydrogen bonding interactions are crucial for the formation of BB protofibril. These findings provide insights into the mechanism that lead to the facial polymorphism of NNQNTF fibrils.

Electronic spectra of uracil, thymine, adenine, guanine, and cytosine in the gas phase and aqueous solution have been studied by extensive time-dependent density functional calculations. Calculations show that the Quantum mechanics/molecular mechanics (QM/MM) geometry optimization based on the molecular dynamics (MD) equilibrated configuration can locate an optimal solvated cluster for the base solvation, and the combined QM/MM and cluster-continuum computational protocol is capable of handling the solvent effect on the excited states of nucleic acid bases and providing realistic absorption spectra in water environment with relatively low computational costs. Generally, the vertical excitation energies in aqueous solution by PCM/TD-X3LYP calculations show excellent agreement with the experimental observations and the maximum deviation is less than 0.2 eV. The present results reveal that the hydrogen bond network around the excited-state base and its dipole moment change may remarkably modify the absorption spectra of nucleic acid bases in aqueous solution.

Physics-based computational approaches to predicting the structure of macromolecules such as proteins are gaining increased use, but there are remaining challenges. In the current work, it is demonstrated that in energy-based prediction methods, the degree of optimization of the sampled structures can influence the prediction results. In particular, discrepancies in the degree of local sampling can bias the predictions in favor of the oversampled structures by shifting the local probability distributions of the minimum sampled energies. In simple systems, it is shown that the magnitude of the errors can be calculated from the energy surface, and for certain model systems, derived analytically. Further, it is shown that for energy wells whose forms differ only by a randomly assigned energy shift, the optimal accuracy of prediction is achieved when the sampling around each structure is equal. Energy correction terms can be used in cases of unequal sampling to reproduce the total probabilities that would occur under equal sampling, but optimal corrections only partially restore the prediction accuracy lost to unequal sampling. For multiwell systems, the determination of the correction terms is a multibody problem; it is shown that the involved cross-correlation multiple integrals can be reduced to simpler integrals. The possible implications of the current analysis for macromolecular structure prediction are discussed.

Protein Lin28 recognizes pre-let-7 microRNAs (miRNAs) through direct interactions between its zinc-knuckle type zinc-finger (ZnF) domains and the terminal loop of pre-let-7, resulting in the inhibition of the synthesis of mature let-7 miRNAs. Despite the physiological importance, the involved conformational changes and energetic factors contributing to the binding affinity and specificity remain unclear. We conducted molecular dynamics (MD) simulations in conjunction with molecular mechanics/generalized born surface area (MM/GBSA) and energy decomposition calculations to investigate the RNA binding-induced conformational changes of the ZnFs and the residual level energetic factors that influence the binding affinity and specificity. We showed that the binding of the RNA results in the inter-domain conformational changes of the two ZnF domains, including changes of the spatial relationships of several nucleobase-binding amino acids. We also observed mutation-induced weakening of the affinity of the Lin28–pre-let-7 binding, which reveals the importance of the stacking interactions between the side-chains of Tyr140, His148, His162 and the bases of nucleic acid G2 and G5 in the specific recognition of pre-let-7 by the ZnFs of Lin28.