The amino acid sequence of a particular target protein, Rv1250, was retrieved from National Center for Biotechnology Information, i.e., NCBI (http://www.ncbi.nlm.nih.gov/) which is a part of the United States National Library of Medicine (NLM), a branch of the National Institutes of Health (NIH). It addresses a sequence of databases related to biomedicine and biotechnology and a significant source for bioinformatics tools and services. For the structure identification, SWISS-MODEL and Phyre2 servers were used which presented only 17.21% of the total predicted structure (Supplementary Figure S1, S2) (template: 6G9X) (Schwede et al., 2003; Kelley et al., 2015; Saini et al., 2021). Therefore, the I-TASSER (Iterative Threading Assembly Refinement) server was used to predict the entire structure of the Rv1250 protein. With the help of I-TASSER tool, five structures were obtained based on the C-Score (confidence score) (http://zhanglab.ccmb.med.umich.edu/I-TASSER/) (Yang and Zhang, 2015). The structure with the least value for C-Score was selected for further docking analysis. These selected tools for structure prediction are automatic homology modeling sever. There is no requirement of the input of a particular template as it automatically investigates a best fitting template of the target and presents results in minutes or hours. Energy minimization of selected template was done by using the SWISS-PDB (Protein Data Bank) Viewer (https://spdbv.vital-it.ch/) which is an automated server that offers an easily accessible interface allowing the investigation of numerous proteins at the same time (Arnold et al., 2006). The proteins can be overlaid in order to infer the structural alignments and relate their active sites or any other relevant parts (Jabeen et al., 2019). The model was validated by using Ramachandran Plot predicted from PROCHECK. Further, the validation of Homology model was done by using Verify-3D (http://servicesn.mbi.ucla.edu/Verify3D/), which ascertains the compatibility of an atomic model with its own amino acid sequence by allocating a structural class based on its environment and position (Luthy et al., 1992; Saini et al., 2019). Additionally, the model was validated by using ProSA integrated webserver (https://prosa.services.came.sbg.ac. at/prosa.php) to approximate deviation of Z-score from the high-resolution structures and to determine the error (Wiederstein and Sippl, 2007; Saini et al., 2021).

The prediction of the active site was done by using Dogsite scorer which is an automated pocket detection tool, moreover it also analyze the drugability of the specific site by taking into consideration its geometric and physiochemical properties (http://dogsite.zbh.uni-hamburg.de) (Malik et al., 2019). The protein model was submitted in the software using this tool and the results were generated. Further validation of the results was done by Castp (http://sts-fw.bioengr.uic.edu/castp/calculation.php). Finally, it was modeled in PyMOL software, a cross-platform molecular graphic tool that has been extensively used for 3D visualization of surfaces, proteins, small molecules, nucleic acids, etc., (Wass et al., 2010; Yuan et al., 2017; Wu et al., 2018). PyMOL provides high-resolution images of macromolecules and it is easy to build a molecule from scratch (Yuan et al., 2017).

The three-Dimensional (3D) structure of the ligand isoniazid was retrieved from the PubChem database (http://www.pubchem.ncbi.nlm.nih.gov/) in a sdf format. This downloaded structure holds the information on the chemical properties and biological activities of small molecules. Further, sdf file was converted into pdb format by Open Babel and energy minimized by Universal Force Field (UFF) in PyRx 0.8 (O’Bayle et al., 2011; Dallankyan & Olson, 2015). The energy minimized isoniazid was considered for molecular docking against Rv1250.

Molecular docking is an extensively used computer simulation protocol to visualize the conformation of a receptor-ligand complex, where the receptor is usually a protein or a nucleic acid molecule and the ligand is either a small molecule or another protein. The docking of a specific ligand to the active site of a precise protein was done by using a server called SWISS-DOCK (http://www.swissdock.ch). There remains an automatic generation of the target protein structure as well as ligand structure for docking. Moreover, docking does not require computational inputs from users as all computations are made by server itself (Grosdidier et al., 2011). Finally, to validate the results of SWISS-DOCK, Autodock vina was used, which is one of the fastest and most widely used open-source programs for molecular docking (Singh et al., 2018; Dalal et al., 2019; Saini et al., 2019; Kesari et al., 2020). The grid box was generated by using Autodock tool and then the site-specific docking was done in the grid box. Hydrogen atoms and Kollman charges (4.03) were added to the protein, while Hydrogen atoms and Gasteiger charges were added on the ligand and both the files were converted in pdbqt format and saved. A potential grid box was generated using AutoGrid4 having spacing of 0.375Å (nearly one-fourth of the length of C-C covalent bond) (Table 1) (Kumari et al., 2021). The box dimensions were as 54 Å × 58 Å × 58Å and center point coordinates were as X = 83.911, Y = 69.166, Z = 77.416. In AutoDock Vina, total nine poses were generated by using the receptor and ligand files together with configuration file encompass grid box properties. An interaction of docking pose with active site residues was observed and the pose with higher binding affinity (−5.3 kcal/mol) was selected (Saini et al., 2019; Kumari et al., 2021).

Molecular dynamics (MD) simulation is a technique that can be used effectively to understand macromolecular structure-to-function relationships (Ahmad & Kesavan, 2021). MD simulation is most widely used to evaluate the stability and to enhance low quality models (Sokkar et al., 2011). Here, it was used to observe the stability of a modelled structure, as well as the stable and static interaction between the ligand and the protein (Kesari et al., 2020; Singh et al., 2020; Saini et al., 2021; Kumari et al., 2022). For MD simulation, firstly GROMACS (Groningen Machine for Chemical Simulations) was used with the GROMOS96 43a1 force field for generating the topology and protein coordinates (Tran et al., 1996; Van Der Spoel et al., 2005). Then, the ligands were processed to create the topology and coordinates files by using PRODRG webserver (http://davapc1.bioch.dundee.ac.uk/cgi-bin/prodrg) (Kumari et al., 2020; 2021). The system was solvated in simple point charge (SPC) in a triclinic box (X = 10.71, Y = 8.73, and Z = 7.52 nm) having a minimum 1 nm distance between protein and edge of the box. Two distinct equilibration phases, i.e., constant number of particles, volume and temperature (NVT) and constant number of particles, pressure, and temperature (NPT) were run for 1 ns (Dalal et al., 2021), while final MD was run for 100 ns with a generation of coordinates at a regular interval of 10 ps. (Singh et al., 2022). The hydrogen bond occupancy and numbers between isoniazid-Rv1250 complex was generated by visual molecular dynamiucs (VMD) (Humphery et al., 1996; Swain et al., 2022).

The Rv1250 protein sequence was obtained from the NCBI. After the sequence retrieval, BLASTp (Basic Local Alignment Search Tool for Protein) was done for to check protein homologues of Rv1250 (Vyas et al., 2012). Among the homoglous structures, Crystal structure of a MFS transporter at 2.54 Å resolution (6G9X) showed a query coverage and identity of 23% and 28.19%, respectively with RV1250 (Binkowski et al., 2003). I-TASSER was used to generate five models of Rv1250 protein with varying C-score values, as shown in Table 2. The model with the highest C-score value, i.e., model-1 (c-score = −2.07) was selected for further performing the docking process (Figure 1).

Model-1 was selected out of 5 models because it had the highest C-score value which signified that the model had a higher confidence (Armando et al., 2018).

This tool also predicted the normalized B-factor, which is a denomination that indicates the degree to which protein residues/atoms have inherent thermal mobility. This value is determined in I-TASSER, by threading the template proteins from the PDB along with the sequence profiles obtained from sequence databases. The reported B-factor profile corresponds to the target protein’s normalized B-factor, designated as B=(B′-u)/s, where B′ is the raw B-factor value, u is the mean and s is the standard deviation of the raw B-factors along the sequence. The resulted graph showed that coils had more mobility than the alpha-helix because the coils had a larger peak.

Ten templates for model building of Rv1250 protein were predicted (Table 3). However the template with a normalized Z-score >1 mean/least Z-score value was used by the I-TASSER for the alignment. The Z-score signifies the indicative of the quality of overall model and is utilized to verify that the input structure belongs to the range of scores usually located in the resident proteins of analogous size (Nawaz et al., 2014).

After building the final protein model, it was prepared for the further docking process (Dar & Mir, 2017). The SWISS PDB viewer performed the energy minimization of the specific protein model (Supplementary Figure S3). In the predicted model, the elevated energies were eliminated using this technique, achieving local minima that were nearby to the native structure. Energy minimization is done as a definite step to resolve minute structural deformations, adverse interactions established during the modeling procedure (Waterhouse et al., 2018).

The Ramachandran plot of Rv1250 generated by PROCHECK showed that 89.9% of the amino acid residues are in the most favoured region, 9.4% in the additional allowed region, 0.6% in the generously allowed region and 0.1% in disallowed region as shown in Figure 2. The results showed that no residues were in the outliners region and similar scores were stated by Dhankhar et al. (2020), Dalal and Kumari (2022).

Furthermore, the 3D verification showed that the predicted model was valid for the further docking process of Rv1250 protein of M. tuberculosis (Abid et al., 2017). Following, results were concluded from the Verify-3D assessment model.

• Residues (82.56%) had an average 3D-1D score ≥0.2

It was observed that the profile score of homology model was 82.56% which indicates that 82.56% of the amino acid residues had an average 3D-1D score of ≥0.2. Further, the model validation was done using ProSA and the z-score was obtained as −3.51 (Figure 4), which indicates that the modeled protein falls in the range of the X-ray elucidated protein structures (Kwofie et al., 2018; Saini et al., 2019; Kumari and Dalal, 2022). It has been also reported that a protein model is considered better with more negative z-score (Yakubu et al., 2017). Therefore, the results demonstrated that the model is reliable and can be considered for the molecular docking and dynamics (Panda et al., 2019).

For the docking process, active site was predicted by the Dogsite scorer and Castp. These servers predicted an active site with sixteen amino acid residues [GLY142, SER275, ALA3232, PRO324, GLY327, MET328, SER333, HIS334, LEU336, CYS337, GLY340, VAL380, PRO383, LEU384, ASP486, PHE487] (Figure 5).

The 3D model with an active site is shown in Figure 6. Once the active site was identified, a suitable ligand, i.e., Isoniazid (molecular formula: C6H7N3O), for binding at this site was selected from the PubChem server, and basic properties and functions of isoniazid were also obtained from the same server. It is evident from the previous research that isoniazid is the best known and well-suited drug to such an infectious disease called TB (Egsmose et al., 1965). Over many decades, isoniazid is considered a highly selective and anti-tuberculosis therapeutic agent (Youatt, 1969; Zhang et al., 1996; Timmins & Deretic, 2006).

Docking of small molecule compounds into a receptor binding site and evaluating the binding likeness of the composite is an imperative section of the structure-based drug design technique (Seeliger & de Groot, 2010). In this study, SWISS-DOCK and PyMOL were used to conduct the protein docking analysis (Dar & Mir, 2017). The protein model and the particular ligand were submitted to the server, and the docked model was received after validation by the SWISS-DOCK software. This software provides an automatic arrangement of protein-ligand structures, the diverse parameter sets and the expedient image and assessment of docking predictions which makes it a choice to be available to extensive viewers (Grosdidier et al., 2011).

The molecular docking results were visualized by the PyMOL (Figure 7). Molecuolar docking performed by using AutoDock vina shows the binding affinity of −5.1 kcal/mol (Table 4). The ligand Isoniazid showed possible hydrogen bonding with ASP331, PHE487, GLY489, ALA490, and ALA491 residues, as displayed in Figure 7. The docking results showed that isoniazid (INH) was the most effective ligand to stop the action of M. tuberculosis by inhibiting the formation of mycolic acid (Takayama et al., 2005; Parulekar et al., 2013), which protects the bacterium from the antibiotic drug. Due to the efficiency of INH-NAD and INH-NAPD adducts to inhibit nucleic acid, and cell wall lipid synthesis, isoniazid provides an effective anti-tuberculosis blend with the efficiency of attacking multiple targets (Timmins & Deretic, 2006).

A molecular dynamic simulation is the tool used to assess protein-ligand system at the atomistic level and articulate on the stability of the protein-ligand complex in the dynamic environment. In the current study, different parameters, for instance, Root Mean Square Deviation (RMSD), Root Mean Square Fluctuation (RMSF), Radius of Gyration (Rg), Solvent Accessible Surface (SASA), and formation of hydrogen bonds were analyzed during the molecular simulation (Table 5).