Essential Dynamics Ensemble Docking for Structure-Based GPCR Drug Discovery

Abstract

The lack of biologically relevant protein structures can hinder rational design of small molecules to target G protein-coupled receptors (GPCRs). While ensemble docking using multiple models of the protein target is a promising technique for structure-based drug discovery, model clustering and selection still need further investigations to achieve both high accuracy and efficiency. In this work, we have developed an original ensemble docking approach, which identifies the most relevant conformations based on the essential dynamics of the protein pocket. This approach is applied to the study of small-molecule antagonists for the PAC1 receptor, a class B GPCR and a regulator of stress. As few as four representative PAC1 models are selected from simulations of a homology model and then used to screen three million compounds from the ZINC database and 23 experimentally validated compounds for PAC1 targeting. Our essential dynamics ensemble docking (EDED) approach can effectively reduce the number of false negatives in virtual screening and improve the accuracy to seek potent compounds. Given the cost and difficulties to determine membrane protein structures for all the relevant states, our methodology can be useful for future discovery of small molecules to target more other GPCRs, either with or without experimental structures.

Keywords: PAC1 receptor; antagonist; computer aided drug design; molecular dynamics; principal component analysis; virtual screening.

FIGURE 1

Cartoon illustrations of the PACAP-bound PAC1R model (PDBID: 6M1I, PAC1very short) and our homology model (template PDBID: 4L6R, PAC1null, simulation snapshot at 500 ns). The PAC1null isoform is more biomedically relevant than the very short isoform. The PACAP peptide is shown as a helix cartoon (pale green); the 21-amino acid ECD insert (see the sequence in Supplementary Figure S1) is shown as a flexible coil (purple). This study focused on docking to the peptide-binding pocket.

FIGURE 2

Four representative PAC1R conformations using in EDED reveal important changes in the binding pocket. The protein is shown in cartoon and the reference ligand is shown as sticks. The histogram of all trajectory frames projected onto the first two principal components of residues within the ligand-binding pocket of PAC1R. Black dots labelled with numbers from 0 to 3 are the representative structures (S0, S1, S2, and S3) determined by the K-means clustering algorithm.

FIGURE 3

Violin plots of the docking score distribution of the top 350 compounds to different receptor models. The dash line shows the –9.0 kcal/mol cutoff used to prioritize compounds for synthesis.

FIGURE 4

Ensemble weighted glide scores ( $\text{Δ}$ G_bind) of 23 experimentally tested compounds. Compounds with strong, modest, and poor ERK inhibitive activity are depicted in green, blue, and red, respectively. Corresponding colored lines represent the average ensemble weighted glide score for that category. A cutoff of –9 kcal/mol was applied for predicted antagonists to be compared to their experimental results showing either strong or medium inhibition (active) or weak inhibition (inactive).

FIGURE 5

Overview of computational workflow for development of PAC1R antagonists. Right Column: selection of input ligands from a structure database (in this example the ZINC15 (Sterling and Irwin, 2015) database). Custom filters were used to select raw structures with desirable properties (molecular weight, logP, etc.). These structures are then prepared using Schrödinger’s ligprep software program. Left column): the PAC1null homology model is constructed from the protein’s sequence, simulated for 500 ns, and the raw coordinates are analyzed. The representative structures are used in ensemble docking. Hit compounds are selected based on visual inspection of the results.