Computational Assessment of Xanthones from African Medicinal Plants as Aldose Reductase Inhibitors

: Diabetes mellitus is a life-threatening non-communicable disease that affects all age groups. Despite the increased attention it has received in recent years, the number of diabetic patients has grown exponentially. These increased cases are attributed to essential enzymes involved in blood glucose regulation. In this study, we attempt to reveal the aldose reductase inhibitory potential of xanthones isolated from African medicinal plants. Ensemble docking, molecular dynamics simulation, density functional theory (DFT), and ADMET methods were employed to identify drug candidates as aldose reductase inhibitors. The ensemble docking results identiﬁed mangostenone B, bangangxanthone A, smeathxanthone B, mangostenone A, and allanxanthone B as potent inhibitors against the aldose reductase enzyme. Molecular dynamics studies showed the xanthones established better binding mode and afﬁnities against the enzyme. Moreover, the electronic properties of the xanthones explained their good pharmacological potentials. Therefore, our ﬁndings suggest that the hit molecules be investigated in vitro and in vivo for drug development against aldose reductase.


Introduction
Aldose reductase (AR) is an enzyme implicated in Type 2 diabetes mellitus (type 2 DM) and plays a vital role in the polyol pathway that catalyzes glucose reduction into sorbitol through the oxidation of NADPH to NADP + [1]. Then, sorbitol dehydrogenase converts the sorbitol to fructose by using NAD as a cofactor. Under normal glucose homeostasis, the pathway represents a minor route of glucose metabolism that operates in parallel with glycolysis, in which the affinity of AR and glucose is low. In the first step of the catalytic mechanism of AR, the anionic pocket of the enzyme, consisting of Asp-43, Tyr-48, Lys-77, His-110, Ser-159, Ans-160, Gln-183, and Tyr-209, forms a network of hydrogen bonds that transfer 4-pro-R hydride from the nicotinamide ring of NADPH to the substrate's carbonyl carbon. Tyr-48 on the active site donates a proton to reduce the carbonyl oxygen atom of the substrate to alcohol during the second step of reduction, while His-110 maintains the substrate's proper orientation in the active site [2]. During hyperglycemic circumstances, highly expressed AR results in two-fold to four-fold accelerated the conversion of glucose to sorbitol. Accumulation of sorbitol in the tissues leads to osmotic swelling, changes in membrane permeability, and oxidative stress [3]. These cell-damaging processes play a significant role in the etiology of diabetic complications such as cataract formation, retinopathy, neuropathy, and nephropathy [4,5].
The over-expression of AR in most diabetic complications leads to a decrease in the histone deacetylase 3 (HDAC3) gene and increased peroxisome proliferator-activated receptor gamma (PPARγ) signaling, also resulting in lipid accumulation in the heart [6]. In addition, the over-expression of AR leads to the repression of other enzymes and biomarkers implicated in diabetes. Inhibitors of AR interact with the elucidated AR amino acid residues and competitively block the active site, preventing substrate access and product formation [2]. Therefore, inhibitors of AR have emerged to be therapeutically significant since they may have the potential to treat long-term diabetic complications. Several efforts and substantial resources have been expended in the development of synthetic aldose reductase (SAR) inhibitors [7]. The developed aldose reductase inhibitors (ARIs) vary structurally, and representative structural classes of ARIs include carboxylic acid derivatives (such as Epalrestat, Alrestatin, Zopalrestat, Zenarestat, Ponalrestat, Lidorestat, and Tolrestat), spirohydantoins and related cyclic amides (such as Sorbinil, Minalrestat, and Fidarestat), and phenolic derivatives (related to Benzopyran-4-one and Chalcone) [8,9]. So far, epalrestat, marketed in Japan and China is the only commercially available SAR. All other inhibitors have failed in clinical trials because of poor pharmacokinetic properties and side effects and even epalrestat has been withdrawn from the market in other countries because of its side effects [4,10]. Moreover, the comparative structural analysis and molecular dynamic (MD) studies of the AR enzyme showed that a single experimental structure is not sufficient to predict all possible binding modes [11]. X-ray studies revealed that there are at least three distinct binding pocket conformations of AR. With this pronounced binding site adaptations, the enzyme can operate on a broad palette of structurally diverse substrates that possess both hydrophobic and hydrophilic groups. These deficiencies result in lower virtual screening hit rates [12]. As a result of the side effects and different structural features needed for the development of potent AR inhibitors, there is a need to explore natural sources to identify safer therapeutic compounds. Plant-derived AR inhibitors are considered suitable alternatives due to their lesser side effects and low cost [13,14].
In Africa, medicinal plants have been extensively used in the treatment of type 2 DM, especially in places where affordability and access to modern healthcare are limited [15]. They are found to be a reservoir of nutrients with remarkable inhibitory and antioxidant activities. In addition, consumers are now mindful of their diets, which shifts recent research trends to natural products as functional food, providing supplementation benefits in diabetes management [16,17]. However, due to cost and limited resources, the isolation of individual molecules required for drug discovery from crude formulations of plants and other natural sources has not been extensively explored in Africa [18]. Xanthones are essential classes of polyphenolic compounds occurring in medicinal plants. Xanthones have been shown in studies to have a variety of pharmacological activities, including antidiabetic, hypolipidemic, anticancer, antioxidant, and antimalarial properties [19][20][21][22][23][24][25]. In addition, they are known to inhibit AR and useful in other diseases caused by accumulation of free radicals in the body [24].
Computational chemistry methods such as molecular docking studies, molecular dynamic simulation, and density functional theory are used to explore the binding mechanisms, electronic properties, and pharmacological potentials of phytoconstituents against targeted diseases [26]. These methods have gained widespread relevance in the drug discovery process because they reduce costs and effort [27][28][29][30].
Using the mentioned in silico approaches, this study provides explicit insights into the top-ranked AR inhibitors' binding affinities, binding modes, and quantum chemical properties, thus identifying new AR inhibitors. It further examines the drug-likeness and ADMET properties of the hit molecules.

Validation of Docking Protocol
The 3D crystal structure of AR (PDB ID: 1EL3) obtained from the protein data bank (http://www.rcsb.org/pdb/ 1 December 2020) was used in the docking study. The AR structure was prepared for docking by removing water molecules, native ligand, co-factor, and ions and then charges and polar hydrogen were added. The structures of 65 xanthones isolated from African medicinal plants were downloaded from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/ 11 December 2020) and energy minimization was performed as described in materials and methods section. In order to validate the docking protocol the native ligand was re-docked at the catalytic site. The result obtained from the re-docking of native ligand of AR enzyme with AutoDockVina [31] showed that the RMSD value estimated between the re-docked and native ligands of AR is 1.53 Å (Figure 1). Furthermore, the interactions established between the native ligand and amino acid residues correlate with those obtained for the re-docked ligand. The re-docked ligand of AR enzyme showed interaction with Trp20, Lys21, Asp43, Val47, Tyr48, Gln49, Val77, Trp79, His110, Trp111, Phe122, Ser159, Asn160, Gln183, Tyr209, Ser210, Pro218, Trp219, Ile260, Cys298, in which, twelve amino acids residues correctly match with those obtained for the native ligand at the active site of the AR enzyme.

Molecular Docking Interaction of Xanthones with Aldose Reductase Enzyme
The estimates of free energies of binding in kcal/mol for xanthones when docked into the active site of AR were used to select the top five hit molecules. The ensemble docking approach [32] was adopted. The AR bare protein was subjected to 100 ns MD simulation with the procedure as described in the methods section. With gmx cluster module the cluster of 21 conformations were obtained by setting a cutoff of 0.08 nm. The structure representing the centroid of each cluster was used in the docking studies. These structures of AR from centroid of each cluster were aligned initially. Each conformation of AR is unique in that it has a different orientation of residues at the binding site. Thus, when we docked the xanthones at the binding site within a grid box large enough to capture the entire binding site, the most favorable xanthone would elicit the better binding free energy for that specific conformation of AR. This molecular docking protocol is anticipated to offer the conformational sampling of the binding site of AR [33]. The results obtained from the ensemble docking studies of xanthones with AR enzyme are given in Table 1. Mangostenone B elicited a significant binding affinity (−11.2 kcal/mol) toward the AR enzyme compared to the binding free energy of co-crystallized ligand (−8.8 kcal/mol). Moreover, an oxygen atom on the benzene ring of the mangostenone B moiety formed one hydrogen bond interaction with Trp20 at 2.75 Å. The molecule established key hydrophobic interactions with Trp20, His110, Trp111, Phe122, Tyr209, Trp219, Cys298, and Leu301. Furthermore, Trp20, His110, Trp111, Phe122, Tyr209, and Trp219 participated in Pi-interaction (Pi-sigma and Pi-Alkyl), while no electrostatic interaction was formed between the amino acid residues at the active site of AR enzyme and the mangostenone B moiety ( Figure 2). Bangangxanthone A elicited binding energy of −10.8 kcal/mol at the enzyme's active site. Two oxygen atoms on the benzene rings of the bangangxanthone A moiety formed two hydrogen bond interactions with two amino acid residues (Cys298 at 3.07 Å and Leu300 at 2.16 Å). Furthermore, several hydrophobic interactions were formed between the hydrogen and oxygen atoms of bangangxanthone A and amino acid residues at the enzyme's active site (Trp20, His110, Phe122, Tyr209, Lys262, Cys298, and Leu300). Furthermore, Pi-interaction (Pi-Pi stacked and Pi-alkyl) was established with five amino acid residues (Trp20, His110, Phe122, Tyr209, and Leu300), while His110 participated in electrostatic interaction with the bangangxanthone A moiety.
Smeathxanthone B showed binding energy of −10.8 kcal/mol after interacting with the active site of the AR enzyme. Trp20, Tyr48, His110, Trp111, Phe122, Tyr209, Lys262 and Cys298 participated in the hydrophobic interactions with the smeathxanthone B moiety. The ligand further established Pi-interactions (Pi-Pi stacked and Pi-alkyl) with Trp20, Tyr48, His110, Trp111, Phe122 and Tyr209, while smeathxanthone A formed no electrostatic interaction with the amino acids at the active site of the AR enzyme.
Mangostenone A showed almost equivalent binding energy (−10.7 kcal/mol) as that of Smeathxanthone B. The oxygen atom of the hydroxyl group on the delocalized benzene ring of mangostenone A formed one hydrogen bond interaction with Trp20 at 2.86 Å. The ligand also established numerous hydrophobic interactions with Trp20, Lys77, His110, Trp111, Phe122, Tyr209, and Trp219. Furthermore, Pi interactions (Pi-alkyl, Pi-Pi T-shaped, and Pi-Pi stacked) were formed with Trp20, His110, Phe122, Tyr209, and Trp219, while one electrostatic interaction was established with His110.
The allanxanthone B elicited binding energy of −10.6 kcal/mol. The keto group of the pyrone ring of allanxanthone B formed one hydrogen bond interaction with Ala299 at 2.10 Å, while the oxygen atoms on its benzene established hydrogen bond interaction with Ala299 at 3.02 Å and Leu300 at 2.52 Å. Furthermore, the amino acid residues Trp20, Trp79, His110, Trp111, Phe122, Trp219, Cys298, Leu300, and Leu301 participated in hydrophobic interactions with the allanxanthone B. No electrostatic interaction was formed while few amino acid residues Trp20, Trp79, His110, Trp111, Phe122, Trp219, Leu300, and Leu301 participated in Pi-interactions between the enzyme and the ligand.
The five-hit molecules formed common hydrophobic interaction with Trp20, His110, Phe122, and Cys298. The hydrogen bond, hydrophobic, electrostatic, and Pi interactions established between the active site of the AR enzyme and the hit molecules are responsible for the high binding affinity and stability elicited by the molecules in the enzyme's binding pocket. Furthermore, AR enzyme is known to possess hydrophobic or specificity pocket, hydrophilic and anion binding pockets. Hence, drug candidates with hydrophilic and hydrophobic groups can exhibit an effective inhibitory property because the hydrophilic functional group can bind to the anion binding pocket thereby hindering proton transfer, while the hydrophobic group binds to the hydrophobic pocket of the enzyme [34]. The top-ranked molecules identified in this study have both hydrophilic and hydrophobic groups such that their hydrophilic group binds to His110 and Tyr48 of the anionic binding pocket, while the hydrophilic group established interaction with the amino acid residues (Trp111, Phe122, Trp219, Cys298, and Leu300) at the specific pocket. These characteristic interactions further establish their tendency to serve as potent AR inhibitors.
Xanthones are a class of polyphenolic compounds with the same skeletal structure and different substituents. The presence of hydroxyl and non-polar alkyl groups (prenylgroups) may contribute to their interaction with the anionic and hydrophobic pocket of the AR enzyme thereby eliciting good binding energy.
The results obtained thus suggest mangostenone B, mangostenoneA, bangangxanthone A, smeathxanthone B, and allanxanthone B are reported as major xanthones in Garcinia mangostana, Garcinia polyantha, Garcinia smeathmannii, and Allanblackia monticola as potent AR inhibitors. In addition, the binding energies elicited by the hit molecules are higher when the results obtained for xanthotoxin, β-sitosterol, gymnemic acid are compared [35,36].

ADMET and Drug-Likeness Analysis
The drug-likeness potential of molecules as possible therapeutic arsenal for the inhibition of AR depends on their absorption, toxicity, and metabolic reaction elicited in diabetic patients. The ADMET predictions were obtained from admetSar-2 online tools (http://lmmd.ecust.edu.cn/admetsar2/ 18 January 2021). Parameters such as human intestinal absorption (HIA) and solubility account for oral bioavailability and dissolution of the ligand for effective metabolic activities. Blood-brain barrier permeability plays a significant role in the passage of the hit molecules to the brain while toxicity assures their safety when developed into drugs to inhibit AR [37,38]. In this study, the five hit molecules selected as inhibitors for AR enzyme elicited ADMET results that potentiate them as good drug candidates (Table 2). The results obtained for the solubility of each ligand are within the acceptable range of −6.5 to 0.5 [39]. Moreover, all the ligands exhibited blood-brain barrier permeability and human intestinal absorption properties. In terms of toxicity, the ADMET studies showed that none of the ligand is carcinogenic while the acute oral toxicity values are within the acceptable range meaning that the ligands are safe and can be developed into potent drugs to inhibit the actions of AR enzyme in type 2 diabetes patients.
Additionally, epalrestat (standard drug) showed similar ADMET characteristics with those obtained for the top-ranked ligands except its inability to cross the blood-brain barrier [40]. Moreover, it obeyed the drug likeness properties. Hence, these molecules can serve as potential drug leads in the inhibition of AR.

Molecular Dynamics Simulation Studies and MM-PBSA Calculations
The ensemble docking studies suggested that the xanthones namely mangostenone B, mangostenone A, mangangxanthone A, smeathxanthone B, and allanxanthone B could potentially inhibit the activity of AR. The scoring function of AutoDockVina is reasonably good in accurately predicting the binding modes and free energy of binding estimates. The solvated environment, the presence of ions, and the protonation states of the residues, however, are not taken into account by the AutoDockVina algorithm. In this case, molecular dynamics provide more insightful evidence for binding site interactions as well as more accurate estimates of binding free energies and binding affinities [41,42]. This is certainly more advantageous than classical docking.
The results of MD simulations of AR and corresponding potential xanthone inhibitors are discussed first. The measurement of protein and ligand RMSD provides good estimates of conformational stability of the system [43]. The deviations from the starting positions of protein or ligand atoms are measured in terms of RMSD. In turn, RMSD measurement provides the protein folding events, the influence of loop flexibility, and binding site adaptation. The RMSD analysis results for the protein showed that the bare protein stabilizes quickly after around 25 ns and thereafter almost remains stable (Supplementary Materials, Figure S1). The complex of mangostenone A gradually raises until around 35 ns and after around 55 ns remains almost steady until the end of simulation with an average RMSD of 0.155 nm. The complex of mangostenone B observed stabilized after 25 ns with average RMSD of 0.15 nm which is well within the acceptable RMSD limit. The complexes with bangangxanthone A observed stabilized throughout the simulation with an average RMSD of 0.15 nm. The complexes of smeathxanthone B and allanxanthone B were also stabilized with average RMSD of 0.15 nm after initial fluctuations which were observed until 25 ns.
The RMSD in the xanthones atoms were also analyzed. Smeathxanthone B atoms have larger deviations compared to other xanthones (Supplementary Materials, Figure S2). The initial 20 ns simulation produced larger deviations in RMSD in the atoms of all the xanthones. Allanxanthone B and mangostenone B atoms also showed slightly higher RMSD with an average of 0.152 nm. The RMSD in mangostenone A atoms was least compared to other xanthones and observed stabilizing at 0.055 nm throughout the simulation except the initial deviations till around 20 ns. The RMSD in the bangangxanthone A remained stable until 60 ns after initial deviations observed till 20 ns. However, thereafter it rises to around 0.125 nm until the end of the simulation. Probably, the initial deviations in RMSD in atoms of all the xanthones under study were due to conformational changes to adapt in the active site of AR. These RMSD values suggest that mangostenone A and bangangxanthone A are more favorable at the active site.
The analysis of radius of gyration (Rg) provides insights into the system's overall compactness when the ligands bind the protein under study [44]. Evidently, all the xanthones stabilize the compact structure of AR almost intact throughout the MDS within the narrow acceptable range of 1.925 to 1.965 nm (Figure 3). This suggests the intact nature of α-helices and other secondary structural features in the AR despite the binding of xanthones. One of the major phenomena of strong binding affinity of ligands is non-bonded interactions such as hydrogen bonds, hydrophobic interactions, and ionic interactions at the binding site [45]. The number of hydrogen bonds formed and the lifetime of such hydrogen bonds formed is in part indicative of the binding affinity. In the case of mangostenone A, one hydrogen bond formed frequently till 75 ns ( Figure 4). The equilibrated system of mangostenone A showed hydrogen bonds with His110 and Trp20. However, both these hydrogen bonds decompose and a new stable hydrogen bond was formed with Ser302. Further, at around 75 ns the hydrogen bond with His110 reformed. No hydrogen bonds were formed after around 75 ns suggesting poor affinity of mangostenone A during the end of simulation. The hydrogen bond analysis was performed setting cutoff of 3.2 Å for bond length range which indicates strong hydrogen bond and cutoff of 120 • for bond angle range. The central carbonyl oxygen atom of xanthone in mangostenone A forms a hydrogen bond with Ser302. However, the most important residue in AR inhibitory mechanism is His110 as it participates in proton transfer during reduction of aldehyde. The oxygen atom of dimethyl chromene ring (the one which has prenyl group not attached) forms a key hydrogen bond with this His110 residue (refer structures of xanthones given in Table 1).
Mangostenone B also formed frequent one hydrogen bond and maximum 3 hydrogen bonds ( Figure 5). The equilibrated system of mangostenone B showed hydrogen bond with Ser302. The initial conformational change in mangostenone B resulted in new hydrogen bonds with residues Trp20 and Val47 at around 24 ns simulation period. All new hydrogen bonds formed with Gly128, Tyr48, and Cys298 at simulation time 52 ns, 69 ns, and 98 ns respectively. However, many frames showed no hydrogen bond formation suggesting less propensity of hydrogen bond formation with mangostenone B. The structure of mangostenone B is very similar to mangostenone A except the attachment of dimethyl chromene ring without prenyl group. However, it results in different hydrogen bond network and it is observed that the key hydrogen bonds are formed with important residues Cys298 and Thr48 which are involved in catalytic mechanism of AR. Further, as the orientation of dimethyl chromene ring is changed in mangostenone B it could not form a hydrogen bond with His110. This result suggests that mangostenone A could offer better AR inhibitory activity than mangostenone B. The visual inspection of frames of mangostenone A and mangostenone B complexes showed the ligands occupied in the binding site. The other non-bonded interactions analyzed in terms of contact analysis are discussed in the following section. Bangangxanthone A showed more frequent hydrogen bond formation and on an average two hydrogen bonds formed throughout the simulation ( Figure 6). In the initial equilibrated trajectory, the hydrogen bonds with binding site residues His110, Trp20, Leu300, and Ala299 were observed. Out of these, the hydrogen bonds with residue His110 were found stable and remained intact till the end of simulation. The hydrogen bond with Ala299 remained intact till the 50 ns simulation; while with Trp20 remained intact till 25 ns. After around 75 ns simulation, a new stable hydrogen bond was found formed with residue Trp111. The structure of bangangxanthone A has more hydroxyl groups and contains a prenylmethyl, methyl chromene ring (not dimethyl chromene as in mangostenone A and mangostenone B). The residues His110 and Trp111 were involved in hydrogen bond formation with the either hydroxyl group on xanthone ring or oxygen atom of chromene ring. These hydrogen bonds are quite stable as can be seen in the hydrogen bond plot. The result suggests that bangangxanthone A could offer more AR inhibitory activity compared to mangostenone A and mangostenone B. In the case of smeathxanthone B one hydrogen bond with residue Trp20 was observed being formed throughout the simulation (Figure 7). Few frames of smeathxanthone B extracted at around 25 ns and 50 ns showed hydrogen bonds with Ala299 and Ser214 respectively. The structure of smeathxanthone B is almost similar to bangangxanthone A except the attachment of chromene ring. Smeathxanthone structure is linear compared to bangangxanthone A. However, in the smeathxanthone B hydrogen bond with His110 or Trp111 were not observed. The hydrogen bond with Trp20 which is also important in catalytic activity of AR was consistently observed throughout the simulation. Possibly, the linear structure of smeathxanthone B offers better AR inhibitory potential than above xanthones. The results of hydrogen bond analysis for allanxanthone B showed one hydrogen bond formed consistently throughout the simulation, while maximum three hydrogen bonds were observed only in few frames ( Figure 8). The visual investigation of frames extracted at different time intervals showed that the residue Ala299 forms a hydrogen bond consistently in all extracted frames. Other than this hydrogen bond, the initial equilibrated trajectory showed hydrogen bonds with Val297 and Leu300. The hydrogen bond with Val297 remained intact till 50 ns, while the hydrogen bond with Leu300 seems to be less stable as it was observed only in initial trajectory. The presence of hydrophobic diprenyl group and dimethyl chromene ring seems to restrict the conformational change in the allanxanthone B structure at the catalytic site. It could form a consistent hydrogen bond with Ala299 and may offer slightly less AR inhibitory potential than smeathxanthone B or bangangxanthone A.
The overall results of hydrogen bond analysis suggests that the residues His110, Trp20, and Ala299 are important in forming a stable hydrogen bond with bangangxanthone A, smeathxanthone B, and allanxanthone B respectively. Consequently, compared to mangostenone A and B, these xanthones may have more favorable binding affinity. Contact analysis between the xanthones under study and catalytic site of AR was performed so as to gauge the distance between the xanthones and the residues of catalytic site during the simulation. The residues Pro218, Trp219, Ala220, Lys221, Val297, and Leu301 were found occupied beyond the average distance of 0.75 nm, while other residues remain within a close average distance of around 0.5 nm in the case of mangostenone B (Supplementary Materials, Figure S3). All the residues as depicted in Figure S3 were found within close contact with average distance of around 0.5 nm in the case of mangostenone A. Around 26 residues were found in close contact with an average distance of around 0.5 nm for the bangangxanthone A and smeathxanthone B, except one residue Val264 in the case of bangangxanthone A. Similarly, 15 residues as depicted in Figure S3 were found within close contact with an average distance of around 0.5 nm for the allanxanthone B. Overall, the results of contact analysis suggests stable contacts of all the xanthones in the catalytic site of AR.
Root mean square fluctuations in atoms of AR residues were also measured (Supplementary Materials, Figure S4). Compared to bare protein, all the xanthone-bound complexes showed fluctuations in the active site residues Trp20, Lys21, Asp43, Val47, Tyr48, Gln49, Val77, Trp79, His110, Trp111, Phe122, Ser159, Asn160, Gln183, Tyr209, Ser210, Pro218, Trp219, Ile260, Cys298. Mangostenone A showed prominent fluctuations in atoms of residues Trp20, Lys21, and Trp79 compared to mangostenone B. Bangangxanthone A also showed fluctuations in these residues compared to other xanthones. These fluctuations are suggestive of binding site adaptability. Principal component analysis (PCA) performed on the MD frames as a part of an essential dynamics analysis helps in understanding the functionally relevant motions in protein and its lowest energy states [46]. In the present work, backbone atoms of AR were considered for PCA. Two principal components (PCs) were extracted from the covariance matrix. The projections of eigenvectors corresponding to each PCs provide the insights of path of motion of protein backbone atoms along with the energy states of AR. The 2D projections of the eigenvectors are shown in (Supplementary Materials, Figure S5). It is clearly evident that the projected motion in bare AR is confined to a narrow range on eigenvector 1 and 2. The distinct projected motion is seen in the case of bangangxanthone A which is split in the range 0 to −2 and 0 to 3 on eigenvector 1. Similarly, in the case of mangostenone A, there are two distinct motions observed which were projected in the range −2 to −4 and −1 to 2 on eigenvector 1. The collective motion in the case of systems with mangostenone B and allanxanthone B are unique, while those in the case of smeathxanthone B are spread over a larger range on both eigenvectors. It can be inferred from this result that the fluctuations in the catalytic site residues in the case of these systems are influenced by the bound xanthone. Further, whether these collective motions in AR when bound to xanthone as captured in eigenvector 1 and eigenvector 2 (or PC1 and PC2) represents the lowest energy states was determined through Gibb's free energy analysis. Lower the value of Gibb's free energy better is the stability of corresponding conformation and the system. The results show that mangostenone A and smeathxanthone B have a greater number of low energy state conformations compared to other xanthones (dark blue regions in Figure 9). Further, mangostenone B and bangangxanthone A have a greater number of low energy conformational states compared to allanxanthone B. The energy range for all the systems is around −3 to 3 kJ/mol. The results of Gibb's free energy estimation suggests better stability of system with mangostenone A, mangostenone B, bangangxanthone A, and smeathxanthone B compared to bare AR. MM-PBSA calculations offer more accurate estimates of binding affinity [47,48] The MM-PBSA calculations were performed on the frames of all the xanthone bound complexes isolated at each 200 ps from the simulation period 50 ns to 100 ns (total 2500 frames). The binding free energy (∆G binding ) for each system was calculated from the van der Waals energy, electrostatic energy, polar solvation energy, and SASA energy (Table 3). In MM-PBSA calculations, the non-bonded interactions such as van der Waals energy and electrostatic energy were measured in terms of Coulomb and Lennard-Jones (LJ) potential functions respectively. These two energy evaluations have the major influence on the binding free energy (∆G binding ) estimates. Smeathxanthone B has the lowest van der Waals energy and thus has the lower ∆G binding energy of −139.426 kJ.mol −1 (Figure 10).  Mangostenone A and allanxanthone B have almost similar and lower van der Waals energies and hence also have almost similar and favorable ∆G binding . The van der Waals energy in the case of bangangxanthone A is lower compared to mangostenone B. However, the polar salvation energy for bangangxanthone A is higher than mangostenone B. Consequently, the ∆G binding for bangangxanthone A and mangostenone B is higher than other three xanthones. Further, as shown in Figure 10 the simulation period from 70 ns and thereafter in the case of mangostenone B and bangangxanthone A shows better ∆G binding almost comparable to other three xanthones. Overall results of MM-PBSA calculations suggest mangostenone A, smeathxanthone B, and allanxanthone B may have the strong binding affinity toward the AR. On the basis of these MM-PBSA calculations, these xanthones could be taken up for further evaluation for the in vitro and in vivo testing.
Further, in order to investigate the key residues at the binding site contributing to the binding free energy, the per-residue binding free energy decomposition studies were performed (Supplementary Materials, Figure S6). The residue Trp20 most significantly contributes to the lowest ∆G binding of smeathxanthone B, while Phe122, Trp11, and Leu300 contribute significantly in the case of mangostenone A. The residue Trp219 alone contributes to in the case of allanxanthone B. The results of per-residue domain decomposition study for mangostenone A and mangostenone B revealed that the contribution of residue Trp111 is important for lower ∆G binding of mangostenone A. In the case of bangangxanthone A, although the residues Trp20, Tyr209, and Leu300 contribute to∆G binding , the acidic residues Asp43 and Glu185 adversely affect the overall ∆G binding .

Global Reactivity Descriptor and Stability Properties of Hit Molecules
The global reactivity descriptor and stability properties of the hit molecules against AR enzyme were calculated and are presented in Table 4. ∆E gap = energy gap; µ = chemical potential; η = hardness; S = softness; χ = electronegativity; ω = electrophilicity index.
The chemical structures of mangostenone B, bangangxanthone A, smeathxanthone B, mangostenone A, and allanxanthone B optimized at B3LYP/6-31G (d,p) level are shown in Figures 11 and 12. The highest occupied molecular orbital (E HOMO ) of a molecule elucidates its ability to donate an electron, while the lowest unoccupied molecular orbital (E LUMO ) providesrelevant information on the electron acceptability of a molecule. Mangostenone B (E HOMO = −5.31 eV) had the highest E HOMO ,thereby suggesting it as the best electron donating moiety among the hit AR inhibitors. Conversely, bangangxanthone A (E LUMO = −1.93 eV) had the lowest E LUMO value, indicating the molecule as the best electron acceptor ( Table 4).
The energy gap (∆E gap ) of a molecule helps in predicting the chemical reactivity, biological activity, and kinetic stability of a molecule [49]. The results obtained showed that mangostenone B exhibit the highest kinetic stability potential. Hence the stability of the hit molecules is in the order of mangostenone B > bangangxanthone A > allanxanthone B > mangostenone A > smeathxanthone B. Also, smeathxanthone B is the highest chemical reactivity potential and is in the order of smeathxanthone B > mangostenone A > allanxanthone B > bangangxanthone A > mangostenone B.
The chemical hardness and softness of a molecule is an essential parameter used to predict its stability and reactivity [50,51]. Hit molecules with lower hardness have a greater potential to exchange charges with the amino acid residues at the active site of the AR enzyme. Hence, it is useful in predicting the bioactivity of a molecule. In our study, Smeathxanthone B had a calculated hardness value of 1.79 eV and a softness value of 0.56 eV. Therefore, the Smeathxanthone B moiety is predicted as the most promising inhibitor against the AR enzyme.
Furthermore, the electronegativity (χ) of a molecule gives relevant information on its ability to attract an electron. It also helps in predicting the possible biological potential of a molecule [52,53]. Hence, the higher the electronegativity value, the greater the tendency to accept an electron. Bangangxanthone A (χ = 3.75 eV) had the highest calculated electronegativity value, thereby making it the best electron-accepting species among the hit molecules.
The chemical potential of a molecule reveals its reactivity and stability tendencies. Higher the chemical potential value, lower is the stability and the higher the reactivity. A highly negative chemical potential translates to easy acceptability of electrons and an uneasy tendency to remove electrons [51]. Our study showed that mangostenone B (−3.35 eV) has a greater reactivity potential and is less stable than mangostenone A (−3.47 eV), allanxanthone B (−3.59 eV), smeathxanthone B (−3.70 eV), and bangangxanthone A (−3.75 eV) (Table 4). Additionally, the electrophilicity index of a molecule is used in evaluating the biological potency of a molecule against a target disease. The higher the electrophilicity index of a molecule, the stronger its electrophilic property. Bangangxanthone A (ω = 3.86 eV) is the highest electrophile when compared to other hit molecules.

Molecular Electrostatic Potential of the Top-Ranked Xanthones
Molecular electrostatic potential (MESP) is a useful technique that has become extremely relevant in explaining the charge distribution pattern of drug leads and bimolecules [54]. It gives important information about the electron-donating and -accepting areas of ligands. Moreover, useful details about the hydrogen bond interactions, chemical reactivity, molecular size, as well as the negative, positive, and neutral electrostatic potential areas of drug like molecules can be obtained from MESP [55]. The MESP illustration diagrams of mangostenone B, bangangxanthone A, smeathxanthone B, mangostenone A, and allanxanthone B are shown in Supplementary Materials, Figure S7. The negative, positive, and neutral electrostatic potential regions on the top-ranked molecules are represented in red, blue, and green colors. The positive site is localized on H49-mangostenone B, H50-bangangxanthone B, H50-smeathxanthone B, H54-mangostenone A, and H51allanxanthone B indicating possible nucleophilic attack on the ligand through this atoms. Moreover, the negative site that is suitable for electrophilic attack was localized on mangostenone B (O3, O4, and O5), bangangxanthone A, and smeathxanthone B (O2, O4, and O5), mangostenone A and allanxanthone B (O2, O3, O4, and O5). These electrophilic and nucleophilic parts of the ligand may participate in hydrogen bonding which can further strengthen their stability in the active site of the AR enzyme and also improve their interaction with the enzyme.

Receptor and Ligand Preparation
The 3D crystal structure of AR (PDB ID: 1EL3) was obtained from the protein data bank (http://www.rcsb.org/pdb/ 1 December 2020). This crystal structure was chosen on the basis of its better resolution (1.70 Å) and as it belongs to homosapiens species. Further, this particular structure has a bound NADP, a benzenesulfonamide ligand and has no missing residues making it a better choice for docking and MD studies. First, the protein was prepared by removing water molecules, native ligand, cofactor, and ions, respectively. Later on, charges and hydrogen were added to the proteins using the MGLTools.
Sixty-five xanthones isolated from African medicinal plants were downloaded from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/ 11 December 2020) in SDF format (see supplementary information). The ligands' energy was minimized under MMFF94x force field using the steepest descent method for 100 steps with a step size of 0.02 using open babel program [56] in PyRx 0.8 software [57].

Validation of Docking Protocol and Molecular Docking Studies
The docking protocol used in this study was validated before virtual screening was performed on the AR enzyme. To validate the docking protocol, residues within 5 Å in the binding site of the enzymes were considered, as shown in Table 5. Later on, water molecules, cofactor, and native ligand were removed from the enzyme's binding site. Finally, the native ligand was removed from the enzyme and re-docked into its binding pockets. The best conformational pose of the re-docked ligand was selected in each case, and their RMSD values were calculated using the PyMol software. An ensemble docking method was employed to obtain the best binding energy for the ligands to perform the docking study. In this method, the bare protein was subjected to 100 ns equilibrated molecular dynamics simulation. The resulting MDS frames were clustered using gromacs clustering method with RMSD cutoff 0.08 nm. The structure representing each cluster's centroid was isolated for docking studies. This process gave 22 conformations of AR protein. These conformational ensembles were aligned with the aid of Chimera program. The proteins were prepared using MGLTools 1.5.2. Polar hydrogen atoms were added to the enzymes, and non-polar hydrogen atoms were merged. Thereafter, the grid box of size 16 Å 3 and center was set at x = 18.551, y = −11.469, z = 17.63 for the AR enzyme, followed by ensemble docking using the AutoDockVina tool [31]. Results with the best binding energy and lowest root mean square deviation (RMSD) value were selected for each ligand. After completing the virtual screening, the hydrogen bond interaction, hydrophobic interaction, Pi-interaction, and electrostatic interaction of the protein-ligand complex of the five-hit molecules were chosen and analyzed with Discovery studio visualizer (2020).

ADME and Drug Likeness Predictions
The absorption, distribution, metabolism, and excretion (ADME) study helps predict the druggability potentials of ligands in the drug discovery process. This study accessed ADMET properties such as solubility, blood-brain barrier, human intestinal absorption, carcinogenicity, and toxicity of the hit molecules using the admetSar-2 online tools (http://lmmd.ecust.edu.cn/ admetsar2/ 3 March 2021). Moreover, drug-likeness properties based on Lipinski rule were studied [58]. The ADMET and drug-likeness predictions are essential in elucidating the drug-likeness properties of the hit molecules as a therapeutic arsenal in the management of type 2 diabetes.

Molecular Dynamics Simulation and MM-PBSA Calculations
The complexes of top five xanthones with AR having lowest binding free energies and the bare AR protein were subjected to the MD simulations using Gromacs2020.4 [59]. The high performance computing (HPC) cluster at Bioinformatics Resources and Applications Facility (BRAF), C-DAC, Pune afforded the MD simulation studies. The topology of the protein was constructed using Chemistry at HARvard Macromolecular Mechanics 36 (CHARMM36) force field, while ligand topologies were constructed using CHARMM36 force field implemented in CGenFF server [60][61][62][63]. The systems were solvated with a simple point charge (SPC) water model and appropriate combination of ions such as sodium and chloride were used to generate the neutral environment [64]. Such solvated and neutral system was subjected to unrestrained energy minimization with steepest descent criteria in order to remove the steric clashes. The resulting system was briefly equilibrated at constant pressure, volume and temperature (300 K) conditions for 100 ps. During NVT equilibration the temperature was held at 300 K using modified Berendsen thermostat (velocity rescaling thermostat) [65] for 100 ps. The NPT equilibration was done where the pressure was hold at 1 bar atmospheric pressure using Berendsenbarostat [66,67]. The equilibrated systems were subjected to 100 ns MD simulation with restrain on covalent bonds using LINear Constraint Solver (LINCS) algorithm [68] and the cut-off value of 12 Å for the long range electrostatics such as Coulomb and Lennard Jones with the Particle Mesh Ewald (PME) method [69]. During production phase MDS, the thermostat and barostat employed were the modified Berendsen thermostat and Parrinello-Rahman barostat respectively [70]. Post MDS, the measurement of deviations in the protein and ligand atoms as RMSD; the extent of fluctuations in the atoms of residues in protein in terms of root mean square fluctuations (RMSF) measurement and radius of gyration (Rg) which is a root mean square distance of a collection of atoms from their common center of mass was undertaken. The hydrogen bond formation was also investigated. Further, the contact analysis using gmxmdmat program was performed to understand the residues within close and within acceptable distance from the xanthones under study. Principal component analysis (PCA) and Gibb's free energy analysis was performed to analyze the most stable conformations and overall stability of AR bound with xanthones. PCA was performed with gmxcovar and gmxanaeig programs. The covariance matrix was constructed during PCA analysis for the backbone atom of AR. The eigenvectors and eigenvalues were obtained after diagonalizing the covariance matrix for each trajectory. Two principal components (PC1 and PC2) were used for the PCA analysis, where the corresponding eigenvectors represented the path of motion. Further, these two PCs were used to analyze the lowest energy state of each complex as Gibb's free energy landscape (Gibb's FEL) with gmx sham program. In the Gibb's FEL plots the deep valleys represent lowest energy states while the boundaries between the deep valleys show intermediate conformations. The most crucial estimate, binding energy of protein-ligand complex (∆G binding ) was estimated through the molecular mechanics-based Poisson Boltzmann surface area continuum solvation calculations (MM-PBSA) [47]. The frames extracted at each 200 ps for the simulation period 50 ns to 100 ns (2500 frames) were used in MM-PBSA calculations. Further, the per-residue domain decomposition studies were carried out on each complex to investigate the contribution of important residues in the binding free energies.

DFT Studies of Top Ranked Aldose Reductase Inhibitors
DFT analysis of the top-ranked xanthones was carried out using the Spartan 14 program with B3LYP functional (Lee-Yang-Parr exchange-correlation functional method) and 6-31G basis set [71]. During the calculations, the values of the frontier orbital energies were computed from the most established conformation of the compounds.

Conclusions
The study explored xanthones isolated from African medicinal plants for their AR inhibitory potential. Molecular docking studies, specifically ensemble docking studies, were performed to identify the potential inhibitors for AR. Five hit molecules were selected based on their binding affinities against the enzyme. The drug-likeness and ADMET studies revealed the selected molecules as suitable inhibitors for the AR enzyme. The phytochemicals were found to stabilize the structure of AR in terms of RMSD, RMSF, Rg, Gibb's free energy and were also found to form key hydrogen interactions at the binding sites. The study finally revealed mangostenone B, bangangxanthone A, smeathxanthone B, mangostenone A, and allanxanthone B as potential AR inhibitors. The MM-PBSA calculations further suggested the strong binding affinities of smeathxanthone B, mangostenone A, and allanxanthone B. These hit molecules could be tested in vitro and in vivo to further establish their potency against the target enzyme. In addition, they can be modified and optimized for the development of improved AR medications with lesser side effects.
Supplementary Materials: The following are available online https://www.mdpi.com/article/10.339 0/computation10090146/s1. Table S1. List of compounds and their PubChem IDs; Figure S1. RMSD in AR backbone atoms when bound to xanthones; Figure S2. RMSD in ligand atoms when bound to AR; Figure S3. Contact map analysis between ligands and residues within 1.5 nm F contact distance during MD simulation; Figure S4. Root mean square fluctuations in AR residues in bare protein and when bound to xanthones under study; Figure S5. 2D projection of eigenvector 1 and eigenvector 2 for all the systems under study; Figure S6. Results of per-residue domain decompostion study; Figure S7 Funding: The authors would like to extend their sincere appreciation to the department of chemistry, Doornfortein campus, University of Johannesburg, South Africa for providing fund for the work.