Valid molecular dynamics simulations of human hemoglobin require a surprisingly large box size

A single mutation in dairy cow-associated H5N1 viruses increases receptor binding breadth

Clade 2.3.4.4b H5N1 is causing an unprecedented outbreak in dairy cows in the United States. To understand if recent H5N1 viruses are changing their receptor use, we screened recombinant hemagglutinin (HA) from historical and recent 2.3.4.4b H5N1 viruses for binding to distinct glycans bearing terminal sialic acids using a glycan microarray. We find that H5 from A/Texas/37/2024, an isolate from the dairy cow outbreak, has increased binding breadth to core glycans bearing terminal α2,3 sialic acids, the avian receptor, compared to historical and recent 2.3.4.4b H5N1 viruses. We do not observe any binding to α2,6 sialic acids, the receptor used by human seasonal influenza viruses. Using molecular dynamics and a cryo-EM structure of A/Texas/37/2024 H5, we show A/Texas/37/2024 H5 is more flexible within the receptor-binding site compared to a 2.3.4.4b H5 from 2022. We identify a single mutation outside of the receptor binding site, T199I, is responsible for increased binding breadth, as it increases receptor binding site flexibility. Together, these data show recent H5N1 viruses are evolving increased receptor binding breadth which could impact the host range and cell types infected with H5N1.

Analysis of the Dynamics of a Complex, Multipathway Reaction: Insulin Dimer Dissociation

The protein hormone insulin forms a homodimer that must dissociate to bind to its receptor. Understanding the kinetics and mechanism of dissociation is essential for the rational design of therapeutic analogs. In addition to its physiological importance, this dissociation process serves as a paradigm for coupled (un)folding and (un)binding. Based on previous free energy simulations, insulin dissociation is thought to involve multiple pathways with comparable free energy barriers. Here, we analyze the mechanism of insulin dimer dissociation using a recently developed computational framework for estimating kinetic statistics from short-trajectory data. These statistics indicate that the likelihood of dissociation (the committor) closely tracks the decrease in the number of (native and nonnative) intermonomer contacts and the increase in the number of water contacts at the dimer interface; the transition state with equal likelihood of association and dissociation corresponds to an encounter complex with relatively few native contacts and many nonnative contacts. We identify four pathways out of the dimer state and quantify their contributions to the rate as well as their exchange by computing reactive fluxes. We show that both the pathways and their extents of exchange can be understood in terms of rotations around three axes of the dimer structure. Our results provide insights into the kinetics of insulin analogs and, more generally, how to characterize complex, multipathway processes.

Structural serology of polyclonal antibody responses to mRNA-1273 and NVX-CoV2373 COVID-19 vaccines

Current COVID-19 vaccines are largely limited in their ability to induce broad, durable immunity against emerging viral variants. Design and development of improved vaccines utilizing existing platforms requires an in-depth understanding of the antigenic and immunogenic properties of available vaccines. Here we examined the antigenicity of two of the original COVID-19 vaccines, mRNA-1273 and NVX-CoV2373, by electron microscopy-based polyclonal epitope mapping (EMPEM) of serum from immunized non-human primates (NHPs) and clinical trial donors. Both vaccines induce diverse polyclonal antibody (pAb) responses to the N-terminal domain (NTD) in addition to the receptor-binding domain (RBD) of the Spike protein, with the NTD supersite being an immunodominant epitope. High-resolution cryo-EMPEM studies revealed extensive pAb responses to and around the supersite with unique angles of approach and engagement. NTD supersite pAbs were also the most susceptible to variant mutations compared to other specificities, indicating that ongoing Spike ectodomain-based vaccine design strategies should consider immuno-masking this site to prevent induction of these strain-specific responses.

CDR L3 Loop Rearrangement Switches Multispecific SPE-7 IgE Antibody From Hapten to Protein Binding

The monoclonal IgE antibody SPE-7 was originally raised against a 2,4-dinitrophenyl (DNP) target. Through its ability to adopt multiple conformations, the antibody is capable of binding to a diverse range of small haptens and large proteins. The present study examines a dataset of experimentally determined crystal structures of the SPE-7 antibody to gain insight into the mechanisms that contribute to its multispecificity. With the emergence of more and more therapeutic antibodies against a huge repertoire of different targets, our research could be of great interest for future drug development. We are able to discriminate between the different paratope-binding states in the conformational ensembles obtained by enhanced sampling molecular dynamics simulations, and to calculate their transition timescales and state probabilities. Furthermore, we describe the key residues responsible for discriminating between the different binding capacities and identify a tryptophan in a central position of the CDR L3 loop as the residue of greatest interest. The overall dynamics of the paratope appear to be mainly influenced by the CDR L3 and CDR L1 loops.

Fit-for-purpose heterodivalent single-domain antibody for gastrointestinal targeting of toxin B from Clostridium difficile

Single-domain antibodies (sdAbs), such as VHHs, are increasingly being developed for gastrointestinal (GI) applications against pathogens to strengthen gut health. However, what constitutes a suitable developability profile for applying these proteins in a gastrointestinal setting remains poorly explored. Here, we describe an in vitro methodology for the identification of sdAb derivatives, more specifically divalent VHH constructs, that display extraordinary developability properties for oral delivery and functionality in the GI environment. We showcase this by developing a heterodivalent VHH construct that cross-inhibits the toxic activity of the glycosyltransferase domains (GTDs) from three different toxinotypes of cytotoxin B (TcdB) from lineages of Clostridium difficile. We show that the VHH construct possesses high stability and binding activity under gastric conditions, in the presence of bile salts, and at high temperatures. We suggest that the incorporation of early developability assessment could significantly aid in the efficient discovery of VHHs and related constructs fit for oral delivery and GI applications.

Optimal Molecular Dynamics System Size for Increased Precision and Efficiency for Epoxy Materials

Molecular dynamics (MD) simulation is an important tool for predicting thermo-mechanical properties of polymer resins at the nanometer length scale, which is particularly important for efficient computationally driven design of advanced composite materials and structures. Because of the statistical nature of modeling amorphous materials on the nanometer length scale, multiple MD models (replicates) are typically built and simulated for statistical sampling of predicted properties. Larger replicates generally provide higher precision in the predictions but result in higher simulation times. Unfortunately, there is insufficient information in the literature to establish guidelines between MD model size and the resulting precision in predicted thermo-mechanical properties. The objective of this study was to determine the optimal MD model size of epoxy resin to balance efficiency and precision. The results show that an MD model size of 15,000 atoms provides for the fastest simulations without sacrificing precision in the prediction of mass density, elastic properties, strength, and thermal properties of epoxy. The results of this study are important for efficient computational process modeling and integrated computational materials engineering (ICME) for the design of next-generation composite materials for demanding applications.

Identification of histidine residues that affect the T/R-state conformations of human hemoglobin using constant pH molecular dynamics simulations

Human hemoglobin (Hb) is a tetrameric protein consisting of two α and two β subunits that can adopt a low-affinity T- and high-affinity R-state conformations. Under physiological pH conditions, histidine (His) residues are the main sites for proton binding or release, and their protonation states can affect the T/R-state conformation of Hb. However, it remains unclear which His residues can effectively affect the Hb conformation. Herein, the impact of the 38 His residues of Hb on its T/R-state conformations was evaluated using constant-pH molecular dynamics (CpHMD) simulations at physiological pH while focusing on the His protonation states. Overall, the protonation states of some His residues were found to be correlated with the Hb conformation state. These residues were mainly located in the proximity of the heme (α87 and β92), and at the α1β2 and α2β1 interfaces (α89 and β97). This correlation may be partly explained by how easily hydrogen bonds can be formed, which depends on the protonation states of the His residues. Taken together, these CpHMD-based findings provide new insights into the identification of titratable His residues α87, α89, β92, and β97 that can affect Hb conformational switching under physiological pH conditions.

Defining the cross-reactivity between peanut allergens Ara h 2 and Ara h 6 using monoclonal antibodies

In peanut allergy, Arachis hypogaea 2 (Ara h 2) and Arachis hypogaea 6 (Ara h 6) are two clinically relevant peanut allergens with known structural and sequence homology and demonstrated cross-reactivity. We have previously utilized X-ray crystallography and epitope binning to define the epitopes on Ara h 2. We aimed to quantitatively characterize the cross-reactivity between Ara h 2 and Ara h 6 on a molecular level using human monoclonal antibodies (mAbs) and structural characterization of allergenic epitopes. We utilized mAbs cloned from Ara h 2 positive single B cells isolated from peanut-allergic, oral immunotherapy-treated patients to quantitatively analyze cross-reactivity between recombinant Ara h 2 (rAra h 2) and Ara h 6 (rAra h 6) proteins using biolayer interferometry and indirect inhibitory ELISA. Molecular dynamics simulations assessed time-dependent motions and interactions in the antibody-antigen complexes. Three epitopes-conformational epitopes 1.1 and 3, and the sequential epitope KRELRNL/KRELMNL-are conserved between Ara h 2 and Ara h 6, while two more conformational and three sequential epitopes are not. Overall, mAb affinity was significantly lower to rAra h 6 than it was to rAra h 2. This difference in affinity was primarily due to increased dissociation of the antibodies from rAra h 6, a phenomenon explained by the higher conformational flexibility of the Ara h 6-antibody complexes in comparison to Ara h 2-antibody complexes. Our results further elucidate the cross-reactivity of peanut 2S albumins on a molecular level and support the clinical immunodominance of Ara h 2.

The Role of Force Fields and Water Models in Protein Folding and Unfolding Dynamics

Protein folding is a fascinating, not fully understood phenomenon in biology. Molecular dynamics (MD) simulations are an invaluable tool to study conformational changes in atomistic detail, including folding and unfolding processes of proteins. However, the accuracy of the conformational ensembles derived from MD simulations inevitably relies on the quality of the underlying force field in combination with the respective water model. Here, we investigate protein folding, unfolding, and misfolding of fast-folding proteins by examining different force fields with their recommended water models, i.e., ff14SB with the TIP3P model and ff19SB with the OPC model. To this end, we generated long conventional MD simulations highlighting the perks and pitfalls of these setups. Using Markov state models, we defined kinetically independent conformational substates and emphasized their distinct characteristics, as well as their corresponding state probabilities. Surprisingly, we found substantial differences in thermodynamics and kinetics of protein folding, depending on the combination of the protein force field and water model, originating primarily from the different water models. These results emphasize the importance of carefully choosing the force field and the respective water model as they determine the accuracy of the observed dynamics of folding events. Thus, the findings support the hypothesis that the water model is at least equally important as the force field and hence needs to be considered in future studies investigating protein dynamics and folding in all areas of biophysics.

Evaluating sheep hemoglobins with MD simulations as an animal model for sickle cell disease

Sickle cell disease (SCD) affects millions worldwide, yet there are few therapeutic options. To develop effective treatments, preclinical models that recapitulate human physiology and SCD pathophysiology are needed. SCD arises from a single Glu-to-Val substitution at position 6 in the β subunit of hemoglobin (Hb), promoting Hb polymerization and subsequent disease. Sheep share important physiological and developmental characteristics with humans, including the same developmental pattern of fetal to adult Hb switching. Herein, we investigated whether introducing the SCD mutation into the sheep β-globin locus would recapitulate SCD's complex pathophysiology by generating high quality SWISS-MODEL sheep Hb structures and performing MD simulations of normal/sickle human (huHbA/huHbS) and sheep (shHbB/shHbS) Hb, establishing how accurately shHbS mimics huHbS behavior. shHbS, like huHbS, remained stable with low RMSD, while huHbA and shHbB had higher and fluctuating RMSD. shHbB and shHbS also behaved identically to huHbA and huHbS with respect to β2-Glu6 and β1-Asp73 (β1-Asn72 in sheep) solvent interactions. These data demonstrate that introducing the single SCD-causing Glu-to-Val substitution into sheep β-globin causes alterations consistent with the Hb polymerization that drives RBC sickling, supporting the development of a SCD sheep model to pave the way for alternative cures for this debilitating, globally impactful disease.

Phage display assisted discovery of a pH-dependent anti-α-cobratoxin antibody from a natural variable domain library

Recycling IgG antibodies bind to their target antigen at physiological pH in the blood stream and release them upon endocytosis when pH levels drop, allowing the IgG antibodies to be recycled into circulation via FcRn-mediated cellular pathways, while the antigens undergo lysosomal degradation. This enables recycling antibodies to achieve comparable therapeutic effect at lower doses than their non-recycling counterparts. The development of such antibodies is typically achieved by histidine doping of their variable regions or by performing in vitro antibody selection campaigns utilizing histidine doped libraries. Both are strategies that may introduce sequence liabilities. Here, we present a methodology that employs a naïve antibody phage display library, consisting of natural variable domains, to discover antibodies that bind α-cobratoxin from the venom of Naja kaouthia in a pH-dependent manner. As a result, an antibody was discovered that exhibits a 7-fold higher off-rate at pH 5.5 than pH 7.4 in bio-layer interferometry experiments. Interestingly, no histidine residues were found in its variable domains, and in addition, the antibody showed pH-dependent binding to a histidine-devoid antigen mutant. As such, the results demonstrate that pH-dependent antigen-antibody binding may not always be driven by histidine residues. By employing molecular dynamics simulations, different protonation states of titratable residues were found, which potentially could be responsible for the observed pH-dependent antigen binding properties of the antibody. Finally, given the typically high diversity of naïve antibody libraries, the methodology presented here can likely be applied to discover recycling antibodies against different targets ab initio without the need for histidine doping.

Structural stability of R-state conformation of carbonmonoxyl sickle and normal hemoglobin dimer

A mutation at the sixth residue, glutamic acid to valine, in beta chain of hemoglobin distorts the entire shape of hemoglobin into a sickle shape. The investigation of the binding mechanisms of different chains of hemoglobin under the mutated condition can give an understanding of the molecular distortion. In this work, we have studied the binding mechanism between two chains in the dimer structure of the R-state conformation of carbonmonoxyl sickle hemoglobin and is compared with that of normal hemoglobin by using molecular dynamics simulations. The binding strength between α-chain (PROA) and β-chain (PROB) in hemoglobin dimer has been analyzed by estimating hydrogen bonds, salt bridges, hydrophobic interactions and non-bonded interactions (electrostatics and van der Waals). The quantitative estimation of aforementioned interactions depicts that the structural stability of normal hemoglobin dimer is found to be greater than that of sickle one. The outcomes of such interactions are also supported by the estimated free energy between the chains in R-state conformation of the dimers. The difference of binding free energy, calculated by utilizing the umbrella sampling technique, is found to be ≈ (0.67 ± 0.06) kcal/mol.Communicated by Ramaswamy H. Sarma.

Structural basis of epitope selectivity and potent protection from malaria by PfCSP antibody L9

A primary objective in malaria vaccine design is the generation of high-quality antibody responses against the circumsporozoite protein of the malaria parasite, Plasmodium falciparum (PfCSP). To enable rational antigen design, we solved a cryo-EM structure of the highly potent anti-PfCSP antibody L9 in complex with recombinant PfCSP. We found that L9 Fab binds multivalently to the minor (NPNV) repeat domain, which is stabilized by a unique set of affinity-matured homotypic, antibody-antibody contacts. Molecular dynamics simulations revealed a critical role of the L9 light chain in integrity of the homotypic interface, which likely impacts PfCSP affinity and protective efficacy. These findings reveal the molecular mechanism of the unique NPNV selectivity of L9 and emphasize the importance of anti-homotypic affinity maturation in protective immunity against P. falciparum.

Structural Characterization of Nanobodies during Germline Maturation

Camelid heavy-chain antibody variable domains (VHH), nanobodies, are the smallest-known functional antibody fragments with high therapeutic potential. In this study, we investigate a VHH binding to hen egg-white lysozyme (HEL). We structurally and dynamically characterized the conformational diversity of four VHH variants to elucidate the antigen-binding process. For two of these antibodies, not only are the dissociation constants known, but also the experimentally determined crystal structures of the VHH in complex with HEL are available. We performed well-tempered metadynamics simulations in combination with molecular dynamics simulations to capture a broad conformational space and to reconstruct the thermodynamics and kinetics of conformational transitions in the antigen-binding site, the paratope. By kinetically characterizing the loop movements of the paratope, we found that, with an increase in affinity, the state populations shift towards the binding competent conformation. The contacts contributing to antigen binding, and those who contribute to the overall stability, show a clear trend towards less variable but more intense contacts. Additionally, these investigated nanobodies clearly follow the conformational selection paradigm, as the binding competent conformation pre-exists within the structural ensembles without the presence of the antigen.

Water dynamics around T0 vs R4 of hemoglobin from local hydrophobicity analysis

The local hydration around tetrameric hemoglobin (Hb) in its T₀ and R₄ conformational substates is analyzed based on molecular dynamics simulations. Analysis of the local hydrophobicity (LH) for all residues at the α₁β₂ and α₂β₁ interfaces, responsible for the quaternary T → R transition, which is encoded in the Monod-Wyman-Changeux model, as well as comparison with earlier computations of the solvent accessible surface area, makes clear that the two quantities measure different aspects of hydration. Local hydrophobicity quantifies the presence and structure of water molecules at the interface, whereas "buried surface" reports on the available space for solvent. For simulations with Hb frozen in its T₀ and R₄ states, the correlation coefficient between LH and buried surface is 0.36 and 0.44, respectively, but it increases considerably if the 95% confidence interval is used. The LH with Hb frozen and flexible changes little for most residues at the interfaces but is significantly altered for a few select ones: Thr41α, Tyr42α, Tyr140α, Trp37β, Glu101β (for T₀) and Thr38α, Tyr42α, Tyr140α (for R₄). The number of water molecules at the interface is found to increase by ∼25% for T₀ → R₄, which is consistent with earlier measurements. Since hydration is found to be essential to protein function, it is clear that hydration also plays an essential role in allostery.

The influence of antibody humanization on shark variable domain (VNAR) binding site ensembles

Sharks and other cartilaginous fish produce new antigen receptor (IgNAR) antibodies, as key part of their humoral immune response and are the phylogenetically oldest living organisms that possess an immunoglobulin (Ig)-based adaptive immune system. IgNAR antibodies are naturally occurring heavy-chain-only antibodies, that recognize antigens with their single domain variable regions (VNARs). In this study, we structurally and biophysically elucidate the effect of antibody humanization of a previously published spiny dogfish VNAR (parent E06), which binds with high affinity to the human serum albumin (HSA). We analyze different humanization variants together with the parental E06 VNAR and the human Vκ1 light chain germline DPK9 antibody to characterize the influence of point mutations in the framework and the antigen binding site on the specificity of VNARs as reported by Kovalenko et al. We find substantially higher flexibility in the humanized variants, reflected in a broader conformational space and a higher conformational entropy, as well as population shifts of the dominant binding site ensembles in solution. A further variant, in which some mutations are reverted, largely restores the conformational stability and the dominant binding minimum of the parent E06. We also identify differences in surface hydrophobicity between the human Vκ1 light chain germline DPK9 antibody, the parent VNAR E06 and the humanized variants. Additional simulations of VNAR-HSA complexes of the parent E06 VNAR and a humanized variant reveal that the parent VNAR features a substantially stronger network of stabilizing interactions. Thus, we conclude that a structural and dynamic understanding of the VNAR binding site upon humanization is a key aspect in antibody humanization.

Atomic-Level View of the Functional Transition in Vertebrate Hemoglobins: The Case of Antarctic Fish Hbs

Tetrameric hemoglobins (Hbs) are prototypal systems for studies aimed at unveiling basic structure-function relationships as well as investigating the molecular/structural basis of adaptation of living organisms to extreme conditions. However, a chronological analysis of decade-long studies conducted on Hbs is illuminating on the difficulties associated with the attempts of gaining functional insights from static structures. Here, we applied molecular dynamics (MD) simulations to explore the functional transition from the T to the R state of the hemoglobin of the Antarctic fish Trematomus bernacchii (HbTb). Our study clearly demonstrates the ability of the MD technique to accurately describe the transition of HbTb from the T to R-like states, as shown by a number of global and local structural indicators. A comparative analysis of the structural states that HbTb assumes in the simulations with those detected in previous MD analyses conducted on HbA (human Hb) highlights interesting analogies (similarity of the transition pathway) and differences (distinct population of intermediate states). In particular, the ability of HbTb to significantly populate intermediate states along the functional pathway explains the observed propensity of this protein to assume these structures in the crystalline state. It also explains some functional data reported on the protein that indicate the occurrence of other functional states in addition to the canonical R and T ones. These findings are in line with the emerging idea that the classical two-state view underlying tetrameric Hb functionality is probably an oversimplification and that other structural states play important roles in these proteins. The ability of MD simulations to accurately describe the functional pathway in tetrameric Hbs suggests that this approach may be effectively applied to unravel the molecular and structural basis of Hbs exhibiting peculiar functional properties as a consequence of the environmental adaptation of the host organism.

Recognition in the Domain of Molecular Chirality: From Noncovalent Interactions to Separation of Enantiomers

It is not a coincidence that both chirality and noncovalent interactions are ubiquitous in nature and synthetic molecular systems. Noncovalent interactivity between chiral molecules underlies enantioselective recognition as a fundamental phenomenon regulating life and human activities. Thus, noncovalent interactions represent the narrative thread of a fascinating story which goes across several disciplines of medical, chemical, physical, biological, and other natural sciences. This review has been conceived with the awareness that a modern attitude toward molecular chirality and its consequences needs to be founded on multidisciplinary approaches to disclose the molecular basis of essential enantioselective phenomena in the domain of chemical, physical, and life sciences. With the primary aim of discussing this topic in an integrated way, a comprehensive pool of rational and systematic multidisciplinary information is provided, which concerns the fundamentals of chirality, a description of noncovalent interactions, and their implications in enantioselective processes occurring in different contexts. A specific focus is devoted to enantioselection in chromatography and electromigration techniques because of their unique feature as "multistep" processes. A second motivation for writing this review is to make a clear statement about the state of the art, the tools we have at our disposal, and what is still missing to fully understand the mechanisms underlying enantioselective recognition.

Implementation of Telescoping Boxes in Adaptive Steered Molecular Dynamics

Long-time dynamical processes, such as those involving protein unfolding and ligand interactions, can be accelerated and realized through steered molecular dynamics (SMD). The challenge has been the extraction of information from such simulations that generalize for complex nonequilibrium processes. The use of Jarzynski's equality opened the possibility of determining the free energy along the steered coordinate, but sampling over the nonequilibrium trajectories is slow to converge. Adaptive steered molecular dynamics (ASMD) and other related techniques have been introduced to overcome this challenge through the use of stages. Here, we take advantage of these stages to address the numerical cost that arises from the required use of very large solvent boxes. We introduce telescoping box schemes within adaptive steered molecular dynamics (ASMD) in which we adjust the solvent box between stages and thereby vary (and optimize) the required number of solvent molecules. We have benchmarked the method on a relatively long α-helical peptide, Ala₃₀, with respect to the potential of mean force and hydrogen bonds. We show that the use of telescoping boxes introduces little numerical error while significantly reducing the computational cost.

Nanobody Paratope Ensembles in Solution Characterized by MD Simulations and NMR

Variable domains of camelid antibodies (so-called nanobodies or VHH) are the smallest antibody fragments that retain complete functionality and therapeutic potential. Understanding of the nanobody-binding interface has become a pre-requisite for rational antibody design and engineering. The nanobody-binding interface consists of up to three hypervariable loops, known as the CDR loops. Here, we structurally and dynamically characterize the conformational diversity of an anti-GFP-binding nanobody by using molecular dynamics simulations in combination with experimentally derived data from nuclear magnetic resonance (NMR) spectroscopy. The NMR data contain both structural and dynamic information resolved at various timescales, which allows an assessment of the quality of protein MD simulations. Thus, in this study, we compared the ensembles for the anti-GFP-binding nanobody obtained from MD simulations with results from NMR. We find excellent agreement of the NOE-derived distance maps obtained from NMR and MD simulations and observe similar conformational spaces for the simulations with and without NOE time-averaged restraints. We also compare the measured and calculated order parameters and find generally good agreement for the motions observed in the ps-ns timescale, in particular for the CDR3 loop. Understanding of the CDR3 loop dynamics is especially critical for nanobodies, as this loop is typically critical for antigen recognition.

Molecular simulations of proteins: From simplified physical interactions to complex biological phenomena

Molecular dynamics simulation is the most popular computational technique for investigating the structural and dynamical behaviour of proteins, in search of the molecular basis of their function. Far from being a completely settled field of research, simulations are still evolving to best capture the essential features of the atomic interactions that govern a protein's inner motions. Modern force fields are becoming increasingly accurate in providing a physical description adequate to this purpose, and allow us to model complex biological systems under fairly realistic conditions. Furthermore, the use of accelerated sampling techniques is improving our access to the observation of progressively larger molecular structures, longer time scales, and more hidden functional events. In this review, the basic principles of molecular dynamics simulations and a number of key applications in the area of protein science are summarized, and some of the most important results are discussed. Examples include the study of the structure, dynamics and binding properties of 'difficult' targets, such as intrinsically disordered proteins and membrane receptors, and the investigation of challenging phenomena like hydration-driven processes and protein aggregation. The findings described provide an overall picture of the current state of this research field, and indicate new perspectives on the road ahead to the upcoming future of molecular simulations.

Solvent Effects on the Menshutkin Reaction

The Menshutkin reaction is a methyl transfer reaction relevant in fields ranging from biochemistry to chemical synthesis. In the present work, the energetics and solvent distributions for NH₃+MeCl and Pyr+MeBr reactions were investigated in explicit solvent (water, methanol, acetonitrile, benzene, cyclohexane) by means of reactive molecular dynamics simulations. For polar solvents (water, methanol, and acetonitrile) and benzene, strong to moderate catalytic effects for both reactions were found, whereas apolar and bulky cyclohexane interacts weakly with the solute and does not show pronounced barrier reduction. The calculated barrier heights for the Pyr+MeBr reaction in acetonitrile and cyclohexane are 23.2 and 28.1 kcal/mol compared with experimentally measured barriers of 22.5 and 27.6 kcal/mol, respectively. The solvent distributions change considerably between reactant and TS but comparatively little between TS and product conformations of the solute. As the system approaches the transition state, correlated solvent motions occur which destabilize the solvent-solvent interactions. This is required for the system to surmount the barrier. Finally, it is found that the average solvent-solvent interaction energies in the reactant, TS, and product state geometries are correlated with changes in the solvent structure around the solute.

Cross-Correlated Motions in Azidolysozyme

The changes in the local and global dynamics of azide-labelled lysozyme compared with that of the wild type protein are quantitatively assessed for all alanine residues along the polypeptide chain. Although attaching -N3 to alanine residues has been considered to be a minimally invasive change in the protein it is found that depending on the location of the alanine residue, the local and global changes in the dynamics differ. For Ala92, the change in the cross-correlated motions are minimal, whereas attaching -N3 to Ala90 leads to pronounced differences in the local and global correlations as quantified by the cross-correlation coefficients of the Cα atoms. We also demonstrate that the spectral region of the asymmetric azide stretch distinguishes between alanine attachment sites, whereas changes in the low frequency, far-infrared region are less characteristic.

Comparing Antibody Interfaces to Inform Rational Design of New Antibody Formats

As the current biotherapeutic market is dominated by antibodies, the design of different antibody formats, like bispecific antibodies and other new formats, represent a key component in advancing antibody therapy. When designing new formats, a targeted modulation of pairing preferences is key. Several existing approaches are successful, but expanding the repertoire of design possibilities would be desirable. Cognate immunoglobulin G antibodies depend on homodimerization of the fragment crystallizable regions of two identical heavy chains. By modifying the dimeric interface of the third constant domain (CH3-CH3), with different mutations on each domain, the engineered Fc fragments form rather heterodimers than homodimers. The first constant domain (CH1-CL) shares a very similar fold and interdomain orientation with the CH3-CH3 dimer. Thus, numerous well-established design efforts for CH3-CH3 interfaces, have also been applied to CH1-CL dimers to reduce the number of mispairings in the Fabs. Given the high structural similarity of the CH3-CH3 and CH1-CL domains we want to identify additional opportunities in comparing the differences and overlapping interaction profiles. Our vision is to facilitate a toolkit that allows for the interchangeable usage of different design tools from crosslinking the knowledge between these two interface types. As a starting point, here, we use classical molecular dynamics simulations to identify differences of the CH3-CH3 and CH1-CL interfaces and already find unexpected features of these interfaces shedding new light on possible design variations. Apart from identifying clear differences between the similar CH3-CH3 and CH1-CL dimers, we structurally characterize the effects of point-mutations in the CH3-CH3 interface on the respective dynamics and interface interaction patterns. Thus, this study has broad implications in the field of antibody engineering as it provides a structural and mechanistical understanding of antibody interfaces and thereby presents a crucial aspect for the design of bispecific antibodies.

Dynamics of camel and human hemoglobin revealed by molecular simulations

Hemoglobin is one of the most widely studied proteins genetically, biochemically, and structurally. It is an oxygen carrying tetrameric protein that imparts the characteristic red color to blood. Each chain of hemoglobin harbors a heme group embedded in a hydrophobic pocket. Several studies have investigated structural variations present in mammalian hemoglobin and their functional implications. However, camel hemoglobin has not been thoroughly explored, especially from a structural perspective. Importantly, very little is known about how the heme group interacts with hemoglobin under varying conditions of osmolarity and temperature. Several experimental studies have indicated that the tense (T) state is more stable than the relaxed (R) state of hemoglobin under normal physiological conditions. Despite the fact that R state is less stable than the T state, no extensive structural dynamics studies have been performed to investigate global quaternary transitions of R state hemoglobin under normal physiological conditions. To evaluate this, several 500 ns all-atom molecular dynamics simulations were performed to get a deeper understanding of how camel hemoglobin behaves under stress, which it is normally exposed to, when compared to human hemoglobin. Notably, camel hemoglobin was more stable under physiological stress when compared to human hemoglobin. Additionally, when compared to camel hemoglobin, cofactor-binding regions of hemoglobin also exhibited more fluctuations in human hemoglobin under the conditions studied. Several differences were observed between the residues of camel and human hemoglobin that interacted with heme. Importantly, distal residues His58 of α hemoglobin and His63 of β hemoglobin formed more sustained interactions, especially at higher temperatures, in camel hemoglobin. These residues are important for oxygen binding to hemoglobin. Thus, this work provides insights into how camel and human hemoglobin differ in their interactions under stress.

The functional role of the hemoglobin-water interface

The interface between hemoglobin (Hb) and its environment, in particular water, is of great physiological relevance. Here, results from in vitro, in vivo, and computational experiments (molecular dynamics simulations) are summarized and put into perspective. One of the main findings from the computations is that the stability of the deoxy, ligand-free T-state (T₀) can be stabilized relative to the deoxy R-state (R₀) only in sufficiently large simulation boxes for the hydrophobic effect to manifest itself. This effect directly influences protein stability and is operative also under physiological conditions. Furthermore, molecular simulations provide a dynamical interpretation of the Perutz model for Hb function. Results from experiments using higher protein concentrations and realistic cellular environments are also discussed. One of the next great challenges for computational studies, which as we show is likely to be taken up in the near future, is to provide a molecular-level understanding of the dynamics of proteins in such crowded environments.

Molecular Dynamics of Hemoglobin Reveals Structural Alterations and Explains the Interactions Driving Sickle Cell Fibrillation

In sickle cell anemia, deoxyhemoglobin deforms RBCs by forming fibrils inside that disintegrate on oxygenation. We studied 100 ns long all-atom molecular dynamics (MD) for sickle and normal hemoglobin fibril models to understand this process, complemented by multiple 1 μs MD for a single tetramer of sickle and normal hemoglobin in deoxy and oxy states. We find that the presence of hydrophobic residues without a bulky side chain at β-6 in hemoglobin is the reason for the stability of the fibrils. Moreover, the free energy landscapes from MD of hemoglobin starting in the tensed (T) state capture the putative transition from T to relaxed (R) state, associated with oxygen binding. The three conformational wells in the landscapes are characterized by the quaternary changes where one αβ dimer rotates with respect to the other. The conformational changes from the oxygenation of sickle hemoglobin hinder the intermolecular contacts necessary for fibril formation.

Quaternary Structure Transitions of Human Hemoglobin: An Atomic-Level View of the Functional Intermediate States

Human hemoglobin (HbA) is one of the prototypal systems used to investigate structure-function relationships in proteins. Indeed, HbA has been used to develop the basic concepts of protein allostery, although the atomic-level mechanism underlying the HbA functionality is still highly debated. This is due to the fact that most of the three-dimensional structural information collected over the decades refers to the endpoints of HbA functional transition with little data available for the intermediate states. Here, we report molecular dynamics (MD) simulations by focusing on the relevance of the intermediate states of the protein functional transition unraveled by the crystallographic studies carried out on vertebrate Hbs. Fully atomistic simulations of the HbA T-state indicate that the protein undergoes a spontaneous transition toward the R-state. The inspection of the trajectory structures indicates that the protein significantly populates the intermediate HL-(C) state previously unraveled by crystallography. In the structural transition, it also assumes the intermediate states crystallographically detected in Antarctic fish Hbs. This finding suggests that HbA and Antarctic fish Hbs, in addition to the endpoints of the transitions, also share a similar deoxygenation pathway despite a distace of hundreds of millions of years in the evolution scale. Finally, using the essential dynamic sampling methodology, we gained some insights into the reverse R to T transition that is not spontaneously observed in classic MD simulations.

Ensembles Are Required to Handle Aleatoric and Parametric Uncertainty in Molecular Dynamics Simulation

Classical molecular dynamics is a computer simulation technique that is in widespread use across many areas of science, from physics and chemistry to materials, biology, and medicine. The method continues to attract criticism due its oft-reported lack of reproducibility which is in part due to a failure to submit it to reliable uncertainty quantification (UQ). Here we show that the uncertainty arises from a combination of (i) the input parameters and (ii) the intrinsic stochasticity of the method controlled by the random seeds. To illustrate the situation, we make a systematic UQ analysis of a widely used molecular dynamics code (NAMD), applied to estimate binding free energy of a ligand-bound to a protein. In particular, we replace the usually fixed input parameters with random variables, systematically distributed about their mean values, and study the resulting distribution of the simulation output. We also perform a sensitivity analysis, which reveals that, out of a total of 175 parameters, just six dominate the variance in the code output. Furthermore, we show that binding energy calculations dampen the input uncertainty, in the sense that the variation around the mean output free energy is less than the variation around the mean of the assumed input distributions, if the output is ensemble-averaged over the random seeds. Without such ensemble averaging, the predicted free energy is five times more uncertain. The distribution of the predicted properties is thus strongly dependent upon the random seed. Owing to this substantial uncertainty, robust statistical measures of uncertainty in molecular dynamics simulation require the use of ensembles in all contexts.

Allosteric regulation in CRISPR/Cas1-Cas2 protospacer acquisition mediated by DNA and Cas2

Cas1 and Cas2 are highly conserved proteins across clustered-regularly-interspaced-short-palindromic-repeat-Cas systems and play a significant role in protospacer acquisition. Based on crystal structure of twofold symmetric Cas1-Cas2 in complex with dual-forked protospacer DNA (psDNA), we conducted all-atom molecular dynamics simulations to study the psDNA binding, recognition, and response to cleavage on the protospacer-adjacent-motif complementary sequence, or PAMc, of Cas1-Cas2. In the simulation, we noticed that two active sites of Cas1 and Cas1’ bind asymmetrically to two identical PAMc on the psDNA captured from the crystal structure. For the modified psDNA containing only one PAMc, as that to be recognized by Cas1-Cas2 in general, our simulations show that the non-PAMc association site of Cas1-Cas2 remains destabilized until after the stably bound PAMc being cleaved at the corresponding association site. Thus, long-range correlation appears to exist upon the PAMc cleavage between the two active sites (∼10 nm apart) on Cas1-Cas2, which can be allosterically mediated by psDNA and Cas2 and Cas2’ in bridging. To substantiate such findings, we conducted repeated runs and further simulated Cas1-Cas2 in complex with synthesized psDNA sequences psL and psH, which have been measured with low and high frequency in acquisition, respectively. Notably, such intersite correlation becomes even more pronounced for the Cas1-Cas2 in complex with psH but remains low for the Cas1-Cas2 in complex with psL. Hence, our studies demonstrate that PAMc recognition and cleavage at one active site of Cas1-Cas2 may allosterically regulate non-PAMc association or even cleavage at the other site, and such regulation can be mediated by noncatalytic Cas2 and DNA protospacer to possibly support the ensued psDNA acquisition.

Pore mutation N617D in the skeletal muscle DHPR blocks Ca2+ influx due to atypical high-affinity Ca2+ binding

Skeletal muscle excitation-contraction (EC) coupling roots in Ca²⁺-influx-independent inter-channel signaling between the sarcolemmal dihydropyridine receptor (DHPR) and the ryanodine receptor (RyR1) in the sarcoplasmic reticulum. Although DHPR Ca²⁺ influx is irrelevant for EC coupling, its putative role in other muscle-physiological and developmental pathways was recently examined using two distinct genetically engineered mouse models carrying Ca²⁺ non-conducting DHPRs: DHPR(N617D) (Dayal et al., 2017) and DHPR(E1014K) (Lee et al., 2015). Surprisingly, despite complete block of DHPR Ca²⁺-conductance, histological, biochemical, and physiological results obtained from these two models were contradictory. Here, we characterize the permeability and selectivity properties and henceforth the mechanism of Ca²⁺ non-conductance of DHPR(N617). Our results reveal that only mutant DHPR(N617D) with atypical high-affinity Ca²⁺ pore-binding is tight for physiologically relevant monovalent cations like Na⁺ and K⁺. Consequently, we propose a molecular model of cooperativity between two ion selectivity rings formed by negatively charged residues in the DHPR pore region.

Uncertainty quantification in classical molecular dynamics

Molecular dynamics simulation is now a widespread approach for understanding complex systems on the atomistic scale. It finds applications from physics and chemistry to engineering, life and medical science. In the last decade, the approach has begun to advance from being a computer-based means of rationalizing experimental observations to producing apparently credible predictions for a number of real-world applications within industrial sectors such as advanced materials and drug discovery. However, key aspects concerning the reproducibility of the method have not kept pace with the speed of its uptake in the scientific community. Here, we present a discussion of uncertainty quantification for molecular dynamics simulation designed to endow the method with better error estimates that will enable it to be used to report actionable results. The approach adopted is a standard one in the field of uncertainty quantification, namely using ensemble methods, in which a sufficiently large number of replicas are run concurrently, from which reliable statistics can be extracted. Indeed, because molecular dynamics is intrinsically chaotic, the need to use ensemble methods is fundamental and holds regardless of the duration of the simulations performed. We discuss the approach and illustrate it in a range of applications from materials science to ligand-protein binding free energy estimation. This article is part of the theme issue 'Reliability and reproducibility in computational science: implementing verification, validation and uncertainty quantification in silico'.

Biomolecular Modeling and Simulation: A Prospering Multidisciplinary Field

We reassess progress in the field of biomolecular modeling and simulation, following up on our perspective published in 2011. By reviewing metrics for the field's productivity and providing examples of success, we underscore the productive phase of the field, whose short-term expectations were overestimated and long-term effects underestimated. Such successes include prediction of structures and mechanisms; generation of new insights into biomolecular activity; and thriving collaborations between modeling and experimentation, including experiments driven by modeling. We also discuss the impact of field exercises and web games on the field's progress. Overall, we note tremendous success by the biomolecular modeling community in utilization of computer power; improvement in force fields; and development and application of new algorithms, notably machine learning and artificial intelligence. The combined advances are enhancing the accuracy andscope of modeling and simulation, establishing an exemplary discipline where experiment and theory or simulations are full partners.

Effect of the Simulation Box Size and Precipitant Concentration on the Behavior of Tetragonal Lysozyme Dimer

Abstract

The 10-nanosecond simulation of a lysozyme dimer, which is a fragment of the tetragonal lysozyme crystal structure, has been carried out by the molecular dynamics method at different simulation box sizes and precipitant concentrations in a solution. The dimer stability has been estimated by calculating the root-mean-square fluctuations of protein atoms. It is shown that the box size does not significantly affect the mobility of protein atoms on a relatively short trajectory, while the effect of the precipitant concentration on this trajectory is noticeable.

Site-selective dynamics of azidolysozyme

The spectroscopic response of and structural dynamics around all azido-modified alanine residues (AlaN₃) in lysozyme are characterized. It is found that AlaN₃ is a positionally sensitive probe for the local dynamics, covering a frequency range of ∼15 cm^-1 for the center frequency of the line shape. This is consistent with findings from selective replacements of amino acids in PDZ2, which reported a frequency span of ∼10 cm^-1 for replacements of Val, Ala, or Glu by azidohomoalanine. For the frequency fluctuation correlation functions, the long-time decay constants τ₂ range from ∼1 to ∼10 ps, which compares with experimentally measured correlation times of 3 ps. Attaching azide to alanine residues can yield dynamics that decays to zero on the few ps time scale (i.e., static component Δ₀ ∼ 0 ps^-1) or to a remaining, static contribution of ∼0.5 ps^-1 (corresponding to 2.5 cm^-1), depending on the local environment on the 10 ps time scale. The magnitude of the static component correlates qualitatively with the degree of hydration of the spectroscopic probe. Although attaching azide to alanine residues is found to be structurally minimally invasive with respect to the overall protein structure, analysis of the local hydrophobicity indicates that the hydration around the modification site differs for modified and unmodified alanine residues, respectively.

Modulation of Nanoparticle Diffusion by Surface Ligand Length and Charge: Analysis with Molecular Dynamics Simulations

To help better interpret experimental measurement of nanoparticle size, it is important to understand how their diffusion depends on the physical and chemical features of surface ligands. In this study, explicit solvent molecular dynamics simulations are used to probe the effect of ligand charge and flexibility on the diffusion of small gold nanoparticles. The results suggest that despite a high bare charge (+18 e), cationic nanoparticles studied here have reduced diffusion constants compared to a hydrophobic gold nanoparticle by merely a modest amount. Increasing the ligand length by 10 CH₂ units also has a limited impact on the diffusion constant. For the three particles studied here, the difference between estimated hydrodynamic radius and radius of gyration is on the order of one solvent layer (3-5 Å), confirming that the significant discrepancies found in the size of similar nanoparticles by recent transmission electron microscopy and dynamic light scattering measurements were due to aggregation under solution conditions. The limited impact of electrostatic friction on the diffusion of highly charged nanoparticles is found to be due to the strong anticorrelation between electrostatic and van der Waals forces between nanoparticle and environment, supporting the generality of recent observation for proteins by Matyushov and co-workers. Including the first shell of solvent molecules as part of the diffusing particle has a minor impact on the total force autocorrelation function but reduces the disparity in relaxation time between the total force and its electrostatic and van der Waals components.

Shark Antibody Variable Domains Rigidify Upon Affinity Maturation-Understanding the Potential of Shark Immunoglobulins as Therapeutics

Sharks and other cartilaginous fish are the phylogenetically oldest living organisms that have antibodies as part of their adaptive immune system. As part of their humoral adaptive immune response, they produce an immunoglobulin, the so-called immunoglobulin new antigen receptor (IgNAR), a heavy-chain only antibody. The variable domain of an IgNAR, also known as V NAR , binds the antigen as an independent soluble domain. In this study, we structurally and dynamically characterized the affinity maturation mechanism of the germline and somatically matured (PBLA8) V NAR to better understand their function and their applicability as therapeutics. We observed a substantial rigidification upon affinity maturation, which is accompanied by a higher number of contacts, thereby contributing to the decrease in flexibility. Considering the static x-ray structures, the observed rigidification is not obvious, as especially the mutated residues undergo conformational changes during the simulation, resulting in an even stronger network of stabilizing interactions. Additionally, the simulations of the V NAR in complex with the hen egg-white lysozyme show that the V NAR antibodies evidently follow the concept of conformational selection, as the binding-competent state already preexisted even without the presence of the antigen. To have a more detailed description of antibody-antigen recognition, we also present here the binding/unbinding mechanism between the hen egg-white lysozyme and both the germline and matured V NAR s. Upon maturation, we observed a substantial increase in the resulting dissociation-free energy barrier. Furthermore, we were able to kinetically and thermodynamically describe the binding process and did not only identify a two-step binding mechanism, but we also found a strong population shift upon affinity maturation toward the native binding pose.

Rapid, accurate, precise and reproducible ligand-protein binding free energy prediction

A central quantity of interest in molecular biology and medicine is the free energy of binding of a molecule to a target biomacromolecule. Until recently, the accurate prediction of binding affinity had been widely regarded as out of reach of theoretical methods owing to the lack of reproducibility of the available methods, not to mention their complexity, computational cost and time-consuming procedures. The lack of reproducibility stems primarily from the chaotic nature of classical molecular dynamics (MD) and the associated extreme sensitivity of trajectories to their initial conditions. Here, we review computational approaches for both relative and absolute binding free energy calculations, and illustrate their application to a diverse set of ligands bound to a range of proteins with immediate relevance in a number of medical domains. We focus on ensemble-based methods which are essential in order to compute statistically robust results, including two we have recently developed, namely thermodynamic integration with enhanced sampling and enhanced sampling of MD with an approximation of continuum solvent. Together, these form a set of rapid, accurate, precise and reproducible free energy methods. They can be used in real-world problems such as hit-to-lead and lead optimization stages in drug discovery, and in personalized medicine. These applications show that individual binding affinities equipped with uncertainty quantification may be computed in a few hours on a massive scale given access to suitable high-end computing resources and workflow automation. A high level of accuracy can be achieved using these approaches.

Surprisingly Fast Interface and Elbow Angle Dynamics of Antigen-Binding Fragments

Fab consist of a heavy and light chain and can be subdivided into a variable (V H and V L ) and a constant region (C H 1 and C L ). The variable region contains the complementarity-determining region (CDR), which is formed by six hypervariable loops, shaping the antigen binding site, the paratope. Apart from the CDR loops, both the elbow angle and the relative interdomain orientations of the V H -V L and the C H 1-C L domains influence the shape of the paratope. Thus, characterization of the interface and elbow angle dynamics is essential to antigen specificity. We studied nine antigen-binding fragments (Fab) to investigate the influence of affinity maturation, antibody humanization, and different light-chain types on the interface and elbow angle dynamics. While the CDR loops reveal conformational transitions in the micro-to-millisecond timescale, both the interface and elbow angle dynamics occur on the low nanosecond timescale. Upon affinity maturation, we observe a substantial rigidification of the V H and V L interdomain and elbow-angle flexibility, reflected in a narrower and more distinct distribution. Antibody humanization describes the process of grafting non-human CDR loops onto a representative human framework. As the antibody framework changes upon humanization, we investigated if both the interface and the elbow angle distributions are changed or shifted. The results clearly showed a substantial shift in the relative V H -V L distributions upon antibody humanization, indicating that different frameworks favor distinct interface orientations. Additionally, the interface and elbow angle dynamics of five antibody fragments with different light-chain types are included, because of their strong differences in elbow angles. For these five examples, we clearly see a high variability and flexibility in both interface and elbow angle dynamics, highlighting the fact that Fab interface orientations and elbow angles interconvert between each other in the low nanosecond timescale. Understanding how the relative interdomain orientations and the elbow angle influence antigen specificity, affinity, and stability has broad implications in the field of antibody modeling and engineering.

Dynamics-Driven Allostery Underlies Ca2+-Mediated Release of SERCA Inhibition by Phospholamban

Sarcoplasmic reticulum (SR) Ca2+-ATPase (SERCA) and phospholamban (PLB) are essential for intracellular Ca2+ transport in myocytes. Ca2+-dependent activation of SERCA-PLB provides a control function that regulates cytosolic and SR Ca2+ levels. Although experimental and computational studies alone have led to a greater insight into SERCA-PLB regulation, the structural mechanisms for Ca2+ binding reversing inhibition of the complex remain poorly understood. Therefore, we have performed atomistic simulations totaling 32.7 μs and cell-based intramolecular fluorescence resonance energy transfer (FRET) experiments to determine structural changes of PLB-bound SERCA in response to binding of a single Ca2+ ion. Complementary MD simulations and FRET experiments showed that open-to-closed transitions in the structure of the headpiece underlie PLB inhibition of SERCA, and binding of a single Ca2+ ion is sufficient to shift the protein population toward a structurally closed structure of the complex. Closure is accompanied by functional interactions between the N-domain β5-β6 loop and the A-domain and the displacement of the catalytic phosphorylation domain toward a competent structure. We propose that reversal of SERCA-PLB inhibition is achieved by stringing together its controlling modules (A-domain and loop Nβ5-β6) with catalytic elements (P-domain) to regulate function during intracellular Ca2+ signaling. We conclude that binding of a single Ca2+ is a critical mediator of allosteric signaling that dictates structural changes and motions that relieve SERCA inhibition by PLB. Understanding allosteric regulation is of paramount importance to guide therapeutic modulation of SERCA and other evolutionarily related ion-motive ATPases.

Antibodies exhibit multiple paratope states influencing VH-VL domain orientations

In the last decades, antibodies have emerged as one of the most important and successful classes of biopharmaceuticals. The highest variability and diversity of an antibody is concentrated on six hypervariable loops, also known as complementarity determining regions (CDRs) shaping the antigen-binding site, the paratope. Whereas it was assumed that certain sequences can only adopt a limited set of backbone conformations, in this study we present a kinetic classification of several paratope states in solution. Using molecular dynamics simulations in combination with experimental structural information we capture the involved conformational transitions between different canonical clusters and additional dominant solution structures occurring in the micro-to-millisecond timescale. Furthermore, we observe a strong correlation of CDR loop movements. Another important aspect when characterizing different paratope states is the relative VH/VL orientation and the influence of the distinct CDR loop states on the VH/VL interface. Conformational rearrangements of the CDR loops do not only have an effect on the relative VH/VL orientations, but also influence in some cases the elbow-angle dynamics and shift the respective distributions. Thus, our results show that antibodies exist as several interconverting paratope states, each contributing to the antibody's properties.

On the importance of statistics in molecular simulations for thermodynamics, kinetics and simulation box size

Computational simulations, akin to wetlab experimentation, are subject to statistical fluctuations. Assessing the magnitude of these fluctuations, that is, assigning uncertainties to the computed results, is of critical importance to drawing statistically reliable conclusions. Here, we use a simulation box size as an independent variable, to demonstrate how crucial it is to gather sufficient amounts of data before drawing any conclusions about the potential thermodynamic and kinetic effects. In various systems, ranging from solvation free energies to protein conformational transition rates, we showcase how the proposed simulation box size effect disappears with increased sampling. This indicates that, if at all, the simulation box size only minimally affects both the thermodynamics and kinetics of the type of biomolecular systems presented in this work.

Water Dynamics Around Proteins: T- and R-States of Hemoglobin and Melittin

The water dynamics, as characterized by the local hydrophobicity (LH), is investigated for tetrameric hemoglobin (Hb) and dimeric melittin. For the T₀ to R₀ transition in Hb, it is found that LH provides additional molecular-level insight into the Perutz mechanism, i.e., the breaking and formation of salt bridges at the α₁/β₂ and α₂/β₁ interface is accompanied by changes in LH. For Hb in cubic water boxes with 90 and 120 Å edge length it is observed that following a decrease in LH as a consequence of reduced water density or change of water orientation at the protein/water interface the α/β interfaces are destabilized; this is a hallmark of the Perutz stereochemical model for the T to R transition in Hb. The present work thus provides a dynamical view of the classical structural model relevant to the molecular foundations of Hb function. For dimeric melittin, earlier results by Cheng and Rossky [ Nature 1998, 392, 696-699] are confirmed and interpreted on the basis of LH from simulations in which the protein structure is frozen. For the flexible melittin dimer, the changes in the local hydration can be as much as 30% greater than for the rigid dimer, reflecting the fact that protein and water dynamics are coupled.

Use of Molecular Dynamics Simulations in Structure-Based Drug Discovery

BACKGROUND

Traditional drug discovery is a lengthy process which involves a huge amount of resources. Modern-day drug discovers various multidisciplinary approaches amongst which, computational ligand and structure-based drug designing methods contribute significantly. Structure-based drug designing techniques require the knowledge of structural information of drug target and drug-target complexes. Proper understanding of drug-target binding requires the flexibility of both ligand and receptor to be incorporated. Molecular docking refers to the static picture of the drug-target complex(es). Molecular dynamics, on the other hand, introduces flexibility to understand the drug binding process.

OBJECTIVE

The aim of the present study is to provide a systematic review on the usage of molecular dynamics simulations to aid the process of structure-based drug design.

METHOD

This review discussed findings from various research articles and review papers on the use of molecular dynamics in drug discovery. All efforts highlight the practical grounds for which molecular dynamics simulations are used in drug designing program. In summary, various aspects of the use of molecular dynamics simulations that underline the basis of studying drug-target complexes were thoroughly explained.

RESULTS

This review is the result of reviewing more than a hundred papers. It summarizes various problems that use molecular dynamics simulations.

CONCLUSION

The findings of this review highlight how molecular dynamics simulations have been successfully implemented to study the structure-function details of specific drug-target complexes. It also identifies the key areas such as stability of drug-target complexes, ligand binding kinetics and identification of allosteric sites which have been elucidated using molecular dynamics simulations.

Common mechanism of thermostability in small α- and β-proteins studied by molecular dynamics

Protein thermostability is important to evolution, diseases, and industrial applications. Proteins use diverse molecular strategies to achieve stability at high temperature, yet reducing the entropy of unfolding seems required. We investigated five small α-proteins and five β-proteins with known, distinct structures and thermostability (T_m ) using multi-seed molecular dynamics simulations at 300, 350, and 400 K. The proteins displayed diverse changes in hydrogen bonding, solvent exposure, and secondary structure with no simple relationship to T_m . Our dynamics were in good agreement with experimental B-factors at 300 K and insensitive to force-field choice. Despite the very distinct structures, the native-state (300 + 350 K) free-energy landscapes (FELs) were significantly broader for the two most thermostable proteins and smallest for the three least stable proteins in both the α- and β-group and with both force fields studied independently (tailed t-test, 95% confidence level). Our results suggest that entropic ensembles stabilize proteins at high temperature due to reduced entropy of unfolding, viz., ΔG = ΔH - TΔS. Supporting this mechanism, the most thermostable proteins were also the least kinetically stable, consistent with broader FELs, typified by villin headpiece and confirmed by specific comparison to a mesophilic ortholog of Thermus thermophilus apo-pyrophosphate phosphohydrolase. We propose that molecular strategies of protein thermostabilization, although diverse, tend to converge toward highest possible entropy in the native state consistent with the functional requirements. We speculate that this tendency may explain why many proteins are not optimally structured and why molten-globule states resemble native proteins so much.

A remark on the efficiency of the double-system/single-box nonequilibrium approach in the SAMPL6 SAMPLing challenge

The alchemical nonequilibrium switching technique was one of several methods in the top tier of performance in the recent SAMPL6 SAMPLing challenge in both accuracy and efficiency. In this paper, in the context of nonequilibrium alchemical switching, we compare the efficiency of the double-system/single-box (DSSB) approach (used in the SAMPL6 challenges) to the standard single-system/double-box method (SSDB). Exploiting the Crooks theorem in a simple but effective test case, we analytically show that the DSSB approach is almost twice as efficient as SSDB for slow near-equilibrium switching but it gives basically no gain over the conventional SSDB approach when the variance of the work distribution exceeds few [Formula: see text], with the potential of producing artifacts and entanglements if not judiciously implemented.

A hallmark of phospholamban functional divergence is located in the N-terminal phosphorylation domain

Sarcoplasmic reticulum Ca²⁺ pump (SERCA) is a critical component of the Ca²⁺ transport machinery in myocytes. There is clear evidence for regulation of SERCA activity by PLB, whose activity is modulated by phosphorylation of its N-terminal domain (residues 1–25), but there is less clear evidence for the role of this domain in PLB’s functional divergence. It is widely accepted that only sarcolipin (SLN), a protein that shares substantial homology with PLB, uncouples SERCA Ca²⁺ transport from ATP hydrolysis by inducing a structural change of its energy-transduction domain; yet, experimental evidence shows that the transmembrane domain of PLB (residues 26–52, PLB_26–52) partially uncouples SERCA in vitro. These apparently conflicting mechanisms suggest that PLB’s uncoupling activity is encoded in its transmembrane domain, and that it is controlled by the N-terminal phosphorylation domain. To test this hypothesis, we performed molecular dynamics simulations (MDS) of the binary complex between PLB_26–52 and SERCA. Comparison between PLB_26–52 and wild-type PLB (PLB_WT) showed no significant changes in the stability and orientation of the transmembrane helix, indicating that PLB_26–52 forms a native-like complex with SERCA. MDS showed that PLB_26–52 produces key intermolecular contacts and structural changes required for inhibition, in agreement with studies showing that PLB_26–52 inhibits SERCA. However, deletion of the N-terminal phosphorylation domain facilitates an order-to-disorder shift in the energy-transduction domain associated with uncoupling of SERCA, albeit weaker than that induced by SLN. This mechanistic evidence reveals that the N-terminal phosphorylation domain of PLB is a primary contributor to the functional divergence among homologous SERCA regulators.

Structure and function of crocodilian hemoglobins and allosteric regulation by chloride, ATP, and CO2

Hemoglobins (Hbs) of crocodilians are reportedly characterized by unique mechanisms of allosteric regulatory control, but there are conflicting reports regarding the importance of different effectors, such as chloride ions, organic phosphates, and CO₂. Progress in understanding the unusual properties of crocodilian Hbs has also been hindered by a dearth of structural information. Here, we present the first comparative analysis of blood properties and Hb structure and function in a phylogenetically diverse set of crocodilian species. We examine mechanisms of allosteric regulation in the Hbs of 13 crocodilian species belonging to the families Crocodylidae and Alligatoridae. We also report new amino acid sequences for the α- and β-globins of these taxa, which, in combination with structural analyses, provide insights into molecular mechanisms of allosteric regulation. All crocodilian Hbs exhibited a remarkably strong sensitivity to CO₂, which would permit effective O₂ unloading to tissues in response to an increase in metabolism during intense activity and diving. Although the Hbs of all crocodilians exhibit similar intrinsic O₂-affinities, there is considerable variation in sensitivity to Cl^- ions and ATP, which appears to be at least partly attributable to variation in the extent of NH₂-terminal acetylation. Whereas chloride appears to be a potent allosteric effector of all crocodile Hbs, ATP has a strong, chloride-independent effect on Hb-O₂ affinity only in caimans. Modeling suggests that allosteric ATP binding has a somewhat different structural basis in crocodilian and mammalian Hbs.

Split the Charge Difference in Two! A Rule of Thumb for Adding Proper Amounts of Ions in MD Simulations

Despite the relevance of properly setting ionic concentrations in Molecular Dynamics (MD) simulations, methods or practical rules to set ionic strength are scarce and rarely documented. Based on a recently proposed thermodynamics method we provide an accurate rule of thumb to define the electrolytic content in simulation boxes. Extending the use of good practices in setting up MD systems is promptly needed to ensure reproducibility and consistency in molecular simulations.