An effective Coarse-grained model for biological simulations: Recent refinements and validations

Unraveling GPCRs Allosteric Modulation. Cannabinoid 1 Receptor as a Case Study

G-protein-coupled receptors (GPCRs) constitute one of the most prominent families of integral membrane receptor proteins that mediate most transmembrane signaling processes. Malfunction of these signal transduction processes is one of the underlying causes of many human pathologies (Parkinson's, Huntington's, heart diseases, etc), provoking that GPCRs are the largest family of druggable proteins. However, these receptors have been targeted traditionally by orthosteric ligands, which usually causes side effects due to the simultaneous targeting of homologous receptor subtypes. Allosteric modulation offers a promising alternative approach to circumvent this problematic and, thus, comprehending its details is a most important task. Here we use the Cannabinoid type-1 receptor (CB1R) in trying to shed light on this issue, focusing on positive allosteric modulation. This is done by using the protein-dipole Langevin-dipole (PDLD) within the linear response approximation (LRA) framework (PDLD/S-2000) along with our coarse-grained (CG) model of membrane proteins to evaluate the dissociation constants (K Bs) and cooperativity factors (αs) for a diverse series of CB1R positive allosteric modulators belonging to the 2-phenylindole structural class, considering CP55940 as an agonist. The agreement with the experimental data evinces that significantly populated allosteric modulator:CB1R and allosteric modulator:CP55940:CB1R complexes have been identified and characterized successfully. Analyzing them, it has been determined that CB1R positive allosteric modulation lies in an outwards displacement of transmembrane α helix (TM) 4 extracellular end and in the regulation of the range of motion of a compound TM7 movement for binary and ternary complexes, respectively. In this respect, we achieved a better comprehension of the molecular architecture of CB1R positive allosteric site, identifying Lys1923.28 and Gly1943.30 as key residues regarding electrostatic interactions inside this cavity, and to rationalize (at both structural and molecular level) the exhibited stereoselectivity in relation to positive allosteric modulation activity by considered CB1R allosteric modulators. Additionally, putative/postulated allosteric binding sites have been screened successfully, identifying the real CB1R positive allosteric site, and most structure-activity relationship (SAR) studies of CB1R 2-phenylindole allosteric modulators have been rationalized. All these findings point out towards the predictive value of the methodology used in the current work, which can be applied to other biophysical systems of interest. The results presented in this study contribute significantly to understand GPCRs allosteric modulation and, hopefully, will encourage a more thorough exploration of the topic.

Acidity Prediction in Arbitrary Solvents: Machine Learning Based on Semiempirical Molecular Orbital Calculation

Due to the nonlinearity of solvent effects, careful solvent selection is essential when using acids in different applications. However, there is a lack of measurements of pKa while systematically changing molecular structures and solvents. Consequently, there was no protocol to predict the acidity in arbitrary environments. This study developed an arbitrary environment pKa prediction protocol by integrating quantum chemical calculations using a polarizable continuum model and machine learning. This protocol constructed models to predict the acidity of biologically relevant molecules in water and candidate superstrong acids in organic solvents. For both systems, the pKa can be predicted with an average absolute error of 1.1 by learning relatively small number of data. This approach can also account for the nonlinear decay of acidity with solvation in different environments (compression effect). The versatility of the protocol extends to its applicability to a wide range of compounds, including those with complex electronic state changes upon proton dissociation, supporting research in diverse fields including, but not limited to, drug discovery and engineering.

Key Interaction Changes Determine the Activation Process of Human Parathyroid Hormone Type 1 Receptor

The parathyroid hormone type 1 receptor (PTH1R) plays a crucial role in modulating various physiological functions and is considered an effective therapeutic target for osteoporosis. However, a lack of detailed molecular and energetic information about PTH1R limits our comprehensive understanding of its activation process. In this study, we performed computational simulations to explore key events in the activation process, such as conformational changes in PTH1R, Gs protein coupling, and the release of guanosine diphosphate (GDP). Our analysis identified kinetic information, including the rate-determining step, transition state, and energy barriers. Free-energy and structural analyses revealed that GDP could be released from the Gs protein when the binding cavity is partially open. Additionally, we predicted important residues, including potential pathogenic mutations, and verified their significance through site-directed mutations. These findings enhance our understanding of class B GPCR activation mechanisms. Furthermore, the methodology employed in this study can be applied to other biophysical systems.

Exploring the Activation Process of the Glycine Receptor

Glycine receptors (GlyR) conduct inhibitory glycinergic neurotransmission in the spinal cord and the brainstem. They play an important role in muscle tone, motor coordination, respiration, and pain perception. However, the mechanism underlying GlyR activation remains unclear. There are five potential glycine binding sites in α1 GlyR, and different binding patterns may cause distinct activation or desensitization behaviors. In this study, we investigated the coupling of protein conformational changes and glycine binding events to elucidate the influence of binding patterns on the activation and desensitization processes of α1 GlyRs. Subsequently, we explored the energetic distinctions between the apical and lateral pathways during α1 GlyR conduction to identify the pivotal factors in the ion conduction pathway preference. Moreover, we predicted the mutational effects of the key residues and verified our predictions using electrophysiological experiments. For the mutants that can be activated by glycine, the predictions of the mutational directions were all correct. The strength of the mutational effects was assessed using Pearson's correlation coefficient, yielding a value of -0.77 between the calculated highest energy barriers and experimental maximum current amplitudes. These findings contribute to our understanding of GlyR activation, identify the key residues of GlyRs, and provide guidance for mechanistic studies on other pLGICs.

Catalytic mechanism study of ATP-citrate lyase during citryl-CoA synthesis process

ATP-citrate lyase (ACLY) is a critical metabolic enzyme and promising target for drug development. The structure determinations of ACLY have revealed its homotetramer states with various subunit symmetries, but catalytic mechanism of ACLY tetramer and the importance of subunit symmetry have not been clarified. Here, we constructed the free energy landscape of ACLY tetramer with arbitrary subunit symmetries and investigated energetic and conformational coupling of subunits during citryl-CoA synthesis process. The optimal conformational pathway indicates that ACLY tetramer encounters three critical conformational barriers and undergoes a loss of rigid-D2 symmetry to gain an energetic advantage. Energetic coupling of conformational changes and biochemical reactions suggests that these biological events are not independent but rather coupled with each other, showing a comparable energy barrier to the experimental data for the rate-limiting step. These findings could contribute to further research on catalytic mechanism, functional modulation, and inhibitor design of ACLY.

Exploring the Activation Mechanism of the GPR183 Receptor

The G protein-coupled receptors (GPCRs) play a pivotal role in numerous biological processes as crucial cell membrane receptors. However, the dynamic mechanisms underlying the activation of GPR183, a specific GPCR, remain largely elusive. To address this, we employed computational simulation techniques to elucidate the activation process and key events associated with GPR183, including conformational changes from inactive to active state, binding interactions with the Gi protein complex, and GDP release. Our findings demonstrate that the association between GPR183 and the Gi protein involves the formation of receptor-specific conformations, the gradual proximity of the Gi protein to the binding pocket, and fine adjustments of the protein conformation, ultimately leading to a stable GPR183-Gi complex characterized by a high energy barrier. The presence of Gi protein partially promotes GPR183 activation, which is consistent with the observation of GPCR constitutive activity test experiments, thus illustrating the reliability of our calculations. Moreover, our study suggests the existence of a stable partially activated state preceding complete activation, providing novel avenues for future investigations. In addition, the relevance of GPR183 for various diseases, such as colitis, the response of eosinophils to Mycobacterium tuberculosis infection, antiviral properties, and pulmonary inflammation, has been emphasized, underscoring its therapeutic potential. Consequently, understanding the activation process of GPR183 through molecular dynamic simulations offers valuable kinetic insights that can aid in the development of targeted therapies.

Mechanistic study of the transmission pattern of the SARS-CoV-2 omicron variant

The omicron variant of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) characterized by 30 mutations in its spike protein, has rapidly spread worldwide since November 2021, significantly exacerbating the ongoing COVID-19 pandemic. In order to investigate the relationship between these mutations and the variant's high transmissibility, we conducted a systematic analysis of the mutational effect on spike-angiotensin-converting enzyme-2 (ACE2) interactions and explored the structural/energy correlation of key mutations, utilizing a reliable coarse-grained model. Our study extended beyond the receptor-binding domain (RBD) of spike trimer through comprehensive modeling of the full-length spike trimer rather than just the RBD. Our free-energy calculation revealed that the enhanced binding affinity between the spike protein and the ACE2 receptor is correlated with the increased structural stability of the isolated spike protein, thus explaining the omicron variant's heightened transmissibility. The conclusion was supported by our experimental analyses involving the expression and purification of the full-length spike trimer. Furthermore, the energy decomposition analysis established those electrostatic interactions make major contributions to this effect. We categorized the mutations into four groups and established an analytical framework that can be employed in studying future mutations. Additionally, our calculations rationalized the reduced affinity of the omicron variant towards most available therapeutic neutralizing antibodies, when compared with the wild type. By providing concrete experimental data and offering a solid explanation, this study contributes to a better understanding of the relationship between theories and observations and lays the foundation for future investigations.

Exploring the activation mechanism of metabotropic glutamate receptor 2

Homomeric dimerization of metabotropic glutamate receptors (mGlus) is essential for the modulation of their functions and represents a promising avenue for the development of novel therapeutic approaches to address central nervous system diseases. Yet, the scarcity of detailed molecular and energetic data on mGlu2 impedes our in-depth comprehension of their activation process. Here, we employ computational simulation methods to elucidate the activation process and key events associated with the mGlu2, including a detailed analysis of its conformational transitions, the binding of agonists, Gi protein coupling, and the guanosine diphosphate (GDP) release. Our results demonstrate that the activation of mGlu2 is a stepwise process and several energy barriers need to be overcome. Moreover, we also identify the rate-determining step of the mGlu2's transition from the agonist-bound state to its active state. From the perspective of free-energy analysis, we find that the conformational dynamics of mGlu2's subunit follow coupled rather than discrete, independent actions. Asymmetric dimerization is critical for receptor activation. Our calculation results are consistent with the observation of cross-linking and fluorescent-labeled blot experiments, thus illustrating the reliability of our calculations. Besides, we also identify potential key residues in the Gi protein binding position on mGlu2, mGlu2 dimer's TM6-TM6 interface, and Gi α5 helix by the change of energy barriers after mutation. The implications of our findings could lead to a more comprehensive grasp of class C G protein-coupled receptor activation.

Energetic and structural insights behind calcium induced conformational transition in calmodulin

Calmodulin (CaM) is a key signaling protein that triggers several cellular and physiological processes inside the cell. Upon binding with calcium ion, CaM undergoes large scale conformational transition from a closed state to an open state that facilitates its interaction with various target protein and regulates their activity. This work explores the origin of the energetic and structural variation of the wild type and mutated CaM and explores the molecular origin for the structural differences between them. We first calculated the sequential calcium binding energy to CaM using the PDLD/S-LRA/β approach. This study  shows a very good correlation with experimental calcium binding energies. Next we calculated the calcium binding energies to the wild type CaM and several mutated CaM systems which were reported experimentally. On the structural aspect, it has been reported experimentally that certain mutation (Q41L-K75I) in calcium bound CaM leads to complete conformational transition from an open to a closed state. By using equilibrium molecular dynamics simulation, free energy calculation and contact frequency map analysis, we have shown that the formation of a cluster of long-range hydrophobic contacts, initiated by the Q41L-K75I CaM variant is the driving force behind its closing motion. This study unravels the energetics and structural aspects behind calcium ion induced conformational changes in wild type CaM and its variant.

Predicting Mutational Effects on Ca2+-Activated Chloride Conduction of TMEM16A Based on a Simulation Study

The dysfunction and defects of ion channels are associated with many human diseases, especially for loss-of-function mutations in ion channels such as cystic fibrosis transmembrane conductance regulator mutations in cystic fibrosis. Understanding ion channels is of great current importance for both medical and fundamental purposes. Such an understanding should include the ability to predict mutational effects and describe functional and mechanistic effects. In this work, we introduce an approach to predict mutational effects based on kinetic information (including reaction barriers and transition state locations) obtained by studying the working mechanism of target proteins. Specifically, we take the Ca2+-activated chloride channel TMEM16A as an example and utilize the computational biology model to predict the mutational effects of key residues. Encouragingly, we verified our predictions through electrophysiological experiments, demonstrating a 94% prediction accuracy regarding mutational directions. The mutational strength assessed by Pearson's correlation coefficient is -0.80 between our calculations and the experimental results. These findings suggest that the proposed methodology is reliable and can provide valuable guidance for revealing functional mechanisms and identifying key residues of the TMEM16A channel. The proposed approach can be extended to a broad scope of biophysical systems.

Interaction of SARS-CoV-2 with host cells and antibodies: experiment and simulation

Severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) is the causative agent of the devastating global COVID-19 pandemic announced by WHO in March 2020. Through unprecedented scientific effort, several vaccines, drugs and antibodies have been developed, saving millions of lives, but the fight against COVID-19 continues as immune escape variants of concern such as Delta and Omicron emerge. To develop more effective treatments and to elucidate the side effects caused by vaccines and therapeutic agents, a deeper understanding of the molecular interactions of SARS-CoV-2 with them and human cells is required. With special interest in computational approaches, we will focus on the structure of SARS-CoV-2 and the interaction of its spike protein with human angiotensin-converting enzyme-2 (ACE2) as a prime entry point of the virus into host cells. In addition, other possible viral receptors will be considered. The fusion of viral and human membranes and the interaction of the spike protein with antibodies and nanobodies will be discussed, as well as the effect of SARS-CoV-2 on protein synthesis in host cells.

Physical pictures of rotation mechanisms of F1- and V1-ATPases: Leading roles of translational, configurational entropy of water

We aim to develop a theory based on a concept other than the chemo-mechanical coupling (transduction of chemical free energy of ATP to mechanical work) for an ATP-driven protein complex. Experimental results conflicting with the chemo-mechanical coupling have recently emerged. We claim that the system comprises not only the protein complex but also the aqueous solution in which the protein complex is immersed and the system performs essentially no mechanical work. We perform statistical-mechanical analyses on V₁-ATPase (the A₃B₃DF complex) for which crystal structures in more different states are experimentally known than for F₁-ATPase (the α₃β₃γ complex). Molecular and atomistic models are employed for water and the structure of V₁-ATPase, respectively. The entropy originating from the translational displacement of water molecules in the system is treated as a pivotal factor. We find that the packing structure of the catalytic dwell state of V₁-ATPase is constructed by the interplay of ATP bindings to two of the A subunits and incorporation of the DF subunit. The packing structure represents the nonuniformity with respect to the closeness of packing of the atoms in constituent proteins and protein interfaces. The physical picture of rotation mechanism of F₁-ATPase recently constructed by Kinoshita is examined, and common points and differences between F₁- and V₁-ATPases are revealed. An ATP hydrolysis cycle comprises binding of ATP to the protein complex, hydrolysis of ATP into ADP and Pi in it, and dissociation of ADP and Pi from it. During each cycle, the chemical compounds bound to the three A or β subunits and the packing structure of the A₃B₃ or α₃β₃ complex are sequentially changed, which induces the unidirectional rotation of the central shaft for retaining the packing structure of the A₃B₃DF or α₃β₃γ complex stabilized for almost maximizing the water entropy. The torque driving the rotation is generated by water with no input of chemical free energy. The presence of ATP is indispensable as a trigger of the torque generation. The ATP hydrolysis or synthesis reaction is tightly coupled to the rotation of the central shaft in the normal or inverse direction through the water-entropy effect.

Exploring the Phospholipid Transport Mechanism of ATP8A1-CDC50

P4-ATPase translocates lipids from the exoplasmic to the cytosolic plasma membrane leaflet to maintain lipid asymmetry distribution in eukaryotic cells. P4-ATPase is associated with severe neurodegenerative and metabolic diseases such as neurological and motor disorders. Thus, it is important to understand its transport mechanism. However, even with progress in X-ray diffraction and cryo-electron microscopy techniques, it is difficult to obtain the dynamic information of the phospholipid transport process in detail. There are still some problems required to be resolved: (1) when does the lipid transport happen? (2) How do the key residues on the transmembrane helices contribute to the free energy of important states? In this work, we explore the phospholipid transport mechanism using a coarse-grained model and binding free energy calculations. We obtained the free energy landscape by coupling the protein conformational changes and the phospholipid transport event, taking ATP8A1-CDC50 (the typical subtype of P4-ATPase) as the research object. According to the results, we found that the phospholipid would bind to the ATP8A1-CDC50 at the early stage when ATP8A1-CDC50 changes from E2P to E2Pi-PL state. We also found that the electrostatic effects play crucial roles in the phospholipid transport process. The information obtained from this work could help us in designing novel drugs for P-type flippase disorders.

A systematic study on the binding affinity of SARS-CoV-2 spike protein to antibodies

The COVID-19 pandemic has caused a worldwide health crisis and economic recession. Effective prevention and treatment methods are urgently required to control the pandemic. However, the emergence of novel SARS-CoV-2 variants challenges the effectiveness of currently available vaccines and therapeutic antibodies. In this study, through the assessment of binding free energies, we analyzed the mutational effects on the binding affinity of the coronavirus spike protein to neutralizing antibodies, patient-derived antibodies, and artificially designed antibody mimics. We designed a scoring method to assess the immune evasion ability of viral variants. We also evaluated the differences between several targeting sites on the spike protein of antibodies. The results presented herein might prove helpful in the development of more effective therapies in the future.

The Nature of Functional Features of Different Classes of G-Protein-Coupled Receptors

G-protein-coupled receptors (GPCRs) are a critical family in the human proteome and are involved in various physiological processes. They are also the most important drug target, with approximately 30% of approved drugs acting on such receptors. The members of the family are divided into six classes based on their structural and functional characteristics. Understanding their structural-functional relationships will benefit us in future drug development. In this article, we investigate the features of protein function, structure, and energy that describe the dynamics of the GPCR activation process between different families. GPCRs straddle the cell membrane and transduce signals from outside the membrane into the cell. During the process, the conformational change in GPCRs that is activated by the binding of signal molecules is essential. During the binding process, different types of signal molecules result in different signal transfer efficiencies. Therefore, the GPCR classes show a variety of structures and activation processes. Based on the experimental crystal structures, we modeled the activation process of the β2 adrenergic receptor (β2AR), glucagon receptor (GCGR), and metabotropic glutamate receptor 2 (mGluR2), which represent class A, B, and C GPCRs, respectively. We calculated their activation free-energy landscapes and analyzed the structure-energy-function relationship. The results show a consistent picture of the activation mechanisms between different types of GPCRs. This could also provide us a way to understand other signal transduction proteins.

Application of Computational Biology and Artificial Intelligence in Drug Design

Traditional drug design requires a great amount of research time and developmental expense. Booming computational approaches, including computational biology, computer-aided drug design, and artificial intelligence, have the potential to expedite the efficiency of drug discovery by minimizing the time and financial cost. In recent years, computational approaches are being widely used to improve the efficacy and effectiveness of drug discovery and pipeline, leading to the approval of plenty of new drugs for marketing. The present review emphasizes on the applications of these indispensable computational approaches in aiding target identification, lead discovery, and lead optimization. Some challenges of using these approaches for drug design are also discussed. Moreover, we propose a methodology for integrating various computational techniques into new drug discovery and design.

Mechanism Study of Proteins under Membrane Environment

Membrane proteins play crucial roles in various physiological processes, including molecule transport across membranes, cell communication, and signal transduction. Approximately 60% of known drug targets are membrane proteins. There is a significant need to deeply understand the working mechanism of membrane proteins in detail, which is a challenging work due to the lack of available membrane structures and their large spatial scale. Membrane proteins carry out vital physiological functions through conformational changes. In the current study, we utilized a coarse-grained (CG) model to investigate three representative membrane protein systems: the TMEM16A channel, the family C GPCRs mGlu2 receptor, and the P4-ATPase phospholipid transporter. We constructed the reaction pathway of conformational changes between the two-end structures. Energy profiles and energy barriers were calculated. These data could provide reasonable explanations for TMEM16A activation, the mGlu2 receptor activation process, and P4-ATPase phospholipid transport. Although they all belong to the members of membrane proteins, they behave differently in terms of energy. Our work investigated the working mechanism of membrane proteins and could give novel insights into other membrane protein systems of interest.

Automated Protein Secondary Structure Assignment from Cα Positions Using Neural Networks

The assignment of secondary structure elements in protein conformations is necessary to interpret a protein model that has been established by computational methods. The process essentially involves labeling the amino acid residues with H (Helix), E (Strand), or C (Coil, also known as Loop). When particular atoms are absent from an input protein structure, the procedure becomes more complicated, especially when only the alpha carbon locations are known. Various techniques have been tested and applied to this problem during the last forty years. The application of machine learning techniques is the most recent trend. This contribution presents the HECA classifier, which uses neural networks to assign protein secondary structure types. The technique exclusively employs Cα coordinates. The Keras (TensorFlow) library was used to implement and train the neural network model. The BioShell toolkit was used to calculate the neural network input features from raw coordinates. The study’s findings show that neural network-based methods may be successfully used to take on structure assignment challenges when only Cα trace is available. Thanks to the careful selection of input features, our approach’s accuracy (above 97%) exceeded that of the existing methods.

Application of Coarse-Grained (CG) Models to Explore Conformational Pathway of Large-Scale Protein Machines

Protein machines are clusters of protein assemblies that function in order to control the transfer of matter and energy in cells. For a specific protein machine, its working mechanisms are not only determined by the static crystal structures, but also related to the conformational transition dynamics and the corresponding energy profiles. With the rapid development of crystallographic techniques, the spatial scale of resolved structures is reaching up to thousands of residues, and the concomitant conformational changes become more and more complicated, posing a great challenge for computational biology research. Previously, a coarse-grained (CG) model aiming at conformational free energy evaluation was developed and showed excellent ability to reproduce the energy profiles by accurate electrostatic interaction calculations. In this study, we extended the application of the CG model to a series of large-scale protein machine systems. The spike protein trimer of SARS-CoV-2, ATP citrate lyase (ACLY) tetramer, and P4-ATPases systems were carefully studied and discussed as examples. It is indicated that the CG model is effective to depict the energy profiles of the conformational pathway between two endpoint structures, especially for large-scale systems. Both the energy change and energy barrier between endpoint structures provide reasonable mechanism explanations for the associated biological processes, including the opening of receptor binding domain (RBD) of spike protein, the phospholipid transportation of P4-ATPase, and the loop translocation of ACLY. Taken together, the CG model provides a suitable alternative in mechanistic studies related to conformational change in large-scale protein machines.

Molecular Modeling is an Enabling Approach to Complement and Enhance Channelopathy Research

Hundreds of human membrane proteins form channels that transport necessary ions and compounds, including drugs and metabolites, yet details of their normal function or how function is altered by genetic variants to cause diseases are often unknown. Without this knowledge, researchers are less equipped to develop approaches to diagnose and treat channelopathies. High-resolution computational approaches such as molecular modeling enable researchers to investigate channelopathy protein function, facilitate detailed hypothesis generation, and produce data that is difficult to gather experimentally. Molecular modeling can be tailored to each physiologic context that a protein may act within, some of which may currently be difficult or impossible to assay experimentally. Because many genomic variants are observed in channelopathy proteins from high-throughput sequencing studies, methods with mechanistic value are needed to interpret their effects. The eminent field of structural bioinformatics integrates techniques from multiple disciplines including molecular modeling, computational chemistry, biophysics, and biochemistry, to develop mechanistic hypotheses and enhance the information available for understanding function. Molecular modeling and simulation access 3D and time-dependent information, not currently predictable from sequence. Thus, molecular modeling is valuable for increasing the resolution with which the natural function of protein channels can be investigated, and for interpreting how genomic variants alter them to produce physiologic changes that manifest as channelopathies. © 2022 American Physiological Society. Compr Physiol 12:3141-3166, 2022.

Predicting Mutational Effects on Receptor Binding of the Spike Protein of SARS-CoV-2 Variants

The pandemic caused by SARS-CoV-2 has cost millions of lives and tremendous social/financial loss. The virus continues to evolve and mutate. In particular, the recently emerged "UK", "South Africa", and Delta variants show higher infectivity and spreading speed. Thus, the relationship between the mutations of certain amino acids and the spreading speed of the virus is a problem of great importance. In this respect, understanding the mutational mechanism is crucial for surveillance and prediction of future mutations as well as antibody/vaccine development. In this work, we used a coarse-grained model (that was used previously in predicting the importance of mutations of N501) to calculate the free energy change of various types of single-site or combined-site mutations. This was done for the UK, South Africa, and Delta mutants. We investigated the underlying mechanisms of the binding affinity changes for mutations at different spike protein domains of SARS-CoV-2 and provided the energy basis for the resistance of the E484 mutant to the antibody m396. Other potential mutation sites were also predicted. Furthermore, the in silico predictions were assessed by functional experiments. The results establish that the faster spreading of recently observed mutants is strongly correlated with the binding-affinity enhancement between virus and human receptor as well as with the reduction of the binding to the m396 antibody. Significantly, the current approach offers a way to predict new variants and to assess the effectiveness of different antibodies toward such variants.

Membrane-binding mechanism of the EEA1 FYVE domain revealed by multi-scale molecular dynamics simulations

Early Endosomal Antigen 1 (EEA1) is a key protein in endosomal trafficking and is implicated in both autoimmune and neurological diseases. The C-terminal FYVE domain of EEA1 binds endosomal membranes, which contain phosphatidylinositol-3-phosphate (PI(3)P). Although it is known that FYVE binds PI(3)P specifically, it has not previously been described of how FYVE attaches and binds to endosomal membranes. In this study, we employed both coarse-grained (CG) and atomistic (AT) molecular dynamics (MD) simulations to determine how FYVE binds to PI(3)P-containing membranes. CG-MD showed that the dominant membrane binding mode resembles the crystal structure of EEA1 FYVE domain in complex with inositol-1,3-diphospate (PDB ID 1JOC). FYVE, which is a homodimer, binds the membrane via a hinge mechanism, where the C-terminus of one monomer first attaches to the membrane, followed by the C-terminus of the other monomer. The estimated total binding energy is ~70 kJ/mol, of which 50-60 kJ/mol stems from specific PI(3)P-interactions. By AT-MD, we could partition the binding mode into two types: (i) adhesion by electrostatic FYVE-PI(3)P interaction, and (ii) insertion of amphipathic loops. The AT simulations also demonstrated flexibility within the FYVE homodimer between the C-terminal heads and coiled-coil stem. This leads to a dynamic model whereby the 200 nm long coiled coil attached to the FYVE domain dimer can amplify local hinge-bending motions such that the Rab5-binding domain at the other end of the coiled coil can explore an area of 0.1 μm2 in the search for a second endosome with which to interact.

Exploring the Activation Process of the β2AR-Gs Complex

G-Protein-coupled receptors (GPCRs) belong to an important family of integral membrane receptor proteins that are essential for a variety of transmembrane signaling process, such as vision, olfaction, and hormone responses. They are also involved in many human diseases (Alzheimer's, heart diseases, etc.) and are therefore common drug targets. Thus, understanding the details of the GPCR activation process is a task of major importance. Various types of crystal structures of GPCRs have been solved either at stable end-point states or at possible intermediate states. However, the detailed mechanism of the activation process is still poorly understood. For example, it is not completely clear when the nucleotide release from the G protein occurs and how the key residues on α5 contribute to the coupling process and further affect the binding specificity. In this work we show by free energy analysis that the guanosine diphosphate (GDP) molecule could be released from the Gs protein when the binding cavity is half open. This occurs during the transition to the Gs open state, which is the rate-determining step in the system conformational change. We also account for the experimentally observed slow-down effects by the change of the reaction barriers after mutations. Furthermore, we identify potential key residues on α5 and validated their significance by site-directed mutagenesis, which illustrates that computational works have predictive value even for complex biophysical systems. The methodology of the current work may be applied to other biophysical systems of interest.

The catalytic dwell in ATPases is not crucial for movement against applied torque

The ATPase-catalysed conversion of ATP to ADP is a fundamental process in biology. During the hydrolysis of ATP, the α₃β₃ domain undergoes conformational changes while the central stalk (γ/D) rotates unidirectionally. Experimental studies have suggested that different catalytic mechanisms operate depending on the type of ATPase, but the structural and energetic basis of these mechanisms remains unclear. In particular, it is not clear how the positions of the catalytic dwells influence the energy transduction. Here we show that the observed dwell positions, unidirectional rotation and movement against the applied torque are reflections of the free-energy surface of the systems. Instructively, we determine that the dwell positions do not substantially affect the stopping torque. Our results suggest that the three resting states and the pathways that connect them should not be treated equally. The current work demonstrates how the free-energy landscape determines the behaviour of different types of ATPases.

Critical Differences between the Binding Features of the Spike Proteins of SARS-CoV-2 and SARS-CoV

The COVID-19 caused by SARS-CoV-2 has spread globally and caused tremendous loss of lives and properties, and it is of utmost urgency to understand its propagation process and to find ways to slow down the epidemic. In this work, we used a coarse-grained model to calculate the binding free energy of SARS-CoV-2 or SARS-CoV to their human receptor ACE2. The investigation of the free energy contribution of the interacting residues indicates that the residues located outside the receptor binding domain are the source of the stronger binding of the novel virus. Thus, the current results suggest that the essential evolution of SARS-CoV-2 happens remotely from the binding domain at the spike protein trimeric body. Such evolution may facilitate the conformational change and the infection process that occurs after the virus is bound to ACE2. By studying the binding pattern between SARS-CoV antibody m396 and SARS-CoV-2, it is found that the remote energetic contribution is missing, which might explain the absence of cross-reactivity of such antibodies.

Combinatorial Approach for Exploring Conformational Space and Activation Barriers in Computer-Aided Enzyme Design

Computer-aided enzyme design is a field of great potential importance for biotechnological applications, medical advances, and a fundamental understanding of enzyme action. However, reaching a predictive ability in this direction is extremely challenging. It requires both the ability to predict quantitatively the activation barriers in cases where the structure and sequence are known and the ability to predict the effect of different mutations. In this work, we propose a protocol for predicting reasonable starting structures of mutants of proteins with known structures and for calculating the activation barriers of the generated mutants. Our approach also allows us to use the predicted structures of the generated mutant to predict structures and activation barriers for subsequent set of mutations. This protocol is used to examine the reliability of the in silico directed evolution of Kemp eliminase and haloalkane dehalogenase. We also used the results of single and double mutations as a base for predicting the effect of transition-state stabilization by multiple concurrent mutations. This strategy seems to be useful in creating an activity funnel that provides a qualitative ranking of the catalytic power of different mutants.

Exploring alternative catalytic mechanisms of the Cas9 HNH domain

Understanding the reaction mechanism of CRISPR-associated protein 9 (Cas9) is crucial for the application of programmable gene editing. Despite the availability of the structures of Cas9 in apo- and substrate-bound forms, the catalytically active structure is still unclear. Our first attempt to explore the catalytic mechanism of Cas9 HNH domain has been based on the reasonable assumption that we are dealing with the same mechanism as endonuclease VII, including the assumption that the catalytic water is in the first shell of the Mg2+ . Trying this mechanism with the cryo-EM structure forced us to induce significant structural change driven by the movement of K848 (or other positively charged residue) close to the active site to facilitate the proton transfer step. In the present study, we explore a second reaction mechanism where the catalytic water is in the second shell of the Mg2+ and assume that the cryo-EM structure by itself is a suitable representation of a catalytic-ready structure. The alternative mechanism indicates that if the active water is from the second shell, then the calculated reaction barrier is lower compared with the corresponding barrier when the water comes from the first shell.

Revisiting the protomotive vectorial motion of F0-ATPase

The elucidation of the detailed mechanism used by F₀ to convert proton gradient to torque and rotational motion presents a major puzzle despite significant biophysical and structural progress. Although the conceptual model has advanced our understanding of the working principles of such systems, it is crucial to explore the actual mechanism using structure-based models that actually reproduce a unidirectional proton-driven rotation. Our previous work used a coarse-grained (CG) model to simulate the action of F₀ However, the simulations were based on a very tentative structural model of the interaction between subunit a and subunit c. Here, we again use a CG model but with a recent cryo-EM structure of cF₁F₀ and also explore the proton path using our water flooding and protein dipole Langevin dipole semimacroscopic formalism with its linear response approximation version (PDLD/S-LRA) approaches. The simulations are done in the combined space defined by the rotational coordinate and the proton transport coordinate. The study reproduced the effect of the protomotive force on the rotation of the F₀ while establishing the electrostatic origin of this effect. Our landscape reproduces the correct unidirectionality of the synthetic direction of the F₀ rotation and shows that it reflects the combined electrostatic coupling between the proton transport path and the c-ring conformational change. This work provides guidance for further studies in other proton-driven mechanochemical systems and should lead (when combined with studies of F₁) to a complete energy transduction picture of the F₀F₁-ATPase system.

A free-energy landscape for the glucagon-like peptide 1 receptor GLP1R

G-protein-coupled receptors (GPCRs) are among the most important receptors in human physiology and pathology. They serve as master regulators of numerous key processes and are involved in as well as cause debilitating diseases. Consequently, GPCRs are among the most attractive targets for drug design and pharmaceutical interventions (>30% of drugs on the market). The glucagon-like peptide 1 (GLP-1) hormone receptor GLP1R is closely involved in insulin secretion by pancreatic β-cells and constitutes a major druggable target for the development of anti-diabetes and obesity agents. GLP1R structure was recently solved, with ligands, allosteric modulators and as part of a complex with its cognate G protein. However, the translation of this structural data into structure/function understanding remains limited. The current study functionally characterizes GLP1R with special emphasis on ligand and cellular partner binding interactions and presents a free-energy landscape as well as a functional model of the activation cycle of GLP1R. Our results should facilitate a deeper understanding of the molecular mechanism underlying GLP1R activation, forming a basis for improved development of targeted therapeutics for diabetes and related disorders.

Exploring the Catalytic Mechanism of Cas9 Using Information Inferred from Endonuclease VII

Elucidating the nature of the gene editing mechanism of CRISPR (Clustered Regularly Interspaced Short Palindromic Repeats) is an important task in view of the role of this breakthrough to the advancement of human medicine. In particular, it is crucial to understand the catalytic mechanism of Cas9 (one of the CRISPR associated proteins) and its role in confirming accurate editing. Thus, we focus in this work on an attempt to analyze the catalytic mechanism of Cas9. Considering the absence of detailed structural information on the active form of Cas9, we use an empirical valence bond (EVB) which is calibrated on the closely related mechanism of T4 endonuclease VII. The calibrated EVB is then used in studying the reaction of Cas9, while trying several structural models. It is found that the catalytic activation requires a large conformational change, where K848 or other positively charged group moves from a relatively large distance toward the scissile phosphate. This conformational change leads to the change in position of the Mg2+ ion and to a major reduction in the activation barrier for the catalytic reaction. Our finding provides an important clue on the nature of the catalytic activation of CAS9 and thus should help in elucidating a key aspect of the gene editing process. For example, the approach used here should be effective in exploring the nature of off target activation and its relationship to the energetics of the unwinding process. This strategy may offer ways to improve the selectivity of Cas9.

Demonstrating aspects of multiscale modeling by studying the permeation pathway of the human ZnT2 zinc transporter

Multiscale modeling provides a very powerful means of studying complex biological systems. An important component of this strategy involves coarse-grained (CG) simplifications of regions of the system, which allow effective exploration of complex systems. Here we studied aspects of CG modeling of the human zinc transporter ZnT2. Zinc is an essential trace element with 10% of the proteins in the human proteome capable of zinc binding. Thus, zinc deficiency or impairment of zinc homeostasis disrupt key cellular functions. Mammalian zinc transport proceeds via two transporter families: ZnT and ZIP; however, little is known about the zinc permeation pathway through these transporters. As a step towards this end, we herein undertook comprehensive computational analyses employing multiscale techniques, focusing on the human zinc transporter ZnT2 and its bacterial homologue, YiiP. Energy calculations revealed a favorable pathway for zinc translocation via alternating access. We then identified key residues presumably involved in the passage of zinc ions through ZnT2 and YiiP, and functionally validated their role in zinc transport using site-directed mutagenesis of ZnT2 residues. Finally, we use a CG Monte Carlo simulation approach to sample the transition between the inward-facing and the outward-facing states. We present our structural models of the inward- and outward-facing conformations of ZnT2 as a blueprint prototype of the transporter conformations, including the putative permeation pathway and participating residues. The insights gained from this study may facilitate the delineation of the pathways of other zinc transporters, laying the foundations for the molecular basis underlying ion permeation. This may possibly facilitate the development of therapeutic interventions in pathological states associated with zinc deficiency and other disorders based on loss-of-function mutations in solute carriers.

Exploring the free-energy landscape of GPCR activation

G-protein-coupled receptors (GPCRs) are a large group of membrane-bound receptor proteins that are involved in a plethora of diverse processes (e.g., vision, hormone response). In mammals, and particularly in humans, GPCRs are involved in many signal transduction pathways and, as such, are heavily studied for their immense pharmaceutical potential. Indeed, a large fraction of drugs target various GPCRs, and drug-development is often aimed at GPCRs. Therefore, understanding the activation of GPCRs is a challenge of major importance both from fundamental and practical considerations. And yet, despite the remarkable progress in structural understanding, we still do not have a translation of the structural information to an energy-based picture. Here we use coarse-grained (CG) modeling to chart the free-energy landscape of the activation process of the β-2 adrenergic receptor (β₂AR) as a representative GPCR. The landscape provides the needed tool for analyzing the processes that lead to activation of the receptor upon binding of the ligand (adrenaline) while limiting constitutive activation. Our results pave the way to better understand the biological mechanisms of action of the β₂AR and GPCRs, from a physical chemistry point of view rather than simply by observing the receptor's behavior physiologically.

Lipid molecules can induce an opening of membrane-facing tunnels in cytochrome P450 1A2

Cytochrome P450 1A2 (P450 1A2, CYP1A2) is a membrane-bound enzyme that oxidizes a broad range of hydrophobic substrates. The structure and dynamics of both the catalytic and trans-membrane (TM) domains of this enzyme in the membrane/water environment were investigated using a multiscale computational approach, including coarse-grained and all-atom molecular dynamics. Starting from the spontaneous self-assembly of the system containing the TM or soluble domain immersed in randomized dilauroyl phosphatidylcholine (DLPC)/water mixture into their respective membrane-bound forms, we reconstituted the membrane-bound structure of the full-length P450 1A2. This structure includes a TM helix that spans the membrane, while being connected to the catalytic domain by a short flexible loop. Furthermore, in this model, the upper part of the TM helix interacts directly with a conserved and highly hydrophobic N-terminal proline-rich segment of the catalytic domain; this segment and the FG loop are immersed in the membrane, whereas the remaining portion of the catalytic domain remains exposed to aqueous solution. The shallow membrane immersion of the catalytic domain induces a depression in the opposite intact layer of the phospholipids. This structural effect may help in stabilizing the position of the TM helix directly beneath the catalytic domain. The partial immersion of the catalytic domain also allows for the enzyme substrates to enter the active site from either aqueous solution or phospholipid environment via several solvent- and membrane-facing tunnels in the full-length P450 1A2. The calculated tunnel dynamics indicated that the opening probability of the membrane-facing tunnels is significantly enhanced when a DLPC molecule spontaneously penetrates into the membrane-facing tunnel 2d. The energetics of the lipid penetration process were assessed by the linear interaction energy (LIE) approximation, and found to be thermodynamically feasible.

The FOF1 ATP synthase: from atomistic three-dimensional structure to the rotary-chemical function

Molecular motors are multi-subunit complexes that are indispensable for accomplishing various tasks of the living cells. One such molecular motor is the F_OF₁ ATP synthase that synthesizes ATP at the expense of the membrane proton gradient. Elucidating the molecular origin of the motor function is challenging despite significant advances in various experimental fields. Currently atomic simulations of whole motor complexes cannot reach to functionally relevant time scales that extend beyond the millisecond regime. Moreover, to reveal the underlying molecular origin of the function, one must model the coupled chemical and conformational events using physically and chemically meaningful multiscaling techniques. In this review, we discuss our approach to model the action of the F₁ and F_O molecular motors, where emphasis is laid on elucidating the molecular origin of the driving force that leads to directional rotation at the expense of ATP hydrolysis or proton gradients. We have used atomic structures of the motors and used hierarchical multiscaling techniques to generate low dimensional functional free energy surfaces of the complete mechano-chemical process. These free energy surfaces were studied further to calculate important characteristics of the motors, such as, rotational torque, temporal dynamics, occurrence of intermittent dwell states, etc. We also studied the result of mutating various parts of the motor domains and our observations correspond very well with the experimental findings. Overall, our studies have generated a cumulative understanding of the motor action, and especially highlight the crucial role of electrostatics in establishing the mechano-chemical coupling.

Simulating the fidelity and the three Mg mechanism of pol η and clarifying the validity of transition state theory in enzyme catalysis

Pol η belongs to the important Y family of DNA polymerases that can catalyze translesion synthesis across sites of damaged DNA. This activity involves the reduced fidelity of Pol η for 8-oxo-7,8-dhyedro-2'-deoxoguanosin(8-oxoG). The fundamental interest in Pol η has grown recently with the demonstration of the importance of a 3rd Mg2+ ion. The current work explores both the fidelity of Pol η and the role of the 3rd metal ion, by using empirical valence bond (EVB) simulations. The simulations reproduce the observed trend in fidelity and shed a new light on the role of the 3rd metal ion. It is found that this ion does not lead to a major catalytic effect, but most probably plays an important role in reducing the product release barrier. Furthermore, it is concluded, in contrast to some implications, that the effect of this metal does not violate transition state theory, and the evaluation of the catalytic effect must conserve the molecular composition upon moving from the reactant to the transition state. Proteins 2017; 85:1446-1453. © 2017 Wiley Periodicals, Inc.

Simulating the dynamics of the mechanochemical cycle of myosin-V

The detailed dynamics of the cycle of myosin-V are explored by simulation approaches, examining the nature of the energy-driven motion. Our study started with Langevin dynamics (LD) simulations on a very coarse landscape with a single rate-limiting barrier and reproduced the stall force and the hand-over-hand dynamics. We then considered a more realistic landscape and used time-dependent Monte Carlo (MC) simulations that allowed trajectories long enough to reproduce the force/velocity characteristic sigmoidal correlation, while also reproducing the hand-over-hand motion. Overall, our study indicated that the notion of a downhill lever-up to lever-down process (popularly known as the powerstroke mechanism) is the result of the energetics of the complete myosin-V cycle and is not the source of directional motion or force generation on its own. The present work further emphasizes the need to use well-defined energy landscapes in studying molecular motors in general and myosin in particular.

Exploring the Dependence of QM/MM Calculations of Enzyme Catalysis on the Size of the QM Region

Although QM/MM calculations are the primary current tool for modeling enzymatic reactions, the reliability of such calculations can be limited by the size of the QM region. Thus, we examine in this work the dependence of QM/MM calculations on the size of the QM region, using the reaction of catechol-O-methyl transferase (COMT) as a test case. Our study focuses on the effect of adding residues to the QM region on the activation free energy, obtained with extensive QM/MM sampling. It is found that the sensitivity of the activation barrier to the size of the QM is rather limited, while the dependence of the reaction free energy is somewhat larger. Of course, the results depend on the inclusion of the first solvation shell in the QM regions. For example, the inclusion of the Mg²⁺ ion can change the activation barrier due to charge transfer effects. However, such effects can easily be included in semiempirical approaches by proper parametrization. Overall, we establish that QM/MM calculations of activation barriers of enzymatic reactions are not highly sensitive to the size of the QM region, beyond the immediate region that describes the reacting atoms.

Torque, chemistry and efficiency in molecular motors: a study of the rotary-chemical coupling in F1-ATPase

Detailed understanding of the action of biological molecular machines must overcome the challenge of gaining a clear knowledge of the corresponding free-energy landscape. An example for this is the elucidation of the nature of converting chemical energy to torque and work in the rotary molecular motor of F₁-ATPase. A major part of the challenge involves understanding the rotary–chemical coupling from a non-phenomenological structure/energy description. Here we focused on using a coarse-grained model of F₁-ATPase to generate a structure-based free-energy landscape of the rotary–chemical process of the whole system. In particular, we concentrated on exploring the possible impact of the position of the catalytic dwell on the efficiency and torque generation of the molecular machine. It was found that the experimentally observed torque can be reproduced with landscapes that have different positions for the catalytic dwell on the rotary–chemical surface. Thus, although the catalysis is undeniably required for torque generation, the experimentally observed position of the catalytic dwell at 80° might not have a clear advantage for the force generation by F₁-ATPase. This further implies that the rotary–chemical couplings in these biological motors are quite robust and their efficiencies do not depend explicitly on the position of the catalytic dwells. Rather, the specific positioning of the dwells with respect to the rotational angle is a characteristic arising due to the structural construct of the molecular machine and might not bear any clear connection to the thermodynamic efficiency for the system.

The control of the discrimination between dNTP and rNTP in DNA and RNA polymerase

Understanding the origin of discrimination between rNTP and dNTP by DNA/RNA polymerases is important both for gaining fundamental knowledge on the corresponding systems and for advancing the design of specific drugs. This work explores the nature of this discrimination by systematic calculations of the transition state (TS) binding energy in RB69 DNA polymerase (gp43) and T7 RNA polymerase. The calculations reproduce the observed trend, in particular when they included the water contribution obtained by the water flooding approach. Our detailed study confirms the idea that the discrimination is due to the steric interaction between the 2'OH and Tyr416 in DNA polymerase, while the electrostatic interaction is the source of the discrimination in RNA polymerase. Proteins 2016; 84:1616-1624. © 2016 Wiley Periodicals, Inc.

Insights into Cotranslational Nascent Protein Behavior from Computer Simulations

Regulation of protein stability and function in vivo begins during protein synthesis, when the ribosome translates a messenger RNA into a nascent polypeptide. Cotranslational processes involving a nascent protein include folding, binding to other macromolecules, enzymatic modification, and secretion through membranes. Experiments have shown that the rate at which the ribosome adds amino acids to the elongating nascent chain influences the efficiency of these processes, with alterations to these rates possibly contributing to diseases, including some types of cancer. In this review, we discuss recent insights into cotranslational processes gained from molecular simulations, how different computational approaches have been combined to understand cotranslational processes at multiple scales, and the new scenarios illuminated by these simulations. We conclude by suggesting interesting questions that computational approaches in this research area can address over the next few years.

Principles and Overview of Sampling Methods for Modeling Macromolecular Structure and Dynamics

Investigation of macromolecular structure and dynamics is fundamental to understanding how macromolecules carry out their functions in the cell. Significant advances have been made toward this end in silico, with a growing number of computational methods proposed yearly to study and simulate various aspects of macromolecular structure and dynamics. This review aims to provide an overview of recent advances, focusing primarily on methods proposed for exploring the structure space of macromolecules in isolation and in assemblies for the purpose of characterizing equilibrium structure and dynamics. In addition to surveying recent applications that showcase current capabilities of computational methods, this review highlights state-of-the-art algorithmic techniques proposed to overcome challenges posed in silico by the disparate spatial and time scales accessed by dynamic macromolecules. This review is not meant to be exhaustive, as such an endeavor is impossible, but rather aims to balance breadth and depth of strategies for modeling macromolecular structure and dynamics for a broad audience of novices and experts.

A Microscopic Capacitor Model of Voltage Coupling in Membrane Proteins: Gating Charge Fluctuations in Ci-VSD

The voltage sensitivity of membrane proteins is reflected in the response of the voltage sensing domains (VSDs) to the changes in membrane potential. This response is implicated in the displacement of positively charged residues, associated with the conformational changes of VSDs. The displaced charges generate nonlinear (i.e., voltage-dependent) capacitance current called the gating current (and its corresponding gating charge), which is a key experimental quantity that characterizes voltage activation in VSMP. However, the relevant theoretical/computational approaches, aimed to correlate the structural information on VSMP to electrophysiological measurements, have been rather limited, posing a broad challenge in computer simulations of VSMP. Concomitant with the development of our coarse-graining (CG) model of voltage coupling, we apply our theoretical framework for the treatments of voltage effects in membrane proteins to modeling the VSMP activation, taking the VSDs (Ci-VSD) derived from the Ciona intestinalis voltage sensitive phosphatase (Ci-VSP) as a model system. Our CG model reproduces the observed gating charge of Ci-VSD activation in several different perspectives. In particular, a new closed-form expression of the gating charge is evaluated in both nonequilibrium and equilibrium ways, while considering the fluctuation-dissipation relation that connects a nonequilibrium measurement of the gating charge to an equilibrium measurement of charge fluctuations (i.e., the voltage-independent linear component of membrane capacitance). In turn, the expression uncovers a novel link that connects an equilibrium measurement of the voltage-independent linear capacitance to a nonequilibrium measurement of the voltage-dependent nonlinear capacitance (whose integral over voltage is equal to the gating charge). In addition, our CG model yields capacitor-like voltage dependent free energy parabolas, resulting in the free energy difference and the free energy barrier for the Ci-VSD activation at "zero" (depolarization) membrane potential. Significantly, the resultant voltage dependent energetics enables a direct evaluation of capacitance-voltage relationship (C-V curve) as well as charge-voltage relationship (Q-V curve) that is in a good agreement with the observed measurement of Ci-VSD voltage activation. Importantly, an extension of our kinetic/thermodynamic model of voltage dependent activation in VSMP allows for novel derivations of voltage-dependent rate constants, whose parameters are expressed by the intrinsic properties of VSMP. These novel closed-form expressions offer a physicochemical foundation for the semiempirical Eyring-type voltage dependent rate equations that have been the cornerstone for the phenomenological (kinetic) descriptions of gating and membrane currents in the mechanistic study of ion channels and transporters. Our extended theoretical framework developed in the present study has potential implications on the roles played by gating charge fluctuations for the spike generations in nerve cells within the framework of the Hodgkin-Huxley-type model.

Force Field Development for Lipid Membrane Simulations

With the rapid development of computer power and wide availability of modelling software computer simulations of realistic models of lipid membranes, including their interactions with various molecular species, polypeptides and membrane proteins have become feasible for many research groups. The crucial issue of the reliability of such simulations is the quality of the force field, and many efforts, especially in the latest several years, have been devoted to parametrization and optimization of the force fields for biomembrane modelling. In this review, we give account of the recent development in this area, covering different classes of force fields, principles of the force field parametrization, comparison of the force fields, and their experimental validation. This article is part of a Special Issue entitled: Biosimulations edited by Ilpo Vattulainen and Tomasz Róg.

Refining the treatment of membrane proteins by coarse-grained models

Obtaining a quantitative description of the membrane proteins stability is crucial for understanding many biological processes. However the advance in this direction has remained a major challenge for both experimental studies and molecular modeling. One of the possible directions is the use of coarse-grained models but such models must be carefully calibrated and validated. Here we use a recent progress in benchmark studies on the energetics of amino acid residue and peptide membrane insertion and membrane protein stability in refining our previously developed coarse-grained model (Vicatos et al., Proteins 2014;82:1168). Our refined model parameters were fitted and/or tested to reproduce water/membrane partitioning energetics of amino acid side chains and a couple of model peptides. This new model provides a reasonable agreement with experiment for absolute folding free energies of several β-barrel membrane proteins as well as effects of point mutations on a relative stability for one of those proteins, OmpLA. The consideration and ranking of different rotameric states for a mutated residue was found to be essential to achieve satisfactory agreement with the reference data.

The ribosome's energy landscape: Recent insights from computation

The ever-increasing capacity of computing resources has extended ribosome calculations from the study of small-scale fluctuations to large-scale barrier-crossing processes. As the field of computational/theoretical biophysics shifts focus to large-scale conformational transitions, there is a growing need for a systematic framework to interpret and analyze ribosome dynamics. To this end, energy landscape principles, largely developed for the study of biomolecular folding, have proven to be invaluable. These tools not only provide a foundation for describing simulations but can be used to reconcile experimental results, as well. In this review, I will discuss recent efforts to employ computational methods to reveal the characteristics of the ribosome's landscape and how these studies can help guide a new generation of experiments that more closely probe the underlying energetics. As a result of these investigations, general principles about ribosome function are beginning to emerge, including that: (1) small-scale fluctuations are the result of structure, rather than detailed energetics, (2) molecular flexibility leads to entropically favored rearrangements, and (3) tRNA dynamics may be accurately described as diffusive movement across an energy landscape.

Equilibrium fluctuation relations for voltage coupling in membrane proteins

A general theoretical framework is developed to account for the effects of an external potential on the energetics of membrane proteins. The framework is based on the free energy relation between two (forward/backward) probability densities, which was recently generalized to non-equilibrium processes, culminating in the work-fluctuation theorem. Starting from the probability densities of the conformational states along the "voltage coupling" reaction coordinate, we investigate several interconnected free energy relations between these two conformational states, considering voltage activation of ion channels. The free energy difference between the two conformational states at zero (depolarization) membrane potential (i.e., known as the chemical component of free energy change in ion channels) is shown to be equivalent to the free energy difference between the two "equilibrium" (resting and activated) conformational states along the one-dimensional voltage couplin reaction coordinate. Furthermore, the requirement that the application of linear response approximation to the free energy functionals of voltage coupling should satisfy the general free energy relations, yields a novel closed-form expression for the gating charge in terms of other basic properties of ion channels. This connection is familiar in statistical mechanics, known as the equilibrium fluctuation-response relation. The theory is illustrated by considering the coupling of a unit charge to the external voltage in the two sites near the surface of membrane, representing the activated and resting states. This is done using a coarse-graining (CG) model of membrane proteins, which includes the membrane, the electrolytes and the electrodes. The CG model yields Marcus-type voltage dependent free energy parabolas for the response of the electrostatic environment (electrolytes etc.) to the transition from the initial to the final configuratinal states, leading to equilibrium free energy difference and free energy barrier that follow the trend of the equilibrium fluctuation relation and the Marcus theory of electron transfer. These energetics also allow for a direct estimation of the voltage dependence of channel activation (Q-V curve), offering a quantitative rationale for a correlation between the voltage dependence parabolas and the Q-V curve, upon site-directed mutagenesis or drug binding. Taken together, by introducing the voltage coupling as the energy gap reaction coordinate, our framework brings new perspectives to the thermodynamic models of voltage activation in voltage-sensitive membrane proteins, offering an a framework for a better understating of the structure-function correlations of voltage gating in ion channels as well as electrogenic phenomena in ion pumps and transporters. Significantly, this formulation also provides a powerful bridge between the CG model of voltage coupling and the conventional macroscopic treatments.

Dissecting the role of the γ-subunit in the rotary-chemical coupling and torque generation of F1-ATPase

Unraveling the molecular nature of the conversion of chemical energy (ATP hydrolysis in the α/β-subunits) to mechanical energy and torque (rotation of the γ-subunit) in F1-ATPase is very challenging. A major part of the challenge involves understanding the rotary-chemical coupling by a nonphenomenological structure-energy description, while accounting for the observed torque generated on the γ-subunit and its change due to mutation of this unit. Here we extend our previous study that used a coarse-grained model of the F1-ATPase to generate a structure-based free energy landscape of the rotary-chemical process. Our quantitative analysis of the landscape reproduced the observed torque for the wild-type enzyme. In doing so, we found that there are several possibilities of torque generation from landscapes with various shapes and demonstrated that a downhill slope along the chemical coordinate could still result in negligible torque, due to ineffective coupling of the chemistry to the γ-subunit rotation. We then explored the relationship between the functionality and the underlying sequence through systematic examination of the effect of various parts of the γ-subunit on free energy surfaces of F1-ATPase. Furthermore, by constructing several types of γ-deletion systems and calculating the corresponding torque generation, we gained previously unknown insights into the molecular nature of the F1-ATPase rotary motor. Significantly, our results are in excellent agreement with recent experimental findings and indicate that the rotary-chemical coupling is primarily established through electrostatic effects, although specific contacts through γ-ionizable residue side chains are not essential for establishing the basic features of the coupling.