Assessing hERG pore models as templates for drug docking using published experimental constraints: the inactivated state in the context of drug block.

Many structurally and therapeutically diverse drugs interact with the human heart K+ channel hERG by binding within the K+ permeation pathway of the open channel, leading to drug-induced 'long QT syndrome'. Drug binding to hERG is often stabilized by inactivation gating. In the absence of a crystal structure, hERG pore homology models have been used to characterize drug interactions. Here we assess potentially inactivated states of the bacterial K+ channel, KcsA, as templates for inactivated state hERG pore models in the context of drug binding using computational docking. Although Flexidock and GOLD docking produced low energy score poses in the models tested, each method selected a MthK K+ channel-based model over models based on the putative inactivated state KcsA structures for each of the 9 drugs tested. The variety of docking poses found indicates that an optimal arrangement for drug binding of aromatic side chains in the hERG pore can be achieved in several different configurations. This plasticity of the drug "binding site" is likely to be a feature of the hERG inactivated state. The results demonstrate that experimental data on specific drug interactions can be used as structural constraints to assess and refine hERG homology models.

Introduction

The effects of most drugs result from binding to target or off-target protein. Characterization of drug binding can provide insight into productive strategies for improving therapeutics and minimizing side effects. Since the physical principles underlying binding are increasingly well understood, drug binding is amenable to computational approaches in which binding sites, drug binding poses, and binding affinities should be accessible to calculation via computational docking.[1,2]

The situation is complicated when an atomic resolution structure of the protein is unavailable. In these cases docking analyses with homology models built on structurally defined templates provide a means of computational assessment of the docking problem. The complexity is increased when the protein can access multiple conformational states and drug binding is state-dependent. The human Ether-à-go-go Related Gene (hERG) product provides one such example. The hERG K+ channel carries the rapid delayed rectifier repolarizing current (IKr) which controls ventricular action repolarization and, thereby, the duration of the QT interval in humans.[3,4] This function is mediated by rapid channel inactivation following channel opening upon membrane depolarization, followed by rapid recovery from inactivation and slow channel closing (deactivation) at repolarizing membrane potentials. HERG is of intense pharmacological interest due to the variety of cardiac and noncardiac drugs that block the channel with potentially fatal consequences.[5,6] The scale of this problem is illustrated by the fact that novel drug candidates are routinely screened against hERG as a key part of the drug development/safety process.[7] Studies on the state-dependence of hERG block indicate that in many cases, especially involving high affinity blockers, the drug binds more strongly as a consequence of inactivation.[8−12] The molecular basis of drug block, including the conformation of the hERG pore in high affinity drug binding states, the nature of the drug binding surface, and the conformations and interactions (poses) of bound drugs, remains poorly defined.

Despite the absence of a hERG crystal structure several structures of K+ channels have been determined, and these provide potential templates for constructing homology models to assess drug binding. On the basis of sequence homologies, especially involving the helix (S6 in hERG) that lines the K+ permeation pathway, most hERG models have been constructed upon the crystal structure templates of KcsA,[13,14] MthK,[15] and KvAP.[16] These models are either unmodified from their templates (e.g., KcsA;[17−19] MthK;[20,21] KvAP[22]), or the model is adjusted after construction based on the expected accessibilities of side chains facing the hERG pore cavity (e.g., KcsA;[23] MthK;[24] KvAP[25]). In some cases models constructed on crystal structure templates have been “relaxed” by molecular dynamics simulation so that the relationship between the hERG model and its template structure is less well-defined (e.g., KcsA;[26] KvAP[27]).

The effects of mutagenesis on drug block in hERG, especially using alanine-scanning mutagenesis,[17] strongly implicate two amino acid residues on the S6 helix, Y652 and F656, whose side chains are expected to be accessible to drug molecules that enter the channel pore from the cytoplasmic side of the membrane when the channel opens.[3−7,9,17] These residues provide 8 aromatic side chains that can make potential hydrophobic, π–π stacking and cation-π interactions with complementary moieties on drug molecules. Alanine replacement of several other amino acids also attenuates drug block in many cases. These include T623 and S624 in a short sequence facing the pore where the pore helix turns into the selectivity filter, V625 within the selectivity filter, G648, V659, and S660 on the S6 helix, and a set of residues including S620 (S620T) and S631 situated above the hERG pore cavity.[17,28−30] These latter residues, and especially N588 (N588K) which lies in a cytoplasmic loop between the top of helix S5 and the pore helix, are unlikely to interact directly with drugs in the open channel pore.[8,11,29] The effects on drug block of mutation of some of these residues (especially N588K[11] and S631A[28] but also S620T[8] and V625A[29,30]) probably arise from attenuation of channel inactivation, with an indirect effect on block of drugs that bind more strongly to the inactivated state.

Crystal structures of putative open inactivated states of the KcsA channel were recently published (PDB:3F5W; PDB:3F7V),[14] and these provide potential templates for hERG pore models. Indeed, one of these structures (3F5W) was used as a template for modeling of the pore region of a full hERG open inactivated state model and docking of hERG blockers.[26] However, there is no a priori reason to expect that the backbone conformation of a putative inactivated state KcsA structure corresponds to the conformation of hERG open-inactivated states that is relevant for drug block. In this paper we address the extent to which hERG homology models can be assessed, and either validated or rejected, in the context of experimental constraints that indirectly define the nature of drug interactions with pore residues. In particular, we focus on hERG pore homology models built onto the putative open-inactivated KcsA structures[14] and compare these with a MthK model that has been used successfully to interpret the effects of pore residue mutation on drug block of hERG.[20,21] We make use of several recent key experimental observations that inform on the essential structural elements of hERG blockers and the likely interactions of these molecules with aromatic residues within the hERG pore.[24,31−33] This analysis gives insight into the nature of the inactivated state of hERG in the context of drug binding.

Methods

Homology Modeling and Model Building

Construction of the hERG pore model built onto the MthK (PDB:1LNQ) template has been described.[20,21] hERG pore models built on putative inactivated state KcsA structure templates (PDB:3F5W and PDB:3F7V)[14] were constructed using the alignment in Figure 1 which defines the parts of the structures (pore helix, selectivity filter, and S6 helix) built into the hERG models. The alignment of sequence corresponding to the S5 helix of hERG onto crystal structure template sequences remains poorly defined[34] and was not included in the models; drugs in the hERG pore cavity are not expected to interact with S5 which does not constitute part of the hERG pore lining. The models were constructed in Insight II (Accelrys, San Diego, CA, USA) by replacing amino acids in the template with the appropriate amino acid in hERG according to the homologies in Figure 1. Side chain rotamers were selected to match the rotamers of the template structure up to the side chain β-atoms. Steric clashes in the models were relieved by reselecting side-chain rotamers having the lowest energies and using additional small manual bond rotations where necessary. In all models K+ ions were inserted into the [1] and [3] positions of the selectivity filter, unless otherwise stated. The models were energy-minimized using 2000 steps of steepest descents using the Discover module of InsightII.

Figure 1

Sequence alignments used to construct hERG pore homology models. Amino acids T632, S624, Y652, and F656 are colored according to the scheme used in structure figures throughout.

In the following text hERG models are defined in terms of the template channel and its PDB accession code (i.e., the hERG pore model built on the MthK, PDB:1LNQ structure template is the “MthK(1LNQ) model”). An overview of the three hERG pore models is shown in Figure 2.

Figure 2

Spatial arrangement of Y652 (pink) and F656 (blue) side chains in hERG pore homology models built onto the crystal structures of MthK (1LNQ) and putative inactivated states of KcsA (3F5W) and (3F7V). Each of the 5 side chain rotamers sampled in GOLD docking runs is shown. This does not represent the full side chain flexibility sampled in GOLD since each bond rotation within each rotamer samples an additional range of angles (of between 10 and 20 o around the specified rotamer torsion angle).[35,36] Flexidock samples a somewhat extended set of side chain conformations since free torsional rotation is allowed for Y652 and F656 side chains.

Computational Docking

The drugs used in this study are shown in Figure 3, and their properties as hERG blockers are compiled in Table 1. The drugs were chosen on the basis of recent structure–activity data that informs on likely arrangements of amino acid side chains in the hERG pore involved in drug binding,[24,31−33] a high representation of drugs that bind preferentially to the inactivated state (Table 1), and a range of IC50 values allowing comparison of hERG block potency with docking parameters that rank drug poses (Table 1). Terfenadine and cisapride were included as examples of “classical” hERG blockers which have been withdrawn from the market due to dangerous side effects resulting from hERG block.[7] Drug molecules were constructed using InsightII. Partial charges were calculated using the Gasteiger–Huckel module in Sybyl 2.0 (Tripos, St. Louis, MO, USA), and the molecules were energy minimized by 5000 steps of steepest descents.

Figure 3

hERG blockers used in docking analysis. For convenience these molecules are referred to as “drugs” throughout although not all of the molecules are prescription drugs. The amino group carrying a positive charge in each drug is indicated with a “+”.

Table 1

Drug Molecules Used in Docking Analysis and Selected Propertiesa

drug	IC50 (ref)	Y652/F656 (ref)	inactivation (ref)
S-bupivacaine	13 μM ([56])	yes ([57])	yes ([58])
cavalli-2	17 nM ([33])	n.a.	n.a.
cavalli-6	2.4 nM ([33])	n.a.	n.a.
cisapride	21 nM ([59])	yes ([17][53])	yes ([48])
chloroquine	2.2 μM ([60])	yes ([62])	yes (weak) ([61])
dofetilide	5 nM ([63])	yes ([52])	yes ([12])
E-4031	13 nM ([64])	yes ([52])	yes ([10])
haloperidol	27 nM ([63])	yes ([47])	yes ([46])
terfenadine	7 nM ([65])	yes ([17], (53))	yes ([12])

The IC50 values are taken from the references cited. All IC50 data were obtained by patch clamp at 34–37 °C of human embryonic kidney (HEK) or Chinese hamster ovary (CHO) cells overexpressing hERG. The table also represents whether IC50 values are reduced by attenuation of inactivation (“inactivation”) and in Y652A and F656A hERG mutants (“Y652/F656”). n.a., data not available.

Flexidock

The Flexidock module of Sybyl 2.0 was run using the default parameter set as previously described.[21] Free side chain flexibility was sampled during docking for the following residues: T623, S624, V625, Y652, F656, and S660. 60,000 generations of the genetic algorithm were used in each run since this number was found to yield energy score convergence. Flexidock outputs tend to be somewhat biased by the starting structure; therefore, each docking analysis was repeated for a total of 240 runs with drug positioned in sets of different starting positions and orientations in the channel. The computational time was around 16 h for each specific drug-model analysis on a single core of a 3.3 GHz Intel ib-2120 processor.

GOLD

In addition to the side chain hydroxyl rotational flexibility utilized in GOLD,[35] side chain flexibility of Y652 and F656 was also sampled. The GOLD rotamer library for these residues[35,36] was supplemented with an additional rotamer for each of Y652 [chi1 −175 (15); chi2 −100 (15)] and F656 [chi1 −60 (15); chi2 80 (15) where the number in brackets specifies the torsion angle range above and below the stated rotamer angle], to incorporate the side chain rotamers observed in a recent simulation that addressed side chain conformers of Y652 and F656 in a hERG model.[37] Due to the large number of rotamers sampled 300,000 steps of the genetic algorithm were used, and 40 docking runs were sampled for each drug. The computational time was around 40 min for each specific drug-model analysis on a single core of a 3.4 GHz Intel i7-3770 processor. Docking series were run twice to sample docking using both Chemscore[38] and ChemPLP[39] scoring functions.

Analysis

GOLD (Chemscore), GOLD (ChemPLP), and Flexidock series were run for each of the 9 drugs in each of the three hERG pore models resulting in a total of 81 (9 × 3 × 3) sets of docking outputs. These outputs were analyzed using Pymol (Schrödinger). For Flexidock runs the 5 lowest energy score poses were inspected to determine the number of interactions between drug and hERG side chains using the criteria compiled in Table 2. A π-stacking interaction was counted if at least two ring carbons on opposite sides of the ring 2-fold axis were within 4.5 Å of partner ring carbons of the adjacent ring (parallel π–π stacking), or if a ring proton was within 4.5 Å of the center of an adjacent ring (T-shaped stacking). A cation-π interaction was counted if the basic N atom of the drug or a proton on a C atom directly bonded to the charged N was within 4 Å of the center of an aromatic ring. Hydrogen bonds were counted if the H atom was within 2.5 Å of the acceptor hydroxyl oxygen. An interaction of the positively charged tertiary nitrogen of the drug with the negative electrostatic field arising from focusing of the pore helix C-terminal dipole charges was counted if the drug nitrogen was within 3 Å of the pore cavity binding site for a K+ ion observed in crystal structures of homologous potassium channels.

Table 2

Distance Criteria Used in Defining Specific Interactions between Drug and hERG Model Pore Residues in Low Energy Score Docked Poses

interaction	distance criteria
parallel π–π	ring centers within 4.5 Å
T-shaped π–π	aromatic H – ring center distance within 4 Å
cation-π	protonated N or adjacent CH within 4 Å of aromatic ring center
H-bond	OH-X or NH-X within 2.5 Å (X is H-bond acceptor)
cation in K+ site	protonated N atom within 3 Å of cavity K+ site

Results

The “Results” are organized as follows: we first assessed Flexidock and GOLD docking outputs for drug-like molecules in the context of K+ channel pores, by docking tetraethylammonium (TEA) and tetrabutylammonium (TBA) into crystal structures of KcsA, in the case of TBA using the structure (KcsA (PDB: 2BOB) from which the molecule was extracted. This analysis provides some validation of these programs for positively charged drug binding to channel pores and yields insight into the nature of interactions that dominate binding poses for the respective programs. In the second section we use GOLD and Flexidock to characterize drug docking into hERG pore homology models and assess their ability to select models that are most likely to correspond to hERG pore conformations compatible with experimental data on drug block. The third section of the Results describes the nature of low energy score docking poses for several drugs in the context of published experimental constraints on binding interactions.

Flexidock and GOLD Docking of Tetraalkylammoniums into KcsA Crystal Structure

Flexidock uses a cut down version of the Tripos force field for optimizing docking poses,[40] whereas GOLD uses empirical fitness functions based largely on optimization of van der Waals surface interactions and hydrogen bonds as selection elements for assessing fitness.[35] These docking protocols are therefore not expected necessarily to give equivalent low energy score docking poses. This is particularly the case for drug block of hERG where the drug binding interactions in the pore cavity may be heterogeneous; in addition a dominant feature of most hERG pore blockers is a secondary, tertiary, or quaternary aliphatic ammonium (Figure 3) that makes a significant contribution to drug binding. The electrostatic component of this contribution is represented in force-field-directed docking but not in empirical protocols unless specifically parametrized. To assess the ability of Flexidock and GOLD to dock molecules of this nature in K+ channel structures we used one of the very limited set of crystal structures of K+ channels with bound positively charged blockers, namely the KcsA crystal structures containing the bound tetraalkylammonium ion, TBA.[41,42]

Figure 4 shows representative low energy score poses for TBA docked into the KcsA crystal structure with potassium ions in the S1 and S3 positions [S1, S3] of the selectivity filter. In both Flexidock and GOLD runs the side chains of I100 and F103, which dominate interactions with TBA in the crystal structure, were allowed full rotamer sampling during docking, and Chemscore was used to rank GOLD poses. Both Flexidock and GOLD position the central ammonium group on the channel axis at a height around 2 Ǻ above the K+ ion binding site effectively overlaying the position of the ammonium nitrogen of TBA in the KcsA crystal structure from which TBA was removed. This position is proposed to be a dehydration transition site for permeant ions.[41] TBA is known to adopt two low energy conformations denoted D2d (4 butyl chains extended in plane) and S4 (tetrahedral).[43] While the electron density for TBA in the closed state KcsA crystal structure is consistent with a planar (D2d) conformation, a significant crystal structure electron density 5–6 Å beneath the ammonium ion position indicates a likely contribution from an S4 conformation similar to that observed for the low energy score TBA Flexidock poses (Figure 4A).[42] In any case both docking programs have localized TBA into the binding site identified in the crystal structure.

Figure 4

Lowest energy score docked outputs for TBA (A, B) and TEA (C, D) using Flexidock (A, C) and GOLD (B, D) docked into KcsA (PDB:2BOB for TBA; PDB:2BOC for TEA). In each case the docked small molecule structure is represented by yellow sticks, and the crystal structure coordinates for TBA are represented by blue sticks. The blue TEA structure in panels C and D was made by editing the TBA coordinates to truncate the butyl chains to ethyl chains. In all runs potassium ions occupied the [1] and [3] positions of the selectivity filter; the [3]K+ ion is shown as a purple sphere.

There is no crystal structure for the bound state for tetraethylammonium (TEA) in the internal site of K+ channels. However TEA is proposed to occupy the same site as TBA; i.e. below the selectivity filter at a position equivalent to the putative dehydration site for permeant ions.[41] Similar to the TBA docking poses, Flexidock localizes TEA near the dehydration site close to that found for TBA (Figure 4C). Since TEA in these poses lies on the channel pore axis and makes few interactions with side chains in the channel pore compared to TBA, the electrostatic contribution to binding, resulting largely from the focused helix dipole charges from the pore helices, dominates binding “energetics” in this model. On the other hand, GOLD maximizes direct interactions between TEA and protein side chains, and a variety of off-axis low energy score poses are obtained (e.g., Figure 4D). These results support the expectation that a force-field based docking analysis is more likely to represent binding modes where electrostatic interactions make a strong contribution to the binding energy.

Assessing hERG Pore Models Built on Putative Inactivated-State KcsA Structure Templates

Successful docking of TBA into its crystal structure binding site, and the representation of electrostatic contributions to binding of TEA using Flexidock, indicates that these programs are likely to produce useful docking output when hERG binding drugs are docked into hERG pore models. We used GOLD and Flexidock to dock each of the 9 drugs in Figure 3 into three hERG pore models. We used a MthK(1LNQ) model which has previously produced docking poses that accord with experimental block of hERG, and specific alanine mutants, for several hERG blockers.[20,21] We assessed two pore models that were constructed on the open and potentially inactivated state KcsA crystal structures[14] (3F5W and 3F7V).

Tables 3 and 4 compile the GOLD docking output, using Chemscore and ChemPLP to rank docking poses. For both Chemscore and ChemPLP the MthK(1LNQ) model outperforms the KcsA(3F5W) and KcsA(3F7V) models for each of the 9 drugs. As might be expected from visual inspection (Figure 2), the arrangement of aromatic side chains in the KcsA(3F7V) model precludes an extensive set of interactions with drug molecules, and this model scores particularly poorly in GOLD runs. However, the results indicate that the MthK(1LNQ) model also affords a more favorable set of interactions with each of the drugs compared to the KcsA(3F5W) model.

Table 3

Best Energy Score for Docking of Drugs into Specified Homology Models Using GOLD and Scoring Docking Poses with ChemPLP

	ChemPLP
drug	MthK(1LNQ)	3F5W	3F7V
S-bupivacaine	73.3	56.8	11.2
cavalli-2	90.1	74.7	26.8
cavalli-6	99.0	87.1	32.1
cisapride	88.4	72.5	18.9
chloroquine	78.1	63.1	21.4
dofetilide	77.3	65.7	22.1
E-4031	79.8	66.5	23.3
haloperidol	83.8	71.9	19.3
terfenadine	97.8	84.4	34.5

Table 4

Best Energy Score for Docking of Drugs into Specified Homology Models Using GOLD and Scoring Docking Poses with Chemscore

	Chemscore
drug	MthK(1LNQ)	3F5W	3F7V
S-bupivacaine	33.5	26.9	–14.4
cavalli-2	44.4	37.8	–4.8
cavalli-6	41.5	35.2	–7.1
cisapride	32.7	30.3	–12.8
chloroquine	34.8	30.6	–10.9
dofetilide	30.5	27.4	–12.8
E-4031	35.2	32.3	–12.2
haloperidol	37.7	33.4	–13.8
terfenidine	51.9	44.1	–1.5

While Flexidock energy scores rank docking poses within a set of runs comprising a particular model and ligand, we have found that comparison of energy scores across a set of different models does not provide reliable ranking of models; this may be because the energy scores contain the internal contribution from the model, and these may outweigh the contributions arising from binding contributions and the internal energy of the bound drug. Instead we assessed the Flexidock output by summing the interactions between drug and hERG model as defined by the criteria in Table 2, for the 5 “best” output docking poses of each docking run according to the Flexidock energy score. Since these interactions are implicitly parametrized within the Tripos force field, this approach allows binding poses to be characterized in terms of a set of specific interactions that can be interpreted in the context of published residue-specific effects of alanine-replacement on drug block. Summations of the interactions for each drug-model set are illustrated in Figure 5. Since alanine-scan effects on drug block, where available, allows independent assessment of the contributions of Y652 and F656 to drug binding we compiled drug interactions involving Y652 and F656 as separate groups. Other interactions comprise hydrogen bonds and the location of the protonated secondary nitrogen atom of the drug in or near the cavity binding site or dehydration site for a K+ ion. This analysis also indicates that each of the 9 drug molecules tested makes more extensive interactions within the MthK(1LNQ) model compared to the hERG models built on the putative inactivated-state KcsA structures.

Figure 5

Summed interactions involved in Flexidock docking of drugs to hERG models built on MthK (1LNQ) and putative inactivated state KcsA (3F5W; 3F7V) crystal structure templates. Interactions from the 5 most favorable energy score output structures were summed. Brown bars: Y652 interactions comprising π-stacking, cation-π and H-bond interactions involving the phenolic hydroxyl group; Orange bars: F656 interactions comprising π-stacking and cation-π interactions; Yellow bars: other interactions comprising the location of the drug protonated amino group in or near the hERG cavity K+ binding site, and hydrogen bond interactions largely involving the side chain hydroxyl group of S624. The totals (green bars) were summed over all specified interactions (as defined in Table 2).

Within the MthK(1LNQ) model both GOLD and Flexidock “scoring” broadly rank computational binding efficacy in line with the efficacy of the drugs as hERG blockers (Figure 6). Assessing docking fitness using both Chemscore and ChemPLP (Tables 3 and 4; Figure 6), the two methanesulfonamide drugs dofetilide and E-4031 “underscore” when compared with their IC50 for drug block. Despite the poor energy scores these molecules appear to make substantial interactions with side chains in the MthK model in GOLD runs (Figures 7 and S1). We emphasize that strong relationships between docking scores and drug block efficacy are not necessarily expected, first because docking methods like GOLD are optimized to dock into crystallographically defined binding sites,[44] and second because IC50 values are nonoptimal measures of drug-receptor interactions (Kd values are preferable). We also emphasize that our analysis does not constitute an assessment of the relative merits of GOLD and Flexidock for docking hERG blockers into homology models;[45] such a comparison would be inappropriate as carried out here, among other reasons because 10–20 times the computational resources were utilized for Flexidock runs compared with GOLD runs. Instead it is the general concordance between the GOLD and Flexidock runs (Tables 3 and 4; Figure 5) that lends confidence that these methods reliably select MthK(1LNQ) over models built on putative inactivated state KcsA structures, as a preferred hERG pore model in the context of drug binding.

Figure 6

Relationship between docking “scores” and hERG channel blocking activity for the drug set described in Table 1 and Figure 3. The CHMPLP score (panel A) is from the “best” pose from GOLD docking. Panel B is the average of the number of interactions as defined in Table 2 and compiled in Figure 5, for the “best” five Flexidock docked poses. The off-diagonal outliers dofetilide (1.) and E-4031 (2.) in panel A were not used in calculating the dotted regression line which has no specific theoretical significance.

Figure 7

Representative low-energy score poses for E-4031 docked into the MthK (1LNQ) hERG pore model using GOLD (CHMPLP scoring) (panel A) and Flexidock (panel B), respectively. E-4031 is represented by yellow bonds. Y652 (mauve), F656 (blue), and T623/S624 (green) side chains are represented as thin sticks. Y652 and F656 side chains that do not make specific drug interactions according to the criteria of Table 2 are omitted. Annotations depict the following interaction types: π–π stacking (black); cation-π (blue); hydrogen bond (green); cation near K+ binding site (blue star). Hydrogen bond annotations in brackets represent intraprotein interactions between Y652 phenolic OH and S624 side chain hydroxyls. The pore subunits of specific Y652 or F656 side chains are numbered in purple in panel B (see text).

Selection of Representative Docking Poses

To what extent can interactions between drugs and models that afford high docking scores be represented as discrete binding poses? This question is relevant to hERG drug docking since there is a wide variety of docking poses in the literature, in many cases describing the interaction of the same drug with different hERG models. Some of the variability is a consequence of the choices of model and docking protocol; the results below, for example, highlight broad differences between Flexidock and GOLD poses. However some of the variability in drug binding poses is likely to be a consequence of the heterogeneous nature of drug binding in the hERG pore (see Discussion). We assessed the extent to which published experimental data can be used to select docking poses for several of the drugs studied. Representative low energy score docking poses for cisapride, dofetilide, terfenadine, and cavalli-2 in the MthK(1LNQ) model are shown in the Supporting Information.

E-4031

Figure 7 illustrates representative low energy poses for E-4031 docked to the hERG MthK(1LNQ) model obtained using GOLD and Flexidock. In each case the drug makes 5 interactions within the hERG pore within the criteria in Table 2. Flexidock poses for E-4031 tend to orient E-4031 vertically with the protonated nitrogen near the cavity K+ site (blue star in Figure 7B). GOLD poses tend to orient E-4031 with the methanesulfonamide directed toward the top of the cavity where it can hydrogen bond with S624 side chains (Figure 7A). In both cases multiple interactions with Y652 and F656 side chains are observed. Two experimental observations indicate that Flexidock-type poses are relevant for hERG block by E-4031 and analogues with high binding affinity. Vilums et al. made a comprehensive study of the effects of varying E-4031 aromatic substituents on hERG block efficacy and found that removing the methanesulfonamide has a negligible effect on hERG block.[32] In fact replacing the methanesulfonamide with a methyl group enhances hERG block by almost 10-fold.[32] This rules out methanesulfonamide hydrogen bonding as a necessary contribution to high affinity hERG block of E-4031 derivatives, although it remains possible that the methanesulfonamide group of E-4031 itself does participate in hydrogen bond interactions. The Flexidock pose in Figure 7B also matches the prediction of Imai et al. who used concatemeric hERG tandem dimers that allow selective alanine replacement of chosen pairs of Y652 or F656 residues.[24] From their experiments E-4031 is proposed to make interactions with adjacent tyrosine side chains (e.g. in subunits 1 and 2) and a single Phe side chain of the following subunit (subunit 3).[24] The pose in Figure 7B accords with this observation. Similar correspondence between predicted hERG subunit interactions from the tandem dimer mutant studies of Imai et al.,[24] and Flexidock poses with the MthK(1LNQ) model, were found for terfenadine and cisapride (Figures S4b and S7 of the Supporting Information).

Cavalli Minimal hERG Blockers

Cavalli et al. recently described an insightful analysis of rationally designed minimal hERG blockers.[33] These high affinity blockers are particularly useful in addressing optimal arrangements of hERG pore aromatic side chains in the context of drug block since they consist simply of 3 benzene (or p-fluorobenzene) rings linked around a protonated secondary nitrogen. High affinity binding must be dominated by aromatic interactions for these molecules. Accordingly, low energy score poses display multiple aromatic interactions for Flexidock docking outputs (Figure 5), and the number of aromatic interactions is greatest in the MthK(1LNQ) model. Figure 8 shows a representative pose for cavalli-2 docked into the MthK(1LNQ) model using Flexidock. In this pose the molecule makes 4 π–π interactions, a cation-π interaction, and a hydrogen bond interaction with a Y652 hydroxyl group. Despite the high affinity nature of drug block by these molecules, multiple interactions of the type illustrated in Figure 8 can be made with the drug bound in a number of different configurations. Figures S5 and S6 in the Supporting Information illustrate some of this diversity in low energy score binding poses for cavalli-2 obtained using both Flexidock and GOLD.

Figure 8

Representative low energy score pose for cavalli-2 docked into the MthK(1LNQ) model using Flexidock. Colors and annotated interactions are the same as described in the legend to Figure 7.

Haloperidol

The inactivation-attenuated mutant S631A reduces hERG drug block by haloperidol.[46] On this basis Durdagi et al. constructed an inactivated state hERG model with a pore domain built onto the KcsA(3F5W) template and used molecular dynamics simulations to identify model states that would facilitate haloperidol binding near S631.[26] However, both GOLD (Tables 3 and 4) and Flexidock (Figure 5) select MthK(1LNQ) over KcsA(3F5W) as a preferred model for haloperidol binding. In each case the drug makes numerous interactions with Y652 and F656. The lowest energy score poses from GOLD (ChemPLP) and Flexidock docking runs are shown in Figure 9, although we emphasize again that both docking methods produce a variety of docked configurations that satisfy multiple interactions with the S6 pore helix aromatic side chains. An analysis of haloperidol block of hERG overexpressed in oocytes indicates that block is strongly attenuated in both Y652A hERG and F656A hERG, supporting the expectation that valid docking poses should show significant interaction with these residues.[47] Since S631 lies far away from the hERG pore cavity and near the selectivity filter close to the extracellular membrane surface,[34] simultaneous interactions of haloperidol with S631 and Y652/F656 are not possible. Consequently, the reduction of haloperidol block in S631A hERG is likely due to the attenuation of inactivation in this mutation and the loss of configurations of S6 residues optimal for drug binding that are a feature of the hERG inactivated state.[48]

Figure 9

Representative low-energy score poses for haloperidol docked into the MthK(1LNQ) hERG pore model using GOLD (ChemPLP scoring) (panel A) and Flexidock (panel B), respectively. Colors and annotated interactions are the same as described in the legend to Figure 7.

Discussion

Since most hERG-blocking drugs access the pore cavity from the intracellular side of the membrane when the channel opens in response to membrane depolarization, open state pore models are likely to best represent the arrangements of key amino side chains that are productive for drug binding. The open state MthK structure (PDB:1LNQ) constitutes a template that has been used successfully to model binding of many drugs with both high and low affinities and variable susceptibilities to Y652A and/or F656A mutation.[20,21,24] Higher resolution structures of the MthK pore excised from the Ca2+-binding domain (PDB: 3LDC; 3LDD)[49] are not good templates for hERG pore states in the context of drug block; the conformation of the S6 helix, and particularly the location of F87 (equivalent to Y652 in hERG), are incompatible with data that support high accessibility of Y652 to drugs (see supplement of ref (21)). A similar situation arises with hERG models built onto the KvAP template. Farid et al. addressed this by manually adjusting the rotational state of the S6 helices in KvAP-based hERG models to recover favorable orientations of the Y652 and F656 side chains.[25] We have found that the MthK(1LNQ) and Farid models are broadly similar in their accordance with alanine-mutation effects on hERG channel block by a number of drugs (see e.g. supplement of ref (21)). Since many drugs, especially high affinity hERG blockers, bind most strongly to the inactivated state of hERG, it is tempting to conclude that these models correspond to structures having S6 helical orientations that are similar to the orientation of S6 in the hERG inactivated state.

How do the putative C-type inactivated state structures of KcsA perform as templates for the hERG pore structure in the context of drug docking? It is straightforward to rule out KcsA:3F7V as a suitable template for the hERG open channel since all of the drugs studied make poor interactions with key residues on the S6 helices. The side chains of Y652 and especially F656 are too widely spaced around the pore helix to allow multiple interactions with aromatic residues on the drugs. Analysis of both Flexidock low energy score output structures (Figure 5) and GOLD scoring (Tables 3 and 4) allows a robust and objective elimination of KcsA:3F7V as a useful template for hERG open state pore models in the context of drug binding.

KcsA:3F5W has previously been used as a structural template for an inactivated state model of hERG.[26] However for each of the 9 drugs tested, both GOLD docking scores (Tables 3 and 4) and analysis of Flexidock poses (Figure 5) selected MthK(1LNQ) over KcsA(3F5W) as a model in which the drugs make a more extensive set of interactions with residues in the hERG pore, especially Y652 and F656. Since the drugs in our set bind more strongly with intact hERG inactivation where this has been determined (Table 1), our analysis supports the conclusion that the MthK(1LNQ) model is more representative of the inactivated state pore structure.

Scoring and Force-Field-Dependence of hERG-Blocking Poses

Although both GOLD and Flexidock uniformly select MthK(1LNQ) as a preferred model for hERG channel block, the low energy score docking poses are somewhat different. Flexidock selects poses in which the drugs are oriented with the positively charged secondary or tertiary amino group oriented toward the top of the pore cavity and either close to the K+ binding site at the focus of the pore helix C-terminal dipole electrostatic field and/or adjacent to a Tyr aromatic ring. GOLD maximizes hydrogen bonds and direct contact between drug and the protein surface lining the pore, and the low energy score poses tend to orient methanesulfonamide groups toward the hydroxyl groups of S624 below the selectivity filter. These differences are not unexpected since Flexidock represents the electrostatic contribution to binding explicitly, and the protonated amino group in all of the drugs studied is expected to make a significant contribution to binding. Recent computational and experimental analysis demonstrates that positively charged amines with multiple aliphatic substituents can interact particularly strongly with aromatic groups.[50,51] Although both GOLD and Flexidock poses are in accord with requirements for both Y652 and F656 as binding determinants for all of the drugs studied, these observations demonstrate that some of the variability in docking poses in hERG homology models in the literature is likely to arise from differences in the docking programs used.

The Role of T623 and S624 in High Affinity hERG Block

Many hERG blockers have significantly reduced binding affinity when tested against T623A and S624A mutants, and this is generally attributed to a role for these residues in hydrogen bond or polar interactions with polar substituents on aromatic rings of drugs in the pore.[29,30] Both dofetilide and E-4031, for example, have considerably reduced block efficacy in T623A and S624A hERG.[52] However, in low energy score docking poses, no hydrogen bond interactions were made with the side chain OH of T623 for any drug, and few hydrogen bonds with the S624 side chain hydroxyl were found in Flexidock poses. For dofetilide and E-4031, enhanced binding affinity is obtained with some analogues in which aromatic ring substituents are removed or replaced with a methyl group,[31,32] demonstrating that interactions of ring substituents with T623 and/or S624 side chains are not critical for high affinity block of hERG by these drugs. Imai et al. suggested that hydrogen bond interactions between the Y652 phenolic OH and T623 side chain hydroxyl might stabilize Y652 side chain rotamers optimized for drug binding,[24] and similar interactions involving the side chain OH of S624 (e.g., Figure 7B) might also play such a role. In addition the side chain hydroxyl groups, especially of S624, might provide a polar environment that favors the location of the positively charged aliphatic amino group near the hERG cavity binding site for a K+ ion,[25] as observed in many of the Flexidock poses.

The Nature of the Inactivated State of hERG

The uniform selection of MthK:1LNQ over KcsA:3F5W as a preferred template for drug docking by GOLD and Flexidock occurs despite differences in binding poses from these docking programs. Additionally, neither docking protocol identifies single binding states for the drugs studied. These observations are consistent with the conclusion that the hERG pore does not contain a conformationally discrete binding site for hERG blockers but that an accessible arrangement of largely aromatic side chains dominates drug binding (see also refs (25)(53)). Successful hERG models converge on structures in which the Y652 and F656 side chains are optimally spaced for multipoint interactions with drug molecules. The interaction of positively charged drugs with multiple aromatic groups in the hERG pore is conceptually similar to the interaction of positively charged, aromatic-rich drugs with multidrug-resistance (MDR) proteins such as P-glycoprotein. Although structurally unrelated to K+ channels, MDR proteins bind positively charged, multiaromatic molecules via clusters of Phe and Tyr residues that project from α-helical segments into a large solvent-accessible cavity (reviewed in ref (54)). Accommodation of diverse drug structures into the MDR binding cavity involves resampling of aromatic side chain rotamers. A similar flexibility of Y652 and F656 rotamers in hERG is probably required to accommodate the diversity of hERG blockers, and rotamer sampling in docking calculations is likely to be important for correctly defining drug poses. The docking analysis also suggests the importance of drug flexibility to optimize multipoint interactions, particularly in Flexidock poses that orient the positively charged aliphatic amine close to the K+ binding site near the top of the pore cavity. This requirement for drug flexibility probably underlies the observations that increasing rigidity around the aliphatic amino group in dofetilide analogues can greatly reduce hERG block.[55]

Conclusions

The results show that experimental data on drug structure, state-dependent block and the effects of channel pore mutation on drug block can be used as constraints to assess pore models for a drug-susceptible channel protein. hERG pore models based on the MthK (PDB:1LNQ) structure broadly accommodate experimental data on drug interactions with the inactivated state that represents a high affinity “receptor” for many drugs. Notably, two versions of the MthK model used in independent studies with different docking methods (this study and that of Imai et al.[24]) produce low energy score poses that accord with the pattern of aromatic side chain-drug interactions in the hERG pore cavity identified using selective mutagenesis in hERG tandem dimers.[24] In the absence of a crystal structure, further refinement of pore conformations of hERG relevant for drug block may be obtained with targeted experimental data, for example to better define the roles of S624, and especially T623, in drug block. These approaches may ultimately lead to hERG pore models that are more usefully predictive of the effects of drugs on hERG channel block in drug development programs.