Leishmania infantum 5’-Methylthioadenosine Phosphorylase presents relevant structural divergence to constitute a potential drug target

Background The 5′-methylthioadenosine phosphorylase (MTAP), an enzyme involved in purine and polyamine metabolism and in the methionine salvage pathway, is considered as a potential drug target against cancer and trypanosomiasis. In fact, Trypanosoma and Leishmania parasites lack de novo purine pathways and rely on purine salvage pathways to meet their requirements. Herein, we propose the first comprehensive bioinformatic and structural characterization of the putative Leishmania infantum MTAP (LiMTAP), using a comparative computational approach. Results Sequence analysis showed that LiMTAP shared higher identity rates with the Trypanosoma brucei (TbMTAP) and the human (huMTAP) homologs as compared to the human purine nucleoside phosphorylase (huPNP). Motifs search using MEME identified more common patterns and higher relatedness of the parasite proteins to the huMTAP than to the huPNP. The 3D structures of LiMTAP and TbMTAP were predicted by homology modeling and compared to the crystal structure of the huMTAP. These models presented conserved secondary structures compared to the huMTAP, with a similar topology corresponding to the Rossmann fold. This confirmed that both LiMTAP and TbMTAP are members of the NP-I family. In comparison to the huMTAP, the 3D model of LiMTAP showed an additional α-helix, at the C terminal extremity. One peptide located in this specific region was used to generate a specific antibody to LiMTAP. In comparison with the active site (AS) of huMTAP, the parasite ASs presented significant differences in the shape and the electrostatic potentials (EPs). Molecular docking of 5′-methylthioadenosine (MTA) and 5′-hydroxyethylthio-adenosine (HETA) on the ASs on the three proteins predicted differential binding modes and interactions when comparing the parasite proteins to the human orthologue. Conclusions This study highlighted significant structural peculiarities, corresponding to functionally relevant sequence divergence in LiMTAP, making of it a potential drug target against Leishmania. Electronic supplementary material The online version of this article (doi: 10.1186/s12900-017-0079-7) contains supplementary material, which is available to authorized users.


Background
Neglected Tropical Diseases (NTDs) correspond to multiple transmissible pathologies that mainly occur in tropical and sub-tropical regions. They affect populations living in poverty with more than a billion people in 149 countries worldwide [1]. Here, we focus on leishmaniases, a group of vector-borne diseases caused by different species of protozoan parasites of the genus Leishmania [2]. Three hundred and 50 million people are at risk of infection and 2 million cases are reported worldwide each year [3]. One to 1.5 million cases of cutaneous leishmaniasis (CL) and 0.2-0.5 million cases of visceral leishmaniasis (VL) are reported annually [3]. VL is mainly caused by Leishmania donovani and Leishmania infantum (L. infantum) species, with an annual death toll of 50,000 cases [3]. Mainstay therapy is based on the use of toxic pentavalent antimonials in long treatment courses [4]. Furthermore, their prolonged use is increasingly inducing parasite drug resistance [5]. Second line drugs, such as pentamidine, miltefosine, and amphotericin B also are toxic, costly or have adverse effects [6]. Therefore, the need for new targets and new drugs is increasingly important, and constitutes research priority.
Search of novel potential drug targets mainly focuses on biochemical and metabolic pathways that show differences between pathogens and their host. Purine salvage, polyamine biosynthesis and thiol metabolism are among the most important metabolic pathways being considered for drug development against diseases caused by Trypanosomatidae parasites [7,8]. Some of the most striking differences between parasites and their mammalian host are found in purine metabolism [9]. In mammals, the de novo and/or the so-called "salvage" pathways ensure the synthesis of the purine nucleotides. To the contrary, most parasites studied rely on the salvage pathways for their purine requirement as they lack the pathways for de novo purine biosynthesis [9]. Therefore, salvage purine metabolism constitutes potentially an excellent target for the rational design of antiparasitic drugs. Among the enzymes involved in purine metabolism, 5′-methylthioadenosine phosphorylase (MTAP) plays a crucial role in purine and polyamine metabolism and in the methionine salvage pathway [10]. The 5′-methylthioadenosine (MTA), natural substrate of MTAPs, is generated during polyamine biosynthesis and is then cleaved to adenine and 5′-methylthioribose-1phosphate [10,11], which are respectively incorporated into the salvage pathways of purine and methionine [12]. MTAP, an entry enzyme to methionine salvage pathway, plays an important role to maintain low intracellular levels of MTA, thus to preserve a proper cellular function. Methionine synthesis, polyamine synthesis, protein trans-methylation and trans-sulfuration pathways are excellent targets for chemotherapeutic intervention against African trypanosomes, which are phylogenetically close to Leishmania parasites [13]. MTAP was described as an interesting chemotherapeutic target in African trypanosomes (Trypanosoma brucei brucei), for which selective transition-state analogues were developed. We cite the 5′-hydroxyethylthio-adenosine (HETA), an analogue of MTA, which is highly metabolized by the Trypanosome MTAP in comparison to the mammalian enzyme [10,14]. Growth inhibition assays showed IC 50 values ≤1 μM for HETA, which was selected among possible candidates for in vivo evaluations. HETA exhibited 70% to 90% cure rates when administered to mice infected with T. brucei brucei [10,14]. Moreover, Leishmania major MTAP (LmMTAP) and Trypanosoma brucei MTAP (TbMTAP) have high druggability indexes (0.8, range: 0 to 1) according to the TDR (Tropical Disease Research) Targets Database (www.tdrtargets.org). This explains our interest to such proteins as promising drug targets against diseases caused by Trypanosomatidae parasites. However so far, no study targeted MTAP in Leishmania.
The enzymes that catalyze the phosphorolytic cleavage of the glycosidic bond in nucleosides are structurally classified in two families called nucleoside phosphorylase-I (NP-I) and nucleoside phosphorylase-II (NP-II) [15]. Members of NP-II family share a common two-domain subunit fold and a dimeric quaternary structure while members of the NP-I family, including MTAP, share a characteristic subunit topology with a trimeric or a hexameric quaternary structure. They accept a range of purine or pyrimidine nucleosides as substrates. Multiple MTAP proteins from archaeal, bacterial and mammalian species were characterized on the enzymatic or structural levels [16][17][18][19][20][21][22][23][24]. It was reported that MTAP functions as a dimer in Mycobacterium smegmatis and Mycobacterium tuberculosis [21,22], as a trimer in human MTAP (huMTAP) [18], in Schistosoma mansoni (SmMTAP) [24] and in a putative bacterial MTAP (PDBids: 4GLF and 4GLJ) [25], and as a hexamer in Sulfolobus solfataricus and Pyrococcus furiosus [19,26]. The first structure of huMTAP was solved at 1.7 Ǻ (PDBid: 1CG6) and the enzyme showed a trimeric quaternary structure very similar to that seen in mammalian purine nucleoside phosphorylase (huPNP) [18]. Since then, other human and bacterial MTAPs were crystallized and deposited in the protein database (www.pdb.org). Despite the existence of multiple quaternary structures within the NP-I family members, the subunit fold is highly conserved [15]. It consists of central β-sheets that form a distorted β-barrel, surrounded by several α-helices, characteristic of the Rossmann fold topology mainly found in proteins that bind nucleotides [15].
The active site (AS) in NP-I family members consists of adjacent phosphate-and nucleoside-binding sites (BSs), mainly constituted by residues from the central β-sheets and the interconnecting loops from one subunit, and residues from an adjacent subunit [15]. While the phosphorolysis reaction is common to all enzymes of the NP-I family, there are significant differences in nucleoside specificity within the family members [15,27].
Characterizing the AS of an enzyme is a key step towards the identification of novel and selective inhibitors that may constitute lead molecules. It is further important to understand and elucidate the binding mode of the enzyme substrate. Such knowledge is valuable for the design of transition-state analogues, which appear to be interesting inhibitors in the case of NP-I family members [28][29][30][31][32][33][34][35]. Indeed the co-crystal structure of the huPNP with acyclovir diphosphate led to the design of a series of 9-substituted 9-deazapurine analogues that showed IC 50 values ranging from 17 to 120 nM [28]. A comparative analysis of the binding mode of MTA versus Methylthio-immucillin-A (MTM), a tight-binding transition state inhibitor of huMTAP, onto the AS of huMTAP revealed differential interactions between the natural substrate and its analogous inhibitor [34].
In this context, we aimed to proceed to a comprehensive comparative bioinformatics and structural characterization of the putative L. infantum MTAP (LiMTAP). Primary sequence alignment (PSA), active site prediction on the PSA, and MEME modeling confirmed closer relationships of LiMTAP and TbMTAP, its orthologue in T. brucei, to huM-TAP than to huPNP. This was further confirmed by EC number predictions using the 3D homology model, thus confirming the annotation of the parasite proteins as 5'-methylthioadenosine phosphorylases. The comparison of the 3D homology models of LiMTAP and TbMTAP proteins to huMTAP was thus conducted showing a global conservation of the 2D topology and 3D structure with notable peculiarities like the occurrence of an extra α helix that brings the N-and C-termini to a close proximity, and strikingly different solvent accessible surface areas (SASAs) and electrostatic potentials (EPs) distributions. Different potential binding sites (BSs) were also identified by SOM-BSfinder, save for a cavity corresponding to the active site (AS) that presented a different shape and a larger size on the parasite proteins. Thus, characterization of the binding mode of MTA and HETA into the AS and docking simulations were performed predicting different binding modes and molecular interactions within the active site cavity. In spite of the relative conservation of the AS residues, in line with different EPs and binding modes, these residues interact differently with MTA in the predicted structural alignments. Such differential interactions also explained the selective inhibition of TbMTAP by HETA. Structural differences between LiMTAP and huMTAP were exploited to produce a polyclonal antibody that is specific to LiMTAP.

Generation of the 2D and 3D structures
The secondary structures of LiMTAP, TbMTAP and huMTAP (PDBid: 1CG6) were detected using the STRIDE program [38] under the Pro-origami web server (http://munk.csse.unimelb.edu.au/pro-origami/; [39]). Advanced options were used to exclude 3 10 helices and pi helices. Topology diagrams were then generated and manually labeled in order to indicate residue numbers for each secondary structure element.
The 3D structure modeling of the parasite MTAP proteins was performed through I-TASSER server [40,41] (http://zhanglab.ccmb.med.umich.edu/I-TASSER). Five models were generated for each submitted protein sequence, ranked according to their C-scores as an estimate of their quality [42]. The model with the best (highest) C-score was retained and further refined through the ModRefiner server [43]. Predictions on functional annotations of each protein based on proteins structurally related to the predicted 3D model were provided by I-TASSER. This included Enzyme Commission (EC) numbers and Gene Ontology (GO) terms predictions using COACH [44].

Surface mapping and active site predictions
The electrostatic potential (EP) on protein surfaces was calculated using the Poisson-Boltzman equation using the Adaptive Poisson-Boltzmann Solver plug-in (APBS) implemented in PyMOL (The PyMOL Molecular Graphics System, Version 1.8 Schrödinger, LLC). Calculations were performed on the refined version of I-TASSER 3D models of LiMTAP and TbMTAP, and on the crystal structure of the huMTAP (PDB ID: 1CG6).
An adapted version of SOM-BSfinder [45], a Self-Organizing Maps-based algorithm [46], was used to identify potential binding sites (BSs) on the parasite proteins. The Enamine Golden Fragments (EGF) collection, containing chemically diverse fragments satisfying the Rule of Three [47] was used as probe library. A 3D Self-Organized Map (SOM) of the atomic coordinates of the docked probes and a Unified Distance Matrix (U-matrix) were generated as previously described [45]. Interaction hot spots revealed through the presence of areas with high neuron densities associated with low values in the U-matrix (U-values) were considered as "high neuron consensus"-areas. Regions with low densities associated with high U-value were considered as barriers separating the docking hot spots. We defined a cutoff (t U ) on the U-values to distinguish potential BSs (consensual binding regions with U-values ≤ t U ) from barriers between BSs (regions with U-values > t U ) as: t u = m u + √ v, where m U and v are respectively the mean and the variance of the U-values.

Molecular docking
MTA, the natural substrate of MTAP and HETA, an inhibitor of the TbMTAP (see Additional file 1: Figure S1), were docked targeting the AS of the huMTAP PDB entry 1CG6 and the two parasite refined models for LiMTAP and TbMTAP, obtained from I-TASSER server. Structure data file (SDF) of the MTA and HETA molecules were downloaded from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/) under the accession numbers 439,176 and CHEMBL191917, respectively. Ligands SDF files and receptors PDB files were converted into the PDBQT format using the Open Babel package [48], as follows: (i) hydrogen atoms were added, (ii) 3D atomic coordinates were generated and (iii) Gasteiger atomic partial charges were calculated. Docking calculations were performed using AutoDock Vina 1.1.2 [49] with its default parameters and up to 20 docking poses for each ligand were asked to be generated. Pairwise atomic Euclidean distance was calculated between each protein and MTA or HETA in their best docking pose, defined as the lowest-energy pose according to the scoring function implemented in AutoDock Vina [49]. Distances lower than or equal to a cutoff of 4Ǻ were considered to determine residues involved in proteinligand interactions. All interactions were examined and irrelevant ones were removed.
To extract endogenous proteins, frozen parasite pellets, kept on ice, were resuspended in 1 mL of lysis buffer (50 mM Tris-HCl (pH 7.4), 0.1 mM disodium EDTA), containing 0.05 mM of Phenylmethanesulfonyl fluoride (PMSF) as inhibitor of proteases. The cells were sonicated (4 × 10 s) to reduce viscosity and were then centrifuged for 15 min at 1600 g, at 4°C. The supernatants were dialyzed during 2 h against a buffer containing 50 mM Tris-HCl (pH 7.4), 0.1 mM disodium EDTA, at 4°C, to eliminate the endogenous phosphate as described previously [35].
Peripheral Blood Mononuclear PBMCs were prepared from heparinized blood, collected from one consented healthy donor (who provided a written consent). The study protocol was approved by the local ethical comittee of the Institut Pasteur de Tunis. The PBMCs were collected by density centrifugation through Lymphocyte Separation medium (Eurobio, France). PBMC were washed two times in 10 ml (1×) PBS at 500 g for 10 min and lysed on ice by sonication (2 × 10 s), in presence of 0.05 mM of PMSF.
Protein concentrations of LV50 and PBMC lysates were determined by the Bicinchoninic acid (BCA) protein assay kit (Sigma, Germany) with bovine serum albumin (BSA) as a standard.

Western blot
Four putative antigenic LiMTAP peptides were predicted using Antigenic, a method described previously [50], implemented in EMBOSS. The peptide (AIVTKPEHIPA ETKQRIAPLVASK), located in the C-terminus of LiM-TAP, was used to produce a polyclonal antibody (Genescript, USA). The specificity of this anti-LiMTAP was assessed by western blot.

Results
Closer relationship of the parasite proteins to the huMTAP and common motif patterns Primary sequence alignment (PSA) of the two parasite proteins, huMTAP and huPNP, was performed using T-Coffee program (Fig. 1). The PSA revealed high divergence on the N-and C-termini of the huPNP compared to the three remaining proteins, along with equivalent identity rates (22-24%). The parasite proteins showed higher identity with the huMTAP (37% and 35%, respectively) than with the huPNP. The Leishmania protein has 60% identity with the Trypanosome protein indicating a higher relatedness between these parasite proteins than with the human ones.
Residues involved in the active site (AS) of the huM-TAP [15] were shaded in pink on the PSA (Fig. 1). Among these 15 residues, five were conserved on the alignment among the four studied sequences and corresponded to A94, F177, M196, T219 and V236 in huM-TAP. Six residues were conserved in the three MTAPs (R60, H61, S178, T197, D220 and D222 in huMTAP). The remaining four residues (T18, T93, V233 and L237) appeared to be specific to the AS of huMTAP. Residues T18 and T93 were replaced by a serine (S) and an asparagine (N) in all three remaining sequences, respectively. Residues V233 and L237 were different from their counterparts in the parasite MTAPs and the huPNP. The V233 was replaced by an alanine in both LiMTAP and TbMTAP (A236 and A240, respectively) and by K254 in the huPNP. The L237 was conserved in the huPNP (L261) but replaced by T240 in LiMTAP and R244 in TbMTAP.
Further analysis of the primary sequences was achieved, on the four protein sequences, through motifs identification using MEME. Eight motifs were asked to be found (Fig. 2), as this is the optimal number of characterized motifs found in PNPs and/or MTAPs [15]. Only two motifs were obtained for the huPNP while four to eight motifs were identified for the huMTAP and the parasite proteins. Interestingly, our analysis identified motifs that embedded many motifs previously described as common and characteristic of members of the NP-I family, and more importantly, they covered all the motifs common to the MTAP proteins [15]. The M1 and M3 motifs were conserved in all 4 proteins (Fig. 2a). M1 included smaller motifs previously described as common to trimeric PNPs and MTAPs. However, at its Cterminal part (AA 207-224 in huMTAP), M1 integrated a motif that was previously defined as specific to MTAP proteins [15]. M3 corresponded to a motif that can be found in all members of the NP-I family [15].
Two motifs, M2 and M6, were conserved among the parasite proteins and huMTAP (Fig. 2a). In fact, M2 Fig. 1 Primary sequence alignment (PSA) of LiMTAP, TbMTAP, huMTAP and huPNP obtained using T-coffee program. The PSA indicates regions of sequence divergence relating to the huPNP or to some regions in the N and C termini of the parasite proteins and in their central part. Pink-shaded residues indicate active site residues of the huMTAP and their corresponding residues on LiMTAP, TbMTAP and huPNP, on the alignment. These residues were described in human MTAP crystals as base-binding sites (cercles), methylthioribose-binding sites (stars) or sulfate/phosphate-binding sites (triangles). The parasite proteins and huMTAP share 11 out of 15 residues. The secondary structures (β-strands in green arrows and α-helices in red lines) are shown above the parasite proteins and huMTAP sequences. They align well to each other confirming correct PSA alignment motif corresponded to a motif already described as specific to MTAP proteins, while M6 motif was not described by Pugmire & Ealick [15]. This could be explained by the small subset of proteins studied herein.
Four motifs M4, M5, M7 and M8 that mapped on variable regions on the PSA (Fig. 1) were specific to the parasite proteins (Fig. 2a). Interestingly, M5 motif contained a sequence (IGIIGGTGL in huMTAP) previously described as a common motif to trimeric PNPs and MTAPs [15] but presented divergent Cand N-termini as compared to the huMTAP. In spite of this common sequence (IGIIGGTGL in huMTAP), the program identified M5 motif as specific to the parasite proteins. Motif 7 was located in the central region (AA 130-148 on LiMTAP) and M4 and M8 motifs covered the C-terminal region on LiMTAP and TbMTAP (Fig. 2a). Overall, the motifs that were specific to the parasite proteins (M4, M5, M7 and M8) had more sequence conservation among the two parasite proteins respectively than the remaining motifs (Fig. 2b). The M8 motif presented the most conserved sequence (Fig. 2). Blast search using the parasite specific MEME motifs confirmed the specificity to the parasite proteins (see Additional file 2: Table S1).
All these results highlighted higher similarities between the two parasite proteins and the huMTAP as compared to huPNP, with a higher relatedness between LiMTAP and TbMTAP than with the human proteins. This was strengthened by the presence of specific motifs on the parasite MTAP sequences and specificities within the active site of the huMTAP that may be further exploited.
3D structure modeling revealed parasite MTAP vs. huMTAP differences A thorough structural analysis was generated for LiMTAP and TbMTAP through the I-TASSER server [44]. Interestingly, templates used for building the parasite protein 3D homology models were the same and corresponded to the putative bacterial MTAP structures having PDB IDs: 4GLF and 4GLJ [25]. Five models were proposed for each protein. For both MTAPs, C-score of the first model among the five generated was significantly higher than those of the remaining models (see Additional file 3: Table  S2). Thus, model 1 was considered as the best one in both cases, and was retained for this study. The 3D model of LiMTAP had a C-score of 1.17, an estimated TM-score of 0.87 ± 0.07 and an estimated root-mean-square deviation Fig. 2 MEME models of the primary sequences of LiMTAP, TbMTAP, huMTAP and huPNP. a MEME generated primary sequence models made of up to 8 motifs on the 4 proteins. The location of these motifs, represented as color-coded boxes, is shown on a graphical illustration of the proteins primary sequence. Numbers above these boxes indicate the first and last amino acids of each MEME motif. The p-value of each MEME modeled protein is indicated. b Logo representation of the sequence composition of the MEME motifs in amino acids RMSD of 3.8 ± 2.6 Ǻ. The 3D model of TbMTAP had a C-score of 0.96, an estimated TM-score of 0.84 ± 0.08 and an estimated RMSD of 4.3 ± 2.9 Ǻ.
The parasite protein models were used as query structures to parse the PDB database for structurally related entries. Resulting protein structures were the same for LiMTAP and TbMTAP models, and corresponded either to MTAPs or to PNPs (Table 1). Their TM-scores exceeded 0.5, reflecting a structural closeness to the target proteins (LiMTAP and TbMTAP) as the TM-score is a length-independent metric that measures the global fold similarity between two proteins, with low sensitivity to local structural variations [51]. Moreover, RMSD of atomic positions were somewhat lower for MTAP entries compared to PNP entries (Table 1). Thus, LiMTAP and TbMTAP proteins presented higher similarity to MTAPs than PNPs.
I-TASSER server also provided EC number [52] predictions for LiMTAP and TbMTAP. It returned five homologous enzymes for each target. In both cases, four enzymes had 2.4.2.28 as EC number, which correspond to MTAPs and one enzyme, with the lowest TM-score, having 2.4.2.1 as EC number, which corresponds to a PNP (see Additional file 4: Table S3). This brings additional confirmation on the annotation of these proteins as 5'-methylthioadenosine phosphorylases.
The LiMTAP and TbMTAP models were further refined using the ModRefiner server. Refined models returned lower RMSD values and higher TM-scores in comparison to initial models (see Additional file 5: Table S4). We performed this step in order to verify that the models that will be used for the next steps present energetically minimized structures. We aligned the refined parasite MTAP structure models on chain A of the PDB entry 1CG6, which corresponds to the huMTAP protein co-crystallized with its natural substrate (MTA) and a sulfate ion (SO 4 2− ) as a co-factor [18]. The three proteins aligned perfectly with high conservation of the secondary structures (see Additional file 6: Figure S2), except for the coil at the Nterminus, the coil linking the α2 helix to the β8 strand and the α7 helix located on the C-terminus of the parasite MTAPs (Fig. 3a). These structural differences also matched high divergence regions on the PSA, i.e., the Nand C-termini and the central region covering residues from 150 to 166 in LiMTAP sequence. This added confidence to our primary sequence alignment at these divergent regions, which were also part of the parasite specific M4, M5, M7 and M8 motifs (Figs. 1 and 2).
The models of LiMTAP and TbMTAP also presented high secondary structure conservation with the huMTAP (Fig. 3b). All three structures of the MTAP proteins presented the same number of β sheets (11), whereas the parasite proteins presented one extra α helix (α7) as compared to the human counterpart, which only counts six α helices (Fig. 3b). Each secondary structure element on the sequences of the parasite proteins matched its counterpart on the human protein sequence, as shown in Fig. 1. Noticeably, the α7 helix appeared at the C-terminal region of the parasite proteins, the most divergent region on the PSA, which corresponds to a C terminal extension on these proteins. On the topology diagrams, this helix brings the C-terminal region of the parasite proteins in close proximity to the N-terminus (Fig. 3b). This is noteworthy, as it may induce significant differences as compared to the huMTAP, provinding a basis for potential specific molecular interactions within a monomer or in the multimeric protein. Moreover, the α7 helix overlapped with MEME M8 motif, identified as specific to LiMTAP and TbMTAP (Fig. 2a). In the same line, the motifs identified as specific to the parasite proteins (M4, M5, M7 and M8) mostly covered α helices ( Table 2) that are exposed on the protein surfaces (Fig. 3a). More importantly, the conserved motifs (M1, M2, M3 and M6) mostly covered β sheets (Table 2), described in the literature as implicated in the active site and/or in trimeric contacts within the subunits of the NP-I family members  [15]. On both models of LiMTAP and TbMTAP, these β sheets formed a central ensemble surrounded by α helices as shown in Fig. 3b. This indicated that both proteins belong to the α/β structural class and to Rossmann fold topology according to the PDB annotations (http://www.rcsb.org). This fold is common to all NP-I family members [15], including the huMTAP (Fig. 3a). All these results demonstrated that our modeled parasite proteins are members of the NP-I family with more closeness to MTAPs than to PNPs. Structural peculiarities correlated with MEME modeling of parasite specific motifs and divergent primary sequence are also highlighted.

Identification of differently sized active sites
The three MTAP proteins (LiMTAP, TbMTAP and huMTAP) were further compared by analysis of their surfaces and electrostatic properties. Solvent accessible   (3) TbMTAP. a The cartoon representation of the 3D models integrates the MEME motifs, which are represented in different colors: yellow for the conserved MEME motifs and red for the parasite specific ones. The protein parts that do not include a MEME motif are colored in grey. The secondary structure elements are labeled on each MTAP model. The N-and C-terminal extremities are marked by N and C, respectively. b Only the topology diagrams of (1) huMTAP and (2) LiMTAP are displayed as the TbMTAP diagram is identical to the LiMTAP one. On these diagrams, the α-helices and β-strands are shown by red cylinders and blue arrows, respectively, and are labeled by numbers starting from the N-terminal extremity. The first and last amino acids of each secondary structure are shown on the corresponding shape, on both huMTAP (1) and LiMTAP (2) diagrams. The connector coils are shown as blue lines. On the parasite protein topology diagrams, the α7 helix brings the C-terminal extremity in a spatially close position to the N-terminus surface areas (SASA) were generated for the three MTAPs and the adaptive Poisson-Boltzmann solver (APBS) was used to calculate their electrostatic potentials (EPs) (Fig. 4a). EP distributions on LiMTAP and TbMTAP surfaces presented equivalent ranges, varying from −42 to +42 and from −45 to +45, respectively (Fig. 4a). However for huMTAP, it presented a wider range varying from −58 to +58 (Fig. 4a). As the EP depends on the nature of the protein residues, these differences may be explained by the fact that the parasite proteins presented more identity to each other than to their human counterpart. However, globally the charge distribution on these proteins was different for each one and unique to it. The protein surfaces also contained different cavities and pockets (Fig. 4b). In order to define and characterize the ASs and other potential BSs on the parasite protein surfaces, we used the SOM-BSfinder algorithm [45]. This probe-mapping method identifies potential BSs as "high neuron consensus"-areas on a 3D map. Four BSs were identified for each target (Fig. 4b). For huMTAP, the first site identified by SOM-BSfinder, BS1, (shown in blue on Fig. 4b) matched the AS defined by the location of the MTA on the co-crystal structure (PDBid: 1CG6) [18]. The first site (BS1) identified on LiMTAP (shown in blue on Fig. 4b) occurred in a cavity having the same location as the AS of the huMTAP, when both proteins were aligned (Fig. 4b, see Additional file 6: Figure S2). Thus, BS1 was considered as the AS of the LiMTAP. In the case of TbMTAP, both BS1 (blue) and BS3 (yellow) (Fig. 4b) occurred at the cavity that corresponded to the AS of huMTAP (Fig. 4b). This suggested that both BSs jointly constituted the AS of TbMTAP, as they occupied two parts of the same cavity. Both proteins, LiMTAP and TbMTAP, presented larger AS cavities than the huMTAP one, which may confer to them specific geometric and physical properties.
These results further point to predicted peculiarities affecting the surface of the proteins and particularly the active sites as a result of global sequence divergence and surface charge distributions.

Molecular docking revealed different ligand binding modes and atomic interactions
Given the differences depicted between the three MTAP proteins at the structural level in general and at the AS in particular, we performed molecular docking of MTA on their ASs in order to predict potential Fig. 4 Solvent-accessible surface area (SASA) of the three MTAP proteins and potential binding sites. a The Electrostatic Potential (EP) projected on proteins SASA is presented in red-white-blue color gradient for negatively, neutral and positively charged regions, respectively. The range of its variation is also presented. b Four binding sites (BSs) were identified by SOM-BSfinder on the SASA of the three MTAPs, each having a different surface color (pink, cyan, or green). They were ranked and colored according to descending cluster size from blue (highest neuron density associated with low U-values) to green, yellow and red (lowest neuron density associated with high U-values). The panels represent: (1) HuMTAP, (2) LiMTAP and (3) TbMTAP impact on the binding modes and interactions with ligands.
Docking of MTA on the huMTAP was performed on the corresponding crystal structure (PDB entry: 1CG6). The lowest-energy pose was very similar to the crystal structure with a good docking score of −7.8. AutoDock vina has been identified by the Community Structure-Activity Resource (CSAR) [53] as one of the most interesting tools in identifying the accurate binding mode of a molecule and predicting the corresponding binding affinity. Thus, we considered the lowest-energy docking poses of MTA on the parasite proteins as the most relevant ones. Docking scores of MTA were of −6.1 and −5.8 on LiMTAP and TbMTAP, respectively. These poses correspond to equivalent binding modes of MTA on LiMTAP and TbMTAP, which differ from the binding mode on huMTAP (Fig. 5a). The purine ring of MTA was oriented to the left on the ASs of the parasite proteins and to the right on the huMTAP (Fig. 5a). In both cases, this part of the MTA molecule was interacting with the negatively charged areas on the corresponding SASA (Fig. 4a).
Residues involved with MTA binding onto the ASs of the three proteins, denoted herein as interacting residues (IRs), were defined as those presenting at least one atom potentially involved in an interaction with at least one ligand atom, at a distance lower than 4 Ǻ. This cutoff was set up according to choices used by other groups to define a potential interaction [18,34]. The nature of the interacting atoms was examined in order to validate the interaction type. We identified structurally equivalent IRs (interacting residues that spatially superimposed) when the protein 3D structure/models were aligned (Fig. 5a). On LiMTAP, ten IRs, G17, H61, H65, I92, N93, A94, M200, M221, M243 and V247, had their structural equivalent on TbMTAP: G18, H62, H66, V93, N94, A95, M204, A225, M247 and V251, respectively (Fig. 5a, Table 3). Six residues on LiMTAP, P69, A94, F181, G199, M200 and T201 had their structural equivalent on huM-TAP: P69, A94, F177, N195, M196 and T197, respectively (Fig. 5a). Noticeably, only the two residues A94 and M200 on LiMTAP structurally aligned simultaneously on huMTAP and TbMTAP (Table 3). Obviously, the huMTAP showed ten specific residues that had no structurally equivalent IRs on the parasite proteins: T18, C95, G96, I194, T219, D220, D222, V231, V233 and V236 (Table 3, Fig. 5a). The T18 residue was described in the literature as interacting with the co-factors SO 4 /PO 4 [15,18], while it appeared herein as interacting with the sulfur atom of MTA at a distance~3.85 Ǻ. On the crystals, residues T219, D220 and D222 interact with the purine ring of MTA, and V233 and V236 are involved in hydrophobic interactions with the methylthio-ribose [15,18]. These MTA interactions were confirmed through our docking results. Residues C95, G96, I194 and V231 were not described in the literature but appeared through our protocol as IRs. In our docking results and on the crystal pose of MTA, these four residues, along with T18, appeared at a distance higher than 3.8Ǻ, which may explain why they were not considered as relevant interactions by   The table illustrates interacting residues (IRs) with MTA at a distance lower than or equal to a cutoff of 4Ǻ. Residues listed in the same column are structurally aligned IRs. As we did not dock the cofactor on the proteins, it was expected not to observe interactions on huMTAP involving the cofactor binding sites (R60, H61, T93). However, three cofactor-binding sites ( a ) were here identified as IR with MTA. ( b ) Corresponds to residues here identified as IRs but not on the crystal structure of huMTAP the crystallographers [15,18]. Herein, we considered all possible interactions as relevant in order to minimize false negatives. On our model, the residues of huMTAP interacting with the purine ring of MTA were F177, T219, D220 and D222, which is consistent with literature [18]. Except for F177, these residues had no structurally equivalent IRs on the parasite proteins (Table 3) although they matched perfectly well on the primary sequence alignment. The residues interacting with the purine ring of MTA were G16, G17, R60, H65, I92, N93, A94, M221, M243 and V247, in LiMTAP and G18, V93, N94, A95, A225, M247, N250, V251 and V254, in TbMTAP (Fig. 5a). Of these residues, I92, M221, M243 and V247 on LiMTAP and their counterparts (V93, A225, M247 and V251) on TbMTAP appeared to establish hydrophobic interactions with the two rings of the purine moiety of MTA (Fig. 5a). This type of hydrophobic pocket/interactions have been depicted with other NP-I family members, namely E. coli Uridine Phosphorylase [54,55]. Interestingly, the residue M196 in huMTAP was interacting with the ribose ring of MTA [18]. Its structural counterparts in LiMTAP (M200) similarly established a hydrogen bond with the 2′-hydroxyl of the ribose (Fig. 5a). However, their trypanosomal counterpart (M204) was involved in hydrophobic interactions with the 5′-methylthioribose part of MTA (Fig. 5a).
In a second step, we docked HETA, a specific inhibitor of TbMTAP [10,35], on all three proteins (Fig. 5b). The lowest-energy docking poses of HETA had scores of −7.5, −6.2 and −5.9 respectively on huMTAP, LiMTAP and TbMTAP. These poses presented similar binding modes to the ones obtained for MTA, with opposite purine orientations of docked HETA on the parasite proteins as compared to huM-TAP (Fig. 5b). For LiMTAP, only residue N68 appeared as specifically interacting with HETA vs. MTA, and V247 as specifically interacting with MTA vs. HETA (Tables 3 and 4). For huMTAP, only residue H65 appeared as specifically interacting with HETA vs. MTA and residues R60, H61, T93 and V233 appeared as specific to the interaction with MTA vs. HETA. For TbMTAP, three residues V96, F185 and T227 appeared as specific to HETA vs. MTA (Tables 3 and 4). They interacted with the ethylthio-ribose part of HETA, which includes an additional ethyl group compared to MTA (see Additional file 1: Figure S1), (Fig. 5b). No such specific interactions could be predicted for LiMTAP and huMTAP (Tables 3 and 4, Fig. 5b). This is consistent with HETA being a specific inhibitor to T. brucei [10,35].
In spite of the active site residues conservation with huMTAP, our results indicated binding modes that are specific to the parasites proteins due to the striking differences predicted on the surface of these proteins.
Moreover, qualitative variations in molecular interactions -with different ligands-are also pointed to within the AS itself as a result of the global sequence divergence. This highlights relevance of natural diversity of LiMTAP and TbMTAP in shaping structural and functional differences.
A peptide unique to LiMTAP could be specifically targeted by a polyclonal antibody Based on differences observed at the primary sequence and the 3D structure levels between LiMTAP and huM-TAP, we selected four peptides on LiMTAP presenting high divergence between both proteins of which 3 corresponded to specific MEME (M4, M5, M8) motifs (Table 5). Peptide number 4 (also corresponding to motif M8) was chosen, according to Genscript recommendations, for having the highest solvent accessibility and antigenicity, to generate an antibody that should be specific to LiMTAP. Peptide 4 included amino acids from 277 to 300 (AIVTKPEHIPAETKQRIAPLVASK), covered the α7 helix on the C-terminus of LiMTAP and comprised six residues identified by I-TASSER as highly exposed, namely E283, H284, A287, E288, Q291 and S299 (Fig. 6a). An antibody was successfully generated against this peptide.
The total Leishmania proteins were extracted, separated on 12% SDS PAGE and transferred to polyvinylidene difluoride membrane and then subjected to western blot analysis using the anti-LiMTAP directed against peptide 4. As a control, we carried out the same experiment on human PBMC lysates. The results showed that anti-LiMTAP antibody recognized only the Leishmania protein extract. No signal was detected with the human PBMC extracts, whereas, they were recognized by the anti-β actin antibody (Fig. 6b). Similar observations were also made using lysates of the THP1 cell line (see Additional file 7: Figure S3).
These results showed that the anti-LiMTAP antibody was specific to the LiMTAP protein and suggested the possibility to exploit the structural differences for the development of specific anti-Leishmania biomolecules.

Discussion
Leishmaniases are neglected tropical diseases having a worldwide distribution and till to date there are no effective vaccines available to prevent them [2]. Mainstay and second line drugs currently used for the treatment of leishmaniasis have serious side effects; resistance to antimony or miltefosine is also increasingly reported [5,56]. Therefore, search for alternative drugs to treat leishmaniasis is a research priority. The metabolic pathways of Leishmania, which are either absent or different from the mammalian host and involved in survival, pathogenesis or drug resistance of the parasite, constitute excellent potential targets for the rational design of Table 4 Residues involved with HETA binding to the three MTAPs The table illustrates interacting residues with HETA at a distance lower than or equal to a cutoff of 4Ǻ. Residues listed in the same column are structurally aligned IRs. Residues were observed that interact in the three MTAPs with both MTA and HETA. Residues specifically involved with HETA binding (vs. MTA) are shown in bold antiparasitic drugs [8,57]. Several reports have shown that targeting the polyamine biosynthesis and purine salvage enzymes against Trypanosomatidae have yielded promising results [10,14,35,58]. Among these enzymes, MTAP plays a crucial role in purine and polyamine metabolism and in the methionine salvage pathway [10]. Selective inhibitors of the T. brucei enzyme have been described, and they showed an in vitro cytotoxicity (IC 50 ) of 10 nM for the 5′deoxy-5′-(hydroxyethylthio)-tubercidin and high cure rates (70% to 90%) of HETA when administrated to mice infected with T. brucei [10,35]. However, to our knowledge, in Leishmania, an organism phylogenetically close to Trypanosoma, the protein MTAP was not considered so far as a potential drug target and was not characterized yet. Given the role of this protein in humans and its putative functional conservation in the studied parasites, it seemed important to characterize the protein in Leishmania infantum, a pathogen that causes visceral leishmaniasis in North Africa, Europe, Asian countries and Latin America. Indeed the knowledge of the 3D structure of a drug target protein is of great importance to conduct structure-based drug discovery. Such characterization allows confirming putative roles and identification of commonalities as well as specific features to an organism, and thus assessing whether natural diversity of such conserved proteins could make them potentially good candidates for drug design. When a protein structure is not or could not be resolved experimentally, homology modeling is one of the most powerful tools to obtain robust models of a protein structure [59,60]. Bioinformatics approaches were used to bring insights into the trypanosomal sirtuin structure and function from L. major, L. infantum, T. brucei and T. cruzi. Structure comparisons with the human protein and molecular docking permitted to highlight specificities that were of interest in predicting specific/selective inhibitors [61]. The Leishmania elongation factor alpha, sharing 82% of identity with its mammalian orthologue, was also successfully used in identifying novel anti-Leishmania molecules in silico [62]. In spite of this high identity rate, selective inhibitors could be identified in silico that targeted a unique structural feature on the parasite protein, resulting from a 12 amino acids long deletion [62]. Herein, we used homology modeling to generate the 3D structure models of LiMTAP and TbMTAP and compare them to other MTAP structures of relevance like in humans. For this purpose, 3D models of the parasite MTAPs were performed by I-TASSER, matching structure predictions with known functional templates [41]. This homology modeling server was already used by other groups for instance to understand functions of human thiol dioxygenase enzymes [63] and to assess stability of the Rabies Virus G protein trimer through molecular dynamics [64].
Through sequence and structure comparison of the putative L. infantum MTAP protein with the human and  Fig. 6 Characterization of the polyclonal antibody directed against an antigenic C-terminal peptide of LiMTAP. a C-terminal location of the surface antigenic peptide unique to LiMTAP. The peptide includes amino acids 277 to 300, colored in yellow or brown. Residues identified as highly exposed through I-TASSER (E283, H284, A287, E288, Q291, S299) are colored in brown. b Specificity of the C-terminal Li-peptide was tested by western blot using respectively the LiMTAP antibody, and the β-actin antibody as a control for human PBMC lysates. The total proteins extracted from Leishmania and human PBMC were resolved on 12% SDS-PAGE gel, transferred to PVDF membrane and then subjected to western blot analysis using anti-LiMTAP T. brucei counterparts, we could depict significant global similarities. Important sequence identity rates, level of active site residues conservation, presence of common MEME motifs to NP-I family members, similar 3D topology, and I-TASSER functional annotations and EC number predictions (2.4.2.28) consolidated the hypothesis that the Leishmania protein is a 5′methylthioadenosine phosphorylase. However, in spite of these commonalities shared by the three proteins, it was possible to identify specific structural features that were congruent with divergence on the primary sequence, itself underlying specific MEME modeled motifs. Relevance of such structural specificities was further confirmed by the identification of a highly antigenic and exposed peptide among four structurally divergent regions/peptides corresponding to one of these specific motifs. Notably, a polyclonal antibody directed against this peptide at the C-terminus proved to recognize specifically the Leishmania protein (and so did not react with human cell extracts). Primary sequence divergence, reflected on surface EP and BS predictions, brought additional information about peculiarities of the parasite protein models. Notably, the predicted active sites (located on similar parts of the proteins) presented different shapes and volumes that prompted looking at the protein-ligand interactions. Indeed, these are of high importance towards a comprehensive study of enzymes in general and drug targets in particular [24,61,65]. Molecular docking of MTA, the natural substrate, and HETA, a well characterized inhibitor in T. brucei, into the AS of the three proteins showed equivalent docking scores between MTA and HETA and lower docking scores (better free energy of binding) on the human target as compared to the parasitic counterparts. This could be due to the fact that docking simulations were performed on a crystal structure of the huMTAP, which corresponds to the biologically optimal conformation for ligand binding. Molecular docking of both MTA and HETA highlighted differential binding modes where the purine ring occupied opposite position in parasite MTAPs to the one in the human protein. The interactions into these pockets defined within a 4 Ǻ range, corresponding to hydrogen bonds or hydrophobic interactions, also revealed qualitative differences at different levels in each comparison. Notably, the docking of MTA molecule predicted specific IRs in huMTAP that had no structural equivalent IRs on the parasite MTAPs. Moreover, the docked MTA into both LiMTAP and TbMTAP involved less hydrophobic IRs with the methylthio part and much more interactions with the purine than those seen in huMTAP. In addition, the interactions with the ethyl group of HETA were characterized by three residues (V96, F185 and T227), only seen in TbMTAP and not in its Leishmania and human counterparts. Noticeably, most residues that are unique to Leishmania or Trypanosoma protein belong to kinetoplastid specific MEME motifs. These were essentially mapping on the surface of the protein and embedding α helices, while the ones shared with the huMTAP were mostly within the central β core. Importantly, our analysis also brought structural explanation to the specific inhibitory effect of TbMTAP by HETA [10], through the presence of specific HETA IRs within the TbMTAP AS. This approach could constitute a basis for the design of non-active mutants and/or the design of transition-state inhibitors [24,35]. In line with this, a recent structural study of the MTAP of Schistosoma mansoni (Sm) highlighted structural features that differentiated this protein from human MTAP, bringing basis for intelligent design of novel SmMTAP inhibitors [24].