Abstract
Leishmania donovani is a parasite that causes visceral leishmaniasis, a severe form of leishmaniasis that affects vital organs. An important target for the treatment of this disease is the protein α-β tubulin, which was modeled in this paper and proposed as a target for the treatment of visceral leishmaniasis. Two classes of compounds were studied, dinitroanilines and oxadiazoles. According to the docking results, dinitroanilines interact better with the L loop domain and oxadiazoles interact better with the colchicine domain.
1. Introduction
Protozoan infections are diseases that affect hundreds of millions of people worldwide. In this scenario, leishmaniasis is a neglected disease caused by protozoan parasites of the genus Leishmania. There are at least 20 species of this genus that give rise to different types of manifestations, ranging from cutaneous injuries to visceral leishmaniasis or Kalazar [1, 2]. The visceral leishmaniasis is the most severe form, because the parasite Leishmania donovani usually migrates to vital organs such as liver, spleen, and bone marrow. The parasite exists in two forms, the mobile promastigotes form, found in the intestine of the female vector, the Sandfly, and the amastigote form, found in the mammalian host, which is the cause of the disease [3]. The pentavalent antimonials are the drugs of first choice for treatment of Kalazar, but they have high toxicity, limited effectiveness, and intolerance [4]. The drugs most used in the treatment of the disease are meglumine antimoniate (Glucantime) and sodium stibogluconate (pentosan), but they are injectables [5]. Currently, an orally administered drug, miltefosine, is known; however, it has teratogenic effects and high cost [4]. Microtubules are crucial to cell division and intracellular organelle transport, hence, recognized as an important target for therapy for several diseases [6]. They are made of repeated α-β tubulin heterodimers, which are associated in parallel forming a wall of microtubules, giving rise to the polymer [7]. Existing drugs used as microtubule inhibitors have high toxicity and their use is limited because of resistance [8]. Therefore, there is need for discovery of new drugs superior to the existing.
In this work, the three-dimensional structure of the Leishmania donovaniα-β tubulin (LdTUB) was built by homology modeling and then applied in molecular docking techniques to study the interaction mode of dinitroaniline and 1,2,4-oxadiazole inhibitors of this protein. The molecular details of the interaction of these inhibitors with the LdTUB are quite unknown. Probably dinitroaniline interacts with the N loop alpha subunit [9] and 1,2,4-oxadiazole interacts with the colchicine domain [10].
2. Experimental
2.1. Homology Modeling
The modeling began with the systematic search in the database Expasy [11] for primary structures of the protein α-β tubulin in LdTUB using the server Protein Knowledgebase (UniProtKB) [12] as a search tool. Thus, we obtained the primary alpha (Q25262_LEIDO, 451 amino acids) and beta (I3W8N7_LEIDO, 442 amino acids) sequences. The crystal structure of pig tubulin heterodimer (protein data bank code: 1TUB) was chosen as the template for the construction of the three-dimensional LdTUB model, as it showed the highest sequence identity therewith [13]. Once the template was chosen, the primary sequence of the target protein was aligned with the template protein using the LALIGN server. The model was constructed on a SWISS-MODEL server and the Swiss-Pdb Viewer program [14, 15]. A trimer was constructed in order to analyze the potential binding sites of the protein.
For validation, we used the program PROCHECK [18, 19], which assesses the stereochemical quality of protein structures and the root mean square deviation (RMSD), which superimposes the structures of proteins and calculates the deviation.
2.2. Inhibitors Studied
Were studied two classes of tubulin inhibitors, described by Bhattacharya et al. [16] (Table 1, Compounds 1–11) and Cottrell et al. [17] (Table 1, Compounds 12–18), which are shown in Table 1. The biological activities of these compounds were reported as a negative logarithm of the concentration capable of inhibiting 50% of the tubulin activity (PIC50 = −logIC50).
2.3. Molecular Docking
The three-dimensional structure of each inhibitor was built using the program PC Spartan Pro [20] and its partial charges were calculated by the semiempirical PM3 method. The compounds were docked into the protein LdTUB using Molegro Virtual Docker (MVD) software [21, 22]. The binding site was limited to a sphere with a 9 Å radius. The MVD allows determining the most probable conformation between protein and ligand by estimating the interaction energy. The best poses are those that have lower energy. The function values used in the anchoring (Docking Scoring Function) are defined by the following equation: where
The term highlights the “piecewise linear potential,” which uses two different sets of parameters as follows: one for approximation to the steric term (van der Waals) between the atoms and the other potential for hydrogen bonding. The second term describes the electrostatic interactions between overloaded atoms. It is a Coulomb potential with a dielectric constant which depends on the distance (). is the internal energy of the ligand:
The first term (double sum) is between all pairs of atoms in the ligand excluding atom pairs that are connected by two bridges. The second term is a torsional energy term, where is the angle of twist ties. The average contribution of the torsion link is used if several torsions are determined. The last term assigns a penalty of 1000 if the distance between two heavy atoms (more than two separate titles) is less than 2.0 Å.
3. Results and Discussion
3.1. Homology Modeling
In the first stage of this work, the primary sequence of the 1TUB dimer was aligned with the sequence of the LdTUB dimer using the LALIGN server. The alignments (Figures 1 and 2) showed that α and β proteins have 85.0% and 79% sequence identity, respectively. These values suggest that 1TUB and LdTUB may share a similar overall structure. According to Ceslovas Venclovas, only if the aligned sequences are over 40–50% identical to each other and have few or no gaps can it be expected that alignments may be accurate in a structural sense [23]. Thus, the alignment was considered satisfactory for modeling the LdTUB tridimensional structure. SWISS-MODEL was used to generate the α and β models and the Maestro (Maestro, version 9.0, Schrödinger, LLC, New York, NY, 2009) program was used for energy minimization with the consistent valence force field (OPLS) until the maximum derivative 1.0 kcal mol−1. After that, the alpha-beta-alpha trimer was constructed and energy was minimized.
3.2. Model Validation
Validation is an essential step that can be implemented at different levels of structural organization. One should evaluate the overall quality of the packing of the protein and the possible errors in structural parameters. The Ramachandran plot (Figures 3 and 4) obtained with the PROCHECK program showed that, in alpha and beta subunits, respectively: (i) 74.4% and 66.5% of the residues are in the most favored regions; (ii) 15.7% and 20.0% of the residues are in additional allowed regions; (iii) 5.7% and 7.0% are of the residues are in generously allowed regions; (iv) 4.2% and 6.5% of the residues in disallowed regions. It can be concluded that LdTUB model is of good quality and is energetically favored and stable. Some amino acid residues are found in disallowed regions; however, it was found that these residues are not part of any protein domain, so it was not necessary to perform any calculations to fit these residues in more detail. The superimposing between the template backbone and target enzymes was available using the Swiss-Pdb Viewer 3.7 program [14, 15].
(a)
(b)
Superimposed tertiary structures of proteins’ (backbone) mold and target found RMSD (root-mean-square deviation) of 0.26 and 0.34 Å for alpha and beta, respectively. It is known from the literature that RMSD values of 0.7 Å are found in identical proteins in different crystal forms, and values in the order of 2 Å are found in homologous proteins with 30% identity [24–27]. Figure 4 shows the 3D LdTUB structure proposed by homology modeling utilizing 1TUB. GDP and GTP cofactors that were added to the model and the complex were slowly relaxed to clean the starting structure.
3.3. Dinitroaniline and Oxadiazole Binding Sites
The colchicine binding sites of LdTUB were calculated and a cavity was observed close to Ile367β, Thr365β, Ile352β, Ser351β, Ser350β, Lys349β, Ile348β, Asn347β, Asn346β, Phe316β, Leu315β, Ala314β, Ser313β, Ala312β, Thr311β, Leu256β, Asn355β, Val354β, Ala353β, Leu352β, Lys351β, Arg250β, Leu249β, Asp248β, Ser247β, Asn246β, Leu245β, Arg240β, Leu239β, Thr236β, Val235β, Fly234β, Val181α, Ala180α, Thr179α, Ser178α, Ala143α, and Asn101α. The vinca binding sites of LdTUB were observed close to Ile355α, Gly354α, Cys353α, Lys352α, Phe351α, Gly350α, Thr349α, Lys336α, Ala333α, Ile332α, Ala331α, Ala330α, Asn329α, Val328α, Lys326α, Pro325α, Tyr319α, Leu317α, Leu248α, Ala247α, Leu224β, Pge221β, Thr220β, Pro219β, Tyr207β, Glu204β, Asn203β, Asp176β, Ser175β, Val174β, and Arg173β. Morrissette et al. analyzed oryzalin (Table 1, Compound 1) resistance in T. gondii [28]. They demonstrated that dinitroanilines bind to α-tubulin. The site identified is located beneath the N loop and is composed of Arg2α, Glu3 α, Val4α, Trp21α, Phe24α, His28α, Ile42α, Asp47α, His61, Pro63, Arg64α, Cys65α, Cys129, Thr239α, Arg243α, and Phe244α residues. Thus, this region was also studied and named N loop domain.
3.4. Docking Studies
Energy score values were calculated to get a better understanding of the differences between the binding modes of each compound in known sites of the tubulin protein and the factors that affect its activity. This procedure was performed for the two classes of compounds studied and similar compounds.
Table 2 presents the experimental values of the inhibition % of leishmanial tubulin assembly at 20 mM dinitroaniline analogues and estimated energy score (MolDockScore) between inhibitor and enzyme obtained with the docking in each tubulin domain. The theoretical results obtained in molecular docking were compared with the experimental results described by the Bhattacharya et al. [16]. Linear regression was used to model the relationship between the independent energy term values and inhibition values by fitting a linear equation to the observed data. The theoretical results obtained by the docking study of the dinitroaniline analogs with colchicine and vinca domains were not well correlated ( and ) with the experimental results. However, the results obtained by the docking of the dinitroaniline analogs with N loop domain were correlated (). Thus, the experimental values agree better with theoretical values from N loop domain results than with the vinca and colchicine. By analyzing the hydrogen bond formed between the eleven compounds (1–11) and the LdTUB active site, we observed (i) compounds 1, 5, 9, and 11 exhibited hydrogen bonds with Ser38; (ii) Compounds 1–4, 6, and 7 interacted with Ile42; (iii) Compounds 1 and 5–7 interacted with His28; (iv) Compounds 1, 2, 4, 6, 7, and 11 interacted with Leu26; (v) Compounds 4 and 7 interacted with Phe24; (vi) Compounds 1, 5, 6, and 9 interacted with Lys40; (vii) all compounds interacted with Cys25. It is evident from Table 2 that Compound 7 interacts more strongly with L loop domain than other compounds. That compound is the most potent because the substituent at the N1 position of the sulfanilamide core is bulkier (Figure 5). Through the docking results, it was observed that the ring at position R2 is directed to Asp39 of the active site, interacting thought cation π interaction.
(a)
(b)
(c)
Table 3 presents the experimental values of the antileishmanial pIC50 and estimated energy between inhibitor and enzyme obtained with the docking in each tubulin domain. The theoretical results obtained in the molecular docking were compared with the experimental results described by Cottrel et al. [17]. Similarly to dinitroaniline inhibitors, linear regression was used for energy term and inhibition value relationships by fitting a linear equation to the observed data. The theoretical results obtained by the docking study of the oxadiazole analogs with N loop and vinca domain were not well correlated ( = 0.06 and = 0.43) with the experimental results. However, the results obtained by the docking of the oxadiazole analogs with colchicine domain were correlated ( = 0.80). This suggests that the experimental values agree better with theoretical values from colchicine domain results than with the vinca and N loop. By analyzing the hydrogen bond formed between the eight compounds (12–18) and the LdTUB active site, we observed (i) Compounds 12, 15, and 18 exhibited hydrogen bonds with Lys251β; (ii) Compounds 12, 13, 14, 16, and 17 interacted with Lys359β; (iii) all compounds interacted with Ser257β. The SCN group of compounds 12–18 is directed to GTP (electrostatic interaction). The more potent oxadiazole is 18 (Table 3). It is the only compound with two substituent groups in the ring.
4. Conclusion
In the present work, the LdTUB active site binding mode interactions of dinitroaniline and oxadiazole analogs were studied. According to the docking results, dinitroanilines interact better with L loop domain and the aromatic ring at the N1 position should be explored. Oxadiazoles interact better with the colchicine domain and those analogs with two substituents in the aromatic ring should be explored.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
The authors thank the Brazilian agencies FAPEMIG, CAPES, and CNPq for funding part of this work.