Development of antidiabetic drugs from benzamide derivatives as glucokinase activator: A computational approach

Hyperglycemia is a condition known for the impairment of insulin secretion and is responsible for diabetes mellitus. Various small molecule inhibitors have been discovered as glucokinase activators. Recent studies on benzamide derivatives showed their importance in the treatment of diabetes as glucokinase activator. The present manuscript showed a computation study on benzamide derivatives to help in the production of potent glucokinase activators. In the present study, pharmacophore development, 3D-QSAR, and docking studies were performed on benzamide derivatives to find out the important features required for the development of a potential glucokinase activator. The generated pharmacophore hypothesis ADRR_1 consisted of essential features required for the activity. The resultant statistical data showed high significant values with R2 > 0.99; 0.98 for the training set and Q2 > 0.52; 0.71 for test set based on atom-based and field-based models, respectively. The potent compound 15b of the series showed a good docking score via binding with different amino acid residues such as (NH…ARG63), (SO2…ARG250, THR65), and π-π staking with (phenyl……TYR214). The virtual screening study used 3563 compounds from ZINC database and screened hit compound ZINC08974524, binds with similar amino acids as shown by compound 15b and crystal ligand with docking scores SP (-11.17 kcal/mol) and XP (-8.43 kcal/mol). Compounds were further evaluated by ADME and MMGBSA parameters. Ligands and ZINC hits showed no violation of Lipinski rules. All the screened compounds showed good synthetic accessibility. The present study may be used by researchers for the development of novel benzamide derivatives as glucokinase activator.


Introduction
Glucokinase is a hexokinase isozyme, consist of 465 amino acids (molecular weight = 50 kD) present in pancreatic b-cells and liver (postprandial state). Glucokinase catalyzes a reaction that involve the transfer of phosphate from ATP to glucose and the generation of glucose 6-phosphate which is the first step in the direction of synthesis of glycogen and glycolysis (Matschinsky 1996;Antoine and Boutin 2009). This reaction is also representing the first ratelimiting step in glucose metabolism. Glucokinase activator (GK-A) worked through two different mechanisms known as lowering the blood glucose level in the liver and increasing insulin secretion in pancreatic b-cells. Therefore, it becomes an interesting target in the present scenario to treat diabetes. Various GK-As have been synthesized, some are under clinical studies and showed promising results to lower blood glucose levels in healthy people and type-2 diabetes mellitus (T2-DM) patients. Glucokinase activator is responsible for several side effects such as hypoglycemia and testicular toxicity (American Diabetes Association 2009; Waring et al. 2011). To eliminate these side effects, frequent dose regimens and dose titration are preferred. There are two different approaches described to treat hypoglycemic state, one is related to the designing of partial activators (Pfefferkorn et al. 2011) and the second approach is to restrict liver-selective glucose activators (Bebernitz et al. 2009;Massa et al. 2011;Pfefferkorn et al. 2012;Park 2012;Park et al. 2013).
In the direction to treat diabetes, various benzamide derivatives have been developed as GK-As (Matschinsky et al. 2010;Mao et al., 2012). Park, et al. identified a novel phenylethyl benzamide GK-A (Park 2012;Park et al. 2013). Furthermore, Park et al. synthesized a series of pyrazole benzamide derivatives as GK-A having 3methylpyridine and 4-phenoxymethyl sulfone groups for the treat- Table 1 In vitro IC 50 and PIC 50 values for GK activity of pyrazole benzamide derivatives (Park 2012;Park et al. 2013  ,5-disubstituted analogues and evaluated them for GK activation activity. These analogues showed considerable antihyperglycemic activity in the animal models (Grewal et al. 2019). In a study, Charaya et al. synthesized thiazole-2-yl benzamide derivatives from benzoic acid and evaluated them for GK activation activity (Charaya et al. 2018). A benzamide derivative PF-04937319 (1) is under phase-1 clinical trial for the treatment of diabetes.
On the basis of previously synthesized compounds (Park 2012;Park et al. 2013), structure-based drug design was performed with pharmacophore development, 3D QSAR, and docking simulations for the determination of potent GK-A by using the PHASE module (Schrodinger). Pharmacophore development determines the important features required for the activity and can be used for 3D QSAR and virtual screening studies. Field and atom-based 3D-QSAR models were developed for the determination of statistical data results of correlation between molecules and their properties. A virtual screening study on the ZINC database generated the potent compounds as GK-A. Docking study revealed the important interactions with amino acids required for activity. On the other hand, the MMGBSA (Molecular Mechanics Generalized Born Surface Area) method was used to predict the binding free energy of the docked molecules. The current findings of the study may be utilized as a guiding tool for the development of novel and effective GK-A.

Software
The 3D-QSAR, pharmacophore modeling and docking were performed by Schrodinger module (Park et al. 2014;Park et al. 2015).

Dataset
A dataset of 43 benzamide derivatives was taken for the development of pharmacophore, 3D-QSAR, virtual screening, and docking studies ( Table 1). The IC 50 values taken from biological activities were converted into pIC 50. Data set were divided into active for 'higher activity' and inactive for 'lower activity' compounds and the remaining were counted in the category of intermediate. The alignment on common scaffold of 43 benzamide derivatives is given in Fig. 1A. The PDB ID-3A0I was used for the docking studies. The crystal ligand of this protein used for comparison of docking outcomes are presented in Fig. 1B (Andreoli et al. 2014).

Preparation of ligands
The molecular structure of benzamide derivatives was converted from 2D to 3D by using LigPrep of Schrodinger software. The OPLS_2005 force field was taken for the preparation of ligands. Further, the prepared molecules were taken for 3D-QSAR and docking simulations (Schrodinger 2017; Ali and Ali 2021).

Pharmacophore development
The pharmacophore model was developed by using PHASE module, where all 43 ligands were aligned on common scaffold and generated conformers using the macromodel search method (Ali and Ali 2021). PHASE module consists six different pharmacophore features including hydrogen bond acceptor (HBA), hydrogen bond donor (HBD), aromatic ring, hydrophobic group, positive ionizable, and negative ionizable groups (Rajeswari et al. 2014). The maximum number of sites was fixed to 5 which further A. Ali Saudi Journal of Biological Sciences 29 (2022) 3313-3325 responsible for the generation of the topmost 20 different hypotheses. For the generation of these hypotheses, ligands were categorized into active, inactive, and intermediate.
The active ligand was set as for those which had pIC 50 value of more than 8 whereas inactive ligands considered for those had pIC 50 value less than 6.6, except these all-other ligands were set as intermediates. The final data consisted of 9 actives and 14 inactive ligands. The hypothesis was generated by using 9 actives and 1 Å box size and 2 Å site distance. The all-generated hypotheses were ranked based on different scores depicted in Table 2 (Crisan et al. 2019;Sakkiah et al. 2014). The hypothesis ADRR_1 (Fig. 2) showed the top scoring features with one HBA, one HBD, and two aromatic rings. The hypothesis determines the essential features required for the binding with receptor for a particular activity. A total of 20 pharmacophore hypotheses were developed which were ranks according to the score hypothesis (Table 2).

3D-QSAR
3D-QSAR models were developed by two techniques known as field-based and atom-based QSAR. The developed models were showed the essential parameters required for activity by correlating the structure features with biological activity. All 43 benzamide derivatives were separated into a training set with 75% and test set with 25% compounds using 4 PLS factors (Kar et al. 2010). Best models were generated by QSAR models (atom and field-based) described in Table 3. Statistics for atom-type fraction and field type fraction were summarized in Tables 4 and 5. In the atom-based model, 35 molecules were selected for the training, and 8 molecules for the test set, whereas in field-based 31 and 12 molecules were taken for the training and the test set, respectively (Table 6). Table 2 The generated pharmacophore hypotheses.

N.
Hypothesis Phasehypo-S Survival-S Inactive-S Site-S Vector-S Volume-S Selectivity-S   The contour maps generated by atom-based and field-based models are presented in Figs. 3 and 4, respectively. These maps are represented as hydrophobic, steric, donor, acceptor, and electrostatic fields (Peng et al. 2017).

Target protein prediction
The potent compounds of the series were selected and submitted to the Swiss-Target-Prediction tool. The said tool predicted var-  Table 6 The IC 50 value (Actual vs predicted) generated by atom-based and field-based 3D-QSAR model using PLS factor 4. ious protein targets, in which GK-A was the most appropriate target for the selected molecules.

Docking analysis
The PDB ID-3A0I consists a three-dimensional structure of GK-A and was downloaded through protein data bank (Tagami et al. 2010). The docking study was performed to determine the binding scores between ligand and receptor by using the Glide module. The SP and XP methodologies were used for all 43 ligands with PDB ID-3A0I (Table 7). The MMGBSA based rescoring technique was used for the prediction of binding free energy calculation between ligands and receptor molecule (Table 8). However, DFT studies can give more accurate assessment of binding/docking (Van Mourik et al., 2014).

ADME prediction studies
The top scored compounds were analyzed by different ADME properties such as drug-likeness, solubility, and pharmacokinetic   studies. The QikProp software was used to calculate ADME properties (Table 9). Further, the Swiss-ADME tool was used to evaluate additional parameters like cytochrome profile of the drug with permeation through different other barriers.

Analysis of pharmacophore modeling
The generated pharmacophore hypotheses with different scores are presented in Table 2. The hypothesis ADRR_1 was chosen as the best hypothesis with phase hypo-score = 1.18, survival score = 5.03, and site score = 0.98. The field-based and atom-based QSAR studies showed reliable statistical parameters with different evaluation factors. The results showed internal validation parameters such as R 2 values 0.99, 0.98; R 2 CV values 0.62, 0.64; Q 2 values 0.52, 0.71; SD values 0.62, 0.19; RMSE values 0.96, 0.66 and F values 299, 271 for atom-based and field-based models, respectively ( Table 3). The scatter plots for both models are shown in Fig. 5.

3D-QSAR
The contour maps showed the correlation between different bioactivities by various substituents on the core moiety (Fig. 3). The correlation between actual and predicted activities of test and training-set compounds for atom-based QSAR is depicted in Fig. 5. The contour maps generated by field-based model include HBA, HBD, steric, electrostatic, and hydrophobic fields (Fig. 4). Different substituent groups on the potent compound (16a) described by various colour responsible for increase or decrease in activity. The correlation between actual and predicted activities of test and training-set compounds for field-based QSAR is depicted in Fig. 5.

Docking and virtual screening studies
The PDB ID-3A0I was taken for docking purposes to evaluate binding interactions of potent ligands and ZINC compounds. The docking scores were compared with the observed activity. Compounds,18g,13b1,19e,19a,16b,16a, and 15b showed binding  interactions with important amino acids required for GK-A. The amino acids bind with both compounds 15b and 18g were ARG250 and THR65, ARG63 as shown in Fig. 6. Compound 16a also showed good binding interactions described in Fig. 7. The docking score of the potent compounds of the series compared with standard drugs such as glucokinase activator 1, piragliatin, ARRY-403, and RO-5305552 are described in Table 8. The binding interactions of other potent compounds 19e and 19 are displayed in Fig. 8.
3.5. ADME properties calculation ADME properties were calculated by the QikProp module, comprise of one another SwissADME tool. These parameters were within the acceptable range for ligands and ZINC hits ( Table 9). The compounds of the series showed drug-likeness properties with no violation of the Lipinski rule. The bioavailability score of compounds showed the value of 0.55.

Discussion
The final hypothesis consisted of one HBD, one HBA, and two aromatic ring structures. The hypothesis (ADRR_1) showed alignment with other molecules of the series and displayed good correlation between structure and bioactivity. The features of the hypothesis were further taken for screening of ZINC compounds from the ZINCPHARMER (http://zincpharmer.csb.pitt.edu/) online tool. The atom based QSAR maps indicated the influence on bioactivity by the addition of substituents on the nucleus. The blue contour maps showed an increase in activity, whereas red maps showed decrease in activity. The compound 16a showed alignment with pharmacophore hypothesis ADRR_1 with different colours on their substituents. The electron-withdrawing group (EWG) substi- Fig. 5. Correlation between test (5 A & C) and training (5B & D) set compounds by using atom and field-based 3D QSAR models.
A. Ali Saudi Journal of Biological Sciences 29 (2022) 3313-3325 tution on benzamide derivatives exhibited increase in the activity, whereas EWG substitution on phenyl ring showed decrease in activity, as represented by red maps (Fig. 3). The HBD group addition showed no changes in activity. The hydrophobic group addition connected to the amide group displayed increase in activity, whereas the introduction of such group at phenyl ring decreased the potency and showed mixed activity throughout the ring. However, the addition of positive and negative ionic groups showed decrease in the activity.
In field based QSAR contour map electrostatic group contains the blue colour at the amide group linked heterocyclic compounds that showed the introduction of electron positive group at the site responsible for increase in the activity. The phenyl ring connected to heterocyclic ring with electron positive group may decrease or increase the activity. The HBA group introduction at phenyl ring connected group may increase the activity. The addition of HBD group at the amide group may increase the activity. The group contains steric field with green color in phenyl ring substituted group may be responsible for increase in activity. The potent compound 15b of the series showed good docking score using interactions with amino acid residues (NH. . .ARG63), (SO 2 . . .ARG250, THR65), and p-p staking with (phenyl. . .. . .TYR214). Compound 16a also showed good binding interactions with amino acid residues such as TRY214, ARG250, THR65, and ARG63 (docking score À11.296 kcal/mol), important for GK-A activity (Fig. 7). The docking scores and the amino acid residues of the potent compounds of the series compared with standard drugs such as glucokinase activator 1 (À3.36 kcal/mol), piragliatin (À8.345 kcal/mol), ARRY-403 (À10.259 kcal/mol), and RO-5305552 (À6.707 kcal/mol) (Table 8). Furthermore, the other potent compounds 19e and 19  A. Ali Saudi Journal of Biological Sciences 29 (2022) 3313-3325 from both the series showed good binding interactions with THR65, ARG63, ARG250, and TYR214 essential amino acids that are important for a GK-A drug (Fig. 8).
The ZINC library was downloaded through the Swiss screening database. Total 3563 compounds were downloaded by ZINC database and further screened by different docking methodologies using Glide module. After applying Lipinski rule, the compounds were filtered through HTVS docking process. The top 50% of the compounds from this process were further taken for SP, and the top 20% compounds were finally taken for XP. Top hit, namely ZINC08974524 ( Fig. 9) showed best docking score in SP (À11.17 kcal/mol), XP (À8.43 kcal/mol), and HTVS (À10.26 kcal/mol). The binding interactions of all the ZINC compounds were similar to crystal ligand interactions. The RMSD value was used as a parameter to check the binding pattern of different com-pounds from crystal ligand. These binding interactions of active compounds and ZINC derivatives displayed similar interaction as shown by crystal ligand of PDB ID-3A0I. The Zinc hit compounds may be used for the in vivo evaluation as GK-A. Furthermore, ZINC screening data suggested that the binding interactions were similar to the compounds taking for 3D QSAR study. So, these compounds may be designed further for the synthesis of potent compounds as GK-A against diabetes.
In ADME analysis, compounds showed solubility range from low to high values. The GI absorption was low for 16a and high for 19a compounds. The SwissTargetPrediction results exhibited that compounds have high target specificity for GK-A. All compounds showed good synthetic accessibility within the range of 2-4. These compounds may be useful for the generation of novel compounds as GK-A. The compounds screened through ZINC data-  base showed good ADME properties (Table 9). However, cytochrome profiling for ZINC hits shown inhibitory activities against CYP3A4, CYP2C19, and CYP2C9. All the compounds of the series and ZINC hits showed noncarcinogenic activities.

SAR optimization
The pharmacophore model, 3D QSAR, virtual screening, and Zinc hit compounds may be used as a basis for the production of novel compounds as GK-A depicted in Fig. 10. The present SAR optimized by the 3D QSAR study revealed that the substitutions on benzamide scaffold with different characteristic features for the development of novel GK-A. The core moiety in the place of thiazole ring can be replaced by some electropositive atoms and hydrophobic groups such as long carbon chain and phenyl ring responsible for increase in activity. The oxygen atom of the benzamide group is more prominent for activity, or if it is replaced by less electronegative groups such as sulphur or nitrogen, then the activity may be diminished. The pyridine ring connected to the benzamide scaffold can be replaced by electropositive groups that results an increase in activity. These features can be used for the further development and synthesis of novel derivatives as GK-A.

Conclusion
The present manuscript revealed a computation study on benzamide derivatives processed in a sequence to produce potent GK-As, where ligands and ZINC hits showed no violation of the Lipinski rules. All screened compounds displayed good synthetic accessibility. Based on 3D-QSAR, pharmacophore development, SwissTargetPrediction, virtual screening, and molecular docking studies, we may design novel molecules with good ADME properties and low toxicity. ADRR_1 is determined as the best pharmacophore in the study. The 3D QSAR study showed the best statistical data by atom-based and field-based models, consecutively. The binding interactions of compounds showed important amino acids required for activity. Fig. 10 showed the importance of different substituents on core moiety which may help for the development of novel compounds as GK-A. Furthermore, the present study may be helpful for researchers as a guiding tool for the development of novel benzamide derivatives as GK-A.

Declaration of Competing Interest
The author declare that there is no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.