Phytocompounds from Himalayan Medicinal Plants as Potential Drugs to Treat Multidrug-Resistant Salmonella typhimurium: An In Silico Approach

Medicinal plants can be used as natural therapeutics to treat diseases in humans. Enteric bacteria possess efflux pumps to remove bile salts from cells to avoid potential membrane damage. Resistance to bile and antibiotics is associated with the survival of Salmonella enterica subspecies enterica serovar Typhimurium (S. typhimurium) within a host. The present study aimed to investigate the binding affinity of major phytocompounds derived from 35 medicinal plants of the North Western Himalayas with the RamR protein (PDB ID 6IE9) of S. typhimurium. Proteins and ligands were prepared using AutoDock software 1.5.6. Molecular docking was performed using AutoDock Vina and MD simulation was performed at 100 ns. Drug likeness and toxicity predictions of hit phytocompounds were evaluated using molinspiration and ProTox II online servers. Moreover, docking, drug likeness, and toxicity results revealed that among all the selected phytocompounds, beta-sitosterol exhibited the most efficacious binding affinity with RamR protein (PDB ID 6IE9) and was nontoxic in nature. MD simulation data revealed that beta-sitosterol in complex with 6IE9 can be used as an antimicrobial. Furthermore, beta-sitosterol is stable in the binding pocket of the target protein; hence, it can be further explored as a drug to inhibit resistance-nodulation-division efflux pumps.


Introduction
Salmonella is a bacterial pathogen that infects the intestinal tract and gallbladder and causes numerous foodborne illnesses in humans. Enteric bacteria, such as Salmonella, tolerate the existence of bile acids for their survival in the gastrointestinal transit and gallbladder [1,2]. Nontyphoidal serovars (NTSs) of Salmonella enterica are the major causes of foodborne illnesses and diarrhea occurring worldwide [3,4]. In S. enterica, the resistancenodulation-division (RND) pump is translated from the acrAB gene regulated by RamA, a transcriptional activator. RamR inhibits the expression of the ramA gene involved in multidrug resistance in Salmonella enterica subspecies enterica serovar Typhimurium (S. typhimurium). S. typhimurium is one of the NTSs causing severe human infections and results in more hospitalizations and mortality worldwide [5]. Additionally, the treatment choices are limited because antibiotics may lead to enhanced shedding of S. typhimurium and its emergence as multidrug-resistant bacteria [6,7]. Presently, it is no longer considered as the first choice of antimicrobial drug due to its resistance. The advent of new resistance mechanisms exists in S. typhimurium, leading to challenges in treating infections. Therefore, alternative therapeutic approaches are required. S. typhimurium comprises of at least nine multidrug efflux systems [8]; among these, the AcrAB-TolC system, containing the AcrB transporter of the RND family, is particularly effective in developing resistance to bile acid [8,9]. RamR is a local transcriptional repressor, which belongs to the TetR family of regulatory proteins [10]; it helps in impairing the ramA gene expression that affects ramA gene transcription resulting in multidrug resistance. Therefore, it is crucial to identify novel pharmacological targets against drug-resistant S. typhimurium.
Medicinal plants have played a pivotal role in treating diseases since the prehistoric period. These plants comprise various phytoconstituents in every part (bark, leaves, flowers, roots, fruits, and seeds), and exhibit high therapeutic value [11].
Herbal medications have recently gained immense interest as they are safe and economic and have been widely used for several years to treat diseases. Bioactive components are secondary metabolites of plants that produce pharmacological and toxicological issues in living organisms. It is difficult to screen each phytoconstituent for toxicity. In drug design, computational techniques play a crucial role in studying the toxicity of chemical and natural compounds as well as their properties [12]. In silico studies, with specific reference to toxicity prediction and molecular docking for each phytochemical in order to determine their therapeutic efficiency, require less time, are economic, and can harm animals [13,14]. Therefore, the present study was designed to investigate the binding affinity of phytocompounds of 35 important medicinal plants of the Northwestern Himalaya with S. typhimurium protein (RamR transcriptional repressor of TetR family) with PDB ID: 6IE9 to prevent inflammatory gastroenteritis.

Ligand Preparation
Seventy major phytocompounds of 35 medicinal plants from Himachal Pradesh, India, were selected for molecular docking analysis. The three-dimensional structures of all the phytocompounds and resistant bile components of S. typhimurium (chenodeoxycholic acid) were downloaded from Pubchem (www.pubchem.com (accessed on 15 May 2021)) in .sdf format, which was finally converted into a PDB file. Each selected ligand (phytocompounds and bile component) was prepared using the open Babel software from the command line on an Ubuntu terminal. Table 1 lists the names of the phytocompounds selected for this study, their plant sources, pharmacological properties, and ethnomedicinal uses.

Protein Preparation
RamR of S. typhimurium [67] was used for molecular docking with major phytocompounds from 35 important medicinal plants (Table 1) found in the northwestern Himalayas of Himachal Pradesh to identify potential inhibitors of S. typhimurium. The 3-D structure protein (PDB ID = 6IE9) was downloaded from the protein databank (http://www.rscb.org/pdb (accessed on 15 May 2021)) as a pentamer, and chain A was extracted for docking using PyMol. Chain A was prepared for docking, and a grid box was set to cover the entire protein (grid box dimensions = 40, 40, 40 Å) and was centered at x, y, z = 11.029, 33.324, 12.359 Å, respectively.

Molecular Docking of Major Phytocompounds of Thirty-Five Medicinal Plants
The AutoDock tool was used to dock the selected ligands to the catalytic triad of proteins, which was further stored as a pdbqt file. Docking was carried out to estimate the population of possible ligand conformations/orientations at the binding site. To align the ligands in the same spatial coordinates, a vina perl script was used [16]. The best conformation was selected with the minimum docked energy after completing the docking search. The pdb complex of protein and ligands was analyzed using Discovery Studio (https://discover.3ds.com/d (accessed on 15 May 2021)) to study the interactions between proteins and ligands. The binding strength of the ligand was calculated as a negative score (kcal/mol).

Drug Likeness Calculations
The drugs were scanned to assess whether the selected phytochemicals met the druglikeness criteria. Lipinski's rule of 5 using Molinspiration (http://www.molinspiration.com (accessed on 15 May 2021)) was used to verify drug likeness attributes, such as the number of hydrogen acceptors <10, number of hydrogen donors <5, molecular weight <500 Da, and partition coefficient log P > 5. The smiles format of all major phytocompounds was uploaded for further screening [68].

ADMET Screening and Toxicity Prediction of Phytocompounds
Absorption, distribution, metabolism, excretion, and toxicity (ADMET) screening was performed to evaluate the absorption, toxicity, and drug-likeness properties of the selected phytocompounds. The 3-D structures of 11 phytocompounds (asiaticoside, beta-sitosterol, bryophyllin A, madecassoside, Mahanimbine, Pennogenin, Rutin, Solasonine, Solamargine, Withaferin A, and Withanone) were saved in smiles format and uploaded on the SWIS-SADME (http://www.swissadme.ch/ (accessed on 15 May 2021)) (Molecular Modeling Group of the SIB (Swiss Institute of Bioinformatics) and PROTOX-II (https://tox-new.charite.de/protox_II/ (accessed on 15 May 2021)) web servers (Charite University of Medicine, Institute for Physiology, Structural Bioinformatics Group, Berlin, Germany) [69][70][71][72] for ADMET screening. SWISSADME is an online tool used to predict ADME and pharmacokinetic and physicochemical features of a molecule, which are the main determinants for clinical trials. Toxicity was evaluated in compounds with LD 50 values ≤50 mg/kg (Class I), >50 mg/kg but <500 mg/kg (Class II), 500 < LD 50 ≤ 5000 mg/kg (Class III), and LD 50 > 5000 mg/kg (Class IV). Classes I, II, and III exhibited less toxicity, whereas Class IV revealed no toxicity [73,74]. Moreover, PROTOX is a rodent oral toxicity server that determines the LD 50 value and the toxicity class of a target molecule [69]. A schematic of the experiment is illustrated in Figure 1. LD50 ≤ 5000 mg/kg (Class III), and LD50 > 5000 mg/kg (Class IV). Classes I, II, and III exhibited less toxicity, whereas Class IV revealed no toxicity [73,74]. Moreover, PROTOX is a rodent oral toxicity server that determines the LD50 value and the toxicity class of a target molecule [69]. A schematic of the experiment is illustrated in Figure 1.

MD Simulation of Protein Ligand Complexes
The ligand-protein complex structure was prepared before MD simulation to remove the structural errors. Extensive 100 ns MD simulation was performed on the Desmond platform to analyze the ligand behavior within the complex [75]. The complex was solvated in a TIP3P (8018 molecules) water model and 0.15 M NaCl (Na: 54.42 mMol and Cl: 49.88 mMol) to mimic a physiological ionic concentration. The molecular mechanics/generalized born surface area (MM/GBSA) was used for binding free energy calculations [75]. MD simulation trajectories were used as inputs to calculate the MM/GBSA of the binding free energies of the ligands and to investigate their binding mechanisms. The thermodynamic binding energy was calculated for every 1000th frame, as the complete MD simulation includes 10,000 frames. The OPLS 2005 force field was used for the MD simulations.

MD Simulation of Protein Ligand Complexes
The ligand-protein complex structure was prepared before MD simulation to remove the structural errors. Extensive 100 ns MD simulation was performed on the Desmond platform to analyze the ligand behavior within the complex [75]. The complex was solvated in a TIP3P (8018 molecules) water model and 0.15 M NaCl (Na: 54.42 mMol and Cl: 49.88 mMol) to mimic a physiological ionic concentration. The molecular mechanics/generalized born surface area (MM/GBSA) was used for binding free energy calculations [75]. MD simulation trajectories were used as inputs to calculate the MM/GBSA of the binding free energies of the ligands and to investigate their binding mechanisms. The thermodynamic binding energy was calculated for every 1000th frame, as the complete MD simulation includes 10,000 frames. The OPLS 2005 force field was used for the MD simulations.     Figures 3 and 4. Furthermore, all the selected phytocompounds were screened for drug likeness and ADME/T.

Drug Likeness Prediction of Active Phytocompounds of Eight Medicinal Plants
Drug likeness was predicted using Molinspiration (server) to study the drug-likeness properties of active phytocompounds, which are based on the Lipinski rule of 5. Lipinski's rule of five was followed by all selected phytocompounds and the bile component, chenodeoxycholic acid, which revealed no violation. Among all the phytocompounds, bryophyllin A, pennogenin, withaferin A, and withanone followed all the rules of drug likeness, whereas beta-sitosterol and mahanimbine exhibited one violation, which was acceptable (Table 3).

Toxicity Prediction of Active Phytocompound and Chenodeoxycholic Acid
The toxicity of phytocompounds was predicted using ProTox-II server and the results are summarized in Table 4. Among all the phytocompounds, only beta-sitosterol exhibited one violation and immune toxicity. Based on the molecular drug likeness and toxicity data, beta-sitosterol was found to be the best phytocompound, which can be used for efflux pump inhibition, and it was further selected for molecular dynamics studies.

MD Simulation of Protein-Ligand Complexes
Molecular dynamics simulation provides insight into the protein-ligand stability and protein structural flexibility of the docked complexes. The root-mean-square deviation (RMSD) plot of beta-sitasterol and 6IE9 complex exhibited significant stability in the protein pocket. These compounds fluctuate within the acceptable range between 3.2 and 5.6 Å, whereas protein Cα RMSD became stable after 25 ns and fluctuates in the range between 4.8 and 6.0 Å ( Figure 5). The trajectory analysis revealed that a sharp change in the ligand RMSD at approximately 45 ns mainly occurred due to the aliphatic chain; this was also observed in the ligand RMSF plot, where fluctuation occurs in atoms 26-30 ( Figures 5 and 6).
Several interactions were responsible for the conformational stability of the compound within the binding pocket, where hydrophobic interactions developed with residues L66, I88, Y92, M126, L130, Val141, Phe155, and L156 (Figure 7). The hydroxyl group revealed hydrogen bonding and water interactions with residues R136, S137, and R131. Moreover, MD simulation supported the docking results, where the compound interacted with residues that were linked with the molecule (Figures 2B and 7). Furthermore, the thermodynamic energy analysis revealed that the average binding free energy was 138.65 ± 19.84 kcal/mol, whereas that of the docked complex was 109.18 kcal/mol (Table 5).
protein structural flexibility of the docked complexes. The root-mean-square deviation (RMSD) plot of beta-sitasterol and 6IE9 complex exhibited significant stability in the protein pocket. These compounds fluctuate within the acceptable range between 3.2 and 5.6 Å, whereas protein Cα RMSD became stable after 25 ns and fluctuates in the range between 4.8 and 6.0 Å ( Figure 5). The trajectory analysis revealed that a sharp change in the ligand RMSD at approximately 45 ns mainly occurred due to the aliphatic chain; this was also observed in the ligand RMSF plot, where fluctuation occurs in atoms 26-30 ( Figures  5 and 6). Several interactions were responsible for the conformational stability of the compound within the binding pocket, where hydrophobic interactions developed with residues L66, I88, Y92, M126, L130, Val141, Phe155, and L156 (Figure 7). The hydroxyl group revealed hydrogen bonding and water interactions with residues R136, S137, and R131. Moreover, MD simulation supported the docking results, where the compound interacted with residues that were linked with the molecule (Figures 2B and 7). Furthermore, the thermodynamic energy analysis revealed that the average binding free energy was 138.65 ± 19.84 kcal/mol, whereas that of the docked complex was 109.18 kcal/mol (Table 5).

Discussion
Poor pharmacological characteristics are the major cause of late-stage failure in drug discovery. Thus, early determination of the inherent medicinal activities of the target compounds is crucial [76]. Moreover, medicinal plant species are abundant in Asia's Himalayan woodlands, and they play a pivotal role in rural livelihoods by producing various valuable food and pharmaceutical commodities [77]. In recent years, the WHO estimated a remarkable increase in the multidrug resistance rate worldwide due to Salmonella strains [78]. Salmonella infections are gaining importance worldwide owing to their socioeconomic impact. Salmonella Typhimurium is one of the most common serovars predominantly associated with clinically reported human salmonellosis in several countries, accounting for at least 15% of infections worldwide [5]. Moreover, at least nine multidrug efflux pumps confer drug resistance in Salmonella; among these, AcrAB is constitutively expressed and is the most potent drug for intrinsic drug resistance [8]. AcrAB is a member of the RND family transporter that cooperates with TolC, an outer membrane component [79]. The AcrAB-TolC system comprising RND transporters can accumulate substrates in the periplasm rather than in the membrane or cytoplasm [80]. A common mechanism of intrinsic resistance to antimicrobial agents in Gram-negative bacteria is represented by the RND family efflux systems, which extrude a broad spectrum of antibiotics and biocides from the periplasm to the exterior of the cell [81].
The present study explored some medicinal plants, including Girardinia diversifolia, which we reported earlier for synergistic and efflux pump inhibitory activity against different strains of S. typhimurium and Staphylococcus aureus [82,83]. In contrast to our study, Mehta et al. [84] reported that methanolic extracts of Pistacia integerrima, Ocimum sanctum, C. asiatica, Momordica charantia, Zingiber officinale, and Withania somnifera exhibited synergistic activity in combination with ciprofloxacin and tetracycline against multidrug resistance. AcrAB-TolC in Salmonella Typhimurium acts as an efflux pump inhibitor. Furthermore, they reported the binding affinity (−8.2 kcal mol −1 ) of lariciresinol with 6EI9 (RamR). Similarly, Luhata et al. [85] reported the antibacterial activity of beta-sitosterol against S. aureus. Sen et al. [86] reported the antibacterial activity of beta-sitosterol against Escherichia coli, Pseudomonas aeruginosa, S. aureus, and Klebsiella pneumoniae. Rolta et al. [87] studied the antibacterial and antifungal activities of phytocompounds of Rheum emodin (emodin, rhein-13c6, and chrysophenodimethy ether) by molecular docking and MD simulations and found that phytocompounds of R. emodin exhibited the best interaction with bacterial and fungal targets. Similarly, Rolta et al. [71] studied the interactions of phytocompounds with the N-protein of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) 2020, and reported that emodin, aloe-emodin, anthrarufin, alizarine, and dantron phytocompounds efficiently inhibit SARS-CoV-2 N-protein. Salaria et al. [88] studied the in vitro and in silico antibacterial and antifungal activities of essential oil, thymol derived from Thymus serpyllum, and validated the docking results via MD simulations. The conformational changes during protein-ligand interactions have been extensively studied via MD simulation methods [89].

Conclusions
The major phytocompounds of 30 fine medicinal plants of the northwestern Himalayas were selected for molecular docking study with the 6EI9 (RamR) target protein of S. typhimurium. Among all the selected phytocompounds, 11 phytocompounds exhibited the best activity compared to the standard drugs. Drug likeness and toxicity data revealed that beta-sitosterol, a major phytocompound of G. diversifolia (Link) Friis, is nontoxic in nature and follows the drug likeness rule. Moreover, MD simulation of betasitosterol in complex with 6EI9 was found to be stable between 0 and 100 ns time period. In this study, we found that beta-sitosterol is a potential plant-based drug for treating S. typhimurium infection. Furthermore, this study needs to be validated through in vitro and in vivo experiments.