Emerging Methods for Ensemble-Based Virtual Screening
Ensemble based virtual screening refers to the use of conformational ensembles from crystal structures, NMR studies or molecular dynamics simulations. It has gained greater acceptance as advances in the theoretical framework, computational algorithms, and software packages enable simulations at longer time scales. Here we focus on the use of computationally generated conformational ensembles and emerging methods that use these ensembles for discovery, such as the Relaxed Complex Scheme or Dynamic Pharmacophore Model. We also discuss the more rigorous physics-based computational techniques such as accelerated molecular dynamics and thermodynamic integration and their applications in improving conformational sampling or the ranking of virtual screening hits. Finally, technological advances that will help make virtual screening tools more accessible to a wider audience in computer aided drug design are discussed.
INTRODUCTION<!>MOVING BEYOND \xe2\x80\x9cLOCK-AND-KEY\xe2\x80\x9d OR SIMPLE \xe2\x80\x9cINDUCED FIT\xe2\x80\x9d<!>CONFORMATIONAL ENSEMBLES<!>Dynamic Pharmacophore Model<!>The Relaxed Complex Scheme<!>ENSEMBLE SELECTION<!>RMSD-BASED CLUSTERING<!>QR-FACTORIZATION<!>Knowledge Gaps<!>Predicting Explicit Water Interactions<!>Promising New Ways to Enhance Sampling of the Receptor Ensemble<!>How Long is Long Enough?<!>Application of more Rigorous Approaches for Predicting Binding Affinities<!>CLOUD COMPUTING, GPU AND HARDWARE ACCELERATION<!>CONCLUSIONS
<p>The computational identification of drug leads out of large compound libraries through receptor-based virtual screening (VS) is a well-established method to predict putative inhibitors for target receptors (reviewed in [1–9]). Despite advances in the underlying algorithms of virtual screening experiments have incrementally improved our ability to discriminate binders from non-binders, a number of obstacles remain. A major outstanding challenge in the practice of virtual screening is the treatment of receptor flexibility, which is especially difficult owing to the many degrees of conformational freedom in target receptors. Steady increases in computational power, coupled with improvements in the underlying algorithms and available structural experimental data, are enabling a new paradigm for virtual screening, wherein computationally predicted ensembles from first-principle simulations are being used in rational drug design efforts. The integration of these more rigorous physics-based methods will have far reaching impact on translational medicine, including the ability to: (1) better understand the structural dynamics of disease-related target receptors, (2) improve our quantitative assessments of ligand-receptor interactions, (3) discover novel modes of ligand binding and inhibition, and (4) develop new therapeutics that are patient-specific and less prone to drug resistance. This review will attempt to summarize the recent computational and theoretical advances that have enabled the development of new methods to integrate larger-scale sampling of receptor space and identify the most biologically relevant structures for drug design. We will attempt to discuss the advantages and drawbacks of current methodologies, but focus on emerging methods that we believe will become more frequently utilized, especially for receptors that exhibit a high degree of flexibility.</p><!><p>A full understanding of molecular recognition presents a problem of intense interest to molecular sciences and the field of computer-aided drug design (CADD) [10–13]. The interactions between ligand molecules and their corresponding receptors are dynamic and complex. Over time, the field has moved away from an original understanding of molecular recognition as a simple, rigid "lock-and-key" mechanism [14], towards an "induced fit" [15]. In the traditional induced-fit theory, the intrinsic plasticity of the receptor is utilized when the protein's structure and function responds to the ligand-binding event through inducible structural motion, and the corresponding rearrangements of the receptor active site upon binding are only capable in the presence of the bound ligand. Many studies describing the induced fit effects in various protein-ligand systems have been presented (for examples, see [16–21]). It is widely accepted that ligands may bind to receptor conformations that occur infrequently in the receptor's dynamics, and that the local motions of active site residues can drastically alter the binding affinity and specificity of ligands to their target. The ability to efficiently sample these rare dynamics and furthermore, to incorporate the resulting conformations into the drug discovery and design protocol is still an active area of investigation. Instead of using static crystal structures in virtual screening, in the traditional lock and key sense, methods that explore large domain motions and significant active site rearrangements in a computationally efficient manner will be required.</p><!><p>An important theoretical notion that presented a modification to traditional induced fit theory is that of ligands binding to pre-existing receptor populations, or conformational ensembles [22, 23]. Borrowed in large part from the more theoretically-leaning field of protein folding, the general idea is that the receptor does not exist in a single native conformation, but more realistically exists in an ensemble of metastates, or energetically low-lying substates along the rugged receptor energy landscape [24, 25]. The ligand is exposed to this ensemble of receptor conformational states and may bind preferentially to one, shifting the equilibrium population towards the favorable bound conformation [26]. These theories are deeply rooted in statistical mechanics and state that the ensemble comprises a statistical distribution of native conformations that essentially move, subject to mechanical forces. Although a review of the elegant protein folding theories is outside the scope of this review (for comprehensive reviews, see [24, 27, 28]), it is important to recognize the influence of these theories on understanding the fundamental nature of protein structure, function and ligand binding [29]. An essential conceptual point that is highly relevant to this review is the recognition that receptors are able to sample productive bound-ligand conformations even while in the apo state. To some extent, this suggests that the "lock and key" is but one of the rare conformations within this unifying scheme, and that conformational selection is an important driving force in ligand binding and recognition.</p><p>Using multiple receptor conformations in CADD and virtual screening ensued not long after the conception of conformational ensembles. These receptor conformations may be obtained from x-ray crystallography, or nuclear magnetic resonance (NMR) experiments, or computationally, from molecular dynamics (MD) or Monte Carlo simulations. In a critical advance, Knegtel et al. [30] presented the first study employing both crystallographically and NMR derived ensembles in an averaged grid-based docking study. Other early and important examples of utilizing multiple experimental structures include FlexE [31], which systematically combined features from multiple crystal structures, and work presented by Osterberg et al., which carefully examined different methods to combine receptor grids stemming from the inclusion of almost two dozen receptor structures of HIV-1 protease [32]. Numerous additional studies continue to enhance the notion that using multiple receptor conformations improves enrichment factors and the ability to predict viable binding modes [22, 33–38].</p><p>Although experimentally derived structures have typically been considered the gold standard for docking and drug design, in many cases, a high correspondence between computationally derived ensembles and experimental structures has been established, though initially limited by insufficient sampling [39]. The progress to date suggests that MD-generated ensembles may closely replicate the structural dynamics of proteins in solution, as shown by NMR experiments [40, 41]. With this in mind, we focus on emerging virtual screening methods that utilize direct sampling of receptor conformational space through MD-based approaches, and note that comprehensive reviews of other methods to include receptor flexibility, including soft docking [42], sidechain rotamer libraries, docking to relevant normal modes [43], induced fit docking [44], and more are presented elsewhere [45–47].</p><!><p>Although the first MD simulation of a protein was performed over 30 years ago [48], the more systematic use of the resulting conformational ensembles in a predictive fashion has only recently been successfully demonstrated. The first experimentally verified study to use multiple computer-generated structures in a systematic ensemble-based approach to discovery was the multiple protein structures based dynamic pharmacophore (DPM) model approach presented by Carlson et al. [49]. In this work, a 500 ps explicitly solvated MD simulation of the apo HIV-1 integrase catalytic domain was performed [50] using the AMBER95 force field [51] and the NWChem v3.2 simulation program [52]. The flexible nature of the integrase system in the active site region precluded resolving the structure completely crystallographically; instead, the authors used a homology model (modeled after the complete structure of Avian Sarcoma Virus integrase) of the flexible region near the active site. Their predicted structure (and the corresponding MD trajectory) was later validated when two additional crystal structures of the full integrase catalytic domain were published, showing a high degree of similarity between the predicted regions and the new crystal information [49]. Eleven conformations from the MD trajectory were used as the representative ensemble to create a dynamic receptor-based pharmacophore model for the integrase system. The bare active sites were flooded with small molecule probes representing different chemical functional groups, which essentially map the most favorable areas in the active site, and are clustered across the ensemble of structures. The resulting DPM was used to search the Available Chemical Database for potential inhibitors, and rerank the known set of active compounds for the integrase. Experimental verification of the predicted set indicated that about one third of the compounds were inhibitory and the reranking with the DPM outperformed any single static pharmacophore.</p><p>Refinements to the DPM method were presented in an application focusing on the flexible HIV-1 protease as a test case [53]. In this work, Meagher and Carlson showed improved performance of the DPM to discriminate binders from non-binders with increasing lengths of MD simulations. This result reinforced not only the importance of including receptor flexibility in CADD efforts, but also the consideration that longer simulation lengths yield improved results for discovery. Presumably, the longer simulation times enhanced the conformational sampling of receptor energy landscape.</p><p>More recently, the DPM was used to discover novel inhibitors of the MDM2-p53 protein [54]. Snapshots extracted from a 2 ns MD simulation of the human MDM3 bound to p53 were extracted every 100 ps. These 21 structures were used to create a 6-site pharmacophore model of the active site. A virtual screening of 35,000 in-house compounds identified 27 compounds, of which 23 were experimentally tested. Four active compounds with unique chemical scaffolds were discovered to inhibit at 50 μM concentration or better. Interestingly, through a close analysis of the MD simulations and crystal structures, the authors had discovered an additional hydrophobic binding area near the known binding cleft [55]. Creating an additional pharmacophore model that utilized this site resulted in the discovery of a fifth compound that inhibited at ~ 20 μM. The discovery of these new compounds is an important step forward in the use of MD-generated structures in virtual screening, and validation of the dynamic pharmacophore method.</p><!><p>Another approach of using MD-generated structural ensembles for drug discovery utilizes a strategy referred to as the relaxed complex scheme (RCS) (Fig. 1). The RCS uses receptor snapshots extracted from MD simulations to search ligand libaries via small molecule docking [56, 57]. It combines the advantages of docking algorithms with dynamic structural information provided by MD simulations, explicitly accounting for the flexibility of both the receptor and docked ligands. This procedure directly applies the concept of conformational ensembles in molecular recognition through docking. Increasingly longer time scale MD simulations enhances our ability to effectively sample the receptor conformational space prior to docking. This scheme has been developed in combination with various MD software packages and AutoDock for the ligand docking [58], although other docking programs can be substituted. The RCS was first applied to the FKBP binding protein [59] and tested using improved re-scoring functions based on MM-PBSA models [56].</p><p>An important example of the method was presented in Schames et al., who performed a 2 ns explicitly solvated MD simulation of HIV-1 integrase in complex with the 5CITEP inhibitor [60]. The simulations revealed an additional cavity adjacent to the integrase active site (i.e. a "trench"), and docking of 5CITEP into this new area proved to be even more favorable than the actual active site based on the AutoDock 3.0 scoring function. Shortly thereafter, Hazuda et al. acknowledged that this new structural understanding, in which essential interactions occur in varying areas of the active site, was invaluable in creating new integrase inhibitors with unique resistance profiles [61]. Just a few years (and rounds of optimization) later, raltegravir was approved for use in humans [62]. This success story provides renewed motivation for the development of computational methodologies that incorporate predictive structural information for flexible receptors systematically.</p><!><p>In order to make these more rigorous yet time-consuming computational methods accessible to a wider audience, various ways to distill the structural ensemble, and thus reduce the requisite computational power, need to be explored. With constant improvements in computing power and the underlying simulation algorithms, we will very soon routinely be sampling into the microsecond timescale [63]. Performing virtual screening experiments on the full set of resulting structures is computationally intractable and likely unnecessary. New metrics to reduce the ensemble systematically – without losing critical structural information – will be important. Several strategies have been developed that select structural information from the resulting ensemble, for both the DPM and the RCS. These approaches, which range from manual selection to more systematic mathematical techniques, both reduce the ensemble to a manageable size and allow the extraction of the most meaningful (i.e. biologically relevant) information.</p><!><p>RMSD-based clustering has been established as an alternate metric that can be used in order to reduce an ensemble to a meaningful set before virtual screening experiments. RMSD-based clustering can be especially useful because it provides information about the most dominant configurations sampled during an MD trajectory. One can argue based on statistics that the most populated clusters of receptor structures would be more frequently encountered by the ligand, energetically stable, and statistically meaningful. Methods remain to be developed that could incorporate both population and energy states, when seeking to create more minimal representative ensembles from MD trajectories.</p><p>Deng et al. presented a refinement of the DPM that combined RMSD-based clustering in a hierarchical fashion with energetic probabilities in order to intelligently select which structures to include in the pharmacophore model [64]. In this study, the authors performed a 1 ns simulation of the HIV integrase in complex with the 5CITEP inhibitor, again using homology modeling to complete part of the active site structure [65]. 1000 snapshots were extracted from the MD trajectory at 1 ps intervals. RMSD-based clustering based on five key active site residues was performed, followed by a relative estimation of the probability of each MD snapshot based on a Boltzmann factor (e−ΔU/RT), wherein the energy difference was based on the difference from the minimum energy structure [66]. The cluster with the largest total sum of the Boltzmann factors was selected as the first cluster. The snapshots from this first cluster were removed from the ensemble and the process repeated to determine the second cluster, and so forth. The final computational ensemble comprised 10 receptor conformations, and collectively represented approximately 50% of the trajectory. LigBuilder [67] was then used to create a DPM based on these 10 representative cluster structures and Catalyst 4.6 [68] was used to search an in-house database of 400 available chemical compounds. A final set of 23 compounds was selected for experimental assays and 9 showed inhibitory capacity at 100 μM or less. Several compounds so identified were missed if only a static pharmacophore model was used.</p><p>Cheng et al. recently extended the RCS for virtual screening by the efficient use of RMSD-based clustering information as well [69] In order to distill the most dominant configurations from the MD simulations, Cheng et al. performed RMSD-based clustering on snapshots extracted every 10 ps from two 40 ns trajectories of avian influenza N1 neuraminidase in the apo form and in complex with the inhibitor oseltamivir [70]. Although the tetramer N1 was used in the simulations, the clustering analyses were carried out on individual monomer protein chains [71]. Therefore, each 40 ns tetramer simulation yielded the equivalent of four times the monomer sampling (160 ns), yielding 16,000 structures for the analyses. To date, this study was based on considerably longer MD trajectories than any other published examples used in CADD (the equivalent of 160 ns for the monomer for each of the apo and holo systems). Therefore an extensive sampling of the receptor configurational space was achieved; in addition, the four monomer copies in the tetramer allows for a multi-copy MD approach, which has been shown to enhance sampling as opposed to a single long trajectory [72]. The RMSD-clustering was performed on a subset of binding-site residues and yielded 10 structures representing the apo ensemble and 5 structures representing the holo ensemble. The central member structure from the three most dominant clusters in the apo and holo sets were screened against the National Cancer Institute Diversity Set 1 (NCIDS1) using AutoDock 4.0 [73]. The final ranking of all ligands from 8 primary screens was determined by taking the weighted average of the docking scores into the full representative ensemble of the holo MD trajectory (Fig. 2). Of the experimentally tested set of 25 compounds, 10 exhibited Ki's of under 500 μM. Of these, 7 compounds were only selected through the screening of the MD-generated structures, which exhibited large structural rearrangements during simulation (Fig. 2) [74]. Eight of the compounds are predicted to at least partially bind to novel binding sites revealed in the simulations. This significant reordering enabled the identification of hits that otherwise would have been missed.</p><p>In a separate system published at almost the same time, Zhong et al. performed an explicitly solvated 5 ns MD simulation of human DNA ligase 1 using the CHARMM27 force field [75] with CMAP corrections [76] and the CHARMM MD engine [77]. Receptor coordinates were extracted every 5 ps (i.e. 1000 conformations total) [78]. They then performed RMSD-based clustering on binding site residues and used representative structures from the four largest clusters in a second round of docking. The top 50,000 compounds from the first phase were docked using DOCK4.0 [79] into the crystal structure as well as representative structures from the five most populated clusters. Compounds were ranked according to the most favorable individual score from the set and 192 compounds tested experimentally. Of the 10 active compounds discovered, nine were chosen based on MD-generated representative clusters.</p><!><p>An alternate technique that has been applied to virtual screening experiments using the RCS is QR-factorization. The technique was originally designed to remove inherent bias in structure databases and distill, from a vast quantity of redundant information, a minimal basis set of protein structures accurately spanning the evolutionary conformation space of a particular protein [80]. In 2008, Amaro et al. incorporated this new method to create a minimal basis set for the configurational space sampled in the MD simulations and used the resulting ensemble in a virtual screening with the RNA editing ligase 1 enzyme in Trypanosoma brucei [81]. In this work, a 20 ns explicitly solvated MD simulation employing the Charmm27 force field [75] and the NAMD2.6 MD program [82] were performed of the ATP-bound enzyme [83]. Snapshots of the receptor only (ATP and all water molecules removed) were extracted every 50 ps, generating a total of 400 receptor conformations. The QR-factorization algorithm reduced the initial set of 400 MD structures to 33. In the first round, virtual screen of the NCIDS1 was performed using the static crystal structure and AutoDock 4.0 [73] to dock and rank approximately 1800 compounds. In order to take receptor flexibility into account and to validate and refine the top hits, the top 2% of the screening hits (corresponding to the top thirty compounds), were redocked into the full and QR-reduced holo MD-ensembles. A comparison of the mean RC energy (i.e. mean energy from docking into all 400 snapshots) and the mean QR-RC energy (i.e. mean energy from docking into the representative 33) for the top thirty compounds indicated that redocking into the QR-reduced set was a more efficient way to capture the effects of receptor flexibility without loss of energetic information. In total, 10 compounds were experimentally tested and 5 were found to inhibit at 10 μM or better. Importantly, two of the best inhibiting compounds would not have been selected for testing based on the crystal structure score alone, therefore, further refinement of the predicted compounds based on the RCS reranking provided an important enrichment of the final ranked set.</p><p>These studies clearly highlight the significance of including MD-generated conformations in the virtual screening protocol for discovering new actives. The identification of entirely novel modes of binding in the hit identification stage provides a key advantage over other methods. This is especially true for very flexible receptors, in which entirely new pockets may be revealed for potential ligand binding. The success of these predictions underscores the importance of incorporating receptor flexibility in docking and scoring methods used in conjunction with virtual screening.</p><!><p>One issue facing many of the more intensive physics-based methods is the knowledge gap between industry and academia, as well as computational and experimental academic groups. Despite an increasing number of academic groups working in CADD, there is still a knowledge gap between industrial and academic groups working in the field. Similarly, computational academic groups do not always publish with experimental validation (in this review, we attempted to only cover computational approaches for which experimental verification was available). The interdisciplinary nature of drug discovery begs for a close interaction between all of these groups. It will be important to leverage knowledge from both camps going forward if we are to maximize progress in improving the "state-of-the-art." In this regard, building bridges that allow the field to utilize data and insight from industry in conjunction with new computational and theoretical methodologies being developed in academia will be especially useful. Additionally, computational and theoretical academicians who are pursing these more rigorous and innovative drug discovery approaches must work more closely with experimentalists to validate their methods. More extensive academic and industrial collaborations in the area of drug discovery for neglected infectious diseases will also be progressively more important given the present climate changes and an increasingly mobile world population.</p><!><p>Another current methodological challenge in virtual screening is the inclusion of explicit ordered water molecules. Water is well-known to play an important role in ligand binding, either through the formation of key lubrication contacts between the receptor and ligand, or by being displaced upon ligand binding [84]. A high percentage of high-resolution crystal structures contain resolved water molecules in or near the active site [85], but the decision as to whether or not to retain them in virtual screening experiments is unclear [86]. A recent study sampling the effects of including different water molecules in docking calculations for 24 receptor targets indicated improved enrichment for approximately half of the receptors studied [87]. A complicating consideration is the potential introduction of hydration artifacts due to the cryogenic temperatures used to trap many receptor structures [88]. MD simulations that include experimentally resolved water molecules can provide easily accessible additional information about these water molecules, such as residence times and dipole order parameters; this information could be used to determine whether to include specific water molecules in virtual screening experiments for a given receptor. More accurate energetic approaches, such as calculating the free energy of binding for individual water molecules can also be considered [89]. Future MD-based ensemble approaches will hopefully make use of the structural water information to find novel candidate compounds in virtual screening strategies.</p><!><p>The ability to sample longer timescales presents yet another challenge for traditional all-atom MD simulations, despite continuous increases in computer power. Sampling routinely into the microsecond timescale and beyond will require the use of new and improved techniques. New and promising emerging approaches that address the general problem of more efficiently sampling receptor conformational space are being pursued; two approaches in particular that we would like to highlight that have benefitted recently from a number of applications on more realistic protein systems are accelerated molecular dynamics (AMD) and generalized Born MD (GB MD). Through the use of a robust potential function, AMD increases the transition of high-energy barriers without any advance knowledge of the underlying energy landscape [90]. This promising method has been shown to sample slow diffusive conformational transitions [91] as well as speed rates of convergence for entropic calculations [92]. Importantly, experimental validation of accurate, enhanced conformational sampling through AMD has been carried out by NMR studies [93]. We anticipate that the structures sampled with this far-reaching technique will facilitate the discovery of new structural understandings that will lead to the discovery of novel therapeutics, and also potentially provide more accurate estimates of the thermodynamics of binding. A second technique to pave the way into the microsecond regime that is ready to be applied to larger biomolecular systems, is generalized Born molecular dynamics (GB MD). GB MD, which employs a continuum representation of the solvent water molecules and salt ions, offers computational efficiency over the more rigorous Poisson-Boltzmann solvers or explicit solvent simulations [94]. The implicit solvent approach essentially reduces solvent friction and has been shown to enhance conformational sampling for a wide variety of peptide, protein, and nucleic acid systems [95–97], and more recently, larger protein systems [98, 99]. Correspondence between GB models and experimentally derived structures has also been established [100, 101], which, similar to AMD, makes it a good choice for extensions of the computational receptor ensemble.</p><!><p>A key question regarding the MD-based approaches for conformational ensemble generation is determining for how long one should simulate. Unfortunately there does not yet appear to be a clear-cut answer to this question. Current successful published studies range from 500 ps to 320 ns in length, which is quite a broad range of timescales. Each of the studies presented some measure of success, mainly the identification of compounds that otherwise would have been missed in static receptor models. One could determine enrichment factors as a function of simulation length for a particular receptor, although it would be difficult to generalize the findings to all possible receptors since different proteins exhibit varied structural dynamics. Moreover, rare events that occur on longer-time scales will only be sampled if the simulations are extended beyond the requisite time point. Including conformational ensembles generated with accelerated sampling approaches, such as AMD or GB MD, may shed light on this aspect. For now, the only clear answer seems to be that including receptor flexibility to some degree, even if only explored for 1 ns of dynamics, is better than a purely rigid model.</p><!><p>Although in virtual screening experiments success is more loosely defined as being able to separate the binders from non-binders, a challenge intimately connected to virtual screening is the development of methods to predict the binding affinities for a series of compounds to a target receptor (reviewed in, among others: [102–107]). Most current methods apply the indirect approach first suggested by Tempe and McCammon, which utilizes a simple thermodynamic cycle for the bound and unbound receptor-ligand complex [108]. The most rigorous treatment of receptor-ligand binding that presently exist are the so-called "computational alchemy" techniques achieved through MD simulations, in which free energy changes are estimated based on coupling-parameter approaches, such as free energy perturbation or thermodynamic integration [109]. These methods describe a higher physical complexity of the binding process and include an extensive sampling of receptor, ligand, and solvent phase spaces 1[07, 110–119]. Importantly, these methods provide more accurate estimates of binding free energies as well as reliable measures of their accuracy. Although they suffer from extraordinary computational costs that essentially prohibit their mainstream application, this too, is changing. More frequently these are being applied to larger ligand sets, and although computational power has not advanced to the point of being able to apply these methodologies to entire ligand libraries, they are being used to handle several dozens of candidate compounds [120]. Additionally, these algorithms now exist as functional modules in the massively parallelizable MD programs Desmond (TI) [121] and NAMD (FEP) [114], allowing the determination of accurate free energies of binding in equivalent wall-clock time of hours-to-days, whereas serial implementations of these algorithms traditionally require(d) weeks-to-months of computing time for a single compound. The utilization of petascale computing infrastructures that will soon become publicly available will continue to drive down the cost of such calculations. Faster and more approximate free-energy based methods, such as molecular mechanics Poisson Boltzmann (generalized Born) surface area (MM-PB(GB)SA) [122, 123], linear interaction energy [124, 125], single-step perturbation [126–129], are less computationally intensive and thus already more easily extendible to larger compound libraries.</p><!><p>The increasing availability of computing power on the order of peta- to exa-FLOPS (floating point operations per second) will make possible longer time course simulations (Fig. 3). Currently, a 100 ns explicit solvent simulation of a 75 kDa protein can be completed in one week on a supercomputer or a local cluster with low latency network. More rigorous physics-based approaches are increasingly possible. The arrival of supercomputers such as the NCSA/IBM Blue Waters with many-core architecture will enable simulations on the order of 1~100 μs on a routine basis. The advance in general purpose graphics processing unit (GPGPU) based processing power, where each GPU contains as many as 512 cores, means that the age of personal supercomputer has dawned upon us. While GPGPU is still challenging to program, it has been used successfully to achieve speedups of tomographic reconstruction codes such as TxBR [130] or molecular dynamics codes such as NAMD [131]. Emerging standards such as OpenCL (open computing language) will likely help make GPGPU or cGPU (integrated CPU/GPU) even more powerful, lower the cost of use, and increase in popularity. What used to take weeks to complete on a TFLOPS (tera-FLOPS) supercomputer could be accomplished in days on a workstation in the near future. More dedicated application and algorithm specific hardware acceleration for molecular dynamics codes has achieved millisecond scale simulations [132]. Although current force fields have yet to be rigorously tested over such longer time scale simulations, pushing time and length scale limits coupled with experimental validation will help identify critical areas for new methodological development.</p><p>Cloud computing refers to the on-demand availability of thousands of processors offered by commercial service providers. Amazon Elastic Cloud 2 (EC2) is offering computation for approximately 10 cents per hour, with enormous impact on virtual screening experiments. Virtual screening services offered as web services are now accessible for interactive use or batch processing of entire library of compounds with transparent access to cloud or cluster resources [133]. The RCS will be available as a computational workflow that will be available to a much wider audience, easily accessible to experimental biologists [134]. Commercial offerings such as the Accelrys Pipeline Pilot tailors to in house chemical compound libraries, whereas the Vision-like environment provides a complete solution to academic researchers in computer aided drug design. Increasingly, both academic and commercial entities will rely on cloud-based virtual screening services that are entirely transparent to the end users. Activities such as the World Community Grid utilize idling computers provides another venue for virtual screening services [135]. From a practical perspective, a brute force approach to cross-dock every resulting structure to every available small molecule compound is counterproductive, regardless of the amount of computational power. Better theories and more efficient computational schemes to allow the selection of conformations most relevant to ligand binding in a predictive manner are still needed. New algorithms such as AutoDock Vina present 10 to 100 fold speed up of AutoDock with equal or better performance [136]. In short, the continued acceleration through algorithmic and hardware improvement, and cloud-based utility computing will make ensemble based virtual screening accessible to a wider audience, and further its validation and optimization.</p><!><p>More rigorous molecular dynamics based approaches to include full receptor flexibility are now being applied in larger-scale virtual screening applications, thanks in part to steady increases in computational power and the development of new and improved methods improving computational efficiency of the underlying strategies. The studies highlighted in this review indicate that these ensemble-based methods can yield material discoveries regardless of the choice of force field or MD program. The general correspondence exhibited among the various strategies is a promising feature and indicative of the quality of many of these methods. MD-generated ensemble-based techniques have the potential to provide breakthrough discoveries for many target receptors, not only through the discovery of potential new ligand binding areas, but also through more accurate estimates of free energies of binding for increasingly larger ligand sets.</p>