The dynamical interplay between a megadalton peptide nanocage and solutes probed by microsecond atomistic MD; implications for design.

Understanding the assembly and dynamics of protein-based supramolecular capsids and cages is of fundamental importance and could lead to applications in synthetic biology and biotechnology. Here we present long and large atomistic molecular dynamics simulations of de novo designed self-assembling protein nanocages (SAGEs) in aqueous media. Microsecond simulations, comprised of ≈42 million atoms for three pre-formed SAGEs of different charges, in the presence of solutes and solvent have been completed. Here, the dynamics, stability and porosity of the peptide networks are explored along with their interactions with ions, small molecules and macromolecular solutes. All assemblies are stable over the µs timescale, and the solutes show a mixture of transport behaviour across or adherence to the fabric of the SAGE particles. Solute proteins largely retained native-like conformation on contact with SAGE. Certain residues of the SAGE peptides are identified as "repeat offenders" for contacting many different solutes, which suggest modifications to reduce non-specific binding. These studies highlight how molecular dynamics can aid the design process of SAGE and similar assemblies for potential applications as diverse as platforms for drug and vaccine delivery and nanoreactors to encapsulate enzyme pathways.

Elfin: An algorithm for the computational design of custom three-dimensional structures from modular repeat protein building blocks.

Computational protein design methods have enabled the design of novel protein structures, but they are often still limited to small proteins and symmetric systems. To expand the size of designable proteins while controlling the overall structure, we developed Elfin, a genetic algorithm for the design of novel proteins with custom shapes using structural building blocks derived from experimentally verified repeat proteins. By combining building blocks with compatible interfaces, it is possible to rapidly build non-symmetric large structures (>1000 amino acids) that match three-dimensional geometric descriptions provided by the user. A run time of about 20min on a laptop computer for a 3000 amino acid structure makes Elfin accessible to users with limited computational resources. Protein structures with controlled geometry will allow the systematic study of the effect of spatial arrangement of enzymes and signaling molecules, and provide new scaffolds for functional nanomaterials.

Atomistic non-adiabatic dynamics of the LH2 complex with a GPU-accelerated ab initio exciton model.

We recently outlined an efficient multi-tiered parallel ab initio excitonic framework that utilizes time dependent density functional theory (TDDFT) to calculate ground and excited state energies and gradients of large supramolecular complexes in atomistic detail - enabling us to undertake non-adiabatic simulations which explicitly account for the coupled anharmonic vibrational motion of all the constituent atoms in a supramolecular system. Here we apply that framework to the 27 coupled bacterio-chlorophyll-a chromophores which make up the LH2 complex, using it to compute an on-the-fly nonadiabatic surface-hopping (SH) trajectory of electronically excited LH2. Part one of this article is focussed on calibrating our ab initio exciton Hamiltonian using two key parameters: a shift δ, which corrects for the error in TDDFT vertical excitation energies; and an effective dielectric constant ε, which describes the average screening of the transition-dipole coupling between chromophores. Using snapshots obtained from equilibrium molecular dynamics simulations (MD) of LH2, we tune the values of both δ and ε through fitting to the thermally broadened experimental absorption spectrum, giving a linear absorption spectrum that agrees reasonably well with experiment. In part two of this article, we construct a time-resolved picture of the coupled vibrational and excitation energy transfer (EET) dynamics in the sub-picosecond regime following photo-excitation. Assuming Franck-Condon excitation of a narrow eigenstate band centred at 800 nm, we use surface hopping to follow a single nonadiabatic dynamics trajectory within the full eigenstate manifold. Consistent with experimental data, this trajectory gives timescales for B800→B850 population transfer (τB800→B850) between 650-1050 fs, and B800 population decay (τ800→) between 10-50 fs. The dynamical picture that emerges is one of rapidly fluctuating LH2 eigenstates that are delocalized over multiple chromophores and undergo frequent crossing on a femtosecond timescale as a result of the atomic vibrations of the constituent chromophores. The eigenstate fluctuations arise from disorder that is driven by vibrational dynamics with multiple characteristic timescales. The scalability of our ab initio excitonic computational framework across massively parallel architectures opens up the possibility of addressing a wide range of questions, including how specific dynamical motions impact both the pathways and efficiency of electronic energy-transfer within large supramolecular systems.

Adaptive free energy sampling in multidimensional collective variable space using boxed molecular dynamics.

The past decade has seen the development of a new class of rare event methods in which molecular configuration space is divided into a set of boundaries/interfaces, and then short trajectories are run between boundaries. For all these methods, an important concern is how to generate boundaries. In this paper, we outline an algorithm for adaptively generating boundaries along a free energy surface in multi-dimensional collective variable (CV) space, building on the boxed molecular dynamics (BXD) rare event algorithm. BXD is a simple technique for accelerating the simulation of rare events and free energy sampling which has proven useful for calculating kinetics and free energy profiles in reactive and non-reactive molecular dynamics (MD) simulations across a range of systems, in both NVT and NVE ensembles. Two key developments outlined in this paper make it possible to automate BXD, and to adaptively map free energy and kinetics in complex systems. First, we have generalized BXD to multidimensional CV space. Using strategies from rigid-body dynamics, we have derived a simple and general velocity-reflection procedure that conserves energy for arbitrary collective variable definitions in multiple dimensions, and show that it is straightforward to apply BXD to sampling in multidimensional CV space so long as the Cartesian gradients ∇CV are available. Second, we have modified BXD to undertake on-the-fly statistical analysis during a trajectory, harnessing the information content latent in the dynamics to automatically determine boundary locations. Such automation not only makes BXD considerably easier to use; it also guarantees optimal boundaries, speeding up convergence. We have tested the multidimensional adaptive BXD procedure by calculating the potential of mean force for a chemical reaction recently investigated using both experimental and computational approaches - i.e., F + CD3CN → DF + D2CN in both the gas phase and a strongly coupled explicit CD3CN solvent. The results obtained using multidimensional adaptive BXD agree well with previously published experimental and computational results, providing good evidence for its reliability.

High performance in silico virtual drug screening on many-core processors.

Drug screening is an important part of the drug development pipeline for the pharmaceutical industry. Traditional, lab-based methods are increasingly being augmented with computational methods, ranging from simple molecular similarity searches through more complex pharmacophore matching to more computationally intensive approaches, such as molecular docking. The latter simulates the binding of drug molecules to their targets, typically protein molecules. In this work, we describe BUDE, the Bristol University Docking Engine, which has been ported to the OpenCL industry standard parallel programming language in order to exploit the performance of modern many-core processors. Our highly optimized OpenCL implementation of BUDE sustains 1.43 TFLOP/s on a single Nvidia GTX 680 GPU, or 46% of peak performance. BUDE also exploits OpenCL to deliver effective performance portability across a broad spectrum of different computer architectures from different vendors, including GPUs from Nvidia and AMD, Intel's Xeon Phi and multi-core CPUs with SIMD instruction sets.

Rapid decomposition and visualisation of protein-ligand binding free energies by residue and by water.

Recent advances in computational hardware, software and algorithms enable simulations of protein-ligand complexes to achieve timescales during which complete ligand binding and unbinding pathways can be observed. While observation of such events can promote understanding of binding and unbinding pathways, it does not alone provide information about the molecular drivers for protein-ligand association, nor guidance on how a ligand could be optimised to better bind to the protein. We have developed the waterswap (C. J. Woods et al., J. Chem. Phys., 2011, 134, 054114) absolute binding free energy method that calculates binding affinities by exchanging the ligand with an equivalent volume of water. A significant advantage of this method is that the binding free energy is calculated using a single reaction coordinate from a single simulation. This has enabled the development of new visualisations of binding affinities based on free energy decompositions to per-residue and per-water molecule components. These provide a clear picture of which protein-ligand interactions are strong, and which active site water molecules are stabilised or destabilised upon binding. Optimisation of the algorithms underlying the decomposition enables near-real-time visualisation, allowing these calculations to be used either to provide interactive feedback to a ligand designer, or to provide run-time analysis of protein-ligand molecular dynamics simulations.

A GPU-accelerated immersive audio-visual framework for interaction with molecular dynamics using consumer depth sensors.

With advances in computational power, the rapidly growing role of computational/simulation methodologies in the physical sciences, and the development of new human-computer interaction technologies, the field of interactive molecular dynamics seems destined to expand. In this paper, we describe and benchmark the software algorithms and hardware setup for carrying out interactive molecular dynamics utilizing an array of consumer depth sensors. The system works by interpreting the human form as an energy landscape, and superimposing this landscape on a molecular dynamics simulation to chaperone the motion of the simulated atoms, affecting both graphics and sonified simulation data. GPU acceleration has been key to achieving our target of 60 frames per second (FPS), giving an extremely fluid interactive experience. GPU acceleration has also allowed us to scale the system for use in immersive 360° spaces with an array of up to ten depth sensors, allowing several users to simultaneously chaperone the dynamics. The flexibility of our platform for carrying out molecular dynamics simulations has been considerably enhanced by wrappers that facilitate fast communication with a portable selection of GPU-accelerated molecular force evaluation routines. In this paper, we describe a 360° atmospheric molecular dynamics simulation we have run in a chemistry/physics education context. We also describe initial tests in which users have been able to chaperone the dynamics of 10-alanine peptide embedded in an explicit water solvent. Using this system, both expert and novice users have been able to accelerate peptide rare event dynamics by 3-4 orders of magnitude.

Computational assay of H7N9 influenza neuraminidase reveals R292K mutation reduces drug binding affinity.

The emergence of a novel H7N9 avian influenza that infects humans is a serious cause for concern. Of the genome sequences of H7N9 neuraminidase available, one contains a substitution of arginine to lysine at position 292, suggesting a potential for reduced drug binding efficacy. We have performed molecular dynamics simulations of oseltamivir, zanamivir and peramivir bound to H7N9, H7N9-R292K, and a structurally related H11N9 neuraminidase. They show that H7N9 neuraminidase is structurally homologous to H11N9, binding the drugs in identical modes. The simulations reveal that the R292K mutation disrupts drug binding in H7N9 in a comparable manner to that observed experimentally for H11N9-R292K. Absolute binding free energy calculations with the WaterSwap method confirm a reduction in binding affinity. This indicates that the efficacy of antiviral drugs against H7N9-R292K will be reduced. Simulations can assist in predicting disruption of binding caused by mutations in neuraminidase, thereby providing a computational 'assay.'

Analysis and assay of oseltamivir-resistant mutants of influenza neuraminidase via direct observation of drug unbinding and rebinding in simulation.

The emergence of influenza drug resistance is a major public health concern. The molecular basis of resistance to oseltamivir (Tamiflu) is investigated using a computational assay involving multiple 500 ns unrestrained molecular dynamics (MD) simulations of oseltamivir complexed with mutants of H1N1-2009 influenza neuraminidase. The simulations, accelerated using graphics processors (GPUs), and using a fully explicit model of water, are of sufficient length to observe multiple drug unbinding and rebinding events. Drug unbinding occurs during simulations of known oseltamivir-resistant mutants of neuraminidase. Molecular-level rationalizations of drug resistance are revealed by analysis of these unbinding trajectories, with particular emphasis on the dynamics of the mutant residues. The results indicate that MD simulations can predict weakening of binding associated with drug resistance. In addition, visualization and analysis of binding site water molecules reveal their importance in stabilizing the binding mode of the drug. Drug unbinding is accompanied by conformational changes, driven by the mutant residues, which results in flooding of a key pocket containing tightly bound water molecules. This displaces oseltamivir, allowing the tightly bound water molecules to be released into bulk. In addition to the role of water, analysis of the trajectories reveals novel behavior of the structurally important 150-loop. Motion of the loop, which can move between an open and closed conformation, is intimately associated with drug unbinding and rebinding. Opening of the loop occurs coincidentally with drug unbinding, and interactions between oseltamivir and the loop seem to aid in the repositioning of the drug back into an approximation of its original binding mode on rebinding. The similarity of oseltamivir to a transition state analogue for neuraminidase suggests that the dynamics of the loop could play an important functional role in the enzyme, with loop closing aiding in binding of the substrate and loop opening aiding the release of the product.

A massively multicore parallelization of the Kohn-Sham energy gradients.

In a previous article [Brown et al., J Chem Theory Comput 2009, 4, 1620], we described a quadrature-based formulation of the Kohn-Sham Coulomb problem that allows for efficient parallelization over thousands of small processor cores. Here, we present the analytic gradients of this modified Kohn-Sham scheme, and describe the parallel implementation of the gradients on a numerical accelerator architecture. We demonstrate an order-of-magnitude acceleration for the combined energy and gradient calculation over a conventional single-core implementation.

Massively Multicore Parallelization of Kohn-Sham Theory.

A multicore parallelization of Kohn-Sham density functional theory is described, using an accelerator technology made by ClearSpeed Technology. Efficiently scaling parallelization over 2304 cores is achieved. To deliver this degree of parallelism, the Coulomb problem is reformulated to use Poisson density fitting with numerical quadrature of the required three-index integrals; extensive testing reveals negligible errors from the additional approximations.