Aminopurine derivatives as putative SopE inhibitors

Abstract

Salmonella outer protein E, SopE, is a virulence factor that is secreted from gram-negative Salmonella enterica, which is pathogenic to humans and is responsible for gastroenteritis and typhoid fever. SopE targets the Rho GTPase family proteins of the host cell and manipulates them by mimicking GTPase regulatory protein guanine nucleotide exchange factor. The aim of this work is the investigation of novel inhibitors against SopE. Structure-based pharmacophore modeling was used to identify structural features that would be important for SopE recognition. Glide fragment library was used for four-point pharmacophore hypothesis development. Small-molecule database screening was performed based on a 3D similarity to the best pharmacophore hypothesis. Binding affinity of filtered database molecules to SopE was predicted quantitatively by molecular docking and scoring using Glide software. Top scoring hits were further analyzed and five molecules were proposed as potent and selective SopE inhibitors. Four out of five proposed molecules were found to be aminopurine derivatives.

• GEF inhibitor
• induced fit docking
• molecular docking
• structure-based pharmacophore building
• virtual screening

Introduction

The Ras (Rat sarcoma) superfamily of GTPases consists of a large number of monomeric G-proteins, which are intracellular cell signaling proteins that take part in many diverse cellular functions, including growth, migration, adhesion, cytoskeletal integrity and survival1,2. One subfamily is the Rho family of small (around 21 kDa) GTPases, which regulates actin dynamics and also participates in cell cycle progression, gene transcription, cell proliferation and motility, vesicular trafficking, phagocytosis and cytokinesis3–6. Guanine nucleotide exchange factors (GEFs) reduce the affinity of Rho GTPases to guanosine diphosphate (GDP) by interacting with their switch regions and the nucleotide containing P-loop6. Upon GEF binding, GDP is substituted by GTP, GTPase is activated and the downstream signaling cascade is turned on7.

Salmonella genus are found widespread in cold- and warm-blooded animals and cause illnesses like typhoid fever and gastroenteritis8. Salmonella outer effector proteins, which are secreted into the host cell by type III secretion system such as SopB, SopE and SopE2, modulate different cellular processes and are therefore involved in various stages of bacterial infection. SopE protein initiates bacterial internalization by membrane ruffling8–10. SopE targets the Rho GTPases Rac1 and Cdc42. SopE acts as a GEF and stimulates nucleotide exchange from GDP to GTP, inducing signal transduction by the host Rho GTPases11. Bacterial uptake of SopE induces actin cytoskeleton rearrangement, which enables the entry of other effector proteins into the host cell11.

The structure of SopE consists of six α-helices in two, three-helix bundles that assume a Λ form and a two stranded β-sheet followed by an eight amino acid long loop (Figure 1). This loop, which connects the two legs of the Λ shape, contains the catalytic core of SopE, a GAGA (residues Gly166, Ala167, Gly168, Ala169) motif (central loop in Figure 1). Catalytic core of SopE does not resemble the catalytic domain of other Rho-GEFs; whose structures have either DH-PH or DHR2 domains12. The crystal structure of the SopE–Cdc42 complex provides information about the interaction between these two proteins as well as clues about the mechanism by which the Cdc42 GTPase is activated by the SopE GEF. The structural analysis of the SopE–Cdc42 complex reveals that insertion of the GAGA loop between the nucleotide-binding regions of Cdc42 leads to major conformational changes in the Switch I (residues 28–40) and Switch II (residues 57–74) regions, resulting in reduced GDP affinity11. Therefore, GAGA loop is crucial for bacterial GEF activity of SopE. Similar conformational changes in these switch regions have been observed by other eukaryotic GEF proteins11. Other than GAGA loop, multiple contacts are observed between Switch I region of Rho GTPase and Asp124 residue of SopE. Ile131 residue of SopE was also reported to have hydrophobic contacts with the Switch I region. Other reported Switch I interactions include Gln194 and Lys198 residues of SopE. GAGA loop residues Gly166, Gly168 and Ala169 interact with the Switch II region. Other Switch II contacting SopE residues are Asp103, Gln109, Ser165 and Thr17411.

Figure 1. Crystal structure of the GEF domain of SopE in complex with Cdc42, shown in cartoon representation (PDB id: 1gzs:A,B). Switch regions of Cdc42 are shown in red and purple. Critical GAGA loop residues are shown in yellow.

Aminopurine is an organic fluorescent molecule that has been proven to be important in various fields and therefore may have a potential for new therapeutic purposes. It is extensively used in nucleic acid studies since it can function as both adenine and guanine analog13,14. Aminopurines were also found to be ATP competitive protein kinase inhibitors in many studies15–17. In this work, we employed a virtual screening (VS) approach that combines structure-based pharmacophore building, database filtering and molecular docking to identify small-drug like molecules that would inhibit GEF activity of SopE. Docking results were further analyzed according to their interactions with target, druglikeness and bioavailability. Up to date, no inhibitors for SopE have been reported. Here we computationally identified a number of aminopurine derivatives as novel potential inhibitors of SopE.

Methods

All VS applications were performed using Schrödinger Suite 2011 (LLC, New York, NY) on Linux platform (Linux Professional Institute, Sacramento, CA) using HP xw6600 Workstation (Hewlett-Packard Co, Palo Alto, CA). The following Schrödinger modules were used: Protein Preparation Wizard18, LigPrep19, ConfGen20, Phase21, QikProp22, Glide23–25 and Induced fit26,27. All modules were accessed via Maestro graphical interface28 (all modules designed by LLC, New York, NY). The workflow depicted in Figure 2 is outlined in detail below.

Figure 2. Schematic workflow of VS and docking of SopE.

Database generation

The 3D structures of the small molecules from the big vendors of ZINC database (Table S1) were prepared via LigPrep and ConfGen scripts, and merged using Phase Database Generation module (LLC, New York, NY). The parameters used in the preparation step are the following: (a) possible ionization and tautomer states were generated at pH 7.0 ± 2.0; (b) chiralities were obtained from the 3D geometry of the structures; (c) for each compound, maximum four low energy isomers, one ring conformation and 10 conformations per rotatable bond were generated; (d) 50 steps of energy minimization were carried out for the conformers with the truncated Newton (TCNG) method using the OPLS_2005 force field and a distance dependent dielectric constant of 4; (e) conformers with high-energy ionization/tautomer states were automatically removed and number of conformers was limited to 100 per compound, which is the default value for conformer generation in Phase Database Generation module; (f) a built-in QikProp script predicted ADME and druglikeness properties of the molecules; (g) molecules that violated Lipinski’s rule or had reactive functional groups were removed. Database preparation yielded a total of 25.8 million conformations from the initial 500 000 molecules.

Receptor protein preparation

Coordinates of the target protein SopE (PDB id: 1gzs) were obtained from RCSB Protein Data Bank (NJ)29. One SopE chain from the crystal structure, which contains the SopE–Cdc42 in dimeric form, was used to prepare the protein within the Protein Preparation Wizard. Bond orders were assigned, hydrogens were added to the heavy atoms and all waters were deleted. Reorientation of certain hydroxyl and thiol groups, amide groups of asparagines, glutamines and imidazole ring of histidines, protonation states of histidines, aspartic acids and glutamic acids were optimized at neutral pH. Using force field OPLS_2005 (NY), minimization was carried out setting maximum heavy atom RMSD to 0.30 Å.

Receptor grid generation

After preparation, receptor grids were generated with Glide by specifying the binding site with a 3D cubic box. Enclosing a box of 20 Å side lengths was placed by centering the critical GAGA loop. Potential energies of receptor atoms confined in this box were calculated without scaling their van der Waals radius. Gly166, Gly168 and Ala169, three GAGA loop residues that are known to interact with Cdc42, were identified as hydrogen bond partners for later use as constraints during docking. Rotation of hydroxyl groups was not allowed.

Structure-based pharmacophore generation

The database, which comprised 500 000 molecules, was filtered using a pharmacophore hypothesis in order to optimize the docking and screening time. Due to the fact that there are no previously identified ligands for SopE; a method that couples structure information with pharmacophore building method, which is called E-pharmacophores (E-pharm), was employed. This methodology makes use of the receptor binding site structural information, and generates pharmacophore sites based on the protein residues that can make interactions upon ligand binding. In order to identify such residues, E-pharm procedure required an initial docking stage. A set of 667 unique small (6–37 atoms, MW range 32–226) fragments30, derived from the molecules in the medicinal chemistry literature (by Schrödinger Company), were used as ligands. The interface between SopE and Cdc42 in the complex structure was examined in detail in order to identify key interacting residues. The contribution of the GAGA loop to Cdc42 recognition has been reported in the literature11. Extra precision (XP) docking without constraints was performed using all fragments and the GAGA loop of SopE as a binding site. XP descriptors were generated to obtain atom-level energy terms such as hydrogen bond interaction, electrostatic interaction, hydrophobic enclosure and π–π stacking interaction. E-pharm script mapped XP energetic terms onto each fragment atom. Pharmacophore sites were determined on fragment atoms with Phase, using the built-in set of six features: hydrogen bond acceptor (A), hydrogen bond donor (D), hydrophobic (H), negative ionizable (N), positive ionizable (P) and aromatic ring (R) site. Energetic contribution of each pharmacophore feature was quantized by assembling energies of the atoms that comprise pharmacophore site. Minimum distance between site points was kept as 2.0 Å and minimum distance between site points with same feature was kept as 4.0 Å. Pharmacophore sites were ranked based on their summed energies, and top ranking sites were combined to create a four-point hypothesis. Acceptor-acceptor-donor-negative ionizable (AADN) hypothesis was selected and the database was pre-filtered based on their 3D similarity to the selected pharmacophore hypothesis, allowing the minimum fitness value to be 0.

Glide docking

A total of 9970 molecules obtained by filtering were used for receptor docking. Glide standard precision (SP) docking was performed with 9970 molecules and top 10% hits based on GlideScore were re-docked to SopE in XP mode, keeping all docking parameters as default. At least one hydrogen bond defined on GAGA loop residues was selected to be satisfied. Using Monte Carlo random search algorithm, ligand poses were generated for each input molecule, and binding affinity of these molecules to SopE was predicted with GlideScore. Potential energies of the docked molecules were also predicted with empirical E-model scoring function. Post-docking minimization was performed with OPLS_2005 force field and at most five poses per ligand were saved. Strain energies of bound and free forms of ligands were calculated and hits having more than 4 kcal/mol energy difference between the two forms received a penalty equal to the quarter of their strain energy difference, which is added to docking score.

Induced fit docking

Induced fit docking (IFD), which allows both receptor side chain and ligand flexibility were carried out with proposed molecules with an inhibition potential towards SopE. IFD was performed to determine the interaction potential of proposed molecules to SopE when receptor-side chain flexibility is induced and no constraints are given. Receptor grid parameters were kept same as Glide docking, except no bonding constraints were determined on GAGA loop residues. Other internal settings and parameters in induced fit workflow were kept at their default states and values.

Results and discussion

In order to find potent inhibitors against SopE, structure-based VS was combined with pre-filtering of small molecule database with pharmacophore building. In this work, a hybrid screening method that combines pharmacophore building and molecular docking was performed to improve screening efficiency. Phamacophore modeling was employed to make use of the binding site shape, size and chemical properties in filtering the database. This allowed us to start the workflow with a larger database, which was then reduced by pharmacophoric search, and the likelihood of finding good-binding ligands within a reasonable screening time could be enhanced. Filtered molecules were docked to receptor SopE in multistep Glide docking protocol and potential inhibitors were identified. A simplified workflow of the VS methodology is presented in Figure 2.

Structure-based pharmacophores

Ligand-based pharmacophore building, which facilitates the database filtering and screening process, requires an initial set of known ligands. However, there are no known ligands identified for SopE to date. Instead, knowledge about the receptor itself was utilized for determining the necessary structural and chemical features in a potential ligand for receptor binding. Fragments were docked to the binding site of the target SopE, near the GAGA loop, as explained in “Methods” section. Docking calculations were performed and protein-fragment interactions were characterized based on energetic contributions of each site on the fragments to the docking score. Thus, coordinates and types of pharmacophore sites (e.g. hydrogen donor, acceptor) that show good binding to GAGA loop were determined, and these sites were ranked.

A total of 21 pharmacophore sites that significantly contribute to the docking score were detected. Figure 3(A) shows all fragment-binding modes that were docked to SopE-binding site. Ranking of all possible pharmacophore sites as well as energetic contributions, coordinates and feature types are listed in Table 1 and are shown in Figure 3(B). In order to create a pharmacophore hypothesis, E-pharm provides two options: (a) allowing the user to select the sites that will constitute the hypothesis and (b) clustering based on site–site distances and total number of site points, specified by the user. Phase allows the number of pharmacophore site points to be between three and seven31. This value was selected as four, because when more site points were included, the likelihood of finding a cluster was decreased. In addition, hypothesis with more site points would yield fewer matches in database, and possibly eliminating potential ligands. The 21 pharmacophore sites were clustered and four high ranking site points that are within 4.2 Å–7.5 Å to each other were determined by E-pharm, to create the hypothesis AADN based on site points A708, N1548, A42 and D1091 (Figure 3B). The angle and distances between the site points are given in Table S2 and Figure S1. D995, which has the highest rank, was not picked because it is at least 10 Å away from the four site points. The AADN hypothesis was superimposed on the Cdc42–SopE complex and the residues that would potentially be replaced by the site points which would prevent Cdc42–SopE interaction were identified. A708 was located near Ala169. Likewise, N1548 site was located near Lys198_SopE, possibly preventing its interaction with Asp38_Cdc42. A42 was occupied by Val36_Cdc42 whereas D1091 was in the vicinity of Thr35_Cdc42. Both residues are in Switch I region of Cdc42 and they were known to interact with Asp124_SopE11.

Figure 3. (A) Binding modes of all fragments superimposed at the binding site of SopE in stick representation. (B) Pharmacophore site points as spheres located on the binding site of SopE, determined by E-pharm. Top six sites based on energetic contribution are labeled. All site points are colored and labeled according to their feature type (donor, acceptor, hydrophobic or ring).

Table 1. Pharmacophore features and ranks by E-pharm. Numbers in feature labels are indices which are automatically assigned during feature generation.

Rank	Feature label	Score	X	Y	Z	Type
1	D995	−1.80	−42.79	58.44	53.22	D
2	A708	−1.70	−30.15	57.8	55.76	A
3	N1548	−1.59	−35.87	53.23	58.95	N
4	D1296	−1.52	−31.51	48.41	53.15	D
5	A42	−1.49	−32.33	53.73	56.24	A
6	D1091	−1.29	−34.84	51.76	54.15	D
7	R2124	−0.98	−42.51	61.08	52.86	R
8	A62	−0.98	−34.14	58.86	54.35	A
9	R2224	−0.95	−41.82	59.85	57.27	R
10	R2233	−0.90	−39.8	58.26	53.84	R
11	D939	−0.90	−39.2	58.7	57.22	D
12	R1961	−0.81	−28.68	56.67	57.28	R
13	R2033	−0.77	−35.52	52.49	56.7	R
14	A81	−0.70	−43.29	56.69	55.55	A
15	D930	−0.70	−29.96	57.53	48.45	D
16	D937	−0.70	−31.93	54.98	52.91	D
17	H1450	−0.70	−31.37	55.91	55.68	H
18	D1282	−0.70	−37.06	56.69	52.44	D
19	R2142	−0.53	−29.74	52.3	54.63	R
20	D999	−0.34	−35.34	56.55	56.15	D
21	R2236	−0.20	−33.94	56.87	54.72	R

Molecular docking

Pharmacophore pre-filtering with the four site point E-pharm hypothesis AADN yielded 9970 compounds from the small molecule database generated from ZINC. The molecules obtained by filtering with AADN hypothesis were docked to SopE with Glide software to predict binding affinity of molecules to SopE, and to investigate ligand–protein interactions.

The SopE coordinates used in this study were taken from the crystal structure of the SopE–Cdc42 complex. The structure of SopE in its apo form is not available, but there are two homologous proteins, BopE32 and SopE233, for which structure determination and relaxation measurements by NMR have been performed. These studies show that the GAGA loop is a flexible structural element. This is an expected result since this loop is a connector between the two arms of the Λ shape and comprises small and flexible residues. The Λ shape is suggested to assume open and closed conformations, in which the open conformation can bind Cdc42 but the closed conformation assumed by BopE would be in steric clash with Cdc42. In our work, the flexibility of the GAGA loop was taken into account by choosing to use the open conformation assumed by SopE in the SopE–Cdc42 complex. Even though only one conformation was used, this conformation allowed the maximum space for docking.

The position and orientation of each ligand relative to the receptor protein was determined and scored with internal scoring function GlideScore. Glide SP and XP docking modes were used, successively. At least one hydrogen bond with Gly166, Gly168 and Ala169 was set as constraint. Four hundred and eighty-four compounds satisfied the constraint criterion and were docked in Glide XP mode. GlideScores of these 484 compounds vary between −8.2 and −0.4 kcal/mol in Glide XP mode. Detailed analysis of docking parameters was carried out for hits with GlideScore −6.5 kcal/mol or higher. This threshold yielded 20 hits out of 484 compounds. 2D structures of all compounds are listed in Table S3.

Molecules were sorted according to their GlideScore values in Table 2. Detailed GlideScore components of each molecule are listed in Table S4. Emodel scoring function (Table 2, last column) includes penalty terms for internal strain energy of the generated poses. Therefore, GlideScore was used for comparing poses of different hits, and Emodel was used to decide between different conformers (or poses) of a particular ligand. Five molecules out of final 20 hits were selected as potent inhibitors of SopE based on the GlideScore and Emodel as well as additional criteria, such as absorption, distribution, metabolism and excretion (ADME) considerations and druglikeness of ligands and strain energy differences of ligands.

Table 2. Docking results of top 20 hits in kcal/mol. ZINC IDs and strain corrected GlideScores for proposed inhibitors are indicated by (*).

#	Title	Glide score	Bound energy	Free energy	Strain energy	Strain penalty	Strain GlideScore	Glide Emodel
1	ZINC03860635	−8.171	16.182	5.816	10.365	1.591	−6.580	−34.444
2	ZINC10274669	−7.948	59.025	26.953	32.072	7.018	−0.930	−57.122
3	ZINC00370772*	−7.681	49.776	45.792	3.984	0.000	−7.681*	−37.519
4	ZINC09214236	−7.567	64.928	64.928	0.000	0.000	−7.567	−33.284
5	ZINC03874923	−7.567	11.363	−5.676	17.039	3.260	−4.307	−50.343
6	ZINC17005625	−7.499	68.736	61.171	7.565	0.891	−6.608	−37.686
7	ZINC17020721*	−7.489	62.430	62.415	0.015	0.000	−7.489*	−38.399
8	ZINC19910989*	−7.456	16.478	13.041	3.437	0.000	−7.456*	−60.893
9	ZINC05033974	−7.410	67.590	54.666	12.924	2.231	−5.179	−44.115
10	ZINC09214236*	−7.330	65.940	64.687	1.253	0.000	−7.330*	−36.937
11	ZINC00375097	−7.324	11.928	11.737	0.191	0.000	−7.324	−39.996
12	ZINC05033974*	−7.274	64.699	63.128	1.572	0.000	−7.274*	−41.142
13	ZINC00388310	−6.960	13.084	10.852	2.232	0.000	−6.960	−43.441
14	ZINC00754304	−6.949	45.244	38.999	6.245	0.561	−6.388	−52.592
15	ZINC09214236	−6.863	69.688	64.003	5.685	0.421	−6.442	−32.037
16	ZINC19910989	−6.841	13.163	8.705	4.459	0.115	−6.726	−57.662
17	ZINC00097607	−6.763	19.312	13.621	5.691	0.423	−6.340	−37.690
18	ZINC39345243	−6.749	55.915	46.046	9.869	1.467	−5.282	−51.248
19	ZINC00155974	−6.700	16.180	14.826	1.353	0.000	−6.700	−41.773
20	ZINC17076964	−6.573	40.677	32.146	8.531	1.133	−5.44	−51.235

Pharmacodynamic properties, oral bioavailability and druglikeness as well as other molecular descriptors for the 20 candidate hits were determined by QikProp (LLC, New York, NY; Table S5). ZINC03860635 and ZINC03874923 have poor predicted human oral absorption (HOA) value, at 0.7% and 0.1%, respectively. ZINC00375097 was also observed to lack hydrophobic solvent accessible surface area with 0 Ų. Therefore, these three molecules were eliminated. Highest predicted HOA was observed for ZINC17076964 at 86.4%. ZINC19910989 and ZINC17076964 showed the highest lipophilicity at 3.231 and 3.433, respectively.

Docking allows ligands to be strained, so that they can fit better into the binding site. Since too much strain may indicate false positives, energies of free and docked conformations of hits were calculated and strain penalties were added to docking score. As a result, 11 out of 20 hits received strain penalties. ZINC10274669 had the highest strain penalty with 7 kcal/mol, which is almost as much as its original docking score. ZINC03860635, ZINC03874923 and 9th ranking ZINC05033974 also received relatively high strain penalties, which lowered their GlideScores considerably.

Molecules with acceptable ADME properties and GlideScore above −7 kcal/mol after strain correction were identified as ZINC00370772, ZINC09214236, ZINC17020721, ZINC19910989 and ZINC05033974. ZINC09214236, ZINC19910989 and ZINC05033974 had multiple poses. Selection of poses of a particular molecule was performed mainly based on their Emodel values. If Emodel value difference was not significant, docking scores were also considered. Poses with too much strain were also eliminated. Accordingly, 10th ranking pose of ZINC09214236, 8th ranking pose of ZINC19910989 and 12th ranking pose of ZINC05033974 were selected. Proposed SopE inhibitors are indicated by star in Table 2, and their ADME properties are given in Table 3.

Table 3. ADME and druglikeness properties of the selected ligands by QikProp.

Molecule*	SASA†	FOSA‡	FISA¶	QPlogHERG§	QplogS||	QPlogPo/w#	% Human oral absorption$
ZINC00370772	415.4	108.6	223.4	−5.653	−1.21	−1.102	54.1
ZINC17020721	463.4	139.2	217.8	−6.036	−2.13	−1.049	55.4
ZINC19910989	694.2	97.3	212.7	−5.439	−4.56	3.187	68.5
ZINC09214236	441.4	160.4	197.5	−5.492	−1.77	−0.975	59.2
ZINC05033974	502.7	221.9	201.7	−5.692	−2.65	0.372	55.7

*ZINC IDs.

†Total solvent accessible surface area in Ų using a probe with a 1.4 Å radius (acceptable range: 300–1000).

‡Hydrophobic component of the SASA (acceptable range: 0–750).

¶Hydrophilic component of the SASA (acceptable range: 7–330).

§Predicted IC50 value for blockage of HERG K+ channels (acceptable range: below −5).

||Predicted aqueous solubility, S in mol/dm⁻³ (acceptable range: −2–6.5).

#Predicted octanol/water partition coefficient (acceptable range: −6.5–0.5).

$Predicted human oral absorption on 0 to 100% scale (<25% is poor and >80% is high).

Structures of final 20 conformations were analyzed and clustered by their 2D structural similarities via ChemMine Web Tool (CA) (http://chemmine.ucr.edu/). Hierarchy of clusters based on pairwise compound similarities were defined using the Atom Pair descriptors and Tanimoto coefficient34. The similarity tree output is provided in Figure S2. The results show that these molecules form two clusters, with one cluster of five compounds and the other cluster with the rest of the compounds. The most significant structural dissimilarity between these clusters was that, the molecules in the first cluster contain an aminopurine ring while molecules in the second cluster do not contain an aminopurine ring. It was observed that four top ranking compounds that were proposed as SopE inhibitors belong to the first cluster.

Binding mode and interaction analysis of the selected hits

Superposition of the proposed molecules based on their interaction with SopE showed that their binding modes are similar (Figure 4F). Proposed molecules with aminopurine groups overlaid well and similar spatial arrangement of S1, S2 and S4 was observed. In contrast, aminopurine group of S5 was faced toward GAGA loop (Figure 4F). The five selected compounds satisfied at least two out of four pharmacophore sites with a 3 Å tolerance, suggesting that the pharmacophore building and docking results converged. The superposition of the selected compounds on the pharmacophore AADN hypothesis and the distances between the pharmacophore sites and the respective ligand atoms are given in Table S6. The interactions of the five proposed inhibitors were analyzed to verify that the ligands made contacts with previously identified important residues of SopE. The binding mode of selected molecules along with receptor hydrogen bond interactions are represented in Figure 4. The list of hydrogen bond interactions between molecules and SopE is given in Table 4.

Figure 4. (A)–(E) Glide XP binding modes of proposed molecules S1–S5 at SopE binding site. SopE is shown in ribbon representation. Ligands and SopE residues interacting with the ligand are shown in stick representation, colored by atom type. Hydrogen bonds between ligands and SopE are indicated by black dashed lines. (F) Binding mode overlay of five proposed molecules. SopE binding site is in surface representation and ligands are in stick representation, colored by molecule.

Table 4. Summary of H-bond interactions between and proposed molecules of SopE. “+” sign indicates presence of a hydrogen bond and number of signs indicates the number of hydrogen bonds between ligand and a particular SopE residue.

		Residue
Database title/ligand	2D structure	Ala117	Lys121	Asp124	Gly168	Asn170	Lys198
ZINC00370772/S1		+	+	+	+
ZINC17020721/S2		+			+	+	++
ZINC19910989/S3		+	+		+		++
ZINC09214236/S4				+	+		+
ZINC05033974/S5			+		++		++

All proposed molecules favored the interaction with Gly168 residue of the GAGA loop, which is reported to be an important promoter of nucleotide release. Mutation of Gly168 to alanine was shown to abolish SopE-mediated nucleotide release11. Therefore, Gly168 interactions were interpreted as a positive outcome. It was observed that Gly168 residue preferred interactions with either hydroxyl group or amine group in ring systems. It should be noted that structure-based VS of SopE was carried out with at least hydrogen bond criterion on GAGA loop residues Gly166, Gly168 and Ala169. Asp124 and Lys198 were also reported to be important residues for stabilizing Switch I regions of small GTPases. Lys198 interactions were mostly with lone oxygen or nitrogen atoms of ligands as well as nitrogen atoms. Two interactions between charged –NH groups and Asp124 residue of SopE were observed. Ala117 residue also participated in hydrophobic interactions. Four out of five proposed SopE inhibitors include an aminopurine ring, which may be a functionally important group upon SopE inhibition. These aminopurine rings bind to Lys 121, Asp124, Gly168 or Lys198 residues of SopE. One of the proposed molecules, S1 or ZINC00370772 is used in the therapy of viral diseases. This molecule was denoted as a novel nucleoside analog and reported to inhibit the replication of several DNA and RNA viruses, including vaccinia, herpes simplex (types 1 and 2), measles and vesicular stomatitis35.

Induced fit docking

IFD was carried out with the proposed molecules, in order to study their binding poses, interactions with the receptor and predicted binding affinity in terms of docking score, when both receptor-side chain and ligand atoms were flexible. IFD binding modes along with their receptor interactions are given in Figure 5 and IFD scores are listed in Table 5.

Figure 5. (A)–(E) IFD binding modes of proposed molecules S1–S5 at SopE binding site. SopE is shown in ribbon representation. Ligands and SopE residues interacting with the ligand are shown in stick representation, colored by atom type. Hydrogen bonds between ligands and SopE are indicated by black dashed lines. (F) Binding mode overlay of five proposed molecules. SopE binding site is in surface representation and ligands are in stick representation, colored by the molecule.

Table 5. IFD Score of proposed molecules in kcal/mol. Docking score of the best pose is presented for each ligand.

Database title	Ligand	IFD Score	GlideScore
ZINC00370772	S1	−7.87	−7.681
ZINC17020721	S2	−7.42	−7.489
ZINC19910989	S3	−4.54	−7.456
ZINC09214236	S4	−7.09	−7.330
ZINC05033974	S5	−8.69	−7.274

In IFD, the highest scoring binding mode of S1 was essentially similar to its Glide XP binding mode (Figure 5A). S1 was able to interact with Ala117, Asp124 and Gly168, as in Glide docking. The Lys121 interaction with the aminopurine ring was not observed in IFD. Instead, Lys198 made a hydrogen bond with S1. Highest scoring binding mode of S2 was also similar to its original Glide XP binding mode in which it was able to interact with Asn170 and Lys198. However interactions with Gly168 were not observed. In the second ranking binding mode of this ligand, which has a completely different orientation (IFD score: −7.18 kcal/mol), the aminopurine ring was able to interact with Ser120, Gly168 and Lys198 (Figure 5B). The highest scoring docking mode of S3 generated in IFD protocol had a completely different binding mode than its original binding mode, in which the molecule showed interactions with Lys123, Ser165 and Ala167 (Figure 5C). This ligand also showed the lowest score in IFD (Table 5). Investigation of the other poses of S3 showed that in all orientations, this ligand preferred to interact with Ala167 of the GAGA loop instead of Gly168. Highest scoring pose of S4, which is considerably different from its original binding pose, was able to interact with Ala117, Asp124 and Gly168 (Figure 5D). IFD binding mode of S5 with highest score was similar to its original docking orientation, interacting with Lys121, Gly168 and Lys198 (Figure 5E). Figure 5(F) shows the overlay of IFD binding modes of proposed molecules. All molecules except S3 had similar binding modes. Overall, induced fit results showed that identified interactions between proposed molecules with aminopurine ring (S1, S2, S4, S5) and SopE in Glide docking were also observed when receptor was flexible, while S3 maintained its contact with the GAGA loop by forming interactions with residue Ala167, different from other ligands.

Selectivity

Binding selectivity of proposed molecules to SopE homologs, which are SopE2 (Pdb id: 1e6e) with 70% amino acid sequence identity36–39 and BopE (Pdb id: 2jok) with 27% amino acid sequence identity to SopE40 were investigated. The catalytic GAGA motif of SopE is conserved in SopE2 (residues 166–169) whereas GAGT motif is present in BopE (residues 171–174)32. Glide XP docking with strain correction was performed on these proteins centering GAGA/GAGT motif around a 20 Å grid box. Proposed molecules were docked to SopE2 and BopE active site with the same docking parameters used for SopE, including at least one hydrogen bond constraint on GAGA/GAGT motif. Since GAGA loop is critical for the GEF function of these proteins, hydrogen bond interaction between proposed ligands and GAGA/GAGT loop was considered as an indicator for potential inhibitory action. No ligand pose could be generated for BopE, indicating that the proposed SopE molecules are not able to bind to BopE with the specified constraints. On the other hand, S1, S3 and S5 were successfully docked to SopE2 satisfying the GAGA motif H-bond constraint. The selectivity docking results are listed in Table 6. None of the hits received strain penalty, and the strain adjusted GlideScores are between −4.58 and −2.92, which is considerably lower than SopE docking results (around −7.50 kcal/mol). However, despite their lower docking scores, hits were able to interact with Gly166 and Ala169. These results suggest that the proposed molecules are selective toward SopE in terms of docking score, but also have inhibition potential against SopE2.

Table 6. Docking results of proposed SopE inhibitors to SopE2 in kcal/mol.

Database title	Ligand	GlideScore	Interacting residue on GAGA loop
ZINC00370772	S1	−4.58	Gly166, Ala169
ZINC19910989	S3	−3.81	Gly166, Ala169
ZINC05033974	S5	−2.92	Ala169

Conclusions

Potent inhibitors of SopE were investigated with computer-aided drug discovery tools. Ligand-based and structure-based VSapproaches were used employing the Schrödinger Suite 2011 modules. Database molecules were screened and pre-filtered, according to their 3D similarity to the selected pharmacophore AADN hypothesis. Filtered database molecules were docked to receptors in two successive docking modes, allowing ligand flexibility. Binding poses of input molecules were determined and scored with GlideScore. Top 20 molecules were subjected to a number of post-docking analyses, including strain energy calculation, assessment of ADME and druglikeness properties. Ultimately, five potential SopE inhibitors molecules with GlideScores around −7.5 kcal/mol and acceptable pharmacodynamic properties were proposed.

The aim of this VS study was to find a potential SopE inhibitor that would hinder the interaction between SopE and Cdc42. The inhibitor binding site was selected such that the inhibitor would interact with the GAGA loop residues that contact the switch regions of Cdc42. Detailed analysis of the interactions of the molecules with SopE showed that the molecules made contacts with Gly168 of the GAGA loop, which interacts with the Switch II of the GTPase Cdc42, and does not stimulate nucleotide exchange on Cdc42 when mutated to alanine. Four out of five proposed molecules were observed to have an aminopurine ring, which may be a functional group for SopE inhibition and needs further investigation in biological assays. Ligand–receptor interactions were investigated with more accurate IFD, and molecules with the aminopurine ring showed favorable interactions, which again points out the importance of this group. Additionally, selectivity of these molecules to known SopE homologs SopE2 and BopE were studied. None of the proposed molecules were able to interact with GAGT motif of BopE whereas three molecules (S1, S3 and S5) interacted with GAGA motif of SopE2 with considerably lower docking scores (−4.5 kcal/mol or less) compared to SopE, implying that proposed molecules are selective toward SopE but have a possibility to inhibit homologous SopE2 GEF as well.

Declaration of interest

The authors report no declarations of interest.

This work is financially supported by the Scientific and Technological Research Council of Turkey (TÜBİTAK) through project 111M444. The authors thank Prof. Sumru Ozkirimli for insightful discussions and her valuable contributions.