Protocol for iterative optimization of modified peptides bound to protein targets.

Peptides are commonly used as therapeutic agents. However, they suffer from easy degradation and instability. Replacing natural by non-natural amino acids can avoid these problems, and potentially improve the affinity towards the target protein. Here, we present a computational pipeline to optimize peptides based on adding non-natural amino acids while improving their binding affinity. The workflow is an iterative computational evolution algorithm, inspired by the PARCE protocol, that performs single-point mutations on the peptide sequence using modules from the Rosetta framework. The modifications can be guided based on the structural properties or previous knowledge of the biological system. At each mutation step, the affinity to the protein is estimated by sampling the complex conformations and applying a consensus metric using various open protein-ligand scoring functions. The mutations are accepted based on the score differences, allowing for an iterative optimization of the initial peptide. The sampling/scoring scheme was benchmarked with a set of protein-peptide complexes where experimental affinity values have been reported. In addition, a basic application using a known protein-peptide complex is also provided. The structure- and dynamic-based approach allows users to optimize bound peptides, with the option to personalize the code for further applications. The protocol, called mPARCE, is available at: https://github.com/rochoa85/mPARCE/ .

Ab initio metadynamics determination of temperature-dependent free-energy landscape in ultrasmall silver clusters.

Ab initio metadynamics enables the extraction of free-energy landscapes having the accuracy of first-principles electronic structure methods. We introduce an interface between the PLUMED code that computes free-energy landscapes and enhanced-sampling algorithms and the Atomic Simulation Environment (ASE) module, which includes several ab initio electronic structure codes. The interface is validated with a Lennard-Jones cluster free-energy landscape calculation by averaging multiple short metadynamics trajectories. We use this interface and analysis to estimate the free-energy landscape of Ag5 and Ag6 clusters at 10, 100, and 300 K with the radius of gyration and coordination number as collective variables, finding at most tens of meV in error. Relative free-energy differences between the planar and non-planar isomers of both clusters decrease with temperature in agreement with previously proposed stabilization of non-planar isomers. Interestingly, we find that Ag6 is the smallest silver cluster where entropic effects at room temperature boost the non-planar isomer probability to a competing state. The new ASE-PLUMED interface enables simulating nanosystem electronic properties under more realistic temperature-dependent conditions.

In silico discovery and biological validation of ligands of FAD synthase, a promising new antimicrobial target.

New treatments for diseases caused by antimicrobial-resistant microorganisms can be developed by identifying unexplored therapeutic targets and by designing efficient drug screening protocols. In this study, we have screened a library of compounds to find ligands for the flavin-adenine dinucleotide synthase (FADS) -a potential target for drug design against tuberculosis and pneumonia- by implementing a new and efficient virtual screening protocol. The protocol has been developed for the in silico search of ligands of unexplored therapeutic targets, for which limited information about ligands or ligand-receptor structures is available. It implements an integrative funnel-like strategy with filtering layers that increase in computational accuracy. The protocol starts with a pharmacophore-based virtual screening strategy that uses ligand-free receptor conformations from molecular dynamics (MD) simulations. Then, it performs a molecular docking stage using several docking programs and an exponential consensus ranking strategy. The last filter, samples the conformations of compounds bound to the target using MD simulations. The MD conformations are scored using several traditional scoring functions in combination with a newly-proposed score that takes into account the fluctuations of the molecule with a Morse-based potential. The protocol was optimized and validated using a compound library with known ligands of the Corynebacterium ammoniagenes FADS. Then, it was used to find new FADS ligands from a compound library of 14,000 molecules. A small set of 17 in silico filtered molecules were tested experimentally. We identified five inhibitors of the activity of the flavin adenylyl transferase module of the FADS, and some of them were able to inhibit growth of three bacterial species: C. ammoniagenes, Mycobacterium tuberculosis, and Streptococcus pneumoniae, where the last two are human pathogens. Overall, the results show that the integrative VS protocol is a cost-effective solution for the discovery of ligands of unexplored therapeutic targets.

PepFun: Open Source Protocols for Peptide-Related Computational Analysis.

Peptide research has increased during the last years due to their applications as biomarkers, therapeutic alternatives or as antigenic sub-units in vaccines. The implementation of computational resources have facilitated the identification of novel sequences, the prediction of properties, and the modelling of structures. However, there is still a lack of open source protocols that enable their straightforward analysis. Here, we present PepFun, a compilation of bioinformatics and cheminformatics functionalities that are easy to implement and customize for studying peptides at different levels: sequence, structure and their interactions with proteins. PepFun enables calculating multiple characteristics for massive sets of peptide sequences, and obtaining different structural observables derived from protein-peptide complexes. In addition, random or guided library design of peptide sequences can be customized for screening campaigns. The package has been created under the python language based on built-in functions and methods available in the open source projects BioPython and RDKit. We present two tutorials where we tested peptide binders of the MHC class II and the Granzyme B protease.

An automated protocol for modelling peptide substrates to proteases.

BACKGROUND: Proteases are key drivers in many biological processes, in part due to their specificity towards their substrates. However, depending on the family and molecular function, they can also display substrate promiscuity which can also be essential. Databases compiling specificity matrices derived from experimental assays have provided valuable insights into protease substrate recognition. Despite this, there are still gaps in our knowledge of the structural determinants. Here, we compile a set of protease crystal structures with bound peptide-like ligands to create a protocol for modelling substrates bound to protease structures, and for studying observables associated to the binding recognition. RESULTS: As an application, we modelled a subset of protease-peptide complexes for which experimental cleavage data are available to compare with informational entropies obtained from protease-specificity matrices. The modelled complexes were subjected to conformational sampling using the Backrub method in Rosetta, and multiple observables from the simulations were calculated and compared per peptide position. We found that some of the calculated structural observables, such as the relative accessible surface area and the interaction energy, can help characterize a protease's substrate recognition, giving insights for the potential prediction of novel substrates by combining additional approaches. CONCLUSION: Overall, our approach provides a repository of protease structures with annotated data, and an open source computational protocol to reproduce the modelling and dynamic analysis of the protease-peptide complexes.

Need for Cross-Validation of Single Particle Cryo-EM.

Cross-validation is used to determine the validity of a model on unseen data by assessing if the model is overfitted to noise. It is widely used in many fields, from artificial intelligence to structural biology in X-ray crystallography and nuclear magnetic resonance. Although there are concerns of map overfitting in cryo-electron microscopy (cryo-EM), cross-validation is rarely used. The problem is that establishing a performance metric of the maps over unseen data (given by 2D-projection images) is difficult due to the low signal-to-noise ratios in the individual particles. Here, I present recent advances for cryo-EM map reconstruction. I highlight that the gold-standard procedure can fail to detect map overfitting in certain cases, showing the necessity of assessing the map quality on unbiased data. Finally, I describe the challenges and advantages of developing a robust cross-validation methodology for cryo-EM.

Flexi-pharma: a molecule-ranking strategy for virtual screening using pharmacophores from ligand-free conformational ensembles.

Computer-aided strategies are useful for reducing the costs and increasing the success-rate in drug discovery. Among these strategies, methods based on pharmacophores (an ensemble of electronic and steric features representing the target active site) are efficient to implement over large compound libraries. However, traditional pharmacophore-based methods require knowledge of active compounds or ligand-receptor structures, and only few ones account for target flexibility. Here, we developed a pharmacophore-based virtual screening protocol, Flexi-pharma, that overcomes these limitations. The protocol uses molecular dynamics (MD) simulations to explore receptor flexibility, and performs a pharmacophore-based virtual screening over a set of MD conformations without requiring prior knowledge about known ligands or ligand-receptor structures for building the pharmacophores. The results from the different receptor conformations are combined using a "voting" approach, where a vote is given to each molecule that matches at least one pharmacophore from each MD conformation. Contrarily to other approaches that reduce the pharmacophore ensemble to some representative models and score according to the matching models or molecule conformers, the Flexi-pharma approach takes directly into account the receptor flexibility by scoring in regards to the receptor conformations. We tested the method over twenty systems, finding an enrichment of the dataset for 19 of them. Flexi-pharma is computationally efficient allowing for the screening of thousands of compounds in minutes on a single CPU core. Moreover, the ranking of molecules by vote is a general strategy that can be applied with any pharmacophore-filtering program.

Predicting the Affinity of Peptides to Major Histocompatibility Complex Class II by Scoring Molecular Dynamics Simulations.

Predicting the binding affinity of peptides able to interact with major histocompatibility complex (MHC) molecules is a priority for researchers working in the identification of novel vaccines candidates. Most available approaches are based on the analysis of the sequence of peptides of known experimental affinity. However, for MHC class II receptors, these approaches are not very accurate, due to the intrinsic flexibility of the complex. To overcome these limitations, we propose to estimate the binding affinity of peptides bound to an MHC class II by averaging the score of the configurations from finite-temperature molecular dynamics simulations. The score is estimated for 18 different scoring functions, and we explored the optimal manner for combining them. To test the predictions, we considered eight peptides of known binding affinity. We found that six scoring functions correlate with the experimental ranking of the peptides significantly better than the others. We then assessed a set of techniques for combining the scoring functions by linear regression and logistic regression. We obtained a maximum accuracy of 82% for the predicted sign of the binding affinity using a logistic regression with optimized weights. These results are potentially useful to improve the reliability of in silico protocols to design high-affinity binding peptides for MHC class II receptors.

Molecular free energy profiles from force spectroscopy experiments by inversion of observed committors.

In single-molecule force spectroscopy experiments, a biomolecule is attached to a force probe via polymer linkers and the total extension of the molecule plus apparatus is monitored as a function of time. In a typical unfolding experiment at constant force, the total extension jumps between two values that correspond to the folded and unfolded states of the molecule. For several biomolecular systems, the committor, which is the probability to fold starting from a given extension, has been used to extract the molecular activation barrier (a technique known as "committor inversion"). In this work, we study the influence of the force probe, which is much larger than the molecule being measured, on the activation barrier obtained by committor inversion. We use a two-dimensional framework in which the diffusion coefficient of the molecule and of the pulling device can differ. We systematically study the free energy profile along the total extension obtained from the committor by numerically solving the Onsager equation and using Brownian dynamics simulations. We analyze the dependence of the extracted barrier on the linker stiffness, molecular barrier height, and diffusion anisotropy and, thus, establish the range of validity of committor inversion. Along the way, we showcase the committor of 2-dimensional diffusive models and illustrate how it is affected by barrier asymmetry and diffusion anisotropy.

The Dimer-of-Trimers Assembly Prevents Catalysis at the Transferase Site of Prokaryotic FAD Synthase.

Flavin mononucleotide (FMN) and flavin-adenine dinucleotide (FAD) are essential flavoprotein cofactors. A riboflavin kinase (RFK) activity catalyzes riboflavin phosphorylation to FMN, which can then be transformed into FAD by an FMN:adenylyltransferase (FMNAT) activity. Two enzymes are responsible for each one of these activities in eukaryotes, whereas prokaryotes have a single bifunctional enzyme, FAD synthase (FADS). FADS folds in two independent modules: the C-terminal with RFK activity and the N-terminal with FMNAT activity. Differences in structure and chemistry for the FMNAT catalysis among prokaryotic and eukaryotic enzymes pointed to the FMNAT activity of prokaryotic FADS as a potential antimicrobial target, making the structural model of the bacterial FMNAT module in complex with substrates relevant to understand the FADS catalytic mechanism and to the discovery of antimicrobial drugs. However, such a crystallographic complex remains elusive. Here, we have used molecular docking and molecular dynamics simulations to generate energetically stable interactions of the FMNAT module of FADS from Corynebacterium ammoniagenes with ATP/Mg2+ and FMN in both the monomeric and dimer-of-trimers assemblies reported for this protein. For the monomer, we have identified the residues that accommodate the reactive phosphates in a conformation compatible with catalysis. Interestingly, for the dimer-of-trimers conformation, we have found that the RFK module negatively influences FMN binding at the interacting FMNAT module. These results agree with calorimetric data of purified samples containing nearly 100% monomer or nearly 100% dimer-of-trimers, indicating that FMN binds to the monomer but not to the dimer-of-trimers. Such observations support regulation of flavin homeostasis by quaternary C. ammoniagenes FADS assemblies.

Assessing the capability of in silico mutation protocols for predicting the finite temperature conformation of amino acids.

Mutation protocols are a key tool in computational biophysics for modelling unknown side chain conformations. In particular, these protocols are used to generate the starting structures for molecular dynamics simulations. The accuracy of the initial side chain and backbone placement is crucial to obtain a stable and quickly converging simulation. In this work, we assessed the performance of several mutation protocols in predicting the most probable conformer observed in finite temperature molecular dynamics simulations for a set of protein-peptide crystals differing only by single-point mutations in the peptide sequence. Our results show that several programs which predict well the crystal conformations fail to predict the most probable finite temperature configuration. Methods relying on backbone-dependent rotamer libraries have, in general, a better performance, but even the best protocol fails in predicting approximately 30% of the mutations.

Transition paths in single-molecule force spectroscopy.

In a typical single-molecule force spectroscopy experiment, the ends of the molecule of interest are connected by long polymer linkers to a pair of mesoscopic beads trapped in the focus of two laser beams. At constant force load, the total extension, i.e., the end-to-end distance of the molecule plus linkers, is measured as a function of time. In the simplest systems, the measured extension fluctuates about two values characteristic of folded and unfolded states, with occasional transitions between them. We have recently shown that molecular (un)folding rates can be recovered from such trajectories, with a small linker correction, as long as the characteristic time of the bead fluctuations is shorter than the residence time in the unfolded (folded) state. Here, we show that accurate measurements of the molecular transition path times require an even faster apparatus response. Transition paths, the trajectory segments in which the molecule (un)folds, are properly resolved only if the beads fluctuate more rapidly than the end-to-end distance of the molecule. Therefore, over a wide regime, the measured rates may be meaningful but not the transition path times. Analytic expressions for the measured mean transition path times are obtained for systems diffusing anisotropically on a two-dimensional free energy surface. The transition path times depend on the properties both of the molecule and of the pulling device.

Discovery of antimicrobial compounds targeting bacterial type FAD synthetases.

The increase of bacterial strains resistant to most of the available antibiotics shows a need to explore novel antibacterial targets to discover antimicrobial drugs. Bifunctional bacterial FAD synthetases (FADSs) synthesise the flavin mononucleotide (FMN) and flavin adenine dinucleotide (FAD). These cofactors act in vital processes as part of flavoproteins, making FADS an essential enzyme. Bacterial FADSs are potential antibacterial targets because of differences to mammalian enzymes, particularly at the FAD producing site. We have optimised an activity-based high throughput screening assay targeting Corynebacterium ammoniagenes FADS (CaFADS) that identifies inhibitors of its different activities. We selected the three best high-performing inhibitors of the FMN:adenylyltransferase activity (FMNAT) and studied their inhibition mechanisms and binding properties. The specificity of the CaFADS hits was evaluated by studying also their effect on the Streptococcus pneumoniae FADS activities, envisaging differences that can be used to discover species-specific antibacterial drugs. The antimicrobial effect of these compounds was also evaluated on C. ammoniagenes, S. pneumoniae, and Mycobacterium tuberculosis cultures, finding hits with favourable antimicrobial properties.

On artifacts in single-molecule force spectroscopy.

In typical force spectroscopy experiments, a small biomolecule is attached to a soft polymer linker that is pulled with a relatively large bead or cantilever. At constant force, the total extension stochastically changes between two (or more) values, indicating that the biomolecule undergoes transitions between two (or several) conformational states. In this paper, we consider the influence of the dynamics of the linker and mesoscopic pulling device on the force-dependent rate of the conformational transition extracted from the time dependence of the total extension, and the distribution of rupture forces in force-clamp and force-ramp experiments, respectively. For these different experiments, we derive analytic expressions for the observables that account for the mechanical response and dynamics of the pulling device and linker. Possible artifacts arise when the characteristic times of the pulling device and linker become comparable to, or slower than, the lifetimes of the metastable conformational states, and when the highly anharmonic regime of stretched linkers is probed at high forces. We also revisit the problem of relating force-clamp and force-ramp experiments, and identify a linker and loading rate-dependent correction to the rates extracted from the latter. The theory provides a framework for both the design and the quantitative analysis of force spectroscopy experiments by highlighting, and correcting for, factors that complicate their interpretation.

Kinetic Ductility and Force-Spike Resistance of Proteins from Single-Molecule Force Spectroscopy.

Ductile materials can absorb spikes in mechanical force, whereas brittle ones fail catastrophically. Here we develop a theory to quantify the kinetic ductility of single molecules from force spectroscopy experiments, relating force-spike resistance to the differential responses of the intact protein and the unfolding transition state to an applied mechanical force. We introduce a class of unistable one-dimensional potential surfaces that encompass previous models as special cases and continuously cover the entire range from ductile to brittle. Compact analytic expressions for force-dependent rates and rupture-force distributions allow us to analyze force-clamp and force-ramp pulling experiments. We find that the force-transmitting protein domains of filamin and titin are kinetically ductile when pulled from their two termini, making them resistant to force spikes. For the mechanostable muscle protein titin, a highly ductile model reconciles data over 10 orders of magnitude in force loading rate from experiment and simulation.

The molecular mechanism of substrate engagement and immunosuppressant inhibition of calcineurin.

Ser/thr phosphatases dephosphorylate their targets with high specificity, yet the structural and sequence determinants of phosphosite recognition are poorly understood. Calcineurin (CN) is a conserved Ca(2+)/calmodulin-dependent ser/thr phosphatase and the target of immunosuppressants, FK506 and cyclosporin A (CSA). To investigate CN substrate recognition we used X-ray crystallography, biochemistry, modeling, and in vivo experiments to study A238L, a viral protein inhibitor of CN. We show that A238L competitively inhibits CN by occupying a critical substrate recognition site, while leaving the catalytic center fully accessible. Critically, the 1.7 Å structure of the A238L-CN complex reveals how CN recognizes residues in A238L that are analogous to a substrate motif, "LxVP." The structure enabled modeling of a peptide substrate bound to CN, which predicts substrate interactions beyond the catalytic center. Finally, this study establishes that "LxVP" sequences and immunosuppressants bind to the identical site on CN. Thus, FK506, CSA, and A238L all prevent "LxVP"-mediated substrate recognition by CN, highlighting the importance of this interaction for substrate dephosphorylation. Collectively, this work presents the first integrated structural model for substrate selection and dephosphorylation by CN and lays the groundwork for structure-based development of new CN inhibitors.

Good Rates From Bad Coordinates: The Exponential Average Time-dependent Rate Approach.

Our ability to calculate rate constants of biochemical processes using molecular dynamics simulations is severely limited by the fact that the time scales for reactions, or changes in conformational state, scale exponentially with the relevant free-energy barrier heights. In this work, we improve upon a recently proposed rate estimator that allows us to predict transition times with molecular dynamics simulations biased to rapidly explore one or several collective variables (CVs). This approach relies on the idea that not all bias goes into promoting transitions, and along with the rate, it estimates a concomitant scale factor for the bias termed the "CV biasing efficiency" γ. First, we demonstrate mathematically that our new formulation allows us to derive the commonly used Infrequent Metadynamics (iMetaD) estimator when using a perfect CV, where γ = 1. After testing it on a model potential, we then study the unfolding behavior of a previously well characterized coarse-grained protein, which is sufficiently complex that we can choose many different CVs to bias, but which is sufficiently simple that we are able to compute the unbiased rate directly. For this system, we demonstrate that predictions from our new Exponential Average Time-Dependent Rate (EATR) estimator converge to the true rate constant more rapidly as a function of bias deposition time than does the previous iMetaD approach, even for bias deposition times that are short. We also show that the γ parameter can serve as a good metric for assessing the quality of the biasing coordinate. We demonstrate that these results hold when applying the methods to an atomistic protein folding example. Finally, we demonstrate that our approach works when combining multiple less-than-optimal bias coordinates, and adapt our method to the related "OPES flooding" approach. Overall, our time-dependent rate approach offers a powerful framework for predicting rate constants from biased simulations.

Bayesian analysis of individual electron microscopy images: towards structures of dynamic and heterogeneous biomolecular assemblies.

We develop a method to extract structural information from electron microscopy (EM) images of dynamic and heterogeneous molecular assemblies. To overcome the challenge of disorder in the imaged structures, we analyze each image individually, avoiding information loss through clustering or averaging. The Bayesian inference of EM (BioEM) method uses a likelihood-based probabilistic measure to quantify the consistency between each EM image and given structural models. The likelihood function accounts for uncertainties in the molecular position and orientation, variations in the relative intensities and noise in the experimental images. The BioEM formalism is physically intuitive and mathematically simple. We show that for experimental GroEL images, BioEM correctly identifies structures according to the functional state. The top-ranked structure is the corresponding X-ray crystal structure, followed by an EM structure generated previously from a superset of the EM images used here. To analyze EM images of highly flexible molecules, we propose an ensemble refinement procedure, and validate it with synthetic EM maps of the ESCRT-I-II supercomplex. Both the size of the ensemble and its structural members are identified correctly. BioEM offers an alternative to 3D-reconstruction methods, extracting accurate population distributions for highly flexible structures and their assemblies. We discuss limitations of the method, and possible applications beyond ensemble refinement, including the cross-validation and unbiased post-assessment of model structures, and the structural characterization of systems where traditional approaches fail. Overall, our results suggest that the BioEM framework can be used to analyze EM images of both ordered and disordered molecular systems.

Microscopic theory, analysis, and interpretation of conductance histograms in molecular junctions.

Molecular electronics break-junction experiments are widely used to investigate fundamental physics and chemistry at the nanoscale. Reproducibility in these experiments relies on measuring conductance on thousands of freshly formed molecular junctions, yielding a broad histogram of conductance events. Experiments typically focus on the most probable conductance, while the information content of the conductance histogram has remained unclear. Here we develop a microscopic theory for the conductance histogram by merging the theory of force-spectroscopy with molecular conductance. The procedure yields analytical equations that accurately fit the conductance histogram of a wide range of molecular junctions and augments the information content that can be extracted from them. Our formulation captures contributions to the conductance dispersion due to conductance changes during the mechanical elongation inherent to the experiments. In turn, the histogram shape is determined by the non-equilibrium stochastic features of junction rupture and formation. The microscopic parameters in the theory capture the junction's electromechanical properties and can be isolated from separate conductance and rupture force (or junction-lifetime) measurements. The predicted behavior can be used to test the range of validity of the theory, understand the conductance histograms, design molecular junction experiments with enhanced resolution and molecular devices with more reproducible conductance properties.