# Molecular Dynamics Simulation and Pharmacoinformatic Integrated Analysis of Bioactive Phytochemicals fromAzadirachta indica (Neem) to Treat Diabetes Mellitus

Authors: Asif Abdullah; Partha Biswas; Md. Sahabuddin; Afiya Mubasharah; Dhrubo Ahmed Khan; Akram Hossain; Tanima Roy; Nishat Md. R. Rafi; Dipta Dey; Md. Nazmul Hasan; Shabana Bibi; Mahmoud Moustafa; Ali Shati; Hesham Hassan; Ruchika Garg
Journal: Journal of Chemistry
 (2023)

---

## Abstract

Diabetes mellitus is a chronic hormonal and metabolic disorder in which our body cannot generate necessary insulin or does not act in response to it, accordingly, ensuing in discordantly high blood sugar (glucose) levels. Diabetes mellitus can lead to systemic dysfunction in the multiorgan system, including cardiac dysfunction, severe kidney disease, lowered quality of life, and increased mortality risk from diabetic complications. To uncover possible therapeutic targets to treat diabetes mellitus, the in silico drug design technique is widely used, which connects the ligand molecules with target proteins to construct a protein-ligand network. To identify new therapeutic targets for type 2 diabetes mellitus,Azadirachta indica is subjected to phytochemical screening using in silico molecular docking, pharmacokinetic behavior analysis, and simulation-based molecular dynamic analysis. This study has analyzed around 63 phytochemical compounds, and the initial selection of the compounds was made by analyzing their pharmacokinetic properties by comparing them with Lipinski’s rule of 5. The selected compounds were subjected to molecular docking. The top four ligand compounds were reported along with the control drug nateglinide based on their highest negative molecular binding affinity. The protein-ligand interaction of selected compounds has been analyzed to understand better how compounds interact with the targeted protein structure. The results of the in silico analysis revealed that 7-Deacetyl-7-oxogedunin had the highest negative docking score of −8.9 Kcal/mol and also demonstrated standard stability in a 100 ns molecular dynamic simulation performed with insulin receptor ectodomain. It has been found that these substances may rank among the essential supplementary antidiabetic drugs for treating type 2 diabetes mellitus. It is suggested that more in vivo and in vitro research studies be carried out to support the conclusions drawn from this in silico research strategy.


---

## Body

## 1. Introduction

The ability of the body to exert control over and use sugar (glucose) as fuel is reduced in type 2 diabetes. Due to this chronic (long-term) illness, too much sugar circulates in the blood. The most notable metabolic dysfunctions associated with type 2 diabetes, often known as diabetes mellitus, are insulin resistance and cell dysfunction. Early in the condition, circulating insulin levels are higher to compensate for insulin resistance, but insulin production eventually becomes insufficient, and hyperglycemia develops [1, 2]. One of the human body’s most significant regulators of energy balance is insulin. Insulin control issues disrupt energy homeostasis, which leads to type 2 diabetes mellitus [3, 4]. The human insulin receptor is a tyrosine kinase that is homodimeric and disulfide-linked (αβ)2 [5, 6]. Even though the receptor’s signaling is essential in various diseases, thorough and atomic-level knowledge is still unclear on how insulin builds up at the receptor and triggers signal transduction. One challenge is that the isolated, soluble receptor ectodomain (sIR), which is suitable for structural biology studies, lacks the remarkable affinity and poor cooperativity of insulin binding that the hormone receptor (hIR) offers. The pancreas secretes insulin hormones [7, 8]. The ability of a cell to absorb blood glucose as an energy source is controlled by insulin through its interaction with the insulin receptor protein found on the cell’s surface. In type 2 diabetes, insulin binds to the insulin receptor as it should. Still, the signal is not carried through into cells resulting in the cells not absorbing sufficient glucose, and the elevated blood glucose levels that result over time harm organs [9–11].Many modern pharmaceuticals are directly or indirectly produced from plants, which have long been a reliable supply of medicines. According to ethnobotanical studies, approximately 800 plants may have antidiabetic qualities [12, 13]. Drug development is an expensive and time-consuming procedure that calls for numerous clinical trials. High throughput virtual screening and de novo structure-based rational drug design are two in silico methods for drug discovery that have shown potential [14–16]. Virtual screening is now widely used to uncover novel compounds like drugs. In silico virtual screening has become a valuable and time- and money-saving addition to in vitro screening for discovering and developing novel effective compounds [17–20]. Ligand-based virtual screening and receptor-based virtual screening are the two main screening procedures for finding potential compounds from a chemical database that are likely to interact favorably with the target binding sites [21, 22].The 3d structure of a protein and protein-ligand complex is essential in molecular modeling for lead discovery [23–25]. Quantitative structure-activity relationships (QSAR), biological tests, and pharmacophore analysis can all be used to improve and create new leads [26–28]. Structure-based drug design assists in creating more potent and significant molecules during the drug development process [29–31].The movement of each molecule in a molecular system is predicted by molecular dynamics (MD) simulations, which depend on a general physics model governing interatomic interactions [32–34]. The general concept of MD simulation is simple. Given the coordinates of every atom in a biomolecular system, it is possible to determine the force each atom experiences from all the other atoms in the system [35, 36] (for example, a protein encircled by a lipid bilayer and potentially some water). Thus, it is possible to forecast each atom’s spatial position as a function of time using Newton’s equations of motion. One, in particular, goes across time, estimating the forces acting on each atom and modifying every atom’s position and speed based on all of those forces. A three-dimensional movie that shows how the system was configured at each point during the simulated period can be compared to the resulting trajectory [37, 38].Finding a ligand that delivers a specific signaling profile and binding to the target is a common goal of drug development. This is true, in particular for signaling receptors [39, 40]. It is possible to use a full agonist, which vigorously stimulates the signaling and activation of the receptor, a partial agonist, which revitalizes signaling only slightly, and a neutral antagonist [41, 42], which inhibits the body’s natural agonists from binding but does not signal on its own, or an inverse agonist, which decreases signaling below basal levels [43, 44]. The drug must stabilize receptor conformational states and, consequently, particular binding pocket conformational states to produce a specific signaling profile. An agonist, for instance, keeps active states above inactive states. The binding pocket’s conformational changes have unique signaling patterns. This information might be available through MD simulations [37, 43].Pharmacokinetic studies are an integral aspect of a new medicine development program. They are used to determine the time course of medicine and significant metabolite concentrations in tubes and other natural fluids to learn about immersion, distribution, metabolism, and elimination. Pharmacokinetics (PK) is commonly utilized in preclinical investigations to evaluate toxicological outcomes. In phase I cure forbearance studies, beast chronic toxicity data can be used to influence lozenge selection and escalation approaches. Pharmacokinetically guided cure escalation (PGDE) methods facilitate the use of preclinical pharmacokinetic data for clinical cure finding in phase I research. In clinical trials, PK is critical for cure detection and escalation investigations [44–46].Every stage in the discovery and development of a drug requires consideration of chemical properties like absorption, distribution, metabolism, excretion, and toxicity (ADMET). It is crucial to develop effective compounds with improved ADMET characteristics [47]. The notion of drug-likeness is a helpful guideline during the early phases of medication development [48]. The first and most well-knownrule-based filters were proposed by Lipinski and colleagues in 1997. They included the molecular weight (MW) 500, the octanol/water partition coefficient (A log P) 5, the number of hydrogen bond donors (HBDs), and the number of hydrogen bond acceptors (HBAs), all of which must be greater than or equal to 5 [49–51]. Ghose and colleagues calculated that, of the 6304 compounds in the CMC (comprehensive medicinal chemistry) database, more than 80% met the following requirements: −0.4 A log P 5.6, 160 MW 480, 40 MR (molar refractivity), and 20 to 70 atoms in total [52]. 7 drug development can be sped up by using filters that compare compounds to drugs based on their physicochemical properties. But 14 numerous studies have demonstrated the limitations of drug-likeness rules or filters based on physicochemical characteristics [53, 54].Several studies have shown the efficacy of neem in treating diabetes meatus. However, the phytochemical screening of ligands, specific pharmacokinetic behavior analysis-related study, and molecular dynamics simulation study of ligands derived from neem is still not performed. Therefore, we aim to analyze the following parameters to have a deep view of new oral drug candidates and potential inhibitors against insulin receptor ectodomains derived from the neem tree via in silico studies.


## 2. Materials and Methods

### 2.1. Preparation of Protein

The RCSB protein data bank (PDB) (https://www.rcsb.org/) contains the 3D experimental tertiary structures of the protein insulin receptor ectodomain. Three hundred six amino acids (AA) make up the protein sequence known as the insulin receptor ectodomain protein with the PDB id 1i44 [55]. The proteins’ PDB structures were assembled using the following criteria: water, metal ions, side chains, and cofactors were excluded, as well as crystallographic protein structure and expressed in humans. Then, the Chimera software mixed polar hydrogen atoms with nonpolar hydrogen atoms. Finally, the Gasteiger charges of the system were determined.


### 2.2. Selection of Ligands

The large variety of chemical spaces covered by phytochemicals derived from naturally occurring medicinal plants can be employed in drug development and discovery. The IMPPAT (Indian medicinal plants, phytochemistry, and therapeutics) database, which contains >1742 Indian medicinal plants and >9500 phytochemical substances, is used to improve natural product-based medication discovery [56, 57]. The chemicals from Azadirachta indica (neem tree) have been compiled using the IMPPAT journal after extensive literature reviews. The compounds discovered in the database were made by accurately classifying the atoms in AutoDock 4, combining nonpolar hydrogens, finding aromatic carbons, and creating a “torsion tree.” It has been found that most of the particles in the compound have the same sort of atom as the AD4 atom.


### 2.3. Molecular Docking

These days, molecular docking is an integral part of structural biology and is primarily employed for CADD. The techniques aid in predicting which small molecule will attach to a target macromolecule (such as a protein, enzyme, or drug) most advantageously [58, 59]. The PyRx virtual screening software integrated with the AutoDock Vina platform was used to conduct a molecular docking analysis to assess the molecular binding energy of the requisite protein with the chosen phytochemicals [60, 61].A library of chemicals can be tested against a specific therapeutic target using PyRx, an open-source virtual screening program. It mostly appears in CADD techniques. PyRx becomes a more dependable CADD tool with the addition of AutoDock 4 and AutoDock Vina as docking wizards and a simple user interface. The ideal protein and ligand binding configurations were found using the molecular docking PyRx tools AutoDock Vina wizard [62]. Utilizing the default setup parameters of the PyRx virtual screening tools, the binding energy with the highest negative number (kcal/mol) was selected as a `preferred drug candidate for further experiments. The protein-ligands complex’s binding interaction was eventually visualized using the Ligplot + Version 2.2 tools [63].


### 2.4. MM-GBSA Analysis

MMGBSA is known as the molecular mechanics-generalized born surface area which can be performed to calculate ligand binding free energies and ligand strain energies for a set of ligands and a single receptor [64, 65]. After completing the interaction binding affinity analysis, the MMGBSA was conducted by utilizing the prime model of Schrödinger suite 2020-3 (maestro application, paid version). This study analyzed the relative binding free affinity of the ligand compounds Vilasinin (PubChem-102090424), Nimbidinin (PubChem-101306757), 7-Deacetyl-7-oxogedunin (PubChem-1886), and PubChem-122801 with the 1I44 protein complex.


### 2.5. Pharmacokinetic Property and QSAR Analysis

In CADD, the term “PK” refers to the computer-based timing of drug absorption, distribution, metabolism, and excretion (ADME). The flow of drugs into and out of the body as a function of intensity and duration is described by PK (ADME) features as a whole [66]. In the early stages of CADD, the integrity and effectiveness of substances are supported and defined by their pharmacokinetic properties. The study used the SwissADME (http://www.swissadme.ch) website to look at the early-stage PK characteristics of the medications we chose [67]. The SwissADME server is a free online web-based server widely used to predict the pharmacokinetic properties and drug-likeness characteristics of simple molecules. Thus, the pharmacokinetic property is analyzed using the SwissADME server. The PKCSM website is used to predict toxicity (https://biosig.unimelb.edu.au/pkcsm/) [68]. Finally, the quantitative structure-activity relationship (QSAR) assessment for the top ligands was carried out on the PASS server (https://www.pharmaexpert.ru/passonline) to validate antidiabetic, anticancer, and antiviral activity.


### 2.6. Molecular Dynamic Simulation

To determine the stability of the targeted protein, insulin receptor ectodomain’s (PDB ID : 1I44) interaction with the four possible ligand molecules was chosen and subjected to 100 nanoseconds MD simulations. The “Desmond v3.6 Program” from Schrödinger (https://www.schrodinger.com/ac) was used to model the molecular dynamics of the protein-ligand complex structures [69].To establish the desired this framework, a predetermined TIP3P water strategy was developed to establish a specific volume with periodic orthorhombic coordinates separated by 10 mm. The necessary ions, for example, 0+ and 0.15 M salts, were randomly added to the solvent solution to neutralize the framework electrically. The solvency protein system was constructed utilizing a ligand complex, and the system framework was reduced using the default protocol. This was accomplished totally inside the Desmond module by applying the force field settings OPLS3e [70]. In NPT assemblies, which were held at 101,325 bar (1 atm) pressure and 300 K with 50 PS capture sessions totaling 1.2 kcal/mole energy coming before them, the Nosé–Hoover temperature combination and the isotropic approach were used.The screenshots of the molecular dynamic simulation were produced using Schrödinger’s maestro application, version 9.5. The Simulations interaction diagram derived from the Desmond modules of the Schrödinger suite has been used to analyze the simulation event and assess the dependability of the MD simulation. The stability of the protein-ligand complex structure was evaluated based on the trajectory performance, root mean square fluctuation (RMSF), root mean square deviation (RMSD), solvent-accessible surface area (SASA) value, intramolecular hydrogen bonds, the radius of gyration (Rg) value, protein-ligand contacts (P-L), the polar surface area (PSA), and MolSA. The root means square deviation (RMSD) in molecular dynamics simulations is the arithmetic mean distance generated by a single atom’s motion over the period of a certain amount of time relative to a reference time [71, 72]. The RMSD of the structural atoms of a protein, such as heavy particles, backbone, C, and sidechain, were first measured before the root mean square deviation of protein-suited ligand compounds from all timescales, which were realigned and compared against the reference time (in our study 100 ns). The following equation (equation (1)) can be used to determine the RMSD of an MD simulation concerning the period of x.(1)RMSD=∑i=0Nmi∗Xi−Yi2M.N stands for the number of selected atoms, rʹ for the bit’s position in the system x after the reference system’s point has been superimposed, and j for the reference time. Most commonly, local changes in the conformational shape of proteins have been detected and tracked using the root mean square fluctuation (RMSF) [73]. The RMSF estimation of a molecular dynamic’s simulation of a protein with 2 residues can be obtained using the continuity equation (2).(2)RMSF=1T∑TjTXitj−xj2.


## 2.1. Preparation of Protein

The RCSB protein data bank (PDB) (https://www.rcsb.org/) contains the 3D experimental tertiary structures of the protein insulin receptor ectodomain. Three hundred six amino acids (AA) make up the protein sequence known as the insulin receptor ectodomain protein with the PDB id 1i44 [55]. The proteins’ PDB structures were assembled using the following criteria: water, metal ions, side chains, and cofactors were excluded, as well as crystallographic protein structure and expressed in humans. Then, the Chimera software mixed polar hydrogen atoms with nonpolar hydrogen atoms. Finally, the Gasteiger charges of the system were determined.


## 2.2. Selection of Ligands

The large variety of chemical spaces covered by phytochemicals derived from naturally occurring medicinal plants can be employed in drug development and discovery. The IMPPAT (Indian medicinal plants, phytochemistry, and therapeutics) database, which contains >1742 Indian medicinal plants and >9500 phytochemical substances, is used to improve natural product-based medication discovery [56, 57]. The chemicals from Azadirachta indica (neem tree) have been compiled using the IMPPAT journal after extensive literature reviews. The compounds discovered in the database were made by accurately classifying the atoms in AutoDock 4, combining nonpolar hydrogens, finding aromatic carbons, and creating a “torsion tree.” It has been found that most of the particles in the compound have the same sort of atom as the AD4 atom.


## 2.3. Molecular Docking

These days, molecular docking is an integral part of structural biology and is primarily employed for CADD. The techniques aid in predicting which small molecule will attach to a target macromolecule (such as a protein, enzyme, or drug) most advantageously [58, 59]. The PyRx virtual screening software integrated with the AutoDock Vina platform was used to conduct a molecular docking analysis to assess the molecular binding energy of the requisite protein with the chosen phytochemicals [60, 61].A library of chemicals can be tested against a specific therapeutic target using PyRx, an open-source virtual screening program. It mostly appears in CADD techniques. PyRx becomes a more dependable CADD tool with the addition of AutoDock 4 and AutoDock Vina as docking wizards and a simple user interface. The ideal protein and ligand binding configurations were found using the molecular docking PyRx tools AutoDock Vina wizard [62]. Utilizing the default setup parameters of the PyRx virtual screening tools, the binding energy with the highest negative number (kcal/mol) was selected as a `preferred drug candidate for further experiments. The protein-ligands complex’s binding interaction was eventually visualized using the Ligplot + Version 2.2 tools [63].


## 2.4. MM-GBSA Analysis

MMGBSA is known as the molecular mechanics-generalized born surface area which can be performed to calculate ligand binding free energies and ligand strain energies for a set of ligands and a single receptor [64, 65]. After completing the interaction binding affinity analysis, the MMGBSA was conducted by utilizing the prime model of Schrödinger suite 2020-3 (maestro application, paid version). This study analyzed the relative binding free affinity of the ligand compounds Vilasinin (PubChem-102090424), Nimbidinin (PubChem-101306757), 7-Deacetyl-7-oxogedunin (PubChem-1886), and PubChem-122801 with the 1I44 protein complex.


## 2.5. Pharmacokinetic Property and QSAR Analysis

In CADD, the term “PK” refers to the computer-based timing of drug absorption, distribution, metabolism, and excretion (ADME). The flow of drugs into and out of the body as a function of intensity and duration is described by PK (ADME) features as a whole [66]. In the early stages of CADD, the integrity and effectiveness of substances are supported and defined by their pharmacokinetic properties. The study used the SwissADME (http://www.swissadme.ch) website to look at the early-stage PK characteristics of the medications we chose [67]. The SwissADME server is a free online web-based server widely used to predict the pharmacokinetic properties and drug-likeness characteristics of simple molecules. Thus, the pharmacokinetic property is analyzed using the SwissADME server. The PKCSM website is used to predict toxicity (https://biosig.unimelb.edu.au/pkcsm/) [68]. Finally, the quantitative structure-activity relationship (QSAR) assessment for the top ligands was carried out on the PASS server (https://www.pharmaexpert.ru/passonline) to validate antidiabetic, anticancer, and antiviral activity.


## 2.6. Molecular Dynamic Simulation

To determine the stability of the targeted protein, insulin receptor ectodomain’s (PDB ID : 1I44) interaction with the four possible ligand molecules was chosen and subjected to 100 nanoseconds MD simulations. The “Desmond v3.6 Program” from Schrödinger (https://www.schrodinger.com/ac) was used to model the molecular dynamics of the protein-ligand complex structures [69].To establish the desired this framework, a predetermined TIP3P water strategy was developed to establish a specific volume with periodic orthorhombic coordinates separated by 10 mm. The necessary ions, for example, 0+ and 0.15 M salts, were randomly added to the solvent solution to neutralize the framework electrically. The solvency protein system was constructed utilizing a ligand complex, and the system framework was reduced using the default protocol. This was accomplished totally inside the Desmond module by applying the force field settings OPLS3e [70]. In NPT assemblies, which were held at 101,325 bar (1 atm) pressure and 300 K with 50 PS capture sessions totaling 1.2 kcal/mole energy coming before them, the Nosé–Hoover temperature combination and the isotropic approach were used.The screenshots of the molecular dynamic simulation were produced using Schrödinger’s maestro application, version 9.5. The Simulations interaction diagram derived from the Desmond modules of the Schrödinger suite has been used to analyze the simulation event and assess the dependability of the MD simulation. The stability of the protein-ligand complex structure was evaluated based on the trajectory performance, root mean square fluctuation (RMSF), root mean square deviation (RMSD), solvent-accessible surface area (SASA) value, intramolecular hydrogen bonds, the radius of gyration (Rg) value, protein-ligand contacts (P-L), the polar surface area (PSA), and MolSA. The root means square deviation (RMSD) in molecular dynamics simulations is the arithmetic mean distance generated by a single atom’s motion over the period of a certain amount of time relative to a reference time [71, 72]. The RMSD of the structural atoms of a protein, such as heavy particles, backbone, C, and sidechain, were first measured before the root mean square deviation of protein-suited ligand compounds from all timescales, which were realigned and compared against the reference time (in our study 100 ns). The following equation (equation (1)) can be used to determine the RMSD of an MD simulation concerning the period of x.(1)RMSD=∑i=0Nmi∗Xi−Yi2M.N stands for the number of selected atoms, rʹ for the bit’s position in the system x after the reference system’s point has been superimposed, and j for the reference time. Most commonly, local changes in the conformational shape of proteins have been detected and tracked using the root mean square fluctuation (RMSF) [73]. The RMSF estimation of a molecular dynamic’s simulation of a protein with 2 residues can be obtained using the continuity equation (2).(2)RMSF=1T∑TjTXitj−xj2.


## 3. Result and Findings

### 3.1. Molecular Docking Result Analysis

The process of drug development has become more reliant on molecular docking. Advances in computer hardware and a rise in the number of and ease of access to small molecule and protein structures, have facilitated the development of new techniques, making docking more popular in both commercial and academic contexts since its inception in the 1980s. The methods for enhancing several drug development tasks with docking have evolved throughout time. Docking was initially created and utilized as a stand-alone technique, but it is now routinely combined with other computational techniques in integrated workflows.The desired intermolecular framework within a drug compound and a macromolecule can be developed through molecular docking. A molecular docking study was first conducted to determine the optimal intermolecular interaction between the phytochemical ligand compounds with the target protein by a group of Python tools called PyRx. Using the AutoDock Vina wizard, we performed molecular docking between 63 phytochemicals and the target protein. The phytochemical compound’s molecular docking revealed binding affinities that range from −4.5 to −8.9 kcal/mol.The top 6.35% of 63 phytochemicals (a total of 4) had a better binding affinity than nateglinide and were chosen based on binding affinity. This study utilized Nateglinide as a control ligand because of its previously reported inhibitory effect against the human insulin receptor ectodomain, which binds with a docking score of −7.1 kcal/mol. The compounds with the best docking score and binding affinity are listed in Table1.Table 1
List of the selected best four ligands and nateglinide (control drug) with their PubChem CID chemical name, two-dimensional (2D) chemical structure, and docking score.
Compound CIDCompound nameChemical structureDocking scoreCID-18867-Deacetyl-7-oxogedunin−8.9CID 122801[1S, 2R, 4S, 10R, 11R, 16R)-6-(Furan-3-yl)-1, 7, 11, 15-pentamethyl-5, 14-dioxo-3 oxapentacyclo [8.8.0.02, 4.02, 7.011, 16] octadec-12-en-18-yl] acetate−8.5CID 102090424Vilasinin−8.2CID 101306757Nimbidinin−7.9CID 5311309 (control)Nateglinide−7.1


### 3.2. MM-GBSA Result Analysis

The greater the negative free energies of the binding value, the stronger the binding between the ligand compound with the targeted protein complex (Table2). The findings showed that the Vilasinin (PubChem CID-102090424) with 1I44 protein complex possessed strong binding free energy −53.644. In contrast, PubChem CID-122801, Nimbidinin (Pubchem-101306757), 7-Deacetyl-7-oxogedunin (PubChem CID-1886), and Nateglinide (PubChem CID-5311309 (Control)) with 1I44 protein complex have been shown excellent binding free energy −48.984, −41.589, −41.461, and −40.984 consecutively.Table 2
MMGBSA calculation by maestro application of schrodinger package software.
Protein IDLigand PubChem CIDLigand nameMMGBSA dG bind score (kcal/mol)1I4418867-Deacetyl-7-oxogedunin−41.4611I44122801[1S, 2R, 4S, 10R, 11R, 16R)-6-(Furan-3-yl)-1,7,11,15-pentamethyl-5, 14-dioxo-3 oxapentacyclo [8.8.0.02, 4.02, 7.011, 16] octadec-12-en-18-yl] acetate−48.9841I44102090424Vilasinin−53.6441I44101306757Nimbidinin−41.5891I445311309 (control)Nateglinide−40.984


### 3.3. Analysis of Protein-Ligand Interactions

The interactions between protein and ligand molecules are necessary for most biological processes. All biological processes depend on the ligand-mediated signaling that occurs through molecular complementarity. At the molecular level, biological recognition is produced by this type of protein-ligand chemical interaction. Establishing specialized sites designed to bind small-molecule ligands having standard affinities adjusted to the requirements of the cell is essential for the evolution of protein function. Cooperation between ligands and the target protein is essential for regulating biological processes that clash with one another. The molecular mechanisms governing protein compatibility changes between low and high-affinity states affect cellular function through collaborative protein-ligand interactions.The creation of small-molecule drug compounds for the treatment of disease is made possible by the atomic resolution structures of protein-ligand complexes. The Ligplot + Version 2.2 tools were used to show how the selected four ligands interacted with the target protein in Table3. In terms of the chemical CID 1886 has been formed with two hydrogen bonds, Trp1200 (A) (3.34 Å) and Asp1229 (A) (3.13 Å), and a hydrophobic bond with Ser1204 (A), Glu1207 (A), Pro1231 (A), Gln1211 (A), Leu1213 (A), Phe1221 (A), Leu1228 (A), Try1227 (A), and Tyr1210 (A).Table 3
Molecular interactions between the chosen phytochemicals and the targeted receptor are tabulated in the docking score.
CompoundsDocking score (kcal/mol)Interaction of the hydrogen bondInteraction of the hydrophobic bond7-Deacetyl-7-oxogedunin (CID-1886)−8.9Trp1200 (A) (3.34 Å) and Asp1229 (A) (3.13 Å)Ser1204 (A), Glu1207 (A), Pro1231 (A), Gln1211 (A), Leu1213 (A), Phe1221 (A), Leu1228 (A), Try1227 (A), and Tyr1210 (A)[1S, 2R, 4S, 10R, 11R, 16R)-6-(Furan-3-yl)-1,7,11,15-pentamethyl-5, 14-dioxo-3 oxapentacyclo [8.8.0.02, 4.02, 7.011, 16] octadec-12-en-18-yl] acetate (CID 122801)−8.5Asp1083 (A) (2.91 Å)Val1010 (A), Gly1149 (A), Val1060 (A), Met1070 (A), Lys1030 (A), Met1139 (A), Asn1137 (A), Arg1136 (A), Gln1004 (A), Gly1082 (A), Gly1003 (A), and Leu1002 (A)Vilasinin (CID 102090424)−8.2Lys1030 (A) (3.08 Å)Ala1028 (A), Met1079 (A), Leu1078 (A), Met1039 (A), Val1010 (A), Leu1002 (A), Gly1082 (A), Gly1003 (A), Arg1136 (A), Asp1083 (A), Gln1004 (A), Asn1137 (A), Ser1006 (A), and Gly1005 (A)Nimbidinin (CID 101306757)−7.9Asp1232 (A) (2.82 Å) and Glu1207 (A) (3.30 Å)Pro1104 (A), Thr1203 (A), Asn1233 (A), Met1109 (A), Pro1235 (A), Pro1103 (A), Ser1204 (A), Leu1205 (A), and Pro1231 (A)Nateglinide (CID 5311309)−7.1Met1079 (A) (2.94 Å)Ala1029 (A), Lys1030 (A), Val1060 (A), Met1139 (A), Met1076 (A), Val1010 (A), Leu1002 (A), Asp1083 (A), Gly1083 (A), Gly1082 (A), and Leu1078 (A)In terms of the chemical CID 122801, it has been revealed that when it comes into contact with the appropriate protein PDB-1I44. From the Ligplot + visualizer, it can be visualized that the hydrogen bond has been formed with Asp1083 (A) (2.91 Å) and supported by 12 hydrophobic bonds (Val1010 (A), Gly1149 (A), Val1060 (A), Met1070 (A), Lys1030 (A), Met1139 (A), Asn1137 (A), Arg1136 (A), Gln1004 (A), Gly1082 (A), Gly1003 (A), and Leu1002 (A)). “The insulin Receptor Ectodomain-Vilasinin” complex was supported by one hydrogen bond Lys1030 (A) (3.08 Å) and fourteen hydrophobic bonds (Ala1028 (A), Met1079 (A), Leu1078 (A), Met1039 (A), Val1010 (A), Leu1002 (A), Gly1082 (A), Gly1003 (A), Arg1136 (A), Asp1083 (A), Gln1004 (A), Asn1137 (A), Ser1006 (A), and Gly1005 (A)) (Figure1), whereas Nimbidinin interacts with insulin receptor ectodomain protease was supported by two hydrogen bonds (Asp1232 (A) (2.82 Å) and Glu1207 (A) (3.30 Å) and nine hydrophobic bonds (Pro1104 (A), Thr1203 (A), Asn1233 (A), Met1109 (A), Pro1235 (A), Pro1103 (A), Ser1204 (A), Leu1205 (A), and Pro1231 (A)) (Figure 1). In addition, the control drug nateglinide possesses the interaction with the targeted protein by one hydrogen bond with Met1079 (A) (2.94 Å) and eleven hydrophobic bonds, including Ala1029 (A), Lys1030 (A), Val1060 (A), Met1139 (A), Met1076 (A), Val1010 (A), Leu1002 (A), Asp1083 (A), Gly1083 (A), Gly1082 (A), and Leu1078 (A).Figure 1
Interaction between the insulin receptor ectodomain (PDB 1I44) with compound CID 1886, CID 122801, CID 102090424, and CID 101306757.


### 3.4. Pharmacokinetic Property and QSAR Analysis

The study of a small molecule candidate’s dynamic motions inside the body, and the ADME characteristics of a chemical that resembles a drug, is known as pharmacokinetics. The term “PK” refers to a class of xenobiotic regulatory procedures that is necessary throughout the drug discovery process, especially in preclinical studies. It employs several mathematical equations to create a model and offers data on the xenobiotics in the body over time. In addition, PK traits can aid in understanding and predicting biological functions, such as a compound’s harmful or beneficial effects. The PKCSM server and SwissADME server were utilized in the study of the necessary PK parameters of the four selected potential phytocompounds that were chosen. The ADME features of the selected compounds, such as lipophilicity (dissolve in oils, fats, and nonpolar solvents), drug-likeness, and water solubility, have been found using the server. In the investigation, the PK characteristics of all of the identified compounds were found to be efficacious and druggable. The compounds studied in this paper have no violation according to Lipinski’s rule of 5 and Ghose and colleague’s rule. Pk properties, including several hydrogen bond donors and acceptors, molecular weight, molar refractivity, and Logp value, are noted in Tables 4 and 5.Table 4
The complete pharmacophore and pharmacokinetic profiling of the selected ligands.
Compound IDMW (g/mol)HBAHBDLogPMRIAbTClCaPCID-1886438.52604.1981115.3499.3940.0930.883CID-122801466.57604.6291124.961000.331.154CID-102090424428.569533.6435117.1197.4030.360.702CID-101306757442.54632.8225117.3186.8471.212.094CID 5311309 (control)317.42224.198191.7599.4050.1080.878Units of physiochemical properties: MW, molecular weight (g/mol); HBA, hydrogen bond acceptor; HBD, hydrogen bond donor; LogP, estimated octanol/coefficient of liquid fraction; MR, molecular refractivity. Units of pharmacokinetics properties: IAb, intestinal absorption (% absorbed); TCl, total clearance (log mL/min/kg); Caco2 permeability (log Papp in 10−6 cm/s).Table 5
Table depicting the toxicity properties of selected compounds showing compound CID, AMES toxicity, max. Tolerated dose (human), hERG I inhibitor, hERG II inhibitor, oral rat acute toxicity, (LD50) oral rat chronic toxicity (LOAEL), hepato-toxicity, skin sensitisation, T. pyriformis toxicity, and minnow toxicity.
Compound CIDAMES toxicityMax. tolerated dose (human)hERG I inhibitorhERG II inhibitorOral rat acute toxicity (LD50)Oral rat chronic toxicity (LOAEL)Hepato-toxicitySkin sensitisationT. pyriformis toxicityMinnow toxicityCID 1886No−0.48NoNo2.7910.66NoNo0.304−0.07CID 122801No−0.658NoNo2.6940.534NoNo0.3060.243CID 102090424No−0.69NoNo2.7161.715YesNo0.3591.129CID 101306757No−0.625NoNo3.712.28NoNo0.2872.354CID 5311309 (control)No0.292NoNo2.0221.977NoNo0.2840.796Table4 has depicted important pharmacokinetic parameters to show the credentials of the selected phytocompounds to be selected as a preferable alternative for present-age diabetic medicine. The molecular weight of the following drug compounds is within the range of 350–480 which firmly supports the ADMET rules of both Lipinski’s Rule of 5 and Ghose and Colleagues’ rules. All the compounds possess more hydrogen bond donors than the control drug. In the case of hydrogen bonds acceptors, CID 1886 and CID 122801 don’t possess any hydrogen bond acceptors, whereas CID 102090424 and CID 101306757 contain 3 hydrogen bond acceptors each, and the control drug CID 5311309 got 2 hydrogen bond acceptors. The logP values are reported as 4.1981, 4.6291, 3.6435, 2.82225, and 4.1981 for compound CID 1886, CID122801, CID 102090424, CID 101306757, and CID 5311309 (Control), respectively. Other important pharmacokinetic parameters including molar refractivity, Intestinal absorption, total clearance, and CaCO2 permeability also have been reported within a standard range as shown in Table 4.In Table5, different toxicity-related parameters for the selected phytocompounds have been demonstrated. From Table 5, we can see that no compounds have shown AMES toxicity or working as hERG I or hERG II inhibitors. The maximum tolerated dose in humans for CID 1886 is 0.48; for CID 122801, the value is −0.658. For CID 102090424 compound, the maximum tolerated dose is −0.69. For CID 101306757 compound, the maximum tolerated dose is −0.625. And lastly, for CID 5311309 (control) compound, the maximum tolerated dose is 0.292. From the column where the rat acute toxicity (LD50) has been depicted, we can see the LD50 values are the compounds are 2.791, 2.694, 2.716, 3.71, and 2.022 for CID 1886, CID 122801, CID 102090424, CID 101306757, and CID 5311309 (control) respectively. The oral rat chronic toxicity (LOAEL) values of those compounds are 0.66, 0.534, 1.715, 2.28, and 1.977 followed by the CID sequence mentioned previously. This work also evaluated the hepatotoxicity and skin sensitiveness of the compounds. It has been found that only CID 102090424 has hepato-toxicity, whereas all the other compounds are reported safe from hepatotoxicity as well as skin sensation parameters in the following in silico study. The T. pyriformis toxicity and minnow toxicity also showed the druggability of the selected phytocompounds.Moreover, the quantitative structure-activity relationship (QSAR) assessment showed that these ligand compounds possessed potent antidiabetic, anticancer, and antiviral activity (Table6). Here, the compounds CID-122801, CID-102090424, CID 101306757, and CID 5311309 (control) possessed potent antidiabetic properties, whereas CID-122801, CID-102090424, CID 101306757, and CID 5311309 (control) have reported antiviral activity. All the compounds also showed potential to be selected as an anticancer compound according to the QSAR analysis.Table 6
QSAR-based bioactivity prediction for ligand validation.
Compounds CIDPrediction of activity spectra for substances (Pa = 0.3 to 0.7)AntidiabeticAnticancerAntiviralCID-1886×√√CID-122801√√√CID-102090424√√√CID-101306757√√×CID 5311309 (control)√√√QSAR: quantitative structure-activity relationship; Pa: prediction of activity score as per PASS server. Pa score of 0.3 to 0.7 signifies moderate activity. ×: Pa score less than 0.3;√: Pa score in the range 0.3 to 0.7.


### 3.5. Analysis of MD Simulation Parameters

The conformational stability and several intramolecular interaction parameters of a protein-ligand complex are studied in real time using MD simulation in CADD. The conformational alteration that a complex system experiences when it is placed in an artificial environment can also be ascertained using this method. To comprehend the conformational variations and rearrangements of the protein in complex, a 100 ns MD simulation of the protein in combination with the particular ligand was carried out in this study. Initially, the terminal snapshots from the 100 ns MDS trajectories were used to examine intermolecular behavior.The average change in the RMSD of the protein-ligand interaction is quite acceptable, with a range of 1–4. If the RMSD value is greater than 1–4, a significant conformational shift in the protein structure has occurred. To examine the conformational change of the target protein combination with the two ligand molecules included by the 4-ligand molecule and control drug, a 100 ns MD simulation was performed, and the related RMSD value was determined.The ligand complex vilasinin, which had fewer fluctuations, had an average RMSD that ranged from 1.2 to 2.4. The RMSD value of the compound, which changed within a very small range and dropped outside the permissible range, indicates that the protein-ligand complex structure depicted in Figure2 is conformationally stable.Figure 2
The RMSD values for the 4 ligand molecules and the nateglinide control for the 1i44 protein model. The selected compounds with the CIDs 1886, 122801, 101306757, 102090424, and 5311309 in relationship with the protein were represented by the hues dark blue, orange, gray, yellow, and light blue, respectively, while the control medicine had the CIDs 5311309 in association with the protein.The RMSF can facilitate in identifying and characterizing the adjacent alterations that arise within the configuration of the protein when a particular ligand compound attaches with remnants. The RMSF of the compounds with the targeted protein model was established for analyzing the changes in protein structure resulting from the attachment of specific ligands with a particular residual site, as illustrated in Figure3.Figure 3
The graph displays the RMSF values for the 1i44 protein model using 4 ligand compounds and 1 control drug that were recovered from the C atoms of the complex system. The markers blue, orange, gray, yellow, and light blue are used to depict the selected compounds and had the CIDs 1886, 122801, 101306757, 102090424, and 5311309, respectively.Alpha helices and beta strands, which are among the most rigid secondary structural elements, were shown to have a minimum observation rate between 7 and 275 amino acid residues. Due to the existence of these domains, the C- and N- termini of the protein exhibit the majority of the variation. Therefore, it is discovered that there is a low fluctuation probability for the displacement of a single atom in the four-ligand complex and control medication under examination in the simulation environment.The radius of gyration (RG) of the interaction between protein and ligand can be determined by the positioning of an interaction system’s atoms across its axis. Since it divulges alteration in complex compactness gradually, the computation of the radius of gyration is one of the most crucial signs to look for when predicting the structural functioning of a macromolecule.The stability of 7-Deacetyl-7-oxogedunin, [(1S, 2R, 4S, 10R, 11R, 16R)-6-(Furan-3-yl)-1,7,11,15,15-pentamethyl-5,14-dioxo-3 oxapenta-cyclo [8.8.0.02, 4.02, 7.011, 16]octadec-12-en-18-yl] acetate, Vilasinin, Nimbidinin, and Nateglinide together with the target protein was also explored in terms of the radius of gyration over a 100 ns simulation run, as depicted in Figure4. The mean rg value for 7-Deacetyl-7-oxogedunin is 3.75, and for CID 122801, CID 102090424, CID 5311309, and Nateglinide (control drug), the average rg value is 4. This result shows that the protein’s binding site does not undergo significant structural alterations upon binding the selected ligand compounds.Figure 4
The radius of gyration (Rg) of the protein-ligand interaction was calculated using the 100 ns simulation. The markers blue, orange, gray, yellow, and light blue are used to depict the selected ligand compounds that contain the CID 1886, CID 122801, CID 101306757, CID 102090424, and CID 5311309, respectively, in association with the protein.The confirmation of structural identity and the actions of biological macromolecules are governed by SASA concentration. When amino acid residues are attached to a surface of a protein, it functions as active sites and/or get interconnected with other molecules and ligands, a molecule’s solventlike behavior (hydrophilic or hydrophobic), and the elements of protein-ligand interactions are better known.The SASA for the selected compounds was in the range of 50 to 40, showing standard exposure of the selected compound with an amino acid residue in the complex systems (Figure5).Figure 5
From the 100 ns simulated interaction diagram, the solvent accessible surface area (SASA) of the protein-ligand interaction compounds was estimated. The markers blue, orange, gray, yellow, and light blue are used to depict the protein in association with the chosen ligand molecules CID 1886, CID 122801, CID 101306757, CID 102090424, and CID 5311309, respectively.The MolSA is equivalent to the Van der Waals surface area, which was discovered using a 1.4 probe radius. In our in silico analysis, each ligand had the standard Van der Waals surface area (Figure6). In addition, only molecule of nitrogen and oxygen atoms contribute to PSA. Here, the targeted protein showed a significant PSA value for each ligand molecule (Figure 7).Figure 6
The 100 ns simulated interaction diagram was used to calculate the molecular surface area (MolSA) of the protein-ligand interaction complexes. The markers blue, orange, gray, yellow, and light blue are used to depict the selected ligand compounds CID 1886, CID 122801, CID 101306757, CID 102090424, and CID 5311309, respectively, in association with the protein.Figure 7
From the 100 ns simulated interaction diagram, the polar surface area (PSA) of the protein-ligand interaction compounds was estimated. The markers blue, orange, gray, yellow, and light blue are used to depict the selected ligand compounds CID 1886, CID 122801, CID 101306757, CID 102090424, and CID 5311309, respectively, in association with the protein.Using the simulation interactions diagram, the intricate structure of a protein with the indicated ligands and their intermolecular interactions have been studied during a simulation time of 100 ns (SID).The interactions between the chosen ligand molecules and proteins are influenced by various factors, including the hydrogen bond, noncovalent bond (hydrophobic bond), ionic bond, and water bridge bond. Figure8 evaluates and depicts the chosen ligand compounds and the control drug. To help build a stable binding with the intended protein, all compounds generated many connections throughout the course of the 100 ns simulation via hydrogen, hydrophobic, ionic, and water bridge bonding and maintained these contacts up until the simulation was over.Figure 8
The interactions between proteins and ligands identified during the 100 ns simulation are shown in the stacked bar charts. In this section, the interactions of the chosen ligand molecules are demonstrated. Figures (a)–(e) in the image above represented, respectively, the intermolecular fractions for the compounds CID 1886, CID 122801, CID 101306757, CID 102090424, and CID 5311309.
(a)(b)(c)(d)(e)


## 3.1. Molecular Docking Result Analysis

The process of drug development has become more reliant on molecular docking. Advances in computer hardware and a rise in the number of and ease of access to small molecule and protein structures, have facilitated the development of new techniques, making docking more popular in both commercial and academic contexts since its inception in the 1980s. The methods for enhancing several drug development tasks with docking have evolved throughout time. Docking was initially created and utilized as a stand-alone technique, but it is now routinely combined with other computational techniques in integrated workflows.The desired intermolecular framework within a drug compound and a macromolecule can be developed through molecular docking. A molecular docking study was first conducted to determine the optimal intermolecular interaction between the phytochemical ligand compounds with the target protein by a group of Python tools called PyRx. Using the AutoDock Vina wizard, we performed molecular docking between 63 phytochemicals and the target protein. The phytochemical compound’s molecular docking revealed binding affinities that range from −4.5 to −8.9 kcal/mol.The top 6.35% of 63 phytochemicals (a total of 4) had a better binding affinity than nateglinide and were chosen based on binding affinity. This study utilized Nateglinide as a control ligand because of its previously reported inhibitory effect against the human insulin receptor ectodomain, which binds with a docking score of −7.1 kcal/mol. The compounds with the best docking score and binding affinity are listed in Table1.Table 1
List of the selected best four ligands and nateglinide (control drug) with their PubChem CID chemical name, two-dimensional (2D) chemical structure, and docking score.
Compound CIDCompound nameChemical structureDocking scoreCID-18867-Deacetyl-7-oxogedunin−8.9CID 122801[1S, 2R, 4S, 10R, 11R, 16R)-6-(Furan-3-yl)-1, 7, 11, 15-pentamethyl-5, 14-dioxo-3 oxapentacyclo [8.8.0.02, 4.02, 7.011, 16] octadec-12-en-18-yl] acetate−8.5CID 102090424Vilasinin−8.2CID 101306757Nimbidinin−7.9CID 5311309 (control)Nateglinide−7.1


## 3.2. MM-GBSA Result Analysis

The greater the negative free energies of the binding value, the stronger the binding between the ligand compound with the targeted protein complex (Table2). The findings showed that the Vilasinin (PubChem CID-102090424) with 1I44 protein complex possessed strong binding free energy −53.644. In contrast, PubChem CID-122801, Nimbidinin (Pubchem-101306757), 7-Deacetyl-7-oxogedunin (PubChem CID-1886), and Nateglinide (PubChem CID-5311309 (Control)) with 1I44 protein complex have been shown excellent binding free energy −48.984, −41.589, −41.461, and −40.984 consecutively.Table 2
MMGBSA calculation by maestro application of schrodinger package software.
Protein IDLigand PubChem CIDLigand nameMMGBSA dG bind score (kcal/mol)1I4418867-Deacetyl-7-oxogedunin−41.4611I44122801[1S, 2R, 4S, 10R, 11R, 16R)-6-(Furan-3-yl)-1,7,11,15-pentamethyl-5, 14-dioxo-3 oxapentacyclo [8.8.0.02, 4.02, 7.011, 16] octadec-12-en-18-yl] acetate−48.9841I44102090424Vilasinin−53.6441I44101306757Nimbidinin−41.5891I445311309 (control)Nateglinide−40.984


## 3.3. Analysis of Protein-Ligand Interactions

The interactions between protein and ligand molecules are necessary for most biological processes. All biological processes depend on the ligand-mediated signaling that occurs through molecular complementarity. At the molecular level, biological recognition is produced by this type of protein-ligand chemical interaction. Establishing specialized sites designed to bind small-molecule ligands having standard affinities adjusted to the requirements of the cell is essential for the evolution of protein function. Cooperation between ligands and the target protein is essential for regulating biological processes that clash with one another. The molecular mechanisms governing protein compatibility changes between low and high-affinity states affect cellular function through collaborative protein-ligand interactions.The creation of small-molecule drug compounds for the treatment of disease is made possible by the atomic resolution structures of protein-ligand complexes. The Ligplot + Version 2.2 tools were used to show how the selected four ligands interacted with the target protein in Table3. In terms of the chemical CID 1886 has been formed with two hydrogen bonds, Trp1200 (A) (3.34 Å) and Asp1229 (A) (3.13 Å), and a hydrophobic bond with Ser1204 (A), Glu1207 (A), Pro1231 (A), Gln1211 (A), Leu1213 (A), Phe1221 (A), Leu1228 (A), Try1227 (A), and Tyr1210 (A).Table 3
Molecular interactions between the chosen phytochemicals and the targeted receptor are tabulated in the docking score.
CompoundsDocking score (kcal/mol)Interaction of the hydrogen bondInteraction of the hydrophobic bond7-Deacetyl-7-oxogedunin (CID-1886)−8.9Trp1200 (A) (3.34 Å) and Asp1229 (A) (3.13 Å)Ser1204 (A), Glu1207 (A), Pro1231 (A), Gln1211 (A), Leu1213 (A), Phe1221 (A), Leu1228 (A), Try1227 (A), and Tyr1210 (A)[1S, 2R, 4S, 10R, 11R, 16R)-6-(Furan-3-yl)-1,7,11,15-pentamethyl-5, 14-dioxo-3 oxapentacyclo [8.8.0.02, 4.02, 7.011, 16] octadec-12-en-18-yl] acetate (CID 122801)−8.5Asp1083 (A) (2.91 Å)Val1010 (A), Gly1149 (A), Val1060 (A), Met1070 (A), Lys1030 (A), Met1139 (A), Asn1137 (A), Arg1136 (A), Gln1004 (A), Gly1082 (A), Gly1003 (A), and Leu1002 (A)Vilasinin (CID 102090424)−8.2Lys1030 (A) (3.08 Å)Ala1028 (A), Met1079 (A), Leu1078 (A), Met1039 (A), Val1010 (A), Leu1002 (A), Gly1082 (A), Gly1003 (A), Arg1136 (A), Asp1083 (A), Gln1004 (A), Asn1137 (A), Ser1006 (A), and Gly1005 (A)Nimbidinin (CID 101306757)−7.9Asp1232 (A) (2.82 Å) and Glu1207 (A) (3.30 Å)Pro1104 (A), Thr1203 (A), Asn1233 (A), Met1109 (A), Pro1235 (A), Pro1103 (A), Ser1204 (A), Leu1205 (A), and Pro1231 (A)Nateglinide (CID 5311309)−7.1Met1079 (A) (2.94 Å)Ala1029 (A), Lys1030 (A), Val1060 (A), Met1139 (A), Met1076 (A), Val1010 (A), Leu1002 (A), Asp1083 (A), Gly1083 (A), Gly1082 (A), and Leu1078 (A)In terms of the chemical CID 122801, it has been revealed that when it comes into contact with the appropriate protein PDB-1I44. From the Ligplot + visualizer, it can be visualized that the hydrogen bond has been formed with Asp1083 (A) (2.91 Å) and supported by 12 hydrophobic bonds (Val1010 (A), Gly1149 (A), Val1060 (A), Met1070 (A), Lys1030 (A), Met1139 (A), Asn1137 (A), Arg1136 (A), Gln1004 (A), Gly1082 (A), Gly1003 (A), and Leu1002 (A)). “The insulin Receptor Ectodomain-Vilasinin” complex was supported by one hydrogen bond Lys1030 (A) (3.08 Å) and fourteen hydrophobic bonds (Ala1028 (A), Met1079 (A), Leu1078 (A), Met1039 (A), Val1010 (A), Leu1002 (A), Gly1082 (A), Gly1003 (A), Arg1136 (A), Asp1083 (A), Gln1004 (A), Asn1137 (A), Ser1006 (A), and Gly1005 (A)) (Figure1), whereas Nimbidinin interacts with insulin receptor ectodomain protease was supported by two hydrogen bonds (Asp1232 (A) (2.82 Å) and Glu1207 (A) (3.30 Å) and nine hydrophobic bonds (Pro1104 (A), Thr1203 (A), Asn1233 (A), Met1109 (A), Pro1235 (A), Pro1103 (A), Ser1204 (A), Leu1205 (A), and Pro1231 (A)) (Figure 1). In addition, the control drug nateglinide possesses the interaction with the targeted protein by one hydrogen bond with Met1079 (A) (2.94 Å) and eleven hydrophobic bonds, including Ala1029 (A), Lys1030 (A), Val1060 (A), Met1139 (A), Met1076 (A), Val1010 (A), Leu1002 (A), Asp1083 (A), Gly1083 (A), Gly1082 (A), and Leu1078 (A).Figure 1
Interaction between the insulin receptor ectodomain (PDB 1I44) with compound CID 1886, CID 122801, CID 102090424, and CID 101306757.


## 3.4. Pharmacokinetic Property and QSAR Analysis

The study of a small molecule candidate’s dynamic motions inside the body, and the ADME characteristics of a chemical that resembles a drug, is known as pharmacokinetics. The term “PK” refers to a class of xenobiotic regulatory procedures that is necessary throughout the drug discovery process, especially in preclinical studies. It employs several mathematical equations to create a model and offers data on the xenobiotics in the body over time. In addition, PK traits can aid in understanding and predicting biological functions, such as a compound’s harmful or beneficial effects. The PKCSM server and SwissADME server were utilized in the study of the necessary PK parameters of the four selected potential phytocompounds that were chosen. The ADME features of the selected compounds, such as lipophilicity (dissolve in oils, fats, and nonpolar solvents), drug-likeness, and water solubility, have been found using the server. In the investigation, the PK characteristics of all of the identified compounds were found to be efficacious and druggable. The compounds studied in this paper have no violation according to Lipinski’s rule of 5 and Ghose and colleague’s rule. Pk properties, including several hydrogen bond donors and acceptors, molecular weight, molar refractivity, and Logp value, are noted in Tables 4 and 5.Table 4
The complete pharmacophore and pharmacokinetic profiling of the selected ligands.
Compound IDMW (g/mol)HBAHBDLogPMRIAbTClCaPCID-1886438.52604.1981115.3499.3940.0930.883CID-122801466.57604.6291124.961000.331.154CID-102090424428.569533.6435117.1197.4030.360.702CID-101306757442.54632.8225117.3186.8471.212.094CID 5311309 (control)317.42224.198191.7599.4050.1080.878Units of physiochemical properties: MW, molecular weight (g/mol); HBA, hydrogen bond acceptor; HBD, hydrogen bond donor; LogP, estimated octanol/coefficient of liquid fraction; MR, molecular refractivity. Units of pharmacokinetics properties: IAb, intestinal absorption (% absorbed); TCl, total clearance (log mL/min/kg); Caco2 permeability (log Papp in 10−6 cm/s).Table 5
Table depicting the toxicity properties of selected compounds showing compound CID, AMES toxicity, max. Tolerated dose (human), hERG I inhibitor, hERG II inhibitor, oral rat acute toxicity, (LD50) oral rat chronic toxicity (LOAEL), hepato-toxicity, skin sensitisation, T. pyriformis toxicity, and minnow toxicity.
Compound CIDAMES toxicityMax. tolerated dose (human)hERG I inhibitorhERG II inhibitorOral rat acute toxicity (LD50)Oral rat chronic toxicity (LOAEL)Hepato-toxicitySkin sensitisationT. pyriformis toxicityMinnow toxicityCID 1886No−0.48NoNo2.7910.66NoNo0.304−0.07CID 122801No−0.658NoNo2.6940.534NoNo0.3060.243CID 102090424No−0.69NoNo2.7161.715YesNo0.3591.129CID 101306757No−0.625NoNo3.712.28NoNo0.2872.354CID 5311309 (control)No0.292NoNo2.0221.977NoNo0.2840.796Table4 has depicted important pharmacokinetic parameters to show the credentials of the selected phytocompounds to be selected as a preferable alternative for present-age diabetic medicine. The molecular weight of the following drug compounds is within the range of 350–480 which firmly supports the ADMET rules of both Lipinski’s Rule of 5 and Ghose and Colleagues’ rules. All the compounds possess more hydrogen bond donors than the control drug. In the case of hydrogen bonds acceptors, CID 1886 and CID 122801 don’t possess any hydrogen bond acceptors, whereas CID 102090424 and CID 101306757 contain 3 hydrogen bond acceptors each, and the control drug CID 5311309 got 2 hydrogen bond acceptors. The logP values are reported as 4.1981, 4.6291, 3.6435, 2.82225, and 4.1981 for compound CID 1886, CID122801, CID 102090424, CID 101306757, and CID 5311309 (Control), respectively. Other important pharmacokinetic parameters including molar refractivity, Intestinal absorption, total clearance, and CaCO2 permeability also have been reported within a standard range as shown in Table 4.In Table5, different toxicity-related parameters for the selected phytocompounds have been demonstrated. From Table 5, we can see that no compounds have shown AMES toxicity or working as hERG I or hERG II inhibitors. The maximum tolerated dose in humans for CID 1886 is 0.48; for CID 122801, the value is −0.658. For CID 102090424 compound, the maximum tolerated dose is −0.69. For CID 101306757 compound, the maximum tolerated dose is −0.625. And lastly, for CID 5311309 (control) compound, the maximum tolerated dose is 0.292. From the column where the rat acute toxicity (LD50) has been depicted, we can see the LD50 values are the compounds are 2.791, 2.694, 2.716, 3.71, and 2.022 for CID 1886, CID 122801, CID 102090424, CID 101306757, and CID 5311309 (control) respectively. The oral rat chronic toxicity (LOAEL) values of those compounds are 0.66, 0.534, 1.715, 2.28, and 1.977 followed by the CID sequence mentioned previously. This work also evaluated the hepatotoxicity and skin sensitiveness of the compounds. It has been found that only CID 102090424 has hepato-toxicity, whereas all the other compounds are reported safe from hepatotoxicity as well as skin sensation parameters in the following in silico study. The T. pyriformis toxicity and minnow toxicity also showed the druggability of the selected phytocompounds.Moreover, the quantitative structure-activity relationship (QSAR) assessment showed that these ligand compounds possessed potent antidiabetic, anticancer, and antiviral activity (Table6). Here, the compounds CID-122801, CID-102090424, CID 101306757, and CID 5311309 (control) possessed potent antidiabetic properties, whereas CID-122801, CID-102090424, CID 101306757, and CID 5311309 (control) have reported antiviral activity. All the compounds also showed potential to be selected as an anticancer compound according to the QSAR analysis.Table 6
QSAR-based bioactivity prediction for ligand validation.
Compounds CIDPrediction of activity spectra for substances (Pa = 0.3 to 0.7)AntidiabeticAnticancerAntiviralCID-1886×√√CID-122801√√√CID-102090424√√√CID-101306757√√×CID 5311309 (control)√√√QSAR: quantitative structure-activity relationship; Pa: prediction of activity score as per PASS server. Pa score of 0.3 to 0.7 signifies moderate activity. ×: Pa score less than 0.3;√: Pa score in the range 0.3 to 0.7.


## 3.5. Analysis of MD Simulation Parameters

The conformational stability and several intramolecular interaction parameters of a protein-ligand complex are studied in real time using MD simulation in CADD. The conformational alteration that a complex system experiences when it is placed in an artificial environment can also be ascertained using this method. To comprehend the conformational variations and rearrangements of the protein in complex, a 100 ns MD simulation of the protein in combination with the particular ligand was carried out in this study. Initially, the terminal snapshots from the 100 ns MDS trajectories were used to examine intermolecular behavior.The average change in the RMSD of the protein-ligand interaction is quite acceptable, with a range of 1–4. If the RMSD value is greater than 1–4, a significant conformational shift in the protein structure has occurred. To examine the conformational change of the target protein combination with the two ligand molecules included by the 4-ligand molecule and control drug, a 100 ns MD simulation was performed, and the related RMSD value was determined.The ligand complex vilasinin, which had fewer fluctuations, had an average RMSD that ranged from 1.2 to 2.4. The RMSD value of the compound, which changed within a very small range and dropped outside the permissible range, indicates that the protein-ligand complex structure depicted in Figure2 is conformationally stable.Figure 2
The RMSD values for the 4 ligand molecules and the nateglinide control for the 1i44 protein model. The selected compounds with the CIDs 1886, 122801, 101306757, 102090424, and 5311309 in relationship with the protein were represented by the hues dark blue, orange, gray, yellow, and light blue, respectively, while the control medicine had the CIDs 5311309 in association with the protein.The RMSF can facilitate in identifying and characterizing the adjacent alterations that arise within the configuration of the protein when a particular ligand compound attaches with remnants. The RMSF of the compounds with the targeted protein model was established for analyzing the changes in protein structure resulting from the attachment of specific ligands with a particular residual site, as illustrated in Figure3.Figure 3
The graph displays the RMSF values for the 1i44 protein model using 4 ligand compounds and 1 control drug that were recovered from the C atoms of the complex system. The markers blue, orange, gray, yellow, and light blue are used to depict the selected compounds and had the CIDs 1886, 122801, 101306757, 102090424, and 5311309, respectively.Alpha helices and beta strands, which are among the most rigid secondary structural elements, were shown to have a minimum observation rate between 7 and 275 amino acid residues. Due to the existence of these domains, the C- and N- termini of the protein exhibit the majority of the variation. Therefore, it is discovered that there is a low fluctuation probability for the displacement of a single atom in the four-ligand complex and control medication under examination in the simulation environment.The radius of gyration (RG) of the interaction between protein and ligand can be determined by the positioning of an interaction system’s atoms across its axis. Since it divulges alteration in complex compactness gradually, the computation of the radius of gyration is one of the most crucial signs to look for when predicting the structural functioning of a macromolecule.The stability of 7-Deacetyl-7-oxogedunin, [(1S, 2R, 4S, 10R, 11R, 16R)-6-(Furan-3-yl)-1,7,11,15,15-pentamethyl-5,14-dioxo-3 oxapenta-cyclo [8.8.0.02, 4.02, 7.011, 16]octadec-12-en-18-yl] acetate, Vilasinin, Nimbidinin, and Nateglinide together with the target protein was also explored in terms of the radius of gyration over a 100 ns simulation run, as depicted in Figure4. The mean rg value for 7-Deacetyl-7-oxogedunin is 3.75, and for CID 122801, CID 102090424, CID 5311309, and Nateglinide (control drug), the average rg value is 4. This result shows that the protein’s binding site does not undergo significant structural alterations upon binding the selected ligand compounds.Figure 4
The radius of gyration (Rg) of the protein-ligand interaction was calculated using the 100 ns simulation. The markers blue, orange, gray, yellow, and light blue are used to depict the selected ligand compounds that contain the CID 1886, CID 122801, CID 101306757, CID 102090424, and CID 5311309, respectively, in association with the protein.The confirmation of structural identity and the actions of biological macromolecules are governed by SASA concentration. When amino acid residues are attached to a surface of a protein, it functions as active sites and/or get interconnected with other molecules and ligands, a molecule’s solventlike behavior (hydrophilic or hydrophobic), and the elements of protein-ligand interactions are better known.The SASA for the selected compounds was in the range of 50 to 40, showing standard exposure of the selected compound with an amino acid residue in the complex systems (Figure5).Figure 5
From the 100 ns simulated interaction diagram, the solvent accessible surface area (SASA) of the protein-ligand interaction compounds was estimated. The markers blue, orange, gray, yellow, and light blue are used to depict the protein in association with the chosen ligand molecules CID 1886, CID 122801, CID 101306757, CID 102090424, and CID 5311309, respectively.The MolSA is equivalent to the Van der Waals surface area, which was discovered using a 1.4 probe radius. In our in silico analysis, each ligand had the standard Van der Waals surface area (Figure6). In addition, only molecule of nitrogen and oxygen atoms contribute to PSA. Here, the targeted protein showed a significant PSA value for each ligand molecule (Figure 7).Figure 6
The 100 ns simulated interaction diagram was used to calculate the molecular surface area (MolSA) of the protein-ligand interaction complexes. The markers blue, orange, gray, yellow, and light blue are used to depict the selected ligand compounds CID 1886, CID 122801, CID 101306757, CID 102090424, and CID 5311309, respectively, in association with the protein.Figure 7
From the 100 ns simulated interaction diagram, the polar surface area (PSA) of the protein-ligand interaction compounds was estimated. The markers blue, orange, gray, yellow, and light blue are used to depict the selected ligand compounds CID 1886, CID 122801, CID 101306757, CID 102090424, and CID 5311309, respectively, in association with the protein.Using the simulation interactions diagram, the intricate structure of a protein with the indicated ligands and their intermolecular interactions have been studied during a simulation time of 100 ns (SID).The interactions between the chosen ligand molecules and proteins are influenced by various factors, including the hydrogen bond, noncovalent bond (hydrophobic bond), ionic bond, and water bridge bond. Figure8 evaluates and depicts the chosen ligand compounds and the control drug. To help build a stable binding with the intended protein, all compounds generated many connections throughout the course of the 100 ns simulation via hydrogen, hydrophobic, ionic, and water bridge bonding and maintained these contacts up until the simulation was over.Figure 8
The interactions between proteins and ligands identified during the 100 ns simulation are shown in the stacked bar charts. In this section, the interactions of the chosen ligand molecules are demonstrated. Figures (a)–(e) in the image above represented, respectively, the intermolecular fractions for the compounds CID 1886, CID 122801, CID 101306757, CID 102090424, and CID 5311309.
(a)(b)(c)(d)(e)


## 4. Discussion

A drug candidate must have good activity against the therapeutic target to achieve the appropriate ADMET properties at a therapeutic dose [74]. ADME, linked to pharmacokinetic characteristics, significantly impacts a drug’s activity. The pharmacokinetic parameters are obligatory to create a medication candidate that will successfully complete the required clinical trials and must be adjusted in time for the drug design process [75–77]. According to Lipinski and Ghose’s drug likeliness qualities, the compounds investigated in this investigation did not exhibit severe toxicity or any other violations. The pharmacokinetics and drug-likeness properties of the investigated compounds were calculated. The compound selected after docking with the most promising result shows no violation in both Lipinski’s rule of 5 and Ghose and colleagues’ rule. The selected ligands followed the following parameters: 160 ≤ MW ≤ 480, −0.4 ≤ WLOGP ≤ 5.6, MLOGP ≤ 4.15, N or O ≤ 10, NH or OH ≤ 5, 40 ≤ Molar Refractivity ≤ 130, and 20 ≤ Atoms ≤ 70. No compound has shown any adverse toxic effect during the in-silico toxicity analysis [78–80]. Furthermore, the QSAR analysis depicted that all the compounds possessed effective antidiabetic, anticancer, and antiviral activity.Analysis of biomolecular interactions and the interface between protein architecture and activity can benefit the development of new medications and performance data from dynamic trajectory analysis, also known as the molecular dynamic simulation (MDS) [18, 81]. The Schrödinger package software from the Desmond application was utilized in our study to execute a 100 ns MDS with the selected physiologic and physicochemical parameters. This simulation trajectory of the simulation tool has allowed accurate analysis of the root mean square deviation (RMSD), root mean square fluctuation (RMSF), the radius of gyration (Rg), hydrogen bond number, and solvent accessible surface area (SASA) [82]. The RMSD of the selected 1i44 model protein was used to evaluate the protein structure’s dependability and conformational fluctuations; a lower value indicates the most stable molecules. RMSD values under 1.5 are typically suggestive of more consistency in docking because RMSD values over 1.5 represent the average binding positions. The RMSD values for the protein-ligand interactions in our study were within a typical range, i.e., average mean values of 4 (the lowest value for the selected ligand compounds is approximately 1.0, and maximum values of 4, indicating a more advantageous docking position and no disruption of the protein-ligand structure). In Figure 2, The average protein fluctuations can be measured using the RMSF as a reference point. RMSF graphs demonstrate how the average protein fluctuations relate to changes at the residue level. The RMSF of the target protein is depicted in Figure 3.The number of intramolecular bonds found in the ligand compound Vilasinin, which is more stable in conformation than the other ligand compound, was higher than the total number of intermolecular bonds created between the macromolecule and its ligand compounds (Figure8). For ligand 7-Deacetyl-7-oxogedunin, the SID diagram observed that the interaction fraction value had reached the highest value at ASP 1229, created by hydrogen bonds and water bridges showing that specific interaction has been maintained over the 100 ns simulation time. The IFV for compound having compound CID 122801 has shown everyday interaction by forming hydrogen bonds, hydrophobic bonds, and water bridges with the highest value of 0.6 with Lys 1030. Vilasinin has formed water bridges and hydrogen bond interactions with the highest value of 2 with ASP 1232. The stacked chart bar for Nimbidinin showed the highest IFV with ARG 1136 by forming hydrogen bonds and water bridges having the highest value of 1. The SID diagram also showed that the highest IFV value of the control drug Nateglinide is formed by water bridges and hydrogen bonds with ARG 1136. The Result shows that the molecular stability of the selected phytocompounds is similar to the control drug nateglinide with the insulin receptor ectodomain.The protein-ligandssolvent-accessible surface area (SASA) was also estimated using simulated trajectories to determine the dimensional changes of the drug-like molecules along the simulation trajectory [15, 83].One of the most significant SASA values is related to the unstable structure, which has hydrophobic amino acid residues close to the water molecule [84, 85]. All chosen compounds have conventional SASA values, as shown by the MDS trajectory’s SASA result (Figure 5). The MolSA and PSA validation graph shows that all ligand compounds have potential advantages over the control medication (Figures 6 and 7).In addition, Rg measures the separation between the end of the protein and the centre mass of the protein. As a result, this metric provides additional information about the protein’s folding properties while measuring how compact the protein molecule is [86]. A significant Rg value indicates slack packing, whereas a lower Rg value indicates compact packing [83]. The summary of the Rg values in Figure 4 reveals that the protein-containing compounds all displayed average compactness.


## 5. Conclusion

Four phytocompounds were selected from a total of 63. The following compounds are projected to have the best drug-like properties based on the virtual-mediated ADMET screening employed in this procedure. These ligand molecules, which are similar to the reference medication nateglinide, have the highest affinities for the insulin receptor, according to molecular docking studies. Notably, these drug-like compounds had longer interactions with the target protein’s active site residues. Interestingly, the highest positive values were observed for the compounds in the simulation study that evaluated the quality of these molecules’ interactions with the target receptor. Due to their ability to block the appropriate insulin receptor, these drugs may delay the development of diabetes mellitus. These top four phytocompounds could be used to treat diabetes, but more (in vivo and in vitro) study is needed to confirm their safety attributes.As the next generation of pharmaceuticals develops, the management of diabetes mellitus will advance. Even though these are only a few examples of the variety of related development initiatives and cutting-edge tactics being utilized to treat type 2 diabetes, the following research has discussed several distinct therapeutic approaches for treating DM. Future generations of antidiabetic drugs will need to find solutions to a variety of real-world issues, including overlapping safety concerns, side effects, drug-drug interactions, dose and futility restrictions, and other issues that currently limit their use in certain situations. Thus, all significant problems would be solved, and the chosen phytocompounds would function as efficient complementary treatments for type 2 Diabetes Mellitus.



---

*Source: 4170703-2023-03-03.xml*