Performance of virtual screening against GPCR homology models: Impact of template selection and treatment of binding site plasticity

Abstract

Rational drug design for G protein-coupled receptors (GPCRs) is limited by the small number of available atomic resolution structures. We assessed the use of homology modeling to predict the structures of two therapeutically relevant GPCRs and strategies to improve the performance of virtual screening against modeled binding sites. Homology models of the D2 dopamine (D2R) and serotonin 5-HT2A receptors (5-HT2AR) were generated based on crystal structures of 16 different GPCRs. Comparison of the homology models to D2R and 5-HT2AR crystal structures showed that accurate predictions could be obtained, but not necessarily using the most closely related template. Assessment of virtual screening performance was based on molecular docking of ligands and decoys. The results demonstrated that several templates and multiple models based on each of these must be evaluated to identify the optimal binding site structure. Models based on aminergic GPCRs showed substantial ligand enrichment and there was a trend toward improved virtual screening performance with increasing binding site accuracy. The best models even yielded ligand enrichment comparable to or better than that of the D2R and 5-HT2AR crystal structures. Methods to consider binding site plasticity were explored to further improve predictions. Molecular docking to ensembles of structures did not outperform the best individual binding site models, but could increase the diversity of hits from virtual screens and be advantageous for GPCR targets with few known ligands. Molecular dynamics refinement resulted in moderate improvements of structural accuracy and the virtual screening performance of snapshots was either comparable to or worse than that of the raw homology models. These results provide guidelines for successful application of structure-based ligand discovery using GPCR homology models.

Fig 1. Structural accuracy of homology models based on templates with different levels of sequence identity.

The average RMSD_TMBB (green and blue bars), RMSD_BSBB (dark green and blue bars), and RMSD_BSSC (grey bars) to the crystal structures were calculated for 50 models per template. Bar charts of structural accuracy using templates with different levels of sequence identity was calculated for (A) the D₂R and (B) 5-HT_2AR.

Fig 2. Relation between structural accuracy and sequence identity.

The average RMSD_TMBB (A and D), RMSD_BSBB (B and E), and RMSD_BSSC (C and F) to the crystal structures for 50 models of the D₂R (A-C) and 5-HT_2AR (D-F) based on different templates. The solid line represents a linear regression and R is Pearson’s correlation coefficient.

Fig 3. Ligand enrichment by homology models.

Distribution of aLogAUC values for (A) D₂R and (B) 5-HT_2AR models based on 16 different templates. Distributions are shown using a boxplot representation. Each boxplot describes the results for 50 models obtained based on one template. The box represents the 50th percentile of the data and the black band shows the median value. The lowest and highest aLogAUC values are represented by the whiskers. The gray boxplot corresponds the results of all 600 models based on aminergic receptor templates. Ensemble enrichments are represented by red filled circles.

Fig 4. Relation between ligand enrichment and sequence identity.

The median aLogAUC values of the D₂R (A-B) and 5-HT_2AR (C-D) homology models based on aminergic templates with different TM (A and C) or BS (B and D) sequence identities. The solid line represents a linear regression and R is Pearson’s correlation coefficient.

Fig 5. Relation between ligand enrichment and structural accuracy.

The median aLogAUC and average RMSD_BSSC to the crystal structures for 50 homology models of the D₂R (A) and 5-HT_2AR (B) based on aminergic templates. The solid line represents a linear regression and R is Pearson’s correlation coefficient.

Fig 6. Influence of template on enrichment of ligand chemotypes.

Enrichment (aLogAUC) of eticlopride- (blue bars), doxepin- (yellow bars), piperidine/piperazine-like (red bars), and all (grey bars) D₂R ligands by D₂R homology models. Homology models based on three different templates (D₃R, H₁R, and 5-HT_2CR) and ensemble enrichments of the D₃R model combined with either the H₁R and 5-HT_2CR models were evaluated.

Fig 7. Structural accuracy of MD simulation snapshots.

RMSD distribution of the TM backbone, BS backbone, and BS side chains of MD snapshots to the crystal structure for the D₂R (A-C) and 5-HT_2AR (D-F). Distributions of RMSD values for the three sets of snapshots based on the Rho-based models (MD_Rho/Apo) and homology models based on the most closely related template in apo (MD_Template/Apo) and holo forms (MD_{Template/Holo}) are shown using a boxplot representation. The box represents the 50th percentile of the data and the black band shows the median value. The lowest and highest RMSD values are represented by the whiskers. The horizontal lines show the RMSD values of the homology model used as starting structure.

Fig 8. Ligand enrichment by MD simulation snapshots.

Distributions of aLogAUC values for 50 homology models (HM) and sets of 50 snapshots from three MD simulations per template for the (A) D₂R based on the D₃R template and (B) 5-HT_2AR based on the 5-HT_2CR template. Results for HM, apo (MD_APO), and holo forms (MD_HOLO) are shown using a boxplot representation. The box represents the 50th percentile of the data and the black band shows the median value. The lowest and highest alogAUC values are represented by the whiskers. The ensemble enrichments of each set of 50 models are represented by red filled circles.