In Silico Screening and Analysis of Potential Inhibitors of Arylamine N-Acetyltransferases ( NATs ) from the Traditional Chinese Medicine : A Study Using Free Available Tools

Arylamine N-acetyltransferases (NATs) are cytosolic enzymes, highly polymorphic, present in both eukaryotes and prokaryotes. These enzymes play an important role in the detoxification and activation of xenobiotics as well as in the synthesis of endogenous compounds. Specific NATs have been pointed out in the literature as possible therapeutic targets. In particular, the human NAT1, for the treatment of certain cancers, and the NAT from M. tuberculosis (TBNAT), for the treatment of tuberculosis. This paper describes an in silico approach to prospect and select potentially inhibitors of NAT1 and TBNAT from the Traditional Chinese Medicine (TCM) using free available tools. A library with ligands from TCM was previously screened in order to select only compounds with optimal pharmacological properties. The affinity of the selected ligands with respect to NAT enzymes was then evaluated by virtual screening (VS). Subsequently, the complexes with the best ligands were submitted to molecular dynamics (MD) simulations aiming to obtain better quality information on affinity and selectivity. The results for one specific ligand, ZINC14690579, indicated its potential for affinity and selectivity. ZINC14690579 structure may represent the discovery of a new scaffold for future development of NAT inhibitors.


Introduction
Arylamine N-acetyltransferases (NATs) are cytosolic enzymes that catalyze the transfer of the acetyl group from acetyl coenzyme A (AcCoA) to the free amino group of arylamines and hydrazines [1][2][3][4].NATs are highly polymorphic and present in both eukaryotes and prokaryotes [5,6].These enzymes play an important role in the detoxification and metabolic activation of xenobiotics [4,7,8] as well as in the synthesis of endogenous compounds [9,10].
potential biomarkers as well as potential drugs for cancer treatment [12][13][14][15].In another therapeutic perspective, studies showed that NAT gene of Mycobacterium tuberculosis (M.tuberculosis) is part of a cluster essential for the survival of the mycobacteria in the environment inside macrophage [10,16,17].This points to NAT from M. tuberculosis (TBNAT) as a potential drug target for antimycobacterial therapy [18,19].
The challenge for the design of selective inhibitors for eukaryotic and prokaryotic NATs begins with the structural similarity shared by these enzymes.In fact, the fold of human NAT1 and NAT2 closely resembles the overall structure of TBNAT [11,20,21].NATs overall fold is composed of three independent domains of approximately the same length.The first two domains, an α-helical bundle and a β-barrel are joined by a linker helix to the third domain, an α/β lid.For native enzymes, the residues of the catalytic triad (Cys-His-Asp) can be easily superimposed.The modular structure of NATs makes these enzymes highly flexible, adding another challenge for the in silico discovery of inhibitors.For instance, previous studies have suggested that there may be multiple conformational states in the region that defines the active site of NAT enzymes [22], which in itself determines that the interpretation of data obtained from computational techniques necessarily take into account the conformational dynamics in the ligand-receptor interaction.
There is currently a wide range of free resources [23][24][25][26][27][28] and free software packages [29][30][31][32][33][34][35] dedicated to drug discovery and design.The availability of free resources and tools allows scientists from all over the world, independent of the financial support, to plan and test hypotheses through in silico experiments.This paper describes an in silico approach to prospect and select potentially selective inhibitors of NAT1 and TBNAT from the Traditional Chinese Medicine (TCM) using free available tools and resources.A library with 36,043 ligands from TCM was previously screened in order to select only compounds with optimal ADMET (absorption, distribution, metabolism and excretion and toxicological properties) pharmacological properties.The affinity of the 1,902 selected ligands with respect to NAT1, NAT2 and TBNAT enzymes was then evaluated by virtual screening (VS).Subsequently, the complexes with the best ligands were submitted to molecular dynamics (MD) simulations aiming to obtain better quality information on the affinity of the ligands.In particular, the MM-PBSA method combined with MD simulations was used to incorporate conformational fluctuations and entropic contributions to the affinity estimation [36][37][38].One potential selective inhibitor for NAT1 and TBNAT was chosen as a promising compound for future development.

Virtual screening and docking results
Figure 1 shows the 2D chemical structures of the ten compounds selected from the virtual screening.Since of none of the ten best ligands found for NAT1 presented better binding score when compared to NAT2 (Table S1 in Supplementary Materials), the selection criterion was based on the best binding score with TBNAT and worse binding score with NAT2.Thus, only the ten best compounds for TBNAT were locally re-evaluated for the purpose of verifying the results generated by the virtual screening with AutoDock Vina [33] (Vina) software implemented in the MTiOpenScreen [28] web server (Table 1).worse binding energy scores for NAT1 when compared to NAT2 and TBNAT.So far we do not have an explanation for this since the binding residues in the homology model of human NAT1 are spatially conserved in relation to the two crystallographic structures available up to the time of the conclusion of this study.Specifically, are deposited in the Protein Data Bank (PDB) [39,40] structures with ID codes 2IJA [41], which refers to the structure of the mutant NAT1 F125S, and 2PQT, which refers to a structure with the cysteine residue at position 68 substituted by a non-standard S-(2-anilino-2-oxoethyl)-l-cysteine (Tyx68) residue.After superimposition it is possible to see that the only conformational exception are the residues His107 and Phe125 in the crystallographic structure 2PQT, which are influenced by the Tyx68 side-chain (Figure S1 in Supplementary Materials).As consequence, in the 2PQT structure the side-chain of His107 shows a torsion about 25° in chi2 angle and a displacement of 1.9 Å of Phe125 side-chain.Even though these differences may not have substantial influence in the docking results since there is no way to anticipate with which residues a non-substrate inhibitor will bind at the enzyme recognition site.Furthermore, as shown in section 2.2., it was possible to observe a spatial repositioning of the ligands during molecular dynamics simulations due to rearrangements in segments of the main chain around the NAT1 recognition site.Vina [33] implemented in the MTiOpenSreen [28] web server.Each compound is identified according to the registration number in the ZINC15 [27] database.
In Table 1 are highlighted in bold the compounds that presented better binding energy for TBNAT than for NAT2 in both Vina and AutoDock: ZINC14496160, ZINC14715626, ZINC14690579, ZINC53007738 and ZINC13480341.The complexes with these compounds were then selected for analysis by molecular dynamics simulation.

Molecular dynamics simulations results
Aiming to increase the sampling all MD simulations were performed in triplicate.The stability of each replicate was evaluated from atom-positional root mean-square deviation (RMSD) analysis (Figure S1 in Supplementary Materials).The last 10 ns of each replicate were extracted and subsequently concatenated forming a single trajectory of 30 ns for further analysis.The conformational stability of the concatenated trajectories for each complex was estimated from the analysis of the mean value and standard deviation of the RMSD (Table S2 in Supplementary Materials).Except for TBNAT complexes, protein RMSD average values for concatenated trajectories were found between 0.17 to 0.20 nm with a maximum fluctuation of ± 0.02 nm, suggesting stability within the ensemble of conformations considered.Similarly, the RMSD average values for the ligands were found between 0.03 to 0.12 nm with a maximum fluctuation of ± 0.03 nm.Protein TBNAT-complexes presented greater variations of RMSD but with stabilization around 0.3 nm in the last 10 ns of the MD simulations (Figure S2 in Supplementary Materials).Ligand RMSDs for the TBNAT complexes ranged from 0.03 to 0.10 nm with a maximum fluctuation of ± 0.01 nm, except for TBNAT-ZINC14496160 which showed a RMSD of 0.78 ± 0.03 nm.The higher mean value of RMSD and the greater amplitude of the conformational fluctuations observed for the TBNAT complexes can be attributed to a greater flexibility and to more significant conformational changes in the enzyme during the simulations.In fact, it is possible that the greater flexibility observed for the MD simulations of TBNAT-ZINC14690579 compared to those observed for NAT1 and NAT2 (Figure S2 and Table S2 in Supplementary Materials) are related to the existence of an additional seventeen amino acid insertion spanning residues 167-183 in human NATs [21].This region, absent in prokaryotic NATs, is reported to increase NAT protein stability [42].The average number of hydrogen bonds (Table S2 in Supplementary Materials) ranged from a minimal of 0.44 (for NAT1-ZINC14715626) to a maximum of 2.09 (for NAT2-ZINC14496160), indicating that the affinity of the ligands occurred mainly through hydrophobic contacts.
Binding score energies of NAT1 and TBNAT inhibitor complexes were analyzed from the concatenated production phase of MD simulations by MM-PBSA method .Table 2 shows MM-PBSA binding score energies for the five compounds select from VS and docking.It is possible to observe that the compound with the greatest binding score difference between TBNAT and NAT2 was ZINC14690579 (ΔΔGbind = -12.7 kcal/mol).This data is in agreement to Table 1, where only ZINC14690579 presented difference in binding score (ΔΔGbind(Vina) = -5.5 kcal/mol; ΔΔGbind(AutoDock) = -16.6 kcal/mol) greater than the standard error of each docking program (which is of about 2-3 kcal/mol in prediction of binding free energy in cross-validation studies) [32,33].
In an unexpected manner, the difference in binding score for compound ZINC14690579 (ΔΔGbind = -15.7 kcal/mol) also demonstrating a potential selectivity factor between NAT1 and NAT2.Thus, the interaction of the respective ligand with TBNAT and NAT1 was subjected to a more detailed analysis.From Table 3 it is possible to observe that the component that most contributes to the interaction of ZINC14690579 with TBNAT and NAT1 is related to the non-polar/hydrophobic contributions (ΔEvdw).In this sense, the sum of the energetic terms allows to conclude that the estimated selectivity profile (predicted in function of ΔΔGbinding) for NAT1/NAT2 and TBNAT/NAT2 is mainly related to the smaller contribution of the non-polar/hydrophobic interactions of NAT2 with ZINC14690579.Specifically for NAT1, a greater contribution of the electrostatic contributions (ΔEelec) is observed when compared to NAT2 and TBNAT.As will be seen below, this must be related to the contribution of basic residues around the NAT1 binding site.In order to determine the contribution of NAT1 residues involved in the interaction with the ligand, free energy decomposition was performed.The number of hydrophobic contacts and hydrogen bonds between NAT1 and ZINC14690579 was also calculated (Table 4).All analysis and calculations were realized considering the 30 ns concatenated trajectory.From Table 4 it is possible to observe that the residues that contributed the most to the interaction are non-polar (Val93, Ile106, Phe125, Leu209, Val216 and Phe217), which is in agreement with the value of ΔEvdw (-38.1 kcal/mol) in Table 3.These residues account for 54.95% of the total hydrophobic contacts with the ligand (Table 4).From Table 3 it can be observed that the value of the term relative to the electrostatic interactions (ΔEelec) between NAT1 and ZINC14690579 (-15.2 kcal/mol) is more than double compared to NAT2 (-6.0 kcal/mol) and TBNAT (-5.7 kcal/mol).From free energy decomposition data in Table 4 it is possible to suggest that the magnitude of the electrostatic contribution between protein and ligand in the NAT1-ZINC14690579 complex derives from the contribution of a group of basic residues (Lys99, Lys100, Arg165, Arg167, Lys188, Lys220 and Arg286).It is interesting to note that five of these residues (Lys99, Lys100, Arg165, Arg167 and Lys188) are distant from the ligand (Figure 2, Panels (a) and (b)), but together they form a basic cluster that contributes significantly to the interaction.Considered together, the fact that the main contributions to the ligand-receptor interaction in the NAT1-ZINC14690579 complex occur through non-polar contacts and long-distance electrostatic interactions suggests that the driving force that maintains the ligand in the active site is mainly the result of a sum of non-directional interactions.In fact, the number of ligand-receptor hydrogen bonds calculated for the NAT1-ZINC14690579 concatenated trajectory was found small, 1.74 ± 0.83 (Table S2 in Supplementary Materials).From Table 4 it can be seen that the most important hydrogen bond interaction was with the backbone of Phe217 (42.2% of the total hydrogen bonds).The other three amino acids with which hydrogen bonds were detected (Phe125, Ser215, Val216) accounted for 34% of the total of these interactions.
Redocking of ZINC14690579 in a structure of NAT1-ZINC14690579 complex collected at 20 ns from MD simulation resulted in a binding energy score of -8.1 kcal/mol.Compared with the value obtained for the docking of ZINC14690579 on the NAT1 homology modeling on error.This finding may indicate that the crystallographic structures of NAT1 may not be the most suitable for docking studies with homology models built from these.
Comparison of the homology modeling structure of NAT1 and a 20ns structure from MD simulation shows important conformational differences in the region of NAT1 binding site (Figure 2, Panel(c)).In particular, the main chain of the segment formed by the residues Ser95 to Lys99 moves 3.45 Å towards the ligand.Following this movement, the backbone and side chain of Val93, the residue that most contributes individually to the protein-ligand interaction, moves around 1.50 closer to the ZINC14690579.In the same sense, a shift of about 1.7 Å in the backbone of the α-helical interdomain (helices α8-α10) allows the side chains of Leu209, Ser215, Val216 and Phe217 to maximize the interaction with the ligand. 1 Energy contribution estimated by the MM-PBSA method [38]. 2 Percentage of hydrophobic contacts and hydrogen bonds calculated with the LIGPLOT [43] program. 3Hydrogen bonds considered only for residues that had more than 2% of the total interactions.

Molecular dynamics simulations of TBNAT-ZINC14690579 complex
Decomposition of free energy and the counting of hydrophobic contacts and hydrogen bonds were performed in the same way as for the NAT1-ZINC14690579 complex.From Table 5 it is possible to observe that all the ten residues that contributed to the interaction are non-polar.This is in agreement with the value of ΔEvdw (-29.1 kcal/mol) in Table 3. Already, the lack of polar residues and hydrogen bonds involved in the interaction with the ligand is in agreement with the low contribution of the electrostatic interactions (ΔEelec = -5.7 kcal/mol).These observations suggest that the protein-ligand interactions in the TBNAT-ZINC14690579 complex are mainly non-directional.For instance, one of the main elements involved in the interaction between TBNAT and ZINC14690579 derives from a cluster of non-polar residues located between Leu99-Ala100 and Ala103-Pro107 (Table 5; Figure 3, Panels (a) and (b)).It is interesting to note that the interaction with this cluster of non-polar residues is a result of a backbone shift of the Leu99-Pro107 segment of approximately 6.2 Å toward the ligand (Figure 3, Panel (c)).Other conformational changes that contribute to maximize the interaction with the ligand involve displacement of the side chains of Phe130, Ala201 and Phe204 in 2.47, 8.31 and 2.97 Å respectively.In particular, changes in the position of residues Ala201 and Phe204 is accompanied by a movement of the backbone of the third domain of approximately 7.5 Å compared to the crystallographic structure of the TBNAT.The obtained results point to a series of concerted rearrangements that result in the maximization of protein-ligand interactions.These observations can be interpreted in the way that the conformational adaptability of the NATs allows these enzymes to perform various functions [21,22]. 1 Energy contribution estimated by the MM-PBSA method [38]. 2 Percentage of hydrophobic contacts and hydrogen bonds calculated with the LIGPLOT [43] program. 3Hydrogen bonds considered only for residues that had more than 2% of the total interactions.A potential anti-tuberculosis activity of ZINC14690579 can be extrapolated from evidence of inhibitory activity of the extract from Premna odorata Blanco (Lamiaceae), a medicinal plant traditionally used in Albay Province, in south-eastern Luzon, Philippines, against M. tuberculosis [44].
It was observed that the most active fraction of the extract of Premna odorata Blanco has a ZINC1469057 analog, the diosmetin (Figure 4).Although not conclusive, this data combined with the

Library preparation
The library containing 36,043 compounds from the Traditional Chinese Medicine (TCM) was downloaded in the SDF (Structures Data File) format from the ZINC database [27], catalog TCM Database@Taiwan [45].Before the use, the TCM library was previously formatted and cleaned (to remove inorganics, large compounds, mixtures, empty structures, salts and duplicates) with Bank Formatter and Bank Cleaner tools, respectively, both implemented in the Mobyle [46] portal.This procedure was followed by the addition of hydrogens and attribution of the 3D structure for each compound with the Standardizer[47] software.Geometries of the ligands were then optimized with the AMMOS [48] software implemented in the Mobyle portal.Afterwards, the library was submitted to the FAF-Drugs3 [5,6] server.Filtering rules were chosen in order to detect and to remove potentially non-bioavailable compounds, substructures with potential toxic groups, aggregators and Pan-Assay interference compounds [26].At the end of this stage, 1,902 compounds were left in the library.

Target preparation
NAT2 and TBNAT structural models were obtained from crystallographic structures deposited in the Protein Data Bank [39,40] (PDB) (www.rcsb.org):NAT2, accession code 2PFR (resolution 1.9 Å) [50] and TBNAT, accession code 4BGF (resolution 2.1 Å) [51].Since the two human NAT1 structures deposited in PDB have missing atoms and modified amino acid residues, the structural model of NAT1 was obtained by homology modelling.The 3D coordinates of wide-type NAT1 (amino-acid sequence access code P18440-1 of UniProtKB [52] database) was generated using the SWISS-MODEL Workspace [53].As template, the crystallographic structure of mutant (F125S) NAT1 (PDB accession code 2IJA, resolution 1.7 Å) [41] was used.Solvent, ions and ligands were removed from the target structures for the subsequent steps.

Virtual screening
Virtual screening with the filtered library was performed with the AudoDock Vina [33], version 1.1.2,software implemented in the MTiOpenScreen [28] web server.All screenings were performed using NAT1, NAT2 and TBNAT as targets.Each enzyme was treated as rigid and each ligand in the MTiOpenScreen returns only the 1000 best ligands for each target.As for NAT1 no compound had more favorable binding score than for NAT2, the ten best compounds were selected considering only TBNAT, i.e., the ten ligands with most favourable binding score for TBNAT and the worst binding score for NAT2.In order to confirm the results of the virtual screening, a second round of docking simulations were locally performed with the ten best compounds selected from MTiOpenScreen.In this round the dockings were performed with AutoDock Vina [33] implemented in the PyRx [35] graphic interface.Except where otherwise indicated, the standard parameters for the programs were used.

Docking with AutoDock 4.2
Autodock [32], version 4.2, a docking program that uses a different approach to the scoring function in comparison to AutoDock Vina, was used to compare the results obtained from the virtual screening.The docking protocol was the same as described for AutoDock Vina, except for the dimensions of the grid box and the use of the Lamarckian Genetic Algorithm.Targets and ligands were prepared for docking simulations with the AutoDockTools (ADT) [55] interface, version 1.5.6.
A cubic box of 40 × 40 × 40 points with a spacing of 0.35 Å between the grid points was generated for each protein target.The grid boxes were centered as for virtual screening protocol.Global search Lamarckian genetic algorithm (LGA) [56] and local search (LS) pseudo-Solis and Wets [57] methods were applied in the docking search.Each single docking simulation consisted of 10 independent runs.The initial population was 150, the maximum number of generations was 27,000 and the maximum number of energy evaluations was 2.5 × 10 6 .Default values were selected for other parameters.

Ligand setup
The topology for each ligand was generated from the Automated Topology Builder (ATB) [58][59] web server.The conformation with better energy score from the docking simulation with Autodock Vina was taken as input structure.Missing hydrogens were added and the ligand structure was uploaded to the ATB server preserving the original geometry.

Molecular dynamics protocol
Molecular dynamics simulations were carried out using GROMOS96 53a6 [60] force field implemented in GROMACS package [34], version 4.5.6.The protonation states of histidines were determined from H++ [61] web server.All systems were simulated in NPT ensemble and periodic boundary conditions.The dimensions of the central box were chosen so that the minimum distance between of any protein atom to the closest box wall was 12 Å.The simulations were carried out using explicit solvent water molecules described by the simple point charge (SPC) [62] model.
Initially, the protein structure in each system was submitted to a maximum of 500 steps of steepest  All simulations were performed in triplicate, using different initial random velocity distributions, for 20 ns with an integration time step of 2 fs.Each system was heated with gradual increments in the following temperatures: 150 K (100 ps), 200 K (100 ps), and 250 K (100 ps).After, the temperature of each system was adjusted to 309.15 K.The first 10 ns of MD simulations was considered as part of the heating (0.30 ns) and the equilibration (9.70 ns) steps and had not been used in the data analysis.The temperature of each group (protein, ligand, solvent and ions) was independently coupled to a thermal bath with a relaxation time of 0.1 ps using the v-rescale thermostat [63].The pressure in each system was weakly coupled to a pressure bath of 1 atm applying an isotropic scaling and 0.5 ps of relaxation time using the Parrinello-Rahman barostat [64,65].Bond lengths were constrained using the LINCS algorithm [66] with 4th order expansion.Electrostatic interactions among non-ligand atoms were evaluated by the smooth particle mesh Ewald (SPME) [67] method with a charge grid spacing of approximately 1.0 Å.The charge grid was interpolated on a cubic grid with the direct sum tolerance set.Lennard-Jones interactions were evaluated using a 14 Å atom-based cutoff.The pair list was updated every 10 steps.

Analysis of the MD trajectories
The last 10ns of each triplicate was concatenated in a 30 ns trajectory for further analysis.As default, the trjconv module of GROMACS was used to restore, if necessary, protein and ligand inside the simulation box, as well as to make the fit of the conformations of each trajectory on the respective reference crystallographic structure.g_rms and g_hbond modules of GROMACS were used, respectively, to calculate the root-mean-square deviation (RMSD) of protein and ligand and the number of hydrogen bonds between protein and ligand.Biding score energies were estimated using the Molecular Mechanics -Poisson Boltzmann Surface Area (MM-PBSA) method, implemented through the module g_mmpbsa in GROMACS by Kumari and co-workers [38].
MM-PBSA based interaction energies were computed on snapshots collected every 2 ps from the DM trajectories.All calculations with g_mmpbsa were performed using default values.
Protein-ligand hydrophobic contacts and hydrogen bonds were computed on snapshots collected every 100 ps from the DM trajectory using the LIGPLOT [43] program.
Visual inspection of the systems was performed using the Visual Molecular Dynamics (VMD) [68] and Chimera [69] software.

Conclusions
The use of free tools and resources allowed us to conclude the objectives of this study with a potential inhibitor identified.ZINC14690579 showed agreement in all methodologies used, indicating a potential for affinity and selectivity for NAT1 and TBNAT in relation to NAT2.This selectivity is mainly related to the smaller contribution of the non-polar/hydrophobic interactions of ZINC14690579 with NAT2 .The results indicate that the binding of ZINC14690579 to both NAT1 and TBNAT involves a set of concerted conformational changes resulting in accommodation of the ligand at the binding site.These initial observations on the conformational behaviour of NATs are in accordance with the theory that the conformational adaptability of these enzymes allows them to perform various functions, and should be taken into account in the in silico study of new inhibitors.
In this sense, ZINC14690579 may represent the discovery of a new scaffold for future development of NAT inhibitors.

Supplementary Materials:
The following are available online: Figure S1: Comparison between the results for the ten best ligands of NAT1 from virtual screening with MTiOpenScreen and the docking with AutoDock 4.2, Figure S2: Time-dependence of atom-positional root mean-square deviation (RMSD) of Cα atoms for NAT-ZINC14690579 complexes, Table S1: Superposition of NAT1 crystallography and molecular modelling structures, Table S2: Average value of structural parameters of concatenated MD trajectories.

Figure 1 .
Figure1.2D chemical structures of the selected compounds from virtual screening with AutoDock Vina[33] implemented in the MTiOpenSreen[28] web server.Each compound is identified according to the registration number in the ZINC15[27] database.

Figure 2 .
Figure 2. Schematic representation of the interaction of NAT1 with ZINC14690579.(a)Representation of main interacting residues of NAT1 with ZINC14690579.NAT1 is represented as transparent ribbon and sticks with shades of red.The ligand is represented as green sticks.The contribution of each residue to the interaction with the ligand is proportional to the increase in red tonality (see Table4).(b) Representation of NAT1-ZINC14690579 in silico complex.NAT1 is represented as gray ribbon and sticks with shades of red.The color bar indicates the relative interaction energy of each residue with the ligand.(c).Structural changes induced by ligand binding at the recognition site of NAT1 depicted by the superimposition of the crystallographic structure of NAT1 (represented gray ribbon) and the 20ns structure from the MD simulation (represented as blue light ribbon).Ligand structure is represented as gray sticks (from docking) and green stick (from DM simulation).Structure of ZINC14690579 on top of panel (c).

Figure 3 .
Figure 3. Schematic representation of the interaction of TBNAT with ZINC14690579.(a)Representation of main interacting residues of TBNAT with ZINC14690579.TBNAT is represented as transparent ribbon and sticks with shades of red.The ligand is represented as green sticks.The contribution of each residue to the interaction with the ligand is proportional to the increase in red tonality (see Table5).(b) Representation of TBNAT-ZINC14690579 in silico complex.TBNAT is represented as gray ribbon and sticks with shades of red.The color bar indicates the relative interaction energy of each residue with the ligand.(c).Structural changes induced by ligand binding at the recognition site of TBNAT depicted by the superposition of the crystallographic structure of TBNAT (represented gray ribbon) and the 20ns structure from the MD simulation (represented as blue light ribbon).Ligand structure is represented as gray sticks (from docking) and green stick (from DM simulation).Red arrows point to the main conformational changes in the recognition site.Structure of the ligand on top of panel (c).

Preprints
(www.preprints.org)| NOT PEER-REVIEWED | Posted: 30 June 2017 doi:10.20944/preprints201706.0132.v1 in silico observations in this work may indicate that ZINC14690579 scaffold may be potentially effective for inhibition of TBNAT and, consequently, treatment of tuberculosis.

Preprints
(www.preprints.org)| NOT PEER-REVIEWED | Posted: 30 June 2017 doi:10.20944/preprints201706.0132.v1library was treated as flexible.Gasteiger[54] partial charges were calculated after addition of all hydrogens.Nonpolar hydrogens of enzyme and ligand were subsequently merged.A box of approximately 14×14×14 Å was generated for the each target.The grid boxes were centralized on the peptidic O of Ile106 (for NAT1), peptidic O of Val106 (for NAT2) and peptidic O of Thr109 (for TBNAT).

Table 1
also shows the binding score values for the same Preprints (

www.preprints.org) | NOT PEER-REVIEWED | Posted: 30 June 2017 doi:10.20944/preprints201706.0132.v1 ligands
[32] the AutoDock[32]software. Autodock uses a completely different methodological approach for score function calculation in comparison to AutoDock Vina.This allows, taking into account the standard error of each program, a preliminary covalidation of the VS.From VS and docking data it is possible to observe that both AutoDock Vina and AutoDock 4.2 tended to have

Table 2 .
Binding score energies for the five compounds with more selective-inhibitory potential from VS revaluated using MM-PBSA method.The energies are in kcal/mol.

Table 3 .
MM-PBSA calculation for NAT-ZINC14690579 complexes.The energies are in kcal/mol.

Table 4 .
Relative contribution of residues to binding between NAT1 and ZINC14690579

Table 5 .
Relative contribution of residues to binding between TBNAT and ZINC14690579

Preprints (www.preprints.org) | NOT PEER-REVIEWED | Posted: 30 June 2017 doi:10.20944/preprints201706.0132.v1 descent
energy minimization.To relax strong solvent-solvent and solvent-protein non-bonded interactions, 100 ps of MD simulation was performed restraining the protein structure.Initial velocities were assigned according to Maxwell distribution.Counter-ions were added to neutralize each system.Table1summarizes the characteristics of each NAT simulated.

Table 6 .
Characteristics of each NAT simulated.