QSAR Studies, Molecular Docking, Molecular Dynamics, Synthesis, and Biological Evaluation of Novel Quinolinone-Based Thiosemicarbazones against Mycobacterium tuberculosis

In this study, a series of novel quinolinone-based thiosemicarbazones were designed in silico and their activities tested in vitro against Mycobacterium tuberculosis (M. tuberculosis). Quantitative structure-activity relationship (QSAR) studies were performed using quinolinone and thiosemicarbazide as pharmacophoric nuclei; the best model showed statistical parameters of R2 = 0.83; F = 47.96; s = 0.31, and was validated by several different methods. The van der Waals volume, electron density, and electronegativity model results suggested a pivotal role in antituberculosis (anti-TB) activity. Subsequently, from this model a new series of quinolinone-thiosemicarbazone 11a–e was designed and docked against two tuberculosis protein targets: enoyl-acyl carrier protein reductase (InhA) and decaprenylphosphoryl-β-D-ribose-2’-oxidase (DprE1). Molecular dynamics simulation over 200 ns showed a binding energy of −71.3 to −12.7 Kcal/mol, suggesting likely inhibition. In vitro antimycobacterial activity of quinolinone-thiosemicarbazone for 11a–e was evaluated against M. bovis, M. tuberculosis H37Rv, and six different strains of drug-resistant M. tuberculosis. All compounds exhibited good to excellent activity against all the families of M. tuberculosis. Several of the here synthesized compounds were more effective than the standard drugs (isoniazid, oxafloxacin), 11d and 11e being the most active products. The results suggest that these compounds may contribute as lead compounds in the research of new potential antimycobacterial agents.


Introduction
Tuberculosis (TB), also known as phthisis, is a chronic severe infectious disease caused by Mycobacterium tuberculosis that affects humans [1]. TB continues to be a serious global health problem, and before the coronavirus pandemic (COVID- 19) it was the main cause of death from an infectious agent, being more deadly than HIV [2]. TB worldwide mortality data show that more than one million people die from this single pathogen annually: TB causes an estimated 10 million new cases and around 1.8 million deaths per year [3]. Moreover, it has been reported that approximately one-quarter of the world population shares the TB infection in a latent phase, representing a high-risk population susceptible to future active TB infection [4]. High morbidity, mortality, and the emergence of new resistant TB strains make this disease a constant source of biomedical research; therefore, a cure for TB has been listed as one of the main goals of the World Health Organization [5,6]. Control of TB in most of the world relies on diagnostics tests such as microscopic-based diagnostics with very low sensitivity and a very old vaccine (Bacille Calmette-Guérin; BCG) and finally, very old therapeutics [7]. Hence, multidrug-resistant (MDR) TB is a growing risk for the population worldwide, with increasingly favorable conditions for the bacteria including the HIV epidemic, other comorbidities such as type 2 diabetes, and low-quality life condition in lower and middle-income countries [5,8].
On the other hand, TB chemotherapy involves first-line oral drugs such as isoniazid, ethambutol, pyrazinamide, and rifampicin for approximately six months; in the first two months of the treatment, all four chemicals are applied, whereas in the last two months, only two antituberculosis (anti-TB) drugs are used: isoniazid and rifampicin. Antibiotic resistance requires a different treatment with second-generation therapeutics [8]. Unfortunately, new TB strains develop resistance to first and second generation [9]. As a result, morbidity and mortality are higher in immunosuppressed HIV patients [10]. Furthermore, the lack of therapeutic adherence due to the long treatment time causes TB-infected patients to abort the treatment, which could increase the resistance of some TB strains [8,11].
The search for new antimicrobial compounds has become a worldwide goal shared by multidisciplinary teams. Recently, bioinformatics teams have focused on compounds designed in silico, thus allowing the discovery of new molecular targets without the experimental laboratory phase. Compound discovery requires investigating the structural features of a compound and its interaction with the cellular pathways that allow it to induce or block a specific human disorder pathway/pathogen. For example, there is a lack of compounds that could deliver efficient results against the TB epidemic; consequently, research and the generation of new chemicals are paramount in the war against TB [11,12]. Novel research could lead to the discovery of several compounds with a high probability of targeting suitable molecules, avoiding the synthesis of inactive compounds and clinical assays, thus reducing the compound discovery time. The different in silico methods, such as quantitative structure-activity relationship (QSAR) and molecular docking, play an essential role in gaining this prior knowledge in compound discovery endeavors [13,14].
In recent years, several compounds have been designed based on cheminformatic tools such as QSAR (quantitative structure-activity relationship) and molecular docking [15,16]. Moreover, additional works have proved that employing these techniques enhances the probability of finding new anti-TB structures, with more probability of success than using only chemical intuition, and further reducing investigation time and costs [17].
Considering the exceptional anti-TB activity of the scaffolds mentioned above, it seems attractive to link both structures, quinolinone and thiosemicarbazide moieties, and prove their potential anti-TB activity. In this sense, this research reports the synthesis of new quinolinone-thiosemicarbazone hybrids and their anti-TB activity evaluation, designed using the combination of quantitative structure-activity relationship (QSAR), molecular docking, and molecular dynamic (MD) studies. Furthermore, the mathematical models herein guided the design of five new compounds that were synthesized and exhibited potent anti-TB in vitro activity against M. tuberculosis.

QSAR Modeling
The minimum energy structures for molecules depicted in Figure 16 (see Section 3) were calculated using Gaussian 16. The minimum energy structures were confirmed using frequency calculation with no negative frequency for all data. In this regard, molecular descriptors were computed using these minimum energy structures together with the anti-TB activity (pMIC). Several multiparametric models were generated using molecular descriptors as independent variables and pMIC as a dependent variable (response variable). The most significant ones are shown in Equations (1)-(3), respectively.  [4][5][6]_v-c_MAS; 75a = I50_B_AB_nCi_2_SS1_T_LGP [1,2,6]_v-h_MAS; 88a = Q1_F_AB_nCi_2_SS6_P_ NSRW_c_MAS; χ = molecular electronegativity, r 2 = correlation coefficient, F = Fisher coefficient, σ = standard deviation, and Q LoO = predictability coefficient (leave-one-out cross-validation). Data analysis revealed that although all generated models were statistically significant, model 3 presented the best combination between statistics and the number of attributes (only three), allowing better handling and interpretation. Therefore, this model was selected for further analysis and to guide the design of new potential anti-TB compounds. The selected model was used to predict biological activity values for a training set of 23 compounds (71.87%) and a prediction set of 9 compounds (28.21%). Figure 2 shows the excellent correspondence between experimental and predicted values for anti-TB activity. In addition, Table S1 (in Supplementary File) shows the values of MICexp, pMIC, their residuals, the corresponding value obtained from model 3, and the difference between MICexp and Q Loo . A linear trend is observed in both diagrams, indicating an excellent internal prediction of model 3. Moreover, the plot of residues vs MICexp and Q Loo residues vs MICexp demonstrates that residues are uniformly closed to zero, corroborating the absence of consequential errors in the amount obtained. Finally, the Pearson correlation between the attributes for model 3 was calculated and is shown in Table 1. According to these results, there is no collinearity between the molecular descriptors present in Equation (3), demonstrating that the equation is not redundant, with attributes independent of each other [41].

Applicability Domain (AD)
The applicability domain analysis was used to detect outliers on the test set and provide reliable and accurate predictions. In this study, the AD was estimated utilizing a William diagram [42]. In this procedure, the leverage threshold (h* = 0.5217) was considered; therefore, a score greater than 0.5217 was used as the criteria to identify if a compound was outside the AD. According to Figure 3, compound 26 fell outside the application domain; therefore, it was removed from the model for future analysis. Thus, according to this AD, the model obtained can be used as a guide for designing hybrids based on quinolines with anti-TB properties since it has a reliable prediction and high confidence [43].

External Validation of the Model QSAR
In Table 2, the different statistical values of the predictive capacity of QSAR model were obtained using the leave-one-out cross-validation (Q LOO ), leave many out cross-validation (Q LMO ), Y-randomization test, and bootstrap (Q boot ). Meanwhile, a r 2 close to 1 means a good fit for the model; the values of Q LOO , Q LMO , and Q boot indicate that the obtained model has high predictability since they have a validation value greater than 0.5 [44]. In addition, Y-scrambling values, 0.1367 for a(R 2 ) and −0.2774 for a(Q 2 ), suggest the statistical model does not correlate by chance. In addition, Tropsha's test validates [45] and the accuracy of predicting pMIC values for model 3. This model met all the test requirements, suggesting an accurate prediction (Table 3). As it can be noted, all of the standard statistical procedures used in this work validate model 3 as a robust, meaningful, and predictable model; therefore, this will operate as a guide for the design of new molecules with potential anti-TB activity.

Design of the Structures
Model 3 indicates that the anti-TB activity might be related to van der Waals descriptors. The parameter has a positive value in the equation, which means that the more positive this descriptor is, the higher the biological activity. The second coefficient in the equation is related to electron density; this descriptor's negative value and magnitude suggest a slight decrease in activity as electron density increases. Finally, the last descriptor is molecular electronegativity; the sign and magnitude suggest decreased activity as electronegativity increases. Considering the importance and the kind of descriptors involved in the QSAR model 3, it started from the most active structure (26, Figure 16), and some modifications were made to find new potential anti-TB compounds.
An important factor to take into account is the free rotation and electronegativity of the atoms present in the thiocarbamide group; free rotation increases the molecular volume, which, according to Equation (3), favors biological activity; on the other hand, the presence of N and S makes the electron density and, therefore, the molecular electronegativity, sensitive to the substituents in this molecular region, which, according to Equation (3), would decrease the biological activity. In this sense, the modifications were focused on the quinoline ring. For the structural change, the simplest isosteric substituents reported in the literature were used and are shown in Figure 4 [46]: H, CH 3 , Cl, and Br. Substitution by methyl groups at positions 6 and 8 slightly decreased the biological activity, while the Cl and Br substituents increased it considerably. This is because the methyl group, which is an electron-donor group, increased the molecular electronegativity and therefore decreased the biological activity. In contrast, the Cl and Br groups, which are electron-withdrawing groups, increased the molecular volume and reduced the molecular electronegativity; therefore, the biological activity was increased in these two compounds (11d and 11e, Table 4).  Table 4 shows the molecular descriptors for the designed compounds according to the selected QSAR model. The values of some parameters used in Lipinski's [47] rule are also shown. As it can be noted, all of the proposed compounds are potentially active against H37Rv anti-TB strains and meet Lipinski rules.

Molecular Orbitals
In the case of antimicrobial activity, several reports suggest structures with, at least, a frontier orbital located in a specific molecular region displaying a broad spectrum of activity [48][49][50][51][52][53][54][55][56][57], compared to compounds whose frontier orbitals are scattered in the whole structures. This latter characteristic is related to the compound fluorescence ability [58][59][60]. In addition, at least one of the two characteristics pointed out above are presented in the anti-TB compounds already reported [58,[61][62][63].
According to Figure S2, the designed compounds could be divided in two kinds: the first one, corresponding to compounds 11b and 11c, exhibiting the HOMO (highest energy occupied orbital) and LUMO (unoccupied orbital of lower energy) distributed among all the structure; in this sense, these compounds might act as nucleophilic (hydrogen acceptor) or electrophilic (hydrogen donor) with the biological target by using both quinolinone and thiosemicarbazide parts on the molecule. In contrast, in the second kind, the electron-withdrawing substituents shift the HOMO to the thiosemicarbazone moiety. In line with these results, the designed compounds meet some characteristics in agreement with the already reported anti-TB compounds. Moreover, compounds 11d and 11e exhibited the HOMO shape more closely than nitro-substituted 1,3-benzothiazoles (BTZ043), a compound family inhibiting the DprE1 protein. These results confirm that the designed compounds meet the electronic characteristic to become molecules potentially active against M. tuberculosis.

Molecular Docking
A molecular docking procedure was carried out to find additional evidence supporting the results obtained from the QSAR method related to the potential antimycobacterial activity of the designed compounds.
In this work, several proteins reported as targets of quinolines and quinolones were evaluated: DNA-gyrase, enoyl-acyl carrier protein reductase, and ATP-synthase [64][65][66][67]; in Table 5 the target with the best scoring values for the designed compounds are shown. The selection criteria were based on the average scoring values for the designed compounds in each protein jointly in comparison with the inhibitor reported for each protein. The protein decaprenylphosphoryl-β-D-ribose-2'-oxidase (DprE1) was also evaluated, given recent interest in it as a fundamental target for developing new antituberculosis compounds because of the structural similarity between its inhibitors and the compounds designed in this research [68]. Table 5 shows a vina score for the designed compounds against proteins InhA and DprE1, which present the best scoring of the designed compounds. The finding of vina score values lower than 5 suggests the possible formation of complexes between the compounds and the selected proteins. According to the obtained results, the designed compounds present suitable affinity towards both proteins, with compounds 11d and 11e having the best performance against both studied proteins.
As shown in Figure 5, all compounds hit the target at the active site. Furthermore, they overlapped the area of action of the compound reported as an inhibitor in each case, with root mean square deviation (RMSD) values less than 3Å. Finally, to obtain more insight into the role of molecular interaction in the ligand-receptor complexes, the molecular interaction for the co-crystallized protein was compared with the 11e protein complex obtained in each case. The respective 2D diagrams are shown in Figure 6.  According to Figure 6, compound 11e interacts with the DprE1 protein analogous to that of the co-crystalized inhibitor (BTZ043). The analysis of molecular interactions in the co-crystallized protein DprE1 suggests that the compound BTZ043 inhibits it mainly by van der Waals interactions where the hydrogen bonds with amino acids Cys387 and Ser59 stand out; these amino acids were reported as fundamental ones in the DprE1 inhibition [69].
In addition, the presence of a π-sulfur interaction (Lys438) can be noted; a strong hydrophobic interaction (Figure 6a). Similarly, in Figure 6b, compound 11e presents the formation of several van der Waals interactions, such as hydrogen bonds (Arg58, Arg54, and Ala53) and a strong π-sulfur type interaction (Cys129). These results show that compound 11e presents the same interactions as the co-crystalized inhibitor but forms different residues within the catalytic pocket. On the other hand, Figure 6c,d shows the molecular interactions between the InhA-NADP complex and compound 11e with the active site of the InhA protein, respectively. It can be noted that the inhibition of this protein by the Inha-NADP complex takes place mainly through dipolar hydrogen bond-type interactions. However, some hydrophobic interactions are essential, such as the formation of π-π bonds and σ-π with amino acids Phe43 and Ile95, respectively.
Interestingly, compound 11e interacted with the active site of the InhA protein in a similar way to the co-crystallized inhibitor. In this regard, hydrogen bonds with amino acids Ala198, Ser19, Ser20, Thr17, and some significant hydrophobic interactions with amino acids Phe43 (π-π) and Ile95(π-σ), respectively, can be highlighted. In general, compound 11e showed excellent values in its vina score and molecular interactions similar to InhA and DprE1 inhibitors already reported, suggesting its potential activity against tuberculosis. Finally, it is noteworthy that the 2D interaction diagrams contrast with the molecular orbitals, indicating that for compound 11e, the thiosemicarbazone fragment has a fundamental role in forming hydrogen bridges. At the same time, the quinolinone moiety would favor hydrophobic interactions with significant energy.

Molecular Dynamics
To evaluate the stability of the complexes formed from the docking calculation and obtain insights into the different interactions that make possible the stability or the lack of it, a 200 ns molecular dynamics simulation was performed for 11c, 11d, and 11e complexes with enoyl-acyl reductase and DprE1 targets. The potential energy, temperature, pressure, volume, and water box size were evaluated to ensure the dynamics ran correctly and reproduced reliable results. For the ligand-enoyl-acyl reductase complexes, the RMSD of the enzyme and the ligands were assessed during 200 ns of simulation ( Figure 7). The results show that the enoyl-ACP reductase reached an equilibrium after 60 ns of simulation at a RMSD value of around 0.3 nm for the three systems. This result suggests that the enzyme does not experience significant conformational changes when interacting with the ligands, which may indicate the possibility of positive interaction with them, translating into a potential inhibition of the enzyme. Analyzing the RMSD of the ligands, values below 0.7 nm were found, which supports the reliability of the docking results previously found. Furthermore, low RMSD values also suggest that the complex formed presents a good stability over time, which supports the possible inhibition. A close inspection at the final step over 200 ns simulation showed that the three ligands stayed in the active site, overlaying between them and with the ligand present in the experimental crystal structure ( Figure 8). To obtain insights on the strength of the interactions, the number of hydrogen bonds (HBs), and the interaction energies (Coulomb and Lennard-Jones) were obtained for the system during the whole length of the simulation (Figure 9). The number of HBs is a good indicator of how strong the ligand-receptor interaction is. As more HBs are formed, the more stable the complex will be. In this sense, compound 11d is the one that may create the most potent interaction with enoyl-acyl reductase, because it formed up to five HBs. Compound 11c formed up to 4 HBs, 1 HB being the most abundant throughout the simulation and between 2 and 3 HBs in the last 75 ns of simulation (when it reaches equilibrium). For compound 11e, although it also managed to make 4 HBs, only 1 HB was present during most of the simulation. These results agree with the Coulomb energy profile, where compound 11d presented the most negative values (around −40 kcal/mol) compared to compounds 11c and 11e, which gave an average value of −10 kcal/mol. The Lennard-Jones potential was similar for the three molecules with values between −30 kcal/mol and −40 kcal/mol suggesting the mode of interaction between the three ligands and enoyl-acyl reductase was similar, which is expected for an inhibition behavior.
The last point of the simulation was analyzed and compared to the (4S)-isonicotinicacetyl-nicotinamide-adenine dinucleotide interaction in the experimental crystal structure ( Figure 8) to obtain more insight into the interaction between the ligands and the enoyl-ACP reductase. In the co-crystal structure (Figure 10a), the dinucleotide formed five HBs with Ile21, Asp64, Val65, Lys165, and Ile194. Furthermore, it had hydrophobic interactions with Ile95, Ile122, and two π-π interactions with Phe43 and Phe149. On the other hand, compound 11c (Figure 10b) showed 2 HBs (Gly14, and Ser94) and 2 hydrophobic interactions (Ile16 and Phe97); compound 11d ( Figure 10c) only formed 2 HBs (Ile95, Gly96); and compound 11e (Figure 10d) 1 HB (Met98) and 1 hydrophobic (Ile95). The dinucleotide managed to form more interactions than the tested ligands, although this was due to the most significant size of the molecule. On the other hand, similar to enoyl-acyl reductase, the DprE1-ligand complexes were evaluated by employing a molecular dynamics procedure, monitoring the results at a different time, and checking their behavior in the potential energy minimization, water box size, the temperature, the pressure, and the volume of the system during the 200 ns. Therefore, these parameters indicated a correct simulation throughout the calculation time, and the different results obtained were analyzed.
First, the pose after the simulation was assessed to see if the ligands managed to stay in the active site or if the interaction was not strong enough and the ligands left the system. DprE1 complex presents an extensive and solvent-accessible active site pocket that allows the entry of FAD and the ligand (Figure 11a). Interactions in open, active sites are less than those in close pockets, requiring stronger interactions to keep the ligand inside. For the three systems studied, the ligands and FAD managed to stay in the active site after the simulation time, suggesting a strong interaction between them and DprE1. Looking at the final pose of FAD in the three systems, their structure was overlaid with FAD in the experimental form (Figure 11b). This result was expected because cofactors-enzyme complexes are very stable.
Furthermore, due to FAD's large size, which manages to penetrate deep inside the enzyme, and the multiple functional groups that can quickly form hydrogen bonds with the enzyme, a stable and robust interaction was very likely to be developed, constraining its movement. On the other side, comparing the final pose of the ligand in the active site with the experimental one found in 6HEZ, a good overlay was not found, although they managed to stay inside the pocket (Figure 11c). This behavior was probably related to the lack of interactions between the ligand and the enzyme. Still, the few interactions with the enzyme and possible interactions with the cofactor may be strong enough to keep the ligand in the active site inhibiting the enzyme. RMSD of the ligands ( Figure 12) shows that compound 11c had a significant change in its conformation compared to the result of the docking calculation. Still, values did not reach 0.3 nm in any of the systems, suggesting good stability in the enzyme-ligand complex. For ligands 11d and 11e, their RMSD values were less than 0.15 nm during most simulation time. For DprE1 in the three complexes, values were very similar, between 0.3 nm and 0.5 nm, during the simulation. In system 11e, a peak was observed at 130 ns, suggesting a particular distortion that disappeared in the last 50 ns of simulation where the RMSD kept stable; therefore, the ligands studied do not produce critical conformational changes to the enzyme, indicating a similar binding mode to that of the co-crystallized ligand. In addition, short-range Coulomb and Lennard-Jones energies were obtained throughout the simulation ( Figure 13) to complement the stability analysis of the complex studied.  As shown in Figure 14, all energies were negative, indicating a suitable attractive interaction between ligands and the target. Regarding Coulomb energies, compound 11c's energy was 10 kcal/mol lower than that of the other two compounds, suggesting a more robust interaction. For Lennard-Jones potential, 11c and 11d presented a similar profile with values around −35 kcal/mol, while 11e energy stabilized during the second half of the simulation at around −20 kcal/mol. Finally, a hydrogen bond (HB) analysis was performed on the three systems ( Figure 14). The results showed that 11d was the ligand that formed the more significant number of HBs during most of the simulation, followed by 11e, and at the end 11c, which agrees with the most negative values found in Coulomb and Lennard-Jones energies. During the 200 ns simulation, 11d managed to form 4 HBs throughout the simulation, becoming 5 in more than half the time. Furthermore, 11d formed 6 and even 7 HBs during a short period. For 11e, 3 HBs were created during most simulation time, with 6 HBs the highest, but only for a short time. Finally, 11c formed only 1 HB in the simulation, while 2 and 3 HBs were formed after half of the time. Interestingly, although it was the compound with the least number of HBs during the simulation, at one point of the calculation, it created 8 HBs, which was the highest of all the systems.
The 2D representation of the final conformation of the complex was analyzed to see which residues were the most important for the interaction (Figure 15). The interactions between FAD and DprE1 were also plotted and HBs with 14 different amino acids were found, explaining the stability of FAD during the simulation. For BTZ043 (experimental), there was one HB with Lys438 and hydrophobic interaction with Val365 where the trifluoro group is located. Compound 11c presented three interactions; 2 HBs (His315 and Leu317), and 1 hydrophobic interaction with Val365. Compound 11d showed only 2 HBs; with Gln336 and Lys438, while 11e had only 1 HB with Thr118. Along with the HBs, any interaction with FAD may also help the stability and inhibit DprE1. Therefore, distances between the different ligands and FAD were measured. Compound 11d presented the shortest average distance (2.3 Å), followed by 11e (3.7 Å), and 11c (7.5 Å). The distance between ligand BTZ043 and FAD for the experimental crystal structure was 4.5 Å. The ligand-protein interaction was analyzed using the MMPBSA approach (Kumari et al. 2014). In this method, the free binding energies (∆G) are computed considering the van der Waals, electrostatic, and solvent-accessible surface area (SASA) energies (Table 6). According to Table 6, all ∆G values were negative, indicating both favorable and stable interaction between the ligand and the receptor; except for 11d, the more significant contribution to the ∆G were van der Waals interactions, followed by electrostatic energies, suggesting all of the interactions were driven mainly by dipolar interactions. Interestingly, the complex formation for compounds 11d-enoyl-ACP reductase was mediated by electrostatic contribution (−77 kcal/mol). Finally, according to Table 6, the compounds interacted more significantly with enoyl-ACP reductase than with DprE1. However, the suitable binding energies (from −12.7 to 71.3 kcal/mol) motivated the synthesis and evaluation of all the designed compounds. The difference in binding energy values obtained for both proteins suggests selectivity in molecular recognition; thus, the compounds designed in this work would have more affinity for the InhA protein than DprE1. On the other hand, although these results are not absolutely indicative of a possible mechanism of action, they do provide additional support for their synthesis and subsequent in vitro testing.

Chemistry
To experimentally validate the results obtained from the computational design, many authors proceed to perform the synthesis of five structures or less [70][71][72][73]. Therefore, Table 7 shows the route employed to prepare the designed compounds using the QSAR procedure, supported by docking and dynamic molecular results. The 3-formyl-2-oxo-quinoline precursors 15a-e were synthesized as presented in Table 7, following the directions of the reported procedures [74]. Two transformations prepared the quinolinone-thiosemicarbazone hybrids 11a-e. The first consisted in the preparation of the chalcone analog derivatives 16a-e by a Claisen-Schmidt condensation between the respective 3-formyl-2-oxo-quinoline 15a-e and an excess of acetone in the presence of 20% NaOH [75]. The mixture was stirred at room temperature for 15 min until all the aldehyde 15 was consumed (monitored by TLC). Then, the reaction mixture was neutralized with AcOH and the precipitate formed was filtered off and washed with a mixture of water/methanol, giving the appropriate chalcone analog derivative 16a-e in 65-88% yields. The second transformation consisted of preparing thiosemicarbazones 11a-e following the general synthetic procedure (Table 7), where the chalcone analog 16 and thiosemicarbazide in EtOH with the addition of AcOH as the catalyst were heated at 80 • C. The synthesized compounds 11a-e afforded yields ranging from 80 to 92%. The structures of the quinolinone-thiosemicarbazone hybrids 11a-e were characterized by IR, NMR (in DMSO-d 6 ), and mass spectrometric analysis.

Antimycobacterial Activity
All the newly synthesized compounds were tested against two different acid-fast slow-growing mycobacteria, M. tuberculosis H37Rv and M. bovis BCG. Isoniazid was used as a reference control, and results were expressed as minimum inhibitory concentration values (MICs, µM) and are presented in Table 8  It is important to note that the experimental values shown in Table 8 (M. tuberculosis H37Rv) are in line with the values predicted by model 3, with compounds 11d and 11e being the most active. The excellent biological activity of the synthesized compounds is attributed to the incorporation of the thiosemicarbazone system in the quinolinone, since it generates an increase in the molecular volume (parameter 24a, see Table 4), specifically due to the free rotation of the NH 2 group (due to its sp 3 hybridization). Added to the above, N and S atoms provide a decrease in molecular electronegativity (χ), which could cause the narrow MIC value among the synthesized products (similar values for each product, see Table 4). These latter results support the hypothesis that the distribution of HOMO in the thiosemicarbazone part seems to be a structural factor related to the antimycobacterial activity.
Complementing this work, the antimycobacterium biological activity for compounds 11a to 11e was tested against four genetically different lineages of M. tuberculosis multiresistant to isoniazid and rifampicin: a M. tuberculosis orphan strain with no match in spolDB4 [76]; M. tuberculosis Beijing, an emerging pathogen in several areas and a predominant endemic strain in others that is frequently associated with drug resistance [77]; M. tuberculosis LAM 9, a predominant genotypic family belonging to the Latin American and Mediterranean (LAM) lineage [78]; and M. tuberculosis Haarlem, whose genotype is ubiquitous worldwide and represents about 25% of the isolates in Europe, Central America, and the Caribbean, suggesting a link with the post-Columbus European colonization [79]. Haarlem strains are actively transmitted in urban settings in Colombia [79]. Additionally, two reference strains, M. tuberculosis rifampicin resistance strain ATCC 35838 and M. tuberculosis isoniazid resistance strain ATCC 35822, were used as internal quality control (Table 9). All compounds exhibited activity in the qualitative range of good to excellent against the different families of multiresistant M. tuberculosis, with MIC values in the range 0.27-66.58 µM (Table 9). Interestingly, many of these compounds showed better activity than the standard drugs (Oxafloxacin, MIC = 2.76 µM). Furthermore, results in Table 9 show that M. tuberculosis Beijing was the most susceptible to compounds developed in this work, with MIC values of 0.27-1.75 µM. Noteworthily, compounds 11c-e showed the best activity against most of the tested M. tuberculosis multiresistant strains. Compound 11e, having a bromine substituent on the quinolinone ring, showed the highest activity against M. tuberculosis Beijing, being tenfold more potent than Oxafloxacin. On the other hand,

Cytotoxicity
The cytotoxic activity of compounds 11a-e was determined by an MTT assay on Vero cells (ATCC ® CCL-81 TM ) with NaDS as a positive control and dimethylsulfoxide (DMSO) as a negative control. The growth of the Vero cell line was expressed as half maximal inhibitory concentration values (IC 50 ). Therefore, the selectivity index (SI) was calculated as the ratio between IC 50 in Vero cells to the concentration that can kill the mycobacteria. The results are summarized in Table 10, revealing that compounds 11a-e were highly cytotoxic with IC 50 values in the range of 1.91-3.53 µM. Compounds 11b and 11e presented minor toxicity. However, many of these compounds showed a selectivity index of >10 against M. bovis BCG and M. tuberculosis H37Rv strains, 11b and 11e being the most potent, indicating that these compounds selectively inhibit M tuberculosis more than normal cells.  In summary, the results obtained in this work suggest that the combination of QSAR and molecular dynamics guided the synthesis of the different compounds with anti-TB activity; moreover, in some compounds, the bioactivity against other strains was more significant than that of the standard compounds. Finally, all the compounds developed in this work showed a selectivity index against M. tuberculosis H37Rv >10, with compound 11e showing the best performance.

Data Collection
The selection criteria were based on compounds bearing quinoline and/or thiosemicarbazone moiety in their structures; furthermore, the guidelines recommended by Alexander Tropsha for the construction of the virtual library were followed [80][81][82][83]. In addition, we calculated the minimum inhibitory concentration of the biological activity using the same protocol. In this work, quinolines are the central nucleus; therefore, a search for antituberculosis quinolines compounds was performed using databases such as SciFinder and PubChem. Subsequently, we searched for thiocarbamides with biological activity against tuberculosis. However, we encountered two obstacles: most of the thiocarbamides were coupled to nuclei other than quinolines outside the application domain. Some thiocarbamides were reported using antituberculosis activities measured in other than IC 50 or mg/mL.
Under the criteria explained above, 32 compounds, including three coupled thiocarbamides and quinolines, met the selection criteria. Thirty-two compounds were selected from the bibliography [84][85][86][87][88] and their structures are shown in Figure 16.
It is important to note that the mechanism of action of the selected structures was not taken into account for the construction of the data, since it is not a requirement for the elaboration of QSAR models according to the different research related to this methodology [13,[89][90][91][92][93][94][95][96]. Nevertheless, the structural diversity of the compounds studied here, with molecular structures similar to compounds already reported as antituberculosis (fluoroquinolones and quinolones), in addition to the great diversity of isosteric substituents on the pharmacophoric nuclei, guarantee a domain of application suitable for the rational design of new drugs.

Minimum Energy Structures and Molecular Descriptors Calculation
The structures were drawn using Avogadro software ( Figure 16) [97], and a minimum energy conformer research was applied. After finding the minimum energy conformers for each structure, these were changed to gjf format and optimized at B3LYP/6-31(d,p) theory level employing Gaussian for Linux [98]. The minimum energy structures were confirmed using frequency calculation. The combination of density functional theory with the set base 6-31G(d,p) was supported by its extraordinary capacity to reproduce experimental results in molecular reaction parameters and QSAR approximations [99][100][101][102][103][104].

Descriptor Calculations
Once the minimum energy structures were assured, three molecular descriptors were computed for the entire data: electronic [105], thermodynamics [106], and topographic [107]. First, electronic descriptors were calculated using conceptual density functional theory [108].
All the electronic molecular descriptors were based on frontier molecular orbitals taken from the output of the optimized structure. Hence, HOMO (highest occupied molecule orbital) and LUMO (lowest unoccupied molecular orbital), were used to calculate other electronic descriptors such as molecular ionization potential (IP), electronic affinity (EA), electronegativity (χ), electronic potential (µe), electrophilic index (ω), molecular hardness (η), molecular softness (s), using the following Equations: IP and EA are ionization potential and electronic affinity, respectively; in both cases, parameters were calculated by assuming that the electronic geometry of the cation or anion remained equal to that of the neutral species [109]. Thus, both were calculated employing Equations (9) and (10).
Thermodynamics descriptors are those related to Gibbs equations: molar refractivity [109], enthalpy (H), entropy (S), Gibbs energy (G), and molar heat capacity; all of them were calculated for the entire data at 25 • C employing the frequency calculation outputs. On the other hand, the topological molecular descriptors such as accessible area, molecular area, solvent excluded volume, H acceptor atoms, H donor atoms, ovality, partition coefficient, Balaban index, cluster counter, molecular topological index, rotatable number of atoms, polar surface area, radius, attribute form, coefficient form, the sum of degrees, sum of degrees of valence, topological diameter, total connectivity, total connectivity of valence, and Wiener index were calculated employing Chembioffice 18.1 software [109]. Finally, topographic 3D molecular descriptors suggested for cheminformatics studies such as QSAR were calculated using QuBiLs-MIDAS [110]. A total of 3031 3D molecular descriptors were estimated based on many algebraic operations considering linear, bilinear, and quadratic indexes for all molecules in the dataset [110].

QSAR Models
The computed molecular descriptors for the entire data were tabulated and used as the independent variable; the anti-TB activity was used as a dependent variable (pMIC = −log(1/MIC)), where experimental MIC (MICexp) represent the minimum inhibitory concentration (in mol/L). A multiple regression analysis (MRA) was performed using the IBM SPSS Statistics 25 software for Windows to find the most molecular descriptors correlated with the biological activity. Different variable selection methods were used for introduction, stepwise, elimination, backward, and forward selection [110].
In all cases, the obtained models were considered "meaningful" by analyzing the statistic quality; statistical parameters such as the regression coefficient (R), the Fisher value (F), and the standard deviation (s) were computed; as already known, the higher the r 2 values, the better the model. In the same line, F values correlated variance (r 2 ) by the number of degrees of freedom with the unexplained variance (1 − r 2 ) by the number of variables in the model. Thus, the higher the percentage of variance explained by the model, the higher the F and, consequently, the better the model. A Pearson correlation was employed to investigate any correlated molecular descriptors in the models. A Pearson's coefficient higher than 0.7 was considered a "strong association" between two independent variables; consequently, one of the correlated independent variables was ruled out.

Applicability Domain (AD)
For the QSAR building and validation, the applicability domain analysis is a fundamental requirement for determining the model's reliability and discarding any incorrect dependent variable. It can be defined as the theoretical spatial section limited by the structure, biological activity, and the descriptors selected in the regression model [42,111]. The AD analysis was carried out using a William plot procedure [112].

Validation Procedures
The robustness of the generated models was examined by employing several statistical validation procedures such as the leave-one-out cross-validation and leave-many-out crossvalidation methods (Q LOO and Q LMO ); for the external validation bootstrapping [113] and y-scrambling [114] were employed. Finally, the Tropsha test criterium was applied to the best models [115].

Molecular Docking
Molecular docking is a versatile tool that allows finding minimum energies of proteinligand interaction structures [116] and has been used as a standard tool for designing and synthesizing new potential anti-TB compounds [69,89,[117][118][119][120][121].
The structure of proteins with their inhibitors was cleaned by withdrawing any water, inhibitor, and non-native ligands molecules. Then, the proteins cleaned were redocked against their inhibitors to ensure that the docking protocols' output showed similarities to the protein retrieved from the PDB web page. For this aim, the autodock vina tool was used [125]. The proteins' PDB ID, active sites, and inhibitors used in the redocking protocols and grillaBox sizes are shown in Table 11.

Molecular Dynamics
Molecular dynamic simulations were carried out to evaluate the different ligandreceptor complexes' stability over time and to obtain insights into the interactions formed. For the calculation, two enzymes were chosen, decaprenylphosphoryl-β-D-ribose 2'epimerase and enoyl-acyl reductase, while the ligands 11c, 46d, and 46e were used. Furthermore, the most favorable conformation of the ligands in the enzymes resulting from docking calculation was chosen as input for the molecular dynamics simulation. For all the calculations, Gromacs 2019 [126] was employed, AMBER99SB-ILDN [127] force field implemented in Gromacs 2019 was picked to build the topology of the enzymes.
On the other hand, the topology was built using the ACPYPE server [126] and the generalized AMBER force field (GAFF) for the ligands. Before the simulation, the enzymeligand complexes were solvated using a three-point water model (TIP3) employing a cube shape and neutralized using chlorine or sodium atoms. Then, the complex was relaxed and equilibrated using constant NVT (number of particles, volume, temperature) and NPT (number of particles, pressure, and temperature) protocols for 100 ps at 300 K [128]. Once the system was equilibrated, the 200 ns simulation was carried out, setting the temperature at 300 K and the pressure to 1 bar [128].

Molecular Orbitals
In the past twenty years, frontier orbitals have played a pivotal role in understanding the molecular interaction between ligand-receptor and a specific biological activity [129][130][131][132][133][134][135][136][137]. Typically, the shape and electron density of boundary orbitals indicate the most reactive regions of a molecule. However, in other cases, compounds with the same biological activity follow a similar pattern of boundary orbitals [48]. The frontier molecular orbitals for the design structures were retrieved from the optimized output using the B3LYP/6-31G(d,p) level of theory; the isosurfaces were generated using 0.02 eV isovalues.

Chemistry
Reagents and solvents were purchased from Sigma (St. Louis, MO, USA) and Merck (Darmstadt, Germany). Melting points were measured using a Stuart SMP3 melting point apparatus. IR spectra were recorded on a Shimadzu FTIR 8400 ATR spectrophotometer. Mass spectra were run on a SHIMADZU-GCMS 2010-DI-2010 spectrometer (equipped with a direct inlet probe) operating at 70 eV. In addition, 1 H and 13 C NMR spectra were recorded on a Bruker Avance 400 spectrophotometer operating at 400 MHz and 100 MHz, respectively, using DMSO-d 6 as solvent and tetramethylsilane as the internal standard. TLC analyses were performed on 0.2 mm pre-coated aluminum plates of silica gel 60 F 254 (Merck, Darmstadt, Germany) and spots visualized were with ultraviolet irradiation.
General procedure for the synthesis of the (E)-3-(3-oxobut-1-en-1-yl)quinolin-2(1H)ones 16a-e. A mixture of acetone (10 mL, 136 mmol), the corresponding 3-formyl-2-oxoquinoline precursors 15a-e (1.0 mmol) and 20% aq NaOH (5 mL) was stirred at room temperature for 15 min and the evolution of the reaction was checked by TLC. Once the reaction was complete, the resulting precipitate was collected by filtration under vacuum and washed several times with water and ethanol to afford 16a-e. No further purification was required.

Antimycobacterial Activity
Antimicrobial activity, was measured by spot culture growth inhibition assay (SPOTi) [138]. In order to evaluate the minimum inhibitory concentration (MIC) values for the different synthetic compounds 11, these were tested against the laboratory strains: M. tuberculosis H37Rv, M. tuberculosis lineage orphan, M. tuberculosis lineage Beijing, M. tuberculosis lineage LAM 9, M. tuberculosis lineage Haarlem, M. tuberculosis ATCC 35838 (rifampicin resistant) and M. tuberculosis ATCC 35822 (isoniazid resistant). Biological activity was tested in a level 3 biosafety laboratory of the Colombian National Health Institute. M. bovis assays were performed at Universidad del Norte at the biological and chemical research laboratory. Briefly, the compound was dissolved in dimethylsulfoxide (DMSO) in a 200 mg/mL stock solution; from this, a serial dilution was made in a sterile 24-well microplate, allowing for evaluation concentrations of 100-0.01 mg/L. First, 2.0 µL of the dilutions were dispensed in each dish. Next, 2.0 mL of Middlebrook 7H10 medium was added appropriately and supplemented with glycerol and OADC (oleic acid, albumin, dextrose, and catalase). The medium was allowed to cool, and then in the middle of each well 2.0 µL of a dilution of the inoculum (containing about 10 6 CFU/mL) of M. tuberculosis H37Rv or M. bovis BCG was added. The isoniazid drug was used as a positive standard. The plates were covered for 2-3 weeks at 37 • C. After the incubation period, the plates were observed and the MIC was determined as the minimum concentration on which growth was not observed. The experiment was repeated on a different day observing exactly the same results.

Cytotoxicity Assay
Vero cells (ATCC ® CCL-81 ™), were cultured in Dulbecco Modified Eagle Medium (DMEM) supplemented with 10% bovine fetal serum and 1% streptomycin-penicillin and passaged twice before the assay in 21 cm 2 cell culture Petri dishes at 37 • C in 5% CO 2 incubator and 100% humidity. Cells were then cultured in 96-well plates for 24 h before the assay to a cell density of 104 cells per well. A 96-well master plate was used to dilute of the compounds in DMEM medium, by diluting 5 mL of the 100 mg/mL DMSO stock of the compounds into 195 mL of the DMEM medium and performing twofold serial dilution with 100 mL. One hundred microliters of this solution was transferred to each plate containing the Vero cells and incubated for 48 h. Five percent NaDS was used as a positive control and DMSO as negative control under the same dilution conditions. Ten milliliters of a freshly prepared MTT (3-(4,5-dimethyl-2-thiazolyl)-2,5-diphenyl-2H-tetrazoliumbromide) solution in PBS at 5 mg/mL was added to each well, and the plates were further incubated for 1.5 h. The media was then removed and 130 mL of DMSO was added to each well. After 30 min of incubation, the absorbance was read at 530 nm on a microplate reader. The experiment was performed in duplicate on different days with different cell cultures and different stock of the compounds. The IC 50 (in µM) values were determined by interpolation from the mean absorbance data of 100% viability (negative control) and 0% viability (positive control).

Conclusions
Thirty-two compounds with anti-TB activity against the M. tuberculosis H37RV strain were retrieved from the literature to obtain possible quantitative structure-anti-TB activity correlations. These compounds were optimized at the B3LYP/6-31G (d,p) level of theory; the minimum energy structures were used to calculate electronic, thermodynamic, and topographical molecular descriptors. As a result, model 3 was chosen according to statistics and the lowest number of attributes, so it was used as a guiding model for generating new compounds. The van der Waals volume, electron density, and electronegativity model results suggest a pivotal role in anti-TB activity; as such, van der Waals volume increases biological activity and electron density and the electronegativity decrease the anti-TB activity. Five new compounds were designed by taking compound 26 as a base structure and adding substituents according to model 3. The designed compounds were compared electronically with reported anti-TB compounds, such as BTZ043 and isoniazid. It was found that compounds with frontier orbitals distributed on opposite sides of the molecules exhibited higher activity values than those without this behavior. In addition, docked compounds against enoyl-acyl carrier protein reductase (InhA) and decaprenylphosphorylβ-D-ribose-2'-oxidase (DprE1) showed affinity values ranging from −8.5 to −7.1 Kcal/mol; compounds 11d and 11e showed the best behavior towards both molecular targets. Proteincompound molecular dynamics calculations resulted in negative binding energy in the range of −71.3 to −12.7 kcal/mol. In all cases, the compounds bound to each target at the catalytic site reported in the literature, which means that there is a high probability that they were inhibitors of these two proteins. The newly designed compounds were synthesized and evaluated for their anti-TB activity against M. bovis BCG, M. tuberculosis H37Rv, and six strains of M. tuberculosis resistant to rifampicin and isoniazid. Most compounds showed better MIC values than the reference drugs (isoniazid and oxafloxacin). Compounds 11d and 11e, with chlorine and bromine groups on the quinolinone ring, were found to be the most active with MIC values of 0.15 µM and 0.13 µM, respectively, against M. tuberculosis H37Rv and selectivity index of >10. This study suggests that 11d and 11e possess great potential in developing of new anti-TB drugs and can serve as hits for in vivo studies to validate their potential as lead compounds.  Figure S2. The frontier orbital for the designed compounds 11a-e from the QSAR model. HOMO: Highest occupied molecular orbital; LUMO: low unoccupied molecular orbital. BTZ: Benzothiazinone derivatives, Figure S3. 1 H and 13 C NMR spectrum of 11a, Figure S4. 1 H, 13 C NMR and DEPT 135 NMR spectrum of 11b, Figure S5. 1 H, 13 C NMR and DEPT 135 NMR spectrum of 11c, Figure S6. 1 H, 13 C NMR and DEPT 135 NMR spectrum of 11d, Figure S7. 1 H, 13 C NMR and DEPT 135 NMR spectrum of 11e, Table S1: Observed and predictive and experimental activities (and their residuals) of the set of quinoline derivatives.
Author Contributions: J.V., investigation; data curation; formal analysis; writing-original draft; J.R.M., S.A.C. and J.L.P., methodology, writing-review and editing; V.R., G.P., L.V., A.B., A.C. and O.V. planned and performed the antimycobacterial assays and review and editing; B.I., R.A., J.Q. and A.I., methodology, writing-review and editing; E.M. and D.I., conceptualization; formal analysis; investigation; resources; data curation; writing-original draft; writing-review and editing; supervision and project administration. All authors approved the final manuscript. All authors have read and agreed to the published version of the manuscript.

Conflicts of Interest:
The authors declare no conflict of interest.