Biomolecular simulations at the exascale: From drug design to organelles and beyond.

The rapid advancement in computational power available for research offers to bring not only quantitative improvements, but also qualitative changes in the field of biomolecular simulation. Here, we review the state of biomolecular dynamics simulations at the threshold to exascale resources becoming available. Both developments in parallel and distributed computing will be discussed, providing a perspective on the state of the art of both. A main focus will be on obtaining binding and conformational free energies, with an outlook to macromolecular complexes and (sub)cellular assemblies.

Selectivity filter mutations shift ion permeation mechanism in potassium channels.

Potassium (K+) channels combine high conductance with high ion selectivity. To explain this efficiency, two molecular mechanisms have been proposed. The "direct knock-on" mechanism is defined by water-free K+ permeation and formation of direct ion-ion contacts in the highly conserved selectivity filter (SF). The "soft knock-on" mechanism involves co-permeation of water and separation of K+ by water molecules. With the aim to distinguish between these mechanisms, crystal structures of the KcsA channel with mutations in two SF residues-G77 and T75-were published, where the arrangements of K+ ions and water display canonical soft knock-on configurations. These data were interpreted as evidence of the soft knock-on mechanism in wild-type channels. Here, we test this interpretation using molecular dynamics simulations of KcsA and its mutants. We show that while a strictly water-free direct knock-on permeation is observed in the wild type, conformational changes induced by these mutations lead to distinct ion permeation mechanisms, characterized by co-permeation of K+ and water. These mechanisms are characterized by reduced conductance and impaired potassium selectivity, supporting the importance of full dehydration of potassium ions for the hallmark high conductance and selectivity of K+ channels. In general, we present a case where mutations introduced at the critical points of the permeation pathway in an ion channel drastically change its permeation mechanism in a nonintuitive manner.

Current State of Open Source Force Fields in Protein-Ligand Binding Affinity Predictions.

In drug discovery, the in silico prediction of binding affinity is one of the major means to prioritize compounds for synthesis. Alchemical relative binding free energy (RBFE) calculations based on molecular dynamics (MD) simulations are nowadays a popular approach for the accurate affinity ranking of compounds. MD simulations rely on empirical force field parameters, which strongly influence the accuracy of the predicted affinities. Here, we evaluate the ability of six different small-molecule force fields to predict experimental protein-ligand binding affinities in RBFE calculations on a set of 598 ligands and 22 protein targets. The public force fields OpenFF Parsley and Sage, GAFF, and CGenFF show comparable accuracy, while OPLS3e is significantly more accurate. However, a consensus approach using Sage, GAFF, and CGenFF leads to accuracy comparable to OPLS3e. While Parsley and Sage are performing comparably based on aggregated statistics across the whole dataset, there are differences in terms of outliers. Analysis of the force field reveals that improved parameters lead to significant improvement in the accuracy of affinity predictions on subsets of the dataset involving those parameters. Lower accuracy can not only be attributed to the force field parameters but is also dependent on input preparation and sampling convergence of the calculations. Especially large perturbations and nonconverged simulations lead to less accurate predictions. The input structures, Gromacs force field files, as well as the analysis Python notebooks are available on GitHub.

N-Acyl-N-Alkyl Sulfonamide Probes for Ligand-Directed Covalent Labeling of GPCRs: The Adenosine A2B Receptor as Case Study.

Small molecular tool compounds play an essential role in the study of G protein-coupled receptors (GPCRs). However, tool compounds most often occupy the orthosteric binding site, hampering the study of GPCRs upon ligand binding. To overcome this problem, ligand-directed labeling techniques have been developed that leave a reporter group covalently bound to the GPCR, while allowing subsequent orthosteric ligands to bind. In this work, we applied such a labeling strategy to the adenosine A2B receptor (A2BAR). We have synthetically implemented the recently reported N-acyl-N-alkyl sulfonamide (NASA) warhead into a previously developed ligand and show that the binding of the A2BAR is not restricted by NASA incorporation. Furthermore, we have investigated ligand-directed labeling of the A2BAR using SDS-PAGE, flow cytometric, and mass spectrometry techniques. We have found one of the synthesized probes to specifically label the A2BAR, although detection was hindered by nonspecific protein labeling most likely due to the intrinsic reactivity of the NASA warhead. Altogether, this work aids the future development of ligand-directed probes for the detection of GPCRs.

Resolving coupled pH titrations using alchemical free energy calculations.

In a protein, nearby titratable sites can be coupled: the (de)protonation of one may affect the other. The degree of this interaction depends on several factors and can influence the measured p K a . Here, we derive a formalism based on double free energy differences ( Δ Δ G ) for quantifying the individual site p K a values of coupled residues. As Δ Δ G values can be obtained by means of alchemical free energy calculations, the presented approach allows for a convenient estimation of coupled residue p K a s in practice. We demonstrate that our approach and a previously proposed microscopic p K a formalism, can be combined with alchemical free energy calculations to resolve pH-dependent protein p K a values. Toy models and both, regular and constant-pH molecular dynamics simulations, alongside experimental data, are used to validate this approach. Our results highlight the insights gleaned when coupling and microstate probabilities are analyzed and suggest extensions to more complex enzymatic contexts. Furthermore, we find that naïvely computed p K a values that ignore coupling, can be significantly improved when coupling is accounted for, in some cases reducing the error by half. In short, alchemical free energy methods can resolve the p K a values of both uncoupled and coupled residues.

New rabies viral resources for multi-scale neural circuit mapping.

Comparisons and linkage between multiple imaging scales are essential for neural circuit connectomics. Here, we report 20 new recombinant rabies virus (RV) vectors that we have developed for multi-scale and multi-modal neural circuit mapping tools. Our new RV tools for mesoscale imaging express a range of improved fluorescent proteins. Further refinements target specific neuronal subcellular locations of interest. We demonstrate the discovery power of these new tools including the detection of detailed microstructural changes of rabies-labeled neurons in aging and Alzheimer's disease mouse models, live imaging of neuronal activities using calcium indicators, and automated measurement of infected neurons. RVs that encode GFP and ferritin as electron microscopy (EM) and fluorescence microscopy reporters are used for dual EM and mesoscale imaging. These new viral variants significantly expand the scale and power of rabies virus-mediated neural labeling and circuit mapping across multiple imaging scales in health and disease.

Molecular basis for human aquaporin inhibition.

Cancer invasion and metastasis are known to be potentiated by the expression of aquaporins (AQPs). Likewise, the expression levels of AQPs have been shown to be prognostic for survival in patients and have a role in tumor growth, edema, angiogenesis, and tumor cell migration. Thus, AQPs are key players in cancer biology and potential targets for drug development. Here, we present the single-particle cryo-EM structure of human AQP7 at 3.2-Å resolution in complex with the specific inhibitor compound Z433927330. The structure in combination with MD simulations shows that the inhibitor binds to the endofacial side of AQP7. In addition, cancer cells treated with Z433927330 show reduced proliferation. The data presented here serve as a framework for the development of AQP inhibitors.

Guidelines for Free-Energy Calculations Involving Charge Changes.

The Coulomb interactions in molecular simulations are inherently approximated due to the finite size of the molecular box sizes amenable to current-day compute power. Several methods exist for treating long-range electrostatic interactions, yet these approaches are subject to various finite-size-related artifacts. Lattice-sum methods are frequently used to approximate long-range interactions; however, these approaches also suffer from artifacts which become particularly pronounced for free-energy calculations that involve charge changes. The artifacts, however, also affect the sampling when plain simulations are performed, leading to a biased ensemble. Here, we investigate two previously described model systems to determine if artifacts continue to play a role when overall neutral boxes are considered, in the context of both free-energy calculations and sampling. We find that ensuring that no net-charge changes take place, while maintaining a neutral simulation box, may be sufficient provided that the simulation boxes are large enough. Addition of salt to the solution (when appropriate) can further alleviate the remaining artifacts in the sampling or the calculated free-energy differences. We provide practical guidelines to avoid finite-size artifacts.

Predicting adverse outcomes in adults with a community-acquired lower respiratory tract infection: a protocol for the development and validation of two prediction models for (i) all-cause hospitalisation and mortality and (ii) cardiovascular outcomes.

BACKGROUND: Community-acquired lower respiratory tract infections (LRTI) are common in primary care and patients at particular risk of adverse outcomes, e.g., hospitalisation and mortality, are challenging to identify. LRTIs are also linked to an increased incidence of cardiovascular diseases (CVD) following the initial infection, whereas concurrent CVD might negatively impact overall prognosis in LRTI patients. Accurate risk prediction of adverse outcomes in LRTI patients, while considering the interplay with CVD, can aid general practitioners (GP) in the clinical decision-making process, and may allow for early detection of deterioration. This paper therefore presents the design of the development and external validation of two models for predicting individual risk of all-cause hospitalisation or mortality (model 1) and short-term incidence of CVD (model 2) in adults presenting to primary care with LRTI. METHODS: Both models will be developed using linked routine electronic health records (EHR) data from Dutch primary and secondary care, and the mortality registry. Adults aged ≥ 40 years with a GP-diagnosis of LRTI between 2016 and 2019 are eligible for inclusion. Relevant patient demographics, medical history, medication use, presenting signs and symptoms, and vital and laboratory measurements will be considered as candidate predictors. Outcomes of interest include 30-day all-cause hospitalisation or mortality (model 1) and 90-day CVD (model 2). Multivariable elastic net regression techniques will be used for model development. During the modelling process, the incremental predictive value of CVD for hospitalisation or all-cause mortality (model 1) will also be assessed. The models will be validated through internal-external cross-validation and external validation in an equivalent cohort of primary care LRTI patients. DISCUSSION: Implementation of currently available prediction models for primary care LRTI patients is hampered by limited assessment of model performance. While considering the role of CVD in LRTI prognosis, we aim to develop and externally validate two models that predict clinically relevant outcomes to aid GPs in clinical decision-making. Challenges that we anticipate include the possibility of low event rates and common problems related to the use of EHR data, such as candidate predictor measurement and missingness, how best to retrieve information from free text fields, and potential misclassification of outcome events.

Accurately Predicting Protein pKa Values Using Nonequilibrium Alchemy.

The stability, solubility, and function of a protein depend on both its net charge and the protonation states of its individual residues. pKa is a measure of the tendency for a given residue to (de)protonate at a specific pH. Although pKa values can be resolved experimentally, theory and computation provide a compelling alternative. To this end, we assess the applicability of a nonequilibrium (NEQ) alchemical free energy method to the problem of pKa prediction. On a data set of 144 residues that span 13 proteins, we report an average unsigned error of 0.77 ± 0.09, 0.69 ± 0.09, and 0.52 ± 0.04 pK for aspartate, glutamate, and lysine, respectively. This is comparable to current state-of-the-art predictors and the accuracy recently reached using free energy perturbation methods (e.g., FEP+). Moreover, we demonstrate that our open-source, pmx-based approach can accurately resolve the pKa values of coupled residues and observe a substantial performance disparity associated with the lysine partial charges in Amber14SB/Amber99SB*-ILDN, for which an underused fix already exists.

Development of an Affinity-Based Probe to Profile Endogenous Human Adenosine A3 Receptor Expression.

The adenosine A3 receptor (A3AR) is a G protein-coupled receptor (GPCR) that exerts immunomodulatory effects in pathophysiological conditions such as inflammation and cancer. Thus far, studies toward the downstream effects of A3AR activation have yielded contradictory results, thereby motivating the need for further investigations. Various chemical and biological tools have been developed for this purpose, ranging from fluorescent ligands to antibodies. Nevertheless, these probes are limited by their reversible mode of binding, relatively large size, and often low specificity. Therefore, in this work, we have developed a clickable and covalent affinity-based probe (AfBP) to target the human A3AR. Herein, we show validation of the synthesized AfBP in radioligand displacement, SDS-PAGE, and confocal microscopy experiments as well as utilization of the AfBP for the detection of endogenous A3AR expression in flow cytometry experiments. Ultimately, this AfBP will aid future studies toward the expression and function of the A3AR in pathologies.

Structure and dynamics of cholesterol-mediated aquaporin-0 arrays and implications for lipid rafts.

Aquaporin-0 (AQP0) tetramers form square arrays in lens membranes through a yet unknown mechanism, but lens membranes are enriched in sphingomyelin and cholesterol. Here, we determined electron crystallographic structures of AQP0 in sphingomyelin/ cholesterol membranes and performed molecular dynamics (MD) simulations to establish that the observed cholesterol positions represent those seen around an isolated AQP0 tetramer and that the AQP0 tetramer largely defines the location and orientation of most of its associated cholesterol molecules. At a high concentration, cholesterol increases the hydrophobic thickness of the annular lipid shell around AQP0 tetramers, which may thus cluster to mitigate the resulting hydrophobic mismatch. Moreover, neighboring AQP0 tetramers sandwich a cholesterol deep in the center of the membrane. MD simulations show that the association of two AQP0 tetramers is necessary to maintain the deep cholesterol in its position and that the deep cholesterol increases the force required to laterally detach two AQP0 tetramers, not only due to protein-protein contacts but also due to increased lipid-protein complementarity. Since each tetramer interacts with four such 'glue' cholesterols, avidity effects may stabilize larger arrays. The principles proposed to drive AQP0 array formation could also underlie protein clustering in lipid rafts.

Interactions between selectivity filter and pore helix control filter gating in the MthK channel.

K+ channel activity can be limited by C-type inactivation, which is likely initiated in part by dissociation of K+ ions from the selectivity filter and modulated by the side chains that surround it. While crystallographic and computational studies have linked inactivation to a "collapsed" selectivity filter conformation in the KcsA channel, the structural basis for selectivity filter gating in other K+ channels is less clear. Here, we combined electrophysiological recordings with molecular dynamics simulations, to study selectivity filter gating in the model potassium channel MthK and its V55E mutant (analogous to KcsA E71) in the pore-helix. We found that MthK V55E has a lower open probability than the WT channel, due to decreased stability of the open state, as well as a lower unitary conductance. Simulations account for both of these variables on the atomistic scale, showing that ion permeation in V55E is altered by two distinct orientations of the E55 side chain. In the "vertical" orientation, in which E55 forms a hydrogen bond with D64 (as in KcsA WT channels), the filter displays reduced conductance compared to MthK WT. In contrast, in the "horizontal" orientation, K+ conductance is closer to that of MthK WT; although selectivity filter stability is lowered, resulting in more frequent inactivation. Surprisingly, inactivation in MthK WT and V55E is associated with a widening of the selectivity filter, unlike what is observed for KcsA and reminisces recent structures of inactivated channels, suggesting a conserved inactivation pathway across the potassium channel family.

Automated relative binding free energy calculations from SMILES to ΔΔG.

In drug discovery, computational methods are a key part of making informed design decisions and prioritising experiments. In particular, optimizing compound affinity is a central concern during the early stages of development. In the last 10 years, alchemical free energy (FE) calculations have transformed our ability to incorporate accurate in silico potency predictions in design decisions, and represent the 'gold standard' for augmenting experiment-driven drug discovery. However, relative FE calculations are complex to set up, require significant expert intervention to prepare the calculation and analyse the results or are provided only as closed-source software, not allowing for fine-grained control over the underlying settings. In this work, we introduce an end-to-end relative FE workflow based on the non-equilibrium switching approach that facilitates calculation of binding free energies starting from SMILES strings. The workflow is implemented using fully modular steps, allowing various components to be exchanged depending on licence availability. We further investigate the dependence of the calculated free energy accuracy on the initial ligand pose generated by various docking algorithms. We show that both commercial and open-source docking engines can be used to generate poses that lead to good correlation of free energies with experimental reference data.

Central cavity dehydration as a gating mechanism of potassium channels.

The hydrophobic gating model, in which ion permeation is inhibited by the hydrophobicity, rather than a physical occlusion of the nanopore, functions in various ion channels including potassium channels. Available research focused on the energy barriers for ion/water conduction due to the hydrophobicity, whereas how hydrophobic gating affects the function and structure of channels remains unclear. Here, we use potassium channels as examples and conduct molecular dynamics simulations to investigate this problem. Our simulations find channel activities (ion currents) highly correlated with cavity hydration level, implying insufficient hydration as a barrier for ion permeation. Enforced cavity dehydration successfully induces conformational transitions between known channel states, further implying cavity dewetting as a key step in the gating procedure of potassium channels utilizing different activation mechanisms. Our work reveals how the cavity dewetting is coupled to structural changes of potassium channels and how it affects channel activity. The conclusion may also apply to other ion channels.

Ion Conduction Mechanisms in Potassium Channels Revealed by Permeation Cycles.

Potassium channels are responsible for the selective yet efficient permeation of potassium ions across cell membranes. Despite many available high-resolution structures of potassium channels, those conformations inform only on static information on the ion permeation processes. Here, we use molecular dynamics simulations and Markov state models to obtain dynamical details of ion permeation. The permeation cycles, expressed in terms of selectivity filter occupancy and representing ion permeation events, are illustrated. We show that the direct knock-on permeation represents the dominant permeation mechanism over a wide range of potassium concentrations, temperatures, and membrane voltages for the pore of MthK. Direct knock-on is also observed in other potassium channels with a highly conserved selectivity filter, demonstrating the robustness of the permeation mechanism. Lastly, we investigate the charge strength dependence of permeation cycles. Our results shed light on the underlying permeation details, which are valuable in studying conduction mechanisms in potassium channels.

APOBEC3B drives PKR-mediated translation shutdown and protects stress granules in response to viral infection.

Double-stranded RNA produced during viral replication and transcription activates both protein kinase R (PKR) and ribonuclease L (RNase L), which limits viral gene expression and replication through host shutoff of translation. In this study, we find that APOBEC3B forms a complex with PABPC1 to stimulate PKR and counterbalances the PKR-suppressing activity of ADAR1 in response to infection by many types of viruses. This leads to translational blockage and the formation of stress granules. Furthermore, we show that APOBEC3B localizes to stress granules through the interaction with PABPC1. APOBEC3B facilitates the formation of protein-RNA condensates with stress granule assembly factor (G3BP1) by protecting mRNA associated with stress granules from RNAse L-induced RNA cleavage during viral infection. These results not only reveal that APOBEC3B is a key regulator of different steps of the innate immune response throughout viral infection but also highlight an alternative mechanism by which APOBEC3B can impact virus replication without editing viral genomes.

Molecular dynamics simulations reveal the importance of amyloid-beta oligomer ß-sheet edge conformations in membrane permeabilization.

Oligomeric aggregates of the amyloid-beta peptide(1-42) (Aß42) are regarded as a primary cause of cytotoxicity related to membrane damage in Alzheimer's disease. However, a dynamical and structural characterization of pore-forming Aß42 oligomers at atomic detail has not been feasible. Here, we used Aß42 oligomer structures previously determined in a membrane-mimicking environment as putative model systems to study the pore formation process in phospholipid bilayers with all-atom molecular dynamics simulations. Multiple Aß42 oligomer sizes, conformations, and N-terminally truncated isoforms were investigated on the multi-µs time scale. We found that pore formation and ion permeation occur via edge conductivity and exclusively for ß-sandwich structures that feature exposed side-by-side ß-strand pairs formed by residues 9 to 21 of Aß42. The extent of pore formation and ion permeation depends on the insertion depth of hydrophilic residues 13 to 16 (HHQK domain) and thus on subtle differences in the overall stability, orientation, and conformation of the aggregates in the membrane. Additionally, we determined that backbone carbonyl and polar side-chain atoms from the edge strands directly contribute to the coordination sphere of the permeating ions. Furthermore, point mutations that alter the number of favorable side-chain contacts correlate with the ability of the Aß42 oligomer models to facilitate ion permeation in the bilayer center. Our findings suggest that membrane-inserted, layered ß-sheet edges are a key structural motif in pore-forming Aß42 oligomers independent of their size and play a pivotal role in aggregate-induced membrane permeabilization.

Enterovirus-Cardiomyocyte Interactions: Impact of Terminally Deleted Genomic RNAs on Viral and Host Functions.

Group B enteroviruses, including coxsackievirus B3 (CVB3), can persistently infect cardiac tissue and cause dilated cardiomyopathy. Persistence is linked to 5' terminal deletions of viral genomic RNAs that have been detected together with minor populations of full-length genomes in human infections. In this study, we explored the functions and interactions of the different viral RNA forms found in persistently infected patients and their putative role(s) in pathogenesis. Since enterovirus cardiac pathogenesis is linked to the viral proteinase 2A, we investigated the effect of different terminal genomic RNA deletions on 2A activity. We discovered that 5' terminal deletions in CVB3 genomic RNAs decreased the levels of 2A proteinase activity but could not abrogate it. Using newly generated viral reporters encoding nano-luciferase, we found that 5' terminal deletions resulted in decreased levels of viral protein and RNA synthesis in singly transfected cardiomyocyte cultures. Unexpectedly, when full-length and terminally deleted forms were cotransfected into cardiomyocytes, a cooperative interaction was observed, leading to increased viral RNA and protein production. However, when viral infections were carried out in cells harboring 5' terminally deleted CVB3 RNAs, a decrease in infectious particle production was observed. Our results provide a possible explanation for the necessity of full-length viral genomes during persistent infection, as they would stimulate efficient viral replication compared to that of the deleted genomes alone. To avoid high levels of viral particle production that would trigger cellular immune activation and host cell death, the terminally deleted RNA forms act to limit the production of viral particles, possibly as trans-dominant inhibitors. IMPORTANCE Enteroviruses like coxsackievirus B3 are able to initiate acute infections of cardiac tissue and, in some cases, to establish a long-term persistent infection that can lead to serious disease sequelae, including dilated cardiomyopathy. Previous studies have demonstrated the presence of 5' terminally deleted forms of enterovirus RNAs in heart tissues derived from patients with dilated cardiomyopathy. These deleted RNAs are found in association with very low levels of full-length enterovirus genomic RNAs, an interaction that may facilitate continued persistence while limiting virus particle production. Even in the absence of detectable infectious virus particle production, these deleted viral RNA forms express viral proteinases at levels capable of causing viral pathology. Our studies provide mechanistic insights into how full-length and deleted forms of enterovirus RNA cooperate to stimulate viral protein and RNA synthesis without stimulating infectious viral particle production. They also highlight the importance of targeting enteroviral proteinases to inhibit viral replication while at the same time limiting the long-term pathologies they trigger.

The 3D structure of lipidic fibrils of α-synuclein.

α-synuclein misfolding and aggregation into fibrils is a common feature of α-synucleinopathies, such as Parkinson's disease, in which α-synuclein fibrils are a characteristic hallmark of neuronal inclusions called Lewy bodies. Studies on the composition of Lewy bodies extracted postmortem from brain tissue of Parkinson's patients revealed that lipids and membranous organelles are also a significant component. Interactions between α-synuclein and lipids have been previously identified as relevant for Parkinson's disease pathology, however molecular insights into their interactions have remained elusive. Here we present cryo-electron microscopy structures of six α-synuclein fibrils in complex with lipids, revealing specific lipid-fibril interactions. We observe that phospholipids promote an alternative protofilament fold, mediate an unusual arrangement of protofilaments, and fill the central cavities of the fibrils. Together with our previous studies, these structures also indicate a mechanism for fibril-induced lipid extraction, which is likely to be involved in the development of α-synucleinopathies. Specifically, one potential mechanism for the cellular toxicity is the disruption of intracellular vesicles mediated by fibrils and oligomers, and therefore the modulation of these interactions may provide a promising strategy for future therapeutic interventions.