Applications of Molecular Modeling in the Design of New Insect Repellents Targeting the Odorant Binding Protein of Anopheles gambiae

Anhanguera Faculty of Brasília, Qs 01, Taguatinga, 72030-700 Brasília-DF, Brazil Laboratory of Molecular Modeling Applied to the Chemical and Biological Defense, Military Institute of Engineering, Praça General Tibúrcio, 80, 22290-270 Rio de Janeiro-RJ, Brazil Institute of Chemistry, University of São Paulo, Av. Trabalhador São-carlense, 400, CP 780, 13560-970 São Carlos-SP, Brazil National Institute of Industrial Property, Praça Mauá, 7, 20081-240 Rio de Janeiro-RJ, Brazil


Introduction
][3][4][5][6][7][8][9][10] The molecules in the cloud penetrate in their cuticular pores and are captured by the olfactory binding proteins (OBP).2][13][14] Studies by Biessmann et al. 11 have proved the importance of OBP for smell recognition in the search for food, in mating and other activities. 3,12,14Among the 60 types of OBP coded in the A. gambiae genoma, the Anopheles gambiae odorant binding protein (AgOBP 1 ) was found in higher levels in the antennas of females when compared to males and its levels were observed to reduce when the females feed on blood. 14espite not having its mechanism of action well elucidated yet, N,N-diethyl-3-methylbenzamide (DEET) is the main repellent in use today against A. gambiae. 2,3,4,6,15sitsanou et al. 14 suggested that this molecule targets AgOBP 1 , blocking the binding site of the attractive molecules, and published the 3D structure of the AgOBP 1 -DEET complex as the first example of a repellent recognized by an OBP.This achievement should give a significant impulse to the drug design of new repellents as it provides the first known target.
Mosquito repellence is one important pharmacological activity that has been recently reported for the indian clove (Syzygium aromaticum) oil.Products containing the clove oil in concentrations of 10 or 20% plus 10% of makaen (Limonella zanthoxylum) presented repellence against Aedes aegypti and Culex quinquefasciatus. 16There is still no proposed mechanism to explain the repelency observed for the components of the indian clove oil but our group believes that these molecules are also able to bind to the AgOBP 1 binding site in the same way as DEET.So, docking and molecular dynamics (MD) techniques were employed to study the interactions between DEET and AgOBP 1 and to compare them to the interactions of this OBP with attractive molecules such as lactic acid, 1-octen-3-ol and indole 11 and the majoritary components of the indian clove oil: eugenol, eugenyl acetate, α-humulene and β-caryophyllene 3,[17][18][19][20][21] (Figure 1).Our results suggest that the repellent molecules have greater affinities for the binding site of AgOBP 1 than the attractive ones and that eugenyl acetate is able to bind more efficiently to AgOBP 1 than DEET.The main features of the AgOBP 1 binding site essential for the design of more effective repellents were also reported.

Docking studies
To perform the calculations in the present work, it was used the dimmeric 3D structure of AgOBP 1 , complexed with DEET and PEG (a 5 unit polyethylene glycol molecule) with resolution of 1.6 Å, available in the Protein Data Bank (PDB) 14 under the code 3N7H.This crystal structure presents one DEET molecule anchored in each subunit at sites located near the interface between the two monomers.
The 3D structures of each compound in Figure 1 were built using the PC Spartan 08 ® program, 22 and their partial atomic charges were calculated by the RM1 semi-empirical method. 23The compounds were docked inside the AgOBP 1 binding sites using the Molegro Virtual Docker (MVD) ® program. 24ockings were done in the two monomers because studies suggest that DEET binds more favorably in the dimeric form of AgOBP 1 since it has been found that a monomer-dimer equilibrium exists in solution and it is possible that the protein dimer is the molecular target of DEET under physiological conditions. 13,14The binding site was restricted to a sphere with a radius of 6 Å around DEET.The ligands and the side chains of 53 additional amino acid residues on the binding site were set to be flexible throughout the dockings.In order to validate the docking protocol used, it was firstly performed the re-docking of DEET on its crystallographic structure.
Due to the stochastic nature of the docking algorithm, about 10 runs were performed for each compound.In each run, 30 poses (conformation and orientation of the ligand) Figure 1.Structures of the molecules studied.Vol. 24, No. 3, 2013   were returned for the analysis of the overlap with DEET in the binding site of the 3D structure of AgOBP 1 .The best pose of each compound was selected for the further steps of MD simulation.

MD simulations
In order to check the docking results, the best conformation of each AgOBP 1 -ligand complex, obtained from the docking studies, was chosen to run MD simulations.However, before this, it was necessary to parameterize the ligands, so they could be recognized by the force field GROMOS96 from the Gromacs program. 25he parameterization was carried out in the Dundee PRODRG Server 26 and the charge distributions were calculated by the ChelpG method using the Gaussian ® program 27 with the B3LYP functional and the 6-31G(d,p) basis function.9][30][31][32][33][34][35] The AgOBP 1 -ligand complexes were simulated using the GROMACS 4.0 package 25 in cubic boxes of approximately 8,214 nm 3 containing about 17,000 water molecules.These systems were minimized using the force field GROMOS96. 25The minimization algorithms used were steepest descent with position restrained (PR) of the ligands and convergence criterion of 100.00 kcal mol -1 Å -1 , followed by steepest descent without PR, conjugate gradients and finally, quasi Newton Raphson until an energy of 1.00 kcal mol -1 Å -1 .The minimized complexes were then submitted to MD simulations in two steps.Initially, it was performed 500 ps of MD at 300 K with PR for the entire system, except the water molecules, in order to ensure a balance of the solvent molecules around the residues of the protein.Subsequently, 20 ns of MD at 300 K without restriction, using 2 fs of integration time and a cutoff of 10 Å for long distance interactions, were performed.A total of 1,000 conformations were obtained during each simulation.In this step, the lists of pairs (pair lists) were updated every 500 steps, all Arg and Lys residues were assigned with positive charges, and the residues Glu and Asp were assigned with negative charges.
The VMD 36 and Swiss-Pdbviewer 24 programs were used to analyze the structures generated after the optimization and MD steps.Plots of variation of total energy, distance, the variation of the Random mean square deviation (RMSD) and H-bonds formed during the MD simulation were generated with the Origin program. 37Qualitative spatial RMSD pictures were generated with MolMol molecular graphics program 38 and the figures of the frames of MD simulations were generated with PyMOL program. 39

Docking studies
Tsitsanou et al. 14 reported the crystal structure of AgOBP 1 complexed with DEET.Their results suggest that, despite the extensive hydrophobic tunnel connecting the two monomers, DEET has affinity to the marginal region of the tunnel, binding in a site at the interface of the two monomers formed by residues from helix 4, 5 and 6 and some additional residues from the other monomer (Table 1).Our docking studies were directed to this pocket, Tsitsanou et al. 14 reported that DEET performs 57 van der Walls interactions, besides interactions with Gly92, Cys95 and Trp114, mediated by H-bonds with water molecules.
In our studies, the molecules in Figure 1 were docked in the binding sites of DEET reported by Tsitsanou et al. 14 The best poses in each monomer were selected according to the criteria of lower docking energy (intermolecular/ electrostatic and H-bond) and best superposition with the crystallographic structure of DEET, in order to obtain a better comprehension of the differences on the binding modes of each molecule and the factors responsible for the affinity to the enzyme and, consequently, the molecular recognition.
RMSD values of 0.4466 and 0.4903 Å, respectively, were observed in our re-docking results and the overlaps to the 3D structures of DEET are shown in Figure 2.1][42] Figures 3 and 4 illustrate the superpositions of the best poses obtained for the attractive and repellent molecules, respectively.It was possible to observe that all these molecules overlap DEET, competing for its binding site, as they adopt similar conformations and interact with the same residues.These conformations were selected for the additional MD simulation studies.
The docking results (Table 2) suggest that eugenyl acetate presents a greater affinity to AgOBP 1 than DEET, with energies of interaction in each monomer of −111.63 and The other monomer Lys93, Arg94 and Leu96 −117.29 kcal mol -1 , respectively.It also forms H-bonds with Arg94 and a water molecule (Figure 5 and Table 2).DEET presented energies of interaction of −107.66 and −109.51kcal mol -1 , respectively, and established H-bond only with a water molecule (Figure 5 and Table 2).
Eugenol, the main component present in the indian clove (S. aromaticum) oil (and likely the component responsible for its pharmacological activities), presented better docking   results than the other attractive compounds studied (indole, 1-octen-3-ol and lactic acid) with interaction energies of -95.14 and -94.04 kcal mol -1 , respectively (Table 2).
Among the compounds with confirmed or potential repellency capacity, β-caryophyllene and α-humulene presented the worst docking results with no H-bonds neither with amino acid residues nor with water molecules.However, these compounds presented energy values similar to 1-octen-3-ol and better than indole and lactic acid.This suggests that they are probably exploring the non-polar interactions inside the hydrophobic tunnel involving the binding sites of AgOBP 1 .
Among the attractive compounds, 1-octen-3-ol presented the best energy values (-81.95 and -84.17 kcal mol -1 ), interacting with residues His77, Lys93 and Arg94, followed by indole (-73.99 and -75.07 kcal mol -1 ) that interacted with Leu73 (Figure 5).Vol. 24, No. 3, 2013   Lactic acid, despite being an important attractive, presented the worst docking energy results (−48.86 and −48.08 kcal mol -1 ) but was capable of establishing interactions with Leu73, Ala88, Trp114 and a water molecule (Figure 6).Despite having attractive activities allegedly based on interactions with AgOBP 1 , lactic acid, 1-octen-3-ol and indole presented energy values inferior to most of the repellent compounds.This result suggests that AgOBP 1 could have more affinity to the repellent molecules than the attractive ones.That led us to propose that the mechanism of repellency is based on the ability of the repellent to block the binding of the attractive molecules to AgOBP 1 .
Our docking studies also suggest that DEET derivatives with potential repellent activity should present one aromatic ring, like DEET, with longer ramifications in order to explore the possible non polar interactions with the extense hydrophobic tunnel.These ramifications should also have polar groups to explore potential H-bonds with residues of the binding site, improving, therefore, the affinity to AgOBP 1 .

Molecular dynamics studies
After the docking studies, the best poses of the molecules were submited to 20 ns of MD simulations in order to check if these compounds would be estabilized in the pockets, via the interactions predicted by the docking studies, and search for new potential interactions.Besides, it was also searched for more information on the hydrophobic tunnel conecting the two monomers in order to better understand the DEET binding sites of AgOBP 1 and propose structures of DEET analogues as new potential repellents.The plots of total energy (data not shown) during the simulations showed that the energy tends to estabilize in the first picoseconds of simulation for all systems studied.
Temporal RMSD calculations were performed on all the atoms of each complex (protein and ligand) to 1,000 frames at every 20 ps, during the 20,000 ps of simulation.Considering that the complexes could fluctuate in the box, each frame was aligned by the least squares method to the first one for the calculation of the standard deviation.In Figures 7 and 8, it can be observed the equilibrations of the systems formed by AgOBP 1 and the attractive and repellent molecules, respectively.It is possible to observe the equilibrium around the first ps of simulation.This behavior was common to all simulations, with deviations never exceeding 0.30 nm (3.0 Å) and 0.15 nm (1.5 Å) for protein and ligand, respectively, except for eugenol that presented more fluctuations during the first ps of simulations (even superior to AgOBP 1 ) with a RMSD of 0.35 nm until 10,000 ps and, then, 0.45 nm until the end of the simulation time.The temporal RMSD results suggest that the compounds accommodate well inside the binding sites during the MD simulations, showing stabilization of the systems and confirming the results obtained by the total energy calculations previously described.The spatial RMSD for each amino acid was also calculated, as illustrated in Figure 9 for the system AgOBP 1 /eugenyl acetate.It provides a qualitative and quantitative view of all regions of the protein during the dynamics.It can be observed that the regions that most fluctuated throughout the MD simulations (larger values of RMSD and major thickening of the tubes in Figure 9) correspond to the C and N-terminal extremities and the regions of loops.Moreover, the residues in the binding sites, the α-helices and β-sheets presented lower fluctuations, revealing to be the more stable regions, as expected.This behavior was common to all systems simulated.
The docking studies suggest that eugenyl acetate is able to form H-bonds with Arg94 and a water molecule.During the MD simulation, however, this compound H-bonded to Trp114 instead and also formed H-bond with a water molecule, but this H-bond did not remain stable during the 20 ns of simulation.The frames of eugenyl acetate selected during the simulations (Figure 10) show that this compound stays stable in the binding site of AgOBP 1 , interacting with Trp114.Additionally, it was observed that its carbonyl group interacted with Leu73 instead of Arg94.These interactions were observed in both monomers and probably contribute more to the stabilization and permanence of eugenyl acetate in the binding site, justifying the better energy values found in the docking studies.
Regarding DEET, the MD simulations showed the formation of an H-bond with residue Trp114 from 10.0 ns of simulation on.This interaction was not observed in the docking studies.In the same way as eugenyl acetate, an H-bond was also observed with a water molecule during the simulation but this interaction was also not stable throughout the simulation time.The frames of DEET collected during the simulation (data not shown) also showed that this compound remained stable in the binding sites of both monomers.
The MD results observed for eugenyl acetate and DEET were analogues to the results found for lactic acid, in which it was also observed an H-Bond with Trp114 and an unstable H-Bond with a water molecule.However, the docking energy values for lactic acid were much inferior than for DEET and eugenyl acetate (Table 2).Based on these results, it can be suggested that the best affinites of eugenyl acetate and DEET for the binding sites of AgOBP 1 are due to non polar interactions, as suggested before by the docking studies.Regarding eugenol, the MD results showed that this was the compound that forms more H-bonds with AgOBP 1 .It was possible to observe up to 4 H-bonds with the prevalence of 1 to 3, since the beginning of the simulation, involving the residues Leu73 (also suggested in the docking studies), His77 and Lys93.The H-Bond with the water molecule, suggested in the docking studies, was also observed in the MD simulation.
The formation of more H-bonds between eugenol and the binding sites of AgOBP 1 is probably due to the bigger mobility of this compound shown in the plot of temporal RSMD in Figure 7.This can be explained by the approximation to residues Leu73, His77 and Lys93, as can be seen in the plots of the distances between the mass centers of eugenol and these residues in Figure 11.
The attractive molecules 1-octen-3-ol and indole established an average of one H-Bond with AgOBP 1 during the simulated time.According to the docking studies, 1-octen-3-ol interacts with His77, Lys93 and Arg94.The MD simulations confirmed interactions with His77 and Lys93 but interactions with His77 were observed only at the end of the simulated time.Results for indole showed that the interaction with Leu73, suggested in the docking studies, was not observed.However, this molecule was able to interact with Trp114 after 10.0 ns of simulation.Additionally, the MD simulations also showed that 1-octen-3-ol and indole were able to establish interactions with water molecules.

Conclusion
The studies developed in the present work allowed us to identify important features of the binding sites of AgOBP 1 complexed with both attractive and repellent molecules, making it possible to observe the fundamental residues for an effective interaction with this protein.
The analysis of the docking results showed that all the studied compounds compete for the DEET binding site in both monomers, interacting with the same residues.Among the compounds with potential repellency, eugenyl acetate presented the best affinity to AgOBP 1 , even better than the DEET, the compound with the best repellency results reported in the literature. 16Eugenol also looked promissing, considering that it presented the best results among the attractive compounds studied.
The MD simulations corroborated the docking results, showing that these compounds bind independently in both monomers, remaining anchored and estabilized during the 20.0 ns of simulation.
Our studies also reveal that new repellent derivatives from DEET should present longer ramifications on the aromatic ring in order to explore non polar interactions inside the extense hydrophobic tunnel connecting the two monomers.These ramifications should also contain polar groups to explore potential H-bonds with polar residues like Leu73.Another important residue that contributed to the stability of the ligands was Trp114.This residue interacted with all of the studied compounds except eugenol and Vol. 24, No. 3, 2013 1-octen-3-ol.These compounds, however, were able to interact with His77 and Lys93, two other important residues that should be considered further in the drug design of new ligands for AgOBP 1 .

Figure 2 .
Figure 2. Re-docking of DEET in the binding sites of AgOBP 1 .The crystalographic structure of DEET is shown in blue.

Figure 3 .
Figure 3. Superpositions of the best poses of the attractive molecules and DEET inside the monomers (a) and (b) of AgOBP 1 .

Figure 4 .
Figure 4. Superpositions of the best poses of the indian clove components and DEET inside the monomers (a) and (b) of AgOBP 1 .

Figure 5 .
Figure 5. Interactions of eugenyl acetate and eugenol in the binding pockets (a) and (b) of AgOBP 1 .

Figure 6 .
Figure 6.Interactions of the attractive molecules in the binding pockets (a) and (b) of AgOBP 1 .

Figure 7 .
Figure 7. Temporal RMSD for AgOBP 1 and the most promising molecules as repellents.

Figure 8 .
Figure 8. Temporal RMSD for AgOBP 1 and the most promising molecules as attractives.

Figure 9 .
Figure 9. Spatial RMSD for the system AgOBP 1 /eugenyl acetate.Eugenyl acetate was omitted from the figure for better clarity.

Figure 10 .
Figure 10.Frames of eugenyl acetate in the binding sites of AgOBP 1 during the 20.0 ns of MD simulation.

Figure 11 .
Figure 11.Variations in the distances between eugenol and the residues Leu73, His77 and Lys93 during the 20.0 ns of MD simulation.

Table 1 .
Main residues of the binding pocket of DEET inside AgOBP 1

Table 2 .
Docking results for the compounds docked in the binding sites of AgOBP 1