Towards New Antidiabetic Development: Molecular Interactions Unfold Potential Agonistic Mechanisms of Naturally Occurring Terpenoids Against Structural Landscape of Glucagon-like Peptide-1 Receptor

1. Introduction

Diabetes mellitus (DM) is a significant public health concern, with the World Health Organization predicting it will be one of the top ten global leading causes of death by 2030 [1]. According to Hossain et al., [2] there are about 536.6 million people currently living with diabetes mellitus. DM is mainly categorized into Type 1 (T1DM) and Type 2 (T2DM). T1DM results from the destruction of pancreatic beta cells, causing insulin deficiency. T2DM on the other hand is caused by insulin resistance. However, genetic susceptibility, obesity and sedentary lifestyles are known to have strong correlation with type 2 diabetes, commonly found in adults [3]. Both DM types cause hyperglycemia, leading to complications like retinopathy, nephropathy, neuropathy, and cardiovascular diseases [4]. Gestational DM also surfaces during pregnancy due to insulin deficiency or hormonal changes. The focus of treatment, especially in T2DM, is normalizing blood glucose levels using hypoglycemics agents like biguanides, thiazolidinediones, sulfonylureas, and DPP-4 inhibitors, along with injectable options like insulin and incretin mimetics [5,6].

One important incretin for T2DM treatment is glucagon-like peptide-1 (GLP-1) [7]. In a glucose-dependent way, GLP-1 binds to the class B GPCR GLP-1 receptor (GLP-1R), which causes insulin release and glucagon suppression [8]. Additionally, GLP-1R activation slows food intake, improves clinical outcomes, delays gastric emptying, promotes satiety, and reduces weight [8]. GLP-1R increases GLP-1 secretion, raises Ca²⁺ levels, and stimulates cyclic AMP synthesis upon agonist binding [9]. A hydrophilic N-terminal domain connected by loops to a transmembrane C-terminal domain is essential for agonist binding in GLP-1R, which comprises 463 amino acids [10]. According to previous studies, GLP-1R exhibits a relative flexibility in aqueous solutions. A disulphide bond between CYS46–CYS71, CYS62–CYS104, and CYS85–CYS126 is also formed between the six highly conserved cysteine residues in the N-terminal region of this protein. Furthermore, class B GPCRs share a high degree of conservation among the ASP67, TRP72, PRO86, ARG102, GLY108, and TRP110 residues, that are essential for the stability [7,11,12]. It has been demonstrated that TRP72 and TRP110 are crucial for agonist binding in GLP-1R, and the crystal structure of the extracellular domain of GLP-1R demonstrates the key location of these conserved residues. The PRO86 residue is crucial for agonist binding, and interactions between these residues in GLP1R are crucial for stabilizing the N-terminal domain [11]. Additionally, studies have demonstrated that TRP33 is required for non-peptide agonist GLP-1R binding, but it is not essential for natural activation of the receptor [12,13]. Previous work already posited the critical need of TRP33 for the function of small-molecule GLP-1R agonists due to its promotion of cAMP production in GLP-1R-expressing cells [14]. The demonstration that non-peptide agonists of the GLP-1R are not required to mimic the extensive receptor contacts formed by peptides within the transmembrane cavity to promote receptor activation will advance the pursuit of non-peptide agonists for therapeutically important class B receptors, including GLP-1R [15,16].

Although synthetic small molecules for T2DM have advanced, natural products (NPs) still have significant potential [17,18,19,20,21]. Drug discovery is sparked by the antidiabetic properties of several plants. For instance, metformin was derived from glargine, a guanidine molecule found in Galega officinalis [22]. NPs extend the chemical field for drug discovery by offering benefits including scaffold diversity and structural complexity. Many groups have dedicated resources to investigation of effect of natural products on GLP1R [23,24]. Phytoconstituents classes like alkaloids (berberine, capsaicin), flavonoids (quercetin, myricetin, luteolin, genistein), polyphenols (curcumin, resveratrol, EGCG), and ginsenosides have been described as acting via the GLP1 pathway [25]. Most of them boost GLP 1R signaling (by downstream cAMP), inhibits DPP 4, improve insulin sensitivity in muscle and adipose tissue, and preserve β cell mass. The limitations of these study lie in their restricted approach, with many only applying high content screening in their investigations [23,24,25,26,27]. A downstream effect without a clear picture of targets of these compounds will limit further optimization efforts. A full understanding of how these molecules potentially interact with the target in biological systems is an insight that molecular dynamics simulation provides.

In this study, we bridged the gap by a full computational investigation of select natural products from plants with known antimalarial compounds against GLP1R. We generated these compounds (244) from published literatures, prepared and docked them to the GLP1R receptor using Autodock Vina®. Our initial findings from molecular docking reveal terpenoid scaffolds as very active against GLP1R. We then sought to understand the dynamic nature of this interactions by molecular simulation for 250ns. We believe that a full terpenoid structure, with a hydrophilic head at one end, and many aliphatic groups in the opposite end provide the needed length for a robust interactions with GLP1R. The hydrophobic effect also played a great role as terpenoid nucleus engaged in robust π interactions (π–π stacking and hydrophobic π-alkyl ) with adjacent residues. Compound B1245 (lupeol) is golden here as it revealed the most energetically favored interactions with GLP1R. Our study is robust, thorough and provide a clear path from natural compounds associated with antimalarial properties to potentially lead agonist of GLP1R.

2. Methods

2.1. Retrieval and Preparation of Proteins

GLP-1R, a crucial protein target linked to the pathogenesis of T2DM was selected for this investigation [28,29,30]. The Protein Data Bank (PDB) (https://www.rcsb.org/) provided the three-dimensional structure of the protein (PDB ID: 7S15) in complex with its co-crystallized ligand. The crystal structure of GLP-1R protein was prepared for docking using the UCSF Chimera 1.14 tool [31]. The structure was cleared of all ions, water molecules, and non-standard residues. The structure editing wizard was used for energy minimization of the crystal structure at 200 steepest descent steps, 0.02 steepest descent steps size (Å), 10 conjugate gradient steps, and 0.02 conjugate gradient steps size (Ᾰ). Hydrogens were added as well as partial charges using the Gasteiger force field [32].

2.2. Ligand Selection and Preparations

A total of 244 natural compounds selected from different medicinal plants with ethnomedicinal usage in the treatment of diabetes were used in this study (Supplementary Table S1). The Structure Data File (SDF) for these substances was compiled from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/). Marvin Sketch 20.22 was used to build structures that were not available in the PubChem database and were converted to 3D. Danuglipron served as the reference drug. Pfizer developed this oral GLP-1R agonist, which is presently undergoing clinical studies [13]. A 3D structure of danuglipron was prepared as previously mentioned after being obtained from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/compound/134611040). The PDB files of all the ligands were optimized and energy minimized using the UCSF Chimera tool, and Gasteiger force field was used to calculate and assign charges to the ligands [33].

2.3. Molecular Docking

Site-specific docking was carried out after initial mapping of the active sites of the protein using the zone wizard on the selection pane in the USCF Chimera program, where the area of amino acids that was less than 5 Å from the co-crystallized ligand was considered the binding pocket for the protein [34]. The docking campaign was conducted using the Auto-dock Vina of PyRx software. The resulting files were uploaded to the Discovery Studio Client 2021 and Chimera 1.14 workspace in order to perform post-docking analysis and visualization [35].

2.4. Molecular Dynamics Simulations (MDS)

For MD simulation, five (5) ligands showing the most favorable binding affinities towards GLP-1R were eventually selected using the cut-off binding energy of -11.0 kcal/mol. The ligands and receptors for the MD simulations were also prepared using the UCSF Chimera. The five ligands complexed with the GLP1R protein, together with the reference danuglipron, were subjected to MD simulation analyses. The protein is characterized using the AMBER Force Field (FF18SB variant) package in the graphics processing unit (GPU) version of the PMEMD.CUDA engine [36]. The stability, compactness, and structural dynamics of the docked complexes were examined with the aid of MD simulations [32]. ANTECHAMBER was used to add partial charges to the ligands. The bound and unbound bio-complexes were neutralized and solvated by introducing H⁺, Na⁺, and Cl^- counter ions using the AMBER18 LEAP module. Furthermore, all-atom explicit solvation was carried out in an orthorhombic TIP3P box containing 10 Å-sized water molecules. With an applied restricted potential of 500 kcal/mol, the initial minimization of 2500 steps was taken into consideration. Full minimization (5000 steps) is the next step, which was completed by a conjugate algorithm without any constraints. In order to maintain a fixed volume and number of atoms while taking into account a canonical ensemble (NVT), the systems were then gradually heated from 0 to 300 K for 50 ps. The Barendsen-Barostat's one bar pressure was applied using the shake algorithm to lessen the hydrogen bond constraint [37,38]. A Langevin thermostat was used to operate the production MD for 200 ns for all systems with a time scale of 2 fs at 300 K and constant pressure of 1 bar. Every 1 ps, the coordinates of the Apo-GLP1R and ligand-GLP-1R complexes were stored, and the MD parameters were examined using the CJTRAJ and PTRAJ modules in AMBER 18/GPU. Principal component analysis (PCA), radius of gyration (RoG), C-α root mean square of deviation (C-α RMSD), and root mean square of fluctuation (C-α RMSF) are among the MD parameters examined [39]. Chimera and the Discovery Studio 2021 client were used to visualize ligand-receptor interactions and perform post-MD structural analysis [40]. For data analysis, Microcal Origin 6.0 software was used to obtain plots of RMSD, RMSF, RoG, PRED, and PCA [41,42].

2.4.1. Binding Free Energy (BFE) Analysis

Molecular dynamics simulation has been a useful computational technique for elucidating the mechanisms of biological processes. Important details on a ligand's binding interactions with a receptor are provided by the thermodynamic binding free energy (BFE). The BFE of a ligand to a protein can be determined using the Molecular Mechanics Poisson-Boltzmann Surface Area (MM/PBSA) approach [43]. Over 200000 images from the 200 ns trajectories were used to compute the energies. Solvents were totally eliminated in the MM/PBSA procedure, and were replaced by a dielectric continuum. The BFE (ΔG) for the compounds- and danuglipron-bound GLP1R systems was calculated using the following expression:

ΔG_bind = ΔG_complex - ΔG_receptor - ΔG_ligand

ΔG_bind = E_gas + G_sol - TΔS

E_gas = E_int + E_vdw + E_ele

G_sol = G_GB + G_SA

G_SA = γSASA

ΔG_bind = gas-phase summation; E_gas = gas-phase energy; G_sol = free solvation energy; TΔS = entropy; E_int= internal energy; E_ele = coulomb and coulomb; E_vdw = van der Waals energy.

2.4.2. Receptor–Ligand Interactions Systems

Using Discovery Studio Visualizer 2021 Client, the active site residues of GLP-1R that form an interaction network with the ligands were analyzed [44]. The different intermolecular forces that occur between GLP-1R and the ligands were identified. The interacting residues of GLP-1R and the types of chemical bonds within the bound complexes were analyzed using Discovery Studio Visualizer [35].

2.4.3. Per-Residue Energy Decomposition (PRED) Analysis

In relation to GLP-1R inhibition, the PRED calculates the energy contributions of the amino acid residues in the binding site of the protein to the overall stability of the chosen compounds. The energetic contributions of individual interacting residues are enhanced and made visible by PRED analysis. All things considered, it improves our comprehension of the connections between protein structure and function. The decomposition of the remaining energy contribution to BFE in the protein-ligand complex was examined using the MM/PBSA method in AMBER18 GPU [32].

2.5. Drug-likeness Studies

In order to make sure that a possible small molecule drug is safe, effective, and able to reach its intended biological target in the body, drug-likeness studies—often ADME/Tox—are essential to drug design. The open-source SwissADME web tool (http://www.swissadme.ch/) was utilized to screen the pharmacokinetic characteristics of the six drugs employed in MDS studies.

3. Results and Discussion

3.1. Molecular Docking Analysis

Medicinal plants play important roles in human wellness. Secondary metabolites from plants are effective in the treatment of several diseases, including diabetes. In this study, 244 Phyto-compounds from medicinal plants used locally in the treatment of diabetes mellitus (DM) were subjected to in silico screening for their antidiabetic potential. These phytochemicals were investigated for their binding potentials against Glucagon-like peptide-1 receptor (GLP-1R).

Resulting from the virtual screening, an initial selection of very promising compounds was done based on ranking with a cut-off threshold of < -10.0 kcal/mol (Table 1). This was done in order to highlight the phytochemicals that bind to the individual target the strongly and most favorably. The ligand (82L) co-crystallized with the GLP-1R had a docking score of -10.7 kcal/mol; we then established a new cut-off binding energy of -11 kcal/mol for our molecules for the purpose of MD simulation, and the compounds advanced are shown in Table 2. The highest binding affinity of these compounds is -13.3 kcal/mol for A1238 (Friedelan-3-one), followed by -12.4, -12.2, -11.7, -11.5, and -11.0 kcal/mol for A1251 (alpha_amyrin), A1245 (beta-amyrin), B1245 (lupeol), B1235 (pyrophorbide), and A1236 (Alpha sitosterol) respectively. Our goal is in the view to unveil potential small molecule oral agonists, in view of this, B1235 (pyrophorbide; -11.5 kcal/mol) was discontinued because it had a molecular weight more than 500 g/mol. Danuglipron which is a Pfizer's small molecule agonist of GLP-1 receptor used as a positive control had binding affinity of -11.1 kcal/mol. As already mentioned in the methodology section, after our iteration was finished, we had five (5) compounds with good binding affinity to move on to MD simulation (Table 2).

The binding energies with the highest negative values show the most stable interactions between the drugs and GLP-1R [45]. Table 2 shows the range of intermolecular interactions that the ligands formed with GLP-1R. They essentially form hydrogen bonding and hydrophobic bonds. Remarkably, there was no hydrogen bond interaction between the two molecules with the highest binding energies, A1238 and A1251. However, others formed hydrogen connections with SER2, ASP163, SER171, THR172, GLN186, CYS191, GLN199, and PHE342. Each molecule has at least three residues with which it reacts with hydrophobic interactions (Table 2).

Table 1. The dock scores of 21 top-performing compounds against GLP-1R protein targets.

S/N	Code	Ligands	Binding energy (kcal/mol)
1.	A1234	Angustidine	-10.4
2.	A1235	Trans-tiliroside	-10.6
3.	A1236	Alpha sitosterol	-11.0
4.	A1238	Friedelan-3-one	-13.3
5.	A1239	Stigmasterol	-10.5
6.	A1242	7-deacetylgedunin	-10.7
7.	A1245	Beta-Amyrin	-12.2
8.	A1251	alpha_amyrin	-12.4
9.	A1252	Arjunglucoside_I	-10.3
10.	A1253	Arjunglucoside_III	-10.4
11.	A1257	Hymenocardine	-10.2
12.	B1235	Pyropheophorbide	-11.5
13.	B1237	Simiarenone	-10.6
14.	B1245	Lupeol	-11.7
15.	B1250	21-O-ethylstrictosamide	-10.2
16.	B1265	Sericoside	-10.3
17.	B1267	Trachelosperoside_E1	-10.3
18.	B1270	13-hydroxy-15-Oxozoapatlin	-10.3
19.	B1277	Arjunglucoside_II	-10.5
20.	B1281	Betulinic_acid	-10.8
21.	B1297	Terminolic_acid	-10.2
22.	B1304	Co-crystallized ligand (82L)	-10.7

Table 1. The dock scores of 21 top-performing compounds against GLP-1R protein targets.

S/N	Code	Ligands	Binding energy (kcal/mol)
1.	A1234	Angustidine	-10.4
2.	A1235	Trans-tiliroside	-10.6
3.	A1236	Alpha sitosterol	-11.0
4.	A1238	Friedelan-3-one	-13.3
5.	A1239	Stigmasterol	-10.5
6.	A1242	7-deacetylgedunin	-10.7
7.	A1245	Beta-Amyrin	-12.2
8.	A1251	alpha_amyrin	-12.4
9.	A1252	Arjunglucoside_I	-10.3
10.	A1253	Arjunglucoside_III	-10.4
11.	A1257	Hymenocardine	-10.2
12.	B1235	Pyropheophorbide	-11.5
13.	B1237	Simiarenone	-10.6
14.	B1245	Lupeol	-11.7
15.	B1250	21-O-ethylstrictosamide	-10.2
16.	B1265	Sericoside	-10.3
17.	B1267	Trachelosperoside_E1	-10.3
18.	B1270	13-hydroxy-15-Oxozoapatlin	-10.3
19.	B1277	Arjunglucoside_II	-10.5
20.	B1281	Betulinic_acid	-10.8
21.	B1297	Terminolic_acid	-10.2
22.	B1304	Co-crystallized ligand (82L)	-10.7

Figure 1. Structures of the top five potential agonists of GLP-1R.

Figure 1. Structures of the top five potential agonists of GLP-1R.

Table 2. Table of binding affinity of top five selected compounds against the GLP1R protein compared to the positive control.

S/N	Code	Compound	Binding Energy	Hydrogen Bond	Hydrophobic Interactions
				Amino Acid Residues
1.	A1236	Alpha-Sitosterol	–11.0	SER2, THR172	LEU3, TRP4, PRO102, LEU106, LYS162, TRP168, PHE195, PHE342
2.	A1238	Friedelan-3-one	–13.3	Nil	TYR110, TRP168, PHE342, LEU345
3.	A1245	Beta-Amyrin	–12.2	GLN199	ARG264, PHE195, MET198, TRP168, LEU106, PHE346, PHE342, LEU345
4.	A1251	Alpha-Amyrin	–12.4	Nil	LEU3, LEU106, LEU107, TRP168, LEU182, PHE195, MET198, PHE342, PHE346
5.	B1245	Lupeol	–11.7	CYS191, GLN186	LEU3, TRP4, LEU107, TYR110, LYS162, LEU166, TRP179, LEU182, PHE342
6.	PC	Danuglipron	–11.1	GLU103, ASP163	TRP4, ALA174, PHE342

Table 2. Table of binding affinity of top five selected compounds against the GLP1R protein compared to the positive control.

S/N	Code	Compound	Binding Energy	Hydrogen Bond	Hydrophobic Interactions
				Amino Acid Residues
1.	A1236	Alpha-Sitosterol	–11.0	SER2, THR172	LEU3, TRP4, PRO102, LEU106, LYS162, TRP168, PHE195, PHE342
2.	A1238	Friedelan-3-one	–13.3	Nil	TYR110, TRP168, PHE342, LEU345
3.	A1245	Beta-Amyrin	–12.2	GLN199	ARG264, PHE195, MET198, TRP168, LEU106, PHE346, PHE342, LEU345
4.	A1251	Alpha-Amyrin	–12.4	Nil	LEU3, LEU106, LEU107, TRP168, LEU182, PHE195, MET198, PHE342, PHE346
5.	B1245	Lupeol	–11.7	CYS191, GLN186	LEU3, TRP4, LEU107, TYR110, LYS162, LEU166, TRP179, LEU182, PHE342
6.	PC	Danuglipron	–11.1	GLU103, ASP163	TRP4, ALA174, PHE342

3.2. Drug-likeness and ADME/Tox Studies

The physicochemical properties of hit and lead compounds are very important in translational drug discovery campaigns. Its evaluation helps reduce attrition rates and failures in later stages of clinical trials. It is always faster, more reliable and less costly to employ in silico toxicity and biological property prediction before validation with traditional in vitro and in vivo models The five top-performing compounds from this study, α-sitosterol, friedelan-3-one, β-amyrin, α-amyrin and lupeol, as potential agonists of GLP-1R were evaluated for their physicochemical properties and pharmacokinetics profiles, affording the prediction and identification of their drug-like characteristics. The drug-like properties of a molecule must be confirmed by balancing the parts of the molecule contributing to its preferential ADMET behavior against the ones for efficacy to reduce attrition risks. Computational science, owing to its deep predictive capabilities, has simplified the prediction of drug-likeness of molecules and aided rational drug design. In Table 3, the physicochemical parameters of the top five hit compounds are shown. Their molecular weight, number of rotatable bonds, hydrogen bond donors and acceptors, solubility, blood-brain barrier (BBB) permeation, gastrointestinal (GI) permeation, ease of synthesis and violation of Lipinski’s rule of 5 were assessed. The qualification of these compounds as druggable was ascertained by predicting their adherence to the Lipinski RO5 recommendations for orally active drugs. The molecular weight, number of rotatable bonds, hydrogen bond acceptors and donors of the compounds showed no violation of the RO5. However, LogP was the only property where the RO5 was violated with values greater than 5 across the compounds. Lipinski has previously reported that any drug molecule that violates two or more of the rules would be considered not orally active. therefore, the five compounds investigated in this study are tipped as orally bioactive compounds [46].

The topological polar surface area (TPSA) of a compound is a crucial parameter in a drug discovery program as it essentially predicts the ability of drug molecules to permeate the blood-brain barrier (BBB). Chemical compounds must be able to pass through the blood-brain barrier in order to target essential receptors and proteins in the brain microenvironment [47]. Generally, substances with a TPSA of less than 140 Å2 can pass through a variety of biological membranes; however, a TPSA of significantly less than 90 Å2 is necessary to pass through the blood-brain barrier [48]. Every compound in Table 3, showed values much lower than 90 Å2, indicating a high capacity to penetrate into the various parts of the brain. Because the GLP1 receptor utilized in this investigation has been shown to be an agonistic receptor linked to aging, it is crucial to also look into the anti-aging potential of these drugs [49]. Furthermore, the TPSA values of the compounds indicated that they have oral activity. The bioavailability scores of these compounds corroborate this as well. With a bioavailability score of 0.55, they have a 55% chance of having an oral bioavailability of at least 10% at biological pH [50].

3.3. Molecular Dynamics Simulations

The computational approach known as molecular dynamics simulation (MDS) is frequently used to study the properties and behavioral changes of protein-ligand interactions [51]. Real-time conformational changes and stability of a protein-ligand interaction can be critically examined and understood with the aid of MDS [52]. This configuration offers comprehensive understanding of conformational changes and intermolecular interactions in complicated systems in simulated environments [53]. They provide atomic-level information of protein movements and interactions, illuminating the connection between protein structure and biological activity. Additionally, dynamic binding interactions with ligands, substrates, or inhibitors are clarified by MD simulations, which demonstrate how these interactions impact secondary structures and biological function in general [48]. All bound and unbound systems were put through 200ns MD simulations in order to understand the conformational changes and rearrangements of the protein (GLP-1R) under investigation. To gain insight into the structural dynamics and equilibrium upon binding, the C-α radius of gyration (RoG), C-α root mean square deviation (RMSD), C-α root mean square fluctuation (RMSF), binding free energy (BFE), and per-residue energy decomposition (PRED) were all extensively examined. All of the evaluated options behaved satisfactorily, according to an analysis of the simulation trajectories. The BFE of the majority of the compounds is similar to that of the positive control (Table 5).

3.3.1. Root Mean Square Deviations (RMSD)

The fluctuation in the C-α backbone is used to determine the stability of the protein-ligand interactions [54]. Perturbations of the protein-ligand complex for a long time is necessary to test the stability of the proteins. The principles of thermodynamics imply that the lower the RMSD values, the more stable the complex, while the higher values have much larger deviations and point to a protein with an unstable structure [52]. RMSD was carried out for all systems to determine their stability during simulations for 200ns. The average RMSD values for all systems are shown in Table 4. Among all the compounds, A1251 showed poor stability upon binding to GLP-1R. Fluctuation is obvious from the beginning of the simulation, and it peaks at around 75ns. It then lowers but remains highly unstable throughout the simulation time. When the average C-α RMSD values of A1236, A1238, and A1251 were compared to that of the unbound protein (Apo), it was found that their interactions with the protein were somewhat unstable because their average RMSD values were greater than those of the unbound protein (3.45 Å). In contrast, the average RMSD value of B1245 (3.45 Å) was the same with that of the unbound protein. Throughout the simulation period, there was little or no variation (Figure 2B). This implies that B1245 bound with GLP-1R have a comparatively better stability (Table 4). Furthermore, compared to the Apo system, the average RMSD of GLP-1R-danuglipron (control) was shown to be lower. This high stability demonstrates that danuglipron can activate GLP-1R, which in turn stimulates β-cells to produce more insulin [55].

3.3.2. Root Mean Square Fluctuations (RMSF)

The structural dynamics of amino acids, which are the building blocks of proteins, are crucial to the overall function of the protein [56]. By measuring the variation of the target protein's residues, RMSF analysis provides information about its flexibility. Understanding the displacement and stability of each GLP-1R residue in the simulation trajectory is made easier with the use of RMSF [57]. Therefore, the conformational changes in the protein's structure with or without ligand binding might be shown by comparing the average RMSF values of bound systems with those of Apo (unbound system). Higher RMSF values suggest more flexibility in residual movement, whereas lower RMSF values indicate fewer conformational changes during the simulation [32].

Danuglipron and the five compounds bound with GLP-1R were analyzed for RMSF. Between residues 1 to 200, contained in the extracellular domain of N-terminal which is crucial for successful agonist binding [58], the RMSF values of the residues were significantly low for B1245 and danuglipron (PC) similar to the Apo (Figure 2C). B1245-GLP-1R, danuglipron-GLP-1R, and the unbound systems were calculated to have average C-α RMSF values of 1.52, 1.51, and 1.71 Å, respectively, as shown in Table 4. Interestingly among all compounds, only B1245 and danuglipron showed to possess lower RMSF values than the unbound system, an indication of higher stability at the active site of the protein. This could indicate that GLP-1R fluctuations decreased when B1245 was bound (Figure 2C), which is consistent with the observation in the RMSD analysis. This bolsters the robust stability found in the RMSD experiment.

3.3.3. Radius of Gyration (RoG)

RoG is regarded as a fundamental indicator of the total size of a chain molecule. According to da Fonseca et al. [57], it can be used to measure how much a protein's structure changes during molecular dynamics simulations. The degree to which a protein's structure is organized tightly or loosely can influence its physiological properties. By evaluating the protein's flexibility and compactness in a biological setting, the RoG compares the protein's structure over time to the hydrodynamic radius that may be objectively observed [59]. Lower RoG values indicate a more structurally compact protein [32]. Unbound GLP-1R, A1236, A1238, A1245, A1251, B1245, and danuglipron have average RoG values of 28.80, 28.45, 28.84, 28.39, 28.61, 28.33, and 27.54 Å, respectively, according to Table 4. According to the RoG plots, GLP-1R achieved more compactness in danuglipron and B1245 (Figure 3B), indicating that the conformational alterations occurred during simulation did not impact the overall structural stability of the protein. This outcome is consistent with the results of the RMSD and RMSF of danuglipron and B1245.

3.3.4. Principal Component Analysis (PCA)

PCA measures the large-scale average motion of a protein and reveals the structures behind the atomic variations. Therefore, during simulations, PCA is used to highlight the most significant motions in a protein's structure [60]. The system's overall motility is given as the sum of the eigen values. Furthermore, PCA can assess a protein's flexibility under different circumstances [61]. The binding interaction between the top 5 drugs and GLP-1R was evaluated using the first two components of the principal components analysis (PCA) (Figure 3C). In comparison to the unbound protein, compound B1245 and danuglipron showed less dispersions along the axes. This validates the previous results of the RoG, RMSD, and RMSF investigations.

Table 4. Maximum, minimum and mean RMSD, RMSF, and ROG of the Apo and potential GLP-1R-bound Systems.

Complex	RMSD (Å) (Mean±SD)	RMSF (Å) (Mean±SD)	RoG (Å) (Mean±SD)
Apo-A1236	4.26±0.69	2.36±1.21	28.45±0.30
Apo-A1238	4.03±0.69	2.05±0.86	28.84±0.26
Apo-A1245	3.90±0.49	2.17±1.20	28.39±0.38
Apo-A1251	5.94±1.06	2.60±1.17	28.61±0.28
Apo-B1245	3.45±0.43	1.52±0.97	28.33±0.18
Apo-PC	3.33±0.40	1.51±0.59	27.54±0.19
Apo-GLP1R	3.45±0.76	1.71±1.13	28.80±0.26

Table 4. Maximum, minimum and mean RMSD, RMSF, and ROG of the Apo and potential GLP-1R-bound Systems.

Complex	RMSD (Å) (Mean±SD)	RMSF (Å) (Mean±SD)	RoG (Å) (Mean±SD)
Apo-A1236	4.26±0.69	2.36±1.21	28.45±0.30
Apo-A1238	4.03±0.69	2.05±0.86	28.84±0.26
Apo-A1245	3.90±0.49	2.17±1.20	28.39±0.38
Apo-A1251	5.94±1.06	2.60±1.17	28.61±0.28
Apo-B1245	3.45±0.43	1.52±0.97	28.33±0.18
Apo-PC	3.33±0.40	1.51±0.59	27.54±0.19
Apo-GLP1R	3.45±0.76	1.71±1.13	28.80±0.26

3.3.5. Binding Free Energy

Molecular mechanics/Poisson-Boltzmann Surface Area (MM/PBSA) is essential when examining drug-receptor interactions because it aids in determining the binding free energy (BFE) between a ligand and its target receptor [53]. The ligand's strength and favorability of binding to the receptor are determined by the BFE [62]. This is essential for a thorough evaluation of its effectiveness. Along with lowering the time and expenses of such drug discoveries, MM/PBSA can predict and then be used to optimize drug designs more accurately, helping to generate novel medications with desirable pharmacodynamic and pharmacokinetic properties [51]. Using the initial 40 ns of the MD trajectory, the BFE was calculated using MM/PBSA for A1236-GLP1R, A1238-GLP1R, A1245-GLP1R, A1251-GLP1R, B1245-GLP1R, and danuglipron-GLP1R.

The binding free energy of all the ligands and danuglipron showed values significantly higher than the one for the co-crystallize ligand 82L (-39.30 kcal/mol) (Table 5). From the results, the known agonist of GLP-1R, danuglipron, has a BFE of -56.47 kcal/mol. The respective ΔGbind (kcal/mol) values for the other bound-GLP1R systems, A1236-GLP1R, A1238-GLP1R, A1245-GLP1R, A1251-GLP1R, and B1245-GLP1R are –55.61, –53.47, –55.90, –56.99, and –57.23, respectively. Previously, danuglipron has been reported to be a small-molecule oral agonist of GLP-1R, with favorable binding and characteristic interactions at the binding pocket of the target [13] [14]. Interestingly, all the five top performing phytochemicals investigated here showed BFE values that is comparable with the standard agonist, danuglipron. By implication, this terpenoids could be said to have characteristic interactions at the binding pocket similar to danuglipron. Top among these compounds is Lupeol- B1245 with energy of binding –57.23 kcal/mol, suggesting the most favorable binding. It is the bond formation that promotes and stabilizes the association of drugs with receptors, therefore, the nature of the chemical bonds that bind drug molecules to their receptors influence the BFE. Based on the molecular interactions observed in the ligand-receptor complexes in this study, it is presumed that the -OH and C=O functional moieties present in the structures of these terpenoids which are similarly present in danuglipron contributed hydrogen bond interactions leading to the high binding energies recorded for the compounds.

Table 5. Summary of MM/PBSA-based computed thermodynamics binding free energy obtained for danuglipron and the proposed GLP-1R agonists.

Energy (kcal/mol)	Apo-A1236	Apo-A1238	Apo-A1245	Apo-A1251	Apo-B1245	Apo-PC	Apo-82L
ΔEvdW[a]	-63.42 ±3.19	–57.43 ±3.37	–57.72 ±2.57	–63.8 ±3.07	–60.98 ±2.93	–66.27 ±3.14	–48.3 ±5.89
ΔEelec	–8.01 ±3.42	1.03 ±2.55	–4.88 ±1.73	–7.32 ±1.81	–6.48 ±3.55	–116.02 ±7.31	–234.95 ±8.03
ΔGgas	–71.43 ±5.04	–56.39 ±4.88	–62.61 ±2.94	–71.21 ±3.77	–67.46 ±4.54	–182.30 ±7.40	–283.32 ±8.34
EGB	23.38 ±2.66	9.77 ±2.75	13.64 ±1.44	21.55 ±1.99	17.36 ±2.44	134.41 ±7.40	250.31 ±6.87
ESA	–7.55 ±0.27	–6.84 ±0.36	–6.94 ±0.23	–7.33 ±0.30	–7.13 ±0.26	–8.59 ±0.30	–6.28 ±0.76
ΔGsolv	15.82 ±2.58	2.92 ±2.64	6.71 ±1.52	14.22 ±2.05	10.23 ±2.35	125.82 ±7.38	244.02 ±36.88
ΔGbind	–55.61 ±4.11	–53.47 ±4.16	–55.90 ±3.59	–56.99 ±3.53	–57.23 ±3.44	–56.47 ±4.06	–39.30 ±6.34

ΔE_vdw = van der Waals, ΔE_elec = electrostatic, ΔG_bind = calculated total free binding energy, ΔG_gas = gas-phase energy, ΔG_solv = solvation free energy, E_SA = solvation energy, E_GB = Gibb’s energy.

Table 5. Summary of MM/PBSA-based computed thermodynamics binding free energy obtained for danuglipron and the proposed GLP-1R agonists.

Energy (kcal/mol)	Apo-A1236	Apo-A1238	Apo-A1245	Apo-A1251	Apo-B1245	Apo-PC	Apo-82L
ΔEvdW[a]	-63.42 ±3.19	–57.43 ±3.37	–57.72 ±2.57	–63.8 ±3.07	–60.98 ±2.93	–66.27 ±3.14	–48.3 ±5.89
ΔEelec	–8.01 ±3.42	1.03 ±2.55	–4.88 ±1.73	–7.32 ±1.81	–6.48 ±3.55	–116.02 ±7.31	–234.95 ±8.03
ΔGgas	–71.43 ±5.04	–56.39 ±4.88	–62.61 ±2.94	–71.21 ±3.77	–67.46 ±4.54	–182.30 ±7.40	–283.32 ±8.34
EGB	23.38 ±2.66	9.77 ±2.75	13.64 ±1.44	21.55 ±1.99	17.36 ±2.44	134.41 ±7.40	250.31 ±6.87
ESA	–7.55 ±0.27	–6.84 ±0.36	–6.94 ±0.23	–7.33 ±0.30	–7.13 ±0.26	–8.59 ±0.30	–6.28 ±0.76
ΔGsolv	15.82 ±2.58	2.92 ±2.64	6.71 ±1.52	14.22 ±2.05	10.23 ±2.35	125.82 ±7.38	244.02 ±36.88
ΔGbind	–55.61 ±4.11	–53.47 ±4.16	–55.90 ±3.59	–56.99 ±3.53	–57.23 ±3.44	–56.47 ±4.06	–39.30 ±6.34

ΔE_vdw = van der Waals, ΔE_elec = electrostatic, ΔG_bind = calculated total free binding energy, ΔG_gas = gas-phase energy, ΔG_solv = solvation free energy, E_SA = solvation energy, E_GB = Gibb’s energy.

Figure 2. (A) A1236, A1238, A1245, A1251, B1245 and Danuglipron (PC) bound to the GLP-1R binding site, (B) conformational stability, C-α atoms RMSD of Apo and the bound systems, (C) flexibility, C-α atoms RMSF of apo and the bound systems, ran for 200 ns MD simulation.

Figure 2. (A) A1236, A1238, A1245, A1251, B1245 and Danuglipron (PC) bound to the GLP-1R binding site, (B) conformational stability, C-α atoms RMSD of Apo and the bound systems, (C) flexibility, C-α atoms RMSF of apo and the bound systems, ran for 200 ns MD simulation.

Figure 3. Inhibitor-induced Conformational changes to GLP-1R protein when inhibitor bound. (A) The Apo form of the protein at 40 ns shows the active site residues (B) RoG of Apo and A1236-, A1238-, A1245-, A1251-, B1245- and danuglipron- bound systems (C) Principal component analyses of Apo and A1236-, A1238-, A1245-, A1251-, B1245- and danuglipron-bound systems.

Figure 3. Inhibitor-induced Conformational changes to GLP-1R protein when inhibitor bound. (A) The Apo form of the protein at 40 ns shows the active site residues (B) RoG of Apo and A1236-, A1238-, A1245-, A1251-, B1245- and danuglipron- bound systems (C) Principal component analyses of Apo and A1236-, A1238-, A1245-, A1251-, B1245- and danuglipron-bound systems.

3.3.6. Per-Residue Energy Decomposition (PRED)

Some compounds have great stability binding to GLP-1R throughout the 200ns simulation time, according to the RMSD and RMSF studies. Additionally, a high BFE indicated favorable binding conformations and interactions with the protein. Accordingly, PRED analysis was performed in order to assess the energy contribution of the active site residues to the compounds' binding and stability at GLP-1R binding pocket [35]. The MMPBSA method on the AMBER18 GPU is used to determine the energy contribution of the essential residues to the binding free energy [32]. The compounds interacted with some essential amino acids of GLP-1R (Figures 4A–9A). Through van der Waals and electrostatic forces, the majority of the residues contribute significantly to the binding (Figures 4B–9B). The amino acids that contributed binding energies of more than -1.0 kcal/mol are considered to produce the most favorable interactions in the complex. The GLP1R residues that interacted with A1236 include SER2 (-1.4 kcal/mol), LEU3 (-2.4 kcal/mol), TRP4 (-1.2 kcal/mol), PRO102 (-0.6 kcal/mol), LEU106 (-1.1 kcal/mol), LYS162 (-1.8 kcal/mol), TRP168 (-2.1 kcal/mol), THR172 (-1.1 kcal/mol), PHE195 (-1.1 kcal/mol), and PHE342 (-1.7 kcal/mol). All the residues contributed high energy towards the binding of A1236 except PRO102 with low energy. However, the contribution of van der Waals interactions in A1236-GLP1R complex is predominantly from TRP4 (-2.6 kcal/mol), TRP168 (-4.0 kcal/mol), and PHE342 (-2.6 kcal/mol), being the only ones above -2.0 kcal/mol. In A1238-GLP1R bound system, the interacting residues are TYR110 (-1.4 kcal/mol), VAL159 (-1.1 kcal/mol), LYS162 (-1.9 kcal/mol), TRP168 (-1.3 kcal/mol), VAL194 (-1.0 kcal/mol), PHE195 (-1.0 kcal/mol), PHE342 (-1.3 kcal/mol), LEU345 (-1.2 kcal/mol), and LEU349 (-0.5 kcal/mol). Based on this, VAL194, PHE195, and LEU349 showed to have unfavorable energies towards the binding. However, favorable van der Waals interactions are observed in TYR110 (-2.2 kcal/mol), TRP168 (-2.3 kcal/mol), and PHE342 (-2.1 kcal/mol) for the complex. In A1245-GLP1R bound system, LEU106 (-0.6 kcal/mol), TRP168 (-1.2 kcal/mol), PHE195 (-1.1 kcal/mol), MET198 (-1.7 kcal/mol), GLN199 (-1.5 kcal/mol), ARG264 (-0.7 kcal/mol), PHE342 (-0.9 kcal/mol), LEU345 (-1.7 kcal/mol), and LEU346 (-0.9 kcal/mol) are the interacting amino acids with their respective energy contributions. From this outcome, the interactions of LEU106, ARG264, PHE342, and LEU346 are considered to be low. Also in this complex, only TRP168 (-2.5 kcal/mol) is considered to produce a favorable van der Waals interaction. Although, GLN199 (-2.1 kcal/mol) provided a complimentary electrostatic interaction to the binding at the active site in this complex.

The amino acids of GLP1R that interacted with A1251 include LEU3 (-1.1 kcal/mol), LEU106 (-1.1 kcal/mol), LEU107 (-0.1 kcal/mol), TRP168 (-1.9 kcal/mol), LEU182 (-1.1 kcal/mol), PHE195 (-1.6 kcal/mol), MET198 (-1.0 kcal/mol), PHE342 (-2.4 kcal/mol), and PHE346 (-0.3 kcal/mol). Observing the interactions, LEU107, MET198, and PHE346 produced weaker energies leading to interactions which are unfavorable in the A1251-GLP1R complex. But the energies from the van der Waals interactions in TRP168 (-3.2 kcal/mol) and PHE342 (-3.0 kcal/mol) provided strong binding in the complex. In the B1245-GLP1R bound system, the interacting amino acid residues include LEU3 (-1.2 kcal/mol), TRP4 (-1.7 kcal/mol), LEU106 (-1.5 kcal/mol), LEU107 (-1.4 kcal/mol), TYR110 (-1.5 kcal/mol), LEU168 (-2.0 kcal/mol), TRP179 (0.0 kcal/mol), LEU182 (-2.1 kcal/mol), GLN186 (-1.3 kcal/mol), CYS191 (-1.0 kcal/mol), and PHE342 (-0.9 kcal/mol). Compound B1245 showed the highest number of interacting residues when compared to other compounds and the positive control, danuglipron. From the PRED outcome, most of the residues contributed high energy towards the binding, however, TRP179, CYS191, and PHE342 did not show favorable binding due to low energies. Based on its van der Waals interactions, they are predominantly from TRP4 (-3.4 kcal/mol), TYR110 (-2.5 kca/mol), and LEU182 (-2.1 kcal/mol) which additionally provided the complex with stability at the binding site. For the danuglipron complexed with GLP1R, the residues that interacted at the binding site are TRP4 (-1.4 kcal/mol), GLU103 (1.8 kcal/mol), LYS162 (-1.6 kcal/mol), ASP163 (1.6 kcal/mol), TRP168 (-4.0 kcal/mol), ALA174 (-1.4 kcal/mol), GLY181 (-0.5 kcal/mol), and PHE342 (-2.0 kcal/mol) with their individual energies contributed. All the residues contributed high energy towards the binding of danuglipron except GLY181 with low energy. However, the contribution of van der Waals interactions in danuglipron-GLP1R complex is principally from GLU103 (-3.1 kcal/mol), LYS162 (-2.2 kcal/mol), TRP168 (-4.1 kcal/mol), and PHE342 (-3.1 kcal/mol), which contributed to the strong binding observed. Also, danuglipron showed high electrostatic energy with GLU103 (–51 kcal/mol) and ASP163 (–57 kcal/mol) contributing significant hydrogen bonding interactions.

Elucidating the chemical basis for ligand-GLP-1R interactions, the interactions between the five pentacyclic triterpenoids and the GLP-1R is predominantly governed by hydrophobic forces, supplemented in some cases by electrostatic and hydrogen bonding interactions (Figures 4A–8A). This aligns with the chemical nature of triterpenoids largely nonpolar frameworks with scattered polar functional groups. The binding pockets of GLP-1R for all triterpenoids are notably enriched with nonpolar and aromatic residues that favor hydrophobic packing involving aromatic residues including tryptophan (TRP168, TRP4) and Phenylalanine (PHE342, PHE195, PHE346) which played central roles across all the ligands. Their large π-systems support π–π stacking and hydrophobic π-alkyl interactions with the planar, bulky ring systems of the triterpenoids. The indole group of TRP is particularly capable of forming snug fits with large nonpolar surfaces, leading to strong van der Waals stabilization. In addition to the aromatic residues, the aliphatic residues such as leucine (LEU3, LEU345, LEU349, LEU107, LEU106, LEU182), Valine (VAL159), Methionine (MET198), and Proline (PRO102) contribute significantly. Their branched hydrocarbon side chains offer complementary nonpolar surfaces that snugly interact with the hydrophobic triterpenoid skeletons. Even Lysine (LYS162), though polar and basic, has a long aliphatic chain that contributes hydrophobically in its neutral region. Moderate hydrophobicity of residues like Tyrosine (TYR110) also promotes interactions, where both the aromatic ring and the hydroxyl group offer dual possibilities for nonpolar and polar contacts.

Though less dominant than hydrophobic effects, electrostatic interactions and hydrogen bonding also contribute to specificity and stabilization, depending on the ligand’s functional groups. Serine (SER2) and Glutamine (GLN186, GLN199) which are their polar side chains (hydroxyl and amide respectively) form potential hydrogen bonds with hydroxyl or carbonyl groups on the ligands. These residues exhibit negative electronic energies, indicating favorable dipole-dipole or H-bonding interactions. Methionine (MET198), while mostly hydrophobic due to its thioether group, the sulfur atom can participate in weak polar interactions, especially via induced dipoles or dipole–dipole interactions, evident in ligands like A1245 and B1244. Arginine (ARG246), though commonly associated with strong positive charges, its hydrophobic methylene backbone can engage in van der Waals contacts, while its terminal guanidinium group can contribute to electrostatics when matched with a complementary polar site. Lysine (LYS162) demonstrates unfavorable electrostatic interactions with some ligands (e.g., A1238) due to charge repulsion likely between its positive amine group and electron-deficient regions on the ligand (e.g., carbonyl oxygens), reducing binding affinity at this site. Cysteine (CYS191) offers weak polar interactions through its thiol group, which can form hydrogen bonds or polar contacts, albeit with lower contribution compared to polar amide or hydroxyl-bearing residues. Overall, across all five pentacyclic triterpenoid ligands, the dominant interaction mechanism is hydrophobic, driven by the ligands’ extensive nonpolar skeletons and the complementary aliphatic and aromatic residues in GLP-1R’s binding site. Where polar interactions occur, they are generally localized to specific residues Glutamine, Serine, and Methionine whose polar groups align well with the functional groups on the triterpenoids. Conversely, residues like Lysine and Arginine may contribute unfavorably or neutrally unless their side chains are positioned to minimize electrostatic conflict. This balance of hydrophobic and polar interactions with target receptors are often exhibited by triterpenoids and similar lipophilic natural products, facilitating specific binding [63].

The investigated terpenoids showed varying degrees of binding energy as calculated by the MMPBSA, but having the highest number of amino acid residues interacting with B1245 contributed to the highest stable system of GLP1R-B1245 (ΔG_bind = -57.23 kcal/mol) relative to other compounds though insignificant. Our findings showed that most binding pocket residues contributed immensely to the interactions through electrostatic and van der Waals forces. These molecular interactions, which are considered necessary for ligand-receptor conformational stability and agonistic potential of the compounds, are essentially found with the amino acids in N-terminal extracellular domain (ECD) of the protein. This domain has been reported to be crucial in agonist binding to GLP-1R [11,12,13] [58,64]. More importantly, a more recent work reported a cryogenic electron microscope structure of GLP-1R revealing a binding pocket requiring a primate-specific tryptophan 33 residue (same as TRP4 in our current study due to early 29 missing amino acids in the protein used) [14]. The binding of danuglipron with this residue was reported to increase insulin levels in primates. This discovery with TRP4 could be explored for further studies leading to oral small-molecule therapies that target GLP1R. Gratifyingly, our current study showed same interaction with TRP4 in GLP1R-danuglipron complex. Additionally, two among our compounds, alpha-sitosterol and lupeol showed similar interaction with TRP4 in this present work. However, this interaction is energetically more profound with lupeol relative to alpha-sitosterol and danuglipron as a result of its more favorable van der Waals force (Figure 8B).

Figure 4. (A) A 3D visualization of the binding interactions of A1236-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP1R with A1236.

Figure 4. (A) A 3D visualization of the binding interactions of A1236-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP1R with A1236.

Figure 5. (A) A 3D visualization of the binding interactions of A1238-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP1R with A1238.

Figure 5. (A) A 3D visualization of the binding interactions of A1238-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP1R with A1238.

Figure 6. (A) A 3D visualization of the binding interactions of A1245-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP1R with A1245.

Figure 6. (A) A 3D visualization of the binding interactions of A1245-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP1R with A1245.

Figure 7. (A) A 3D visualization of the binding interactions of the A1251-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP1R with A1251.

Figure 7. (A) A 3D visualization of the binding interactions of the A1251-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP1R with A1251.

Figure 8. (A) A 3D visualization of the binding interactions of the B1245-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP1R with B1245.

Figure 8. (A) A 3D visualization of the binding interactions of the B1245-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP1R with B1245.

Figure 9. (A) A 3D visualization of the binding interactions of PC-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP-1R with PC.

Figure 9. (A) A 3D visualization of the binding interactions of PC-GLP1R bound system, (B) Individual energy contributions (PRED) of the interacting residues of GLP-1R with PC.

According to reports, the nineteen (19) medicinal herbs used in this study can be used to treat diabetes mellitus in traditional medicine [65]. This implies that the phytochemicals obtained from any of these plants may include a lead component or compounds that could be used to treat diabetes. In this study, 244 phytochemicals (found in 19 plants) were screened against GLP-1R. The five best-performing compounds were produced by Hymenocardia acida, Spondias mombin, and Alstonia boonei. The potential of Hymenocardia acida in the treatment of diabetes mellitus has been emphasized by many pharmacological investigations. As evidenced by its capacity to lower blood glucose levels and enhance insulin secretion and sensitivity, the plant extract has shown notable antidiabetic effects in a variety of in vitro and in vivo models [66]. The phyto-compounds that demonstrated the most promising effects were lupeol (B1245) and Friedelan-3-one (A1238) from Hymenocardia acida, α-sitosterol (A1236) from Spondias mombin, and beta amyrin (A1245) and alpha amyrin (A1251) from Alstonia boonei. These compounds’ ability to stimulate GLP-1R activities based on their strong binding energies and special molecular interactions in the extracellular and transmembrane domains could open up future studies of terpenoids to provide insulinotropic effects via GLP-1R agonism. Because they have fewer side effects, research into natural herbal phytoconstituents for the treatment of diabetes has gained more attention. One family of phytoconstituents that shows great promise for the treatment of diabetes mellitus is terpenoids. The improvement of insulin sensitivity is one of the main ways that terpenoids mitigate the effects of diabetes. By increasing insulin sensitivity at the target locations, secondary metabolites such as terpenoids can counteract the molecular pathways that lead to the development of insulin resistance [67]. According to earlier research, terpenoids increase insulin sensitivity, glucose uptake, and glycogen synthesis by activating the AMPK pathway, PPAR-γ, incretins, iNOS, and PPAR-γ, and by inhibiting the expression of α-glucosidase, plasma lipase, aldose reductase, tyrosine phosphatase, α-amylase, ghrelin, SGLT1, GLUT2, DPP4, DGAT, HSD, TXINP, and PTB-1B [68,69,70]. Nevertheless, as far as we are aware, no research has been done on terpenoids with GLP-1R agonism characteristics. This observation could prompt further experimental works involving this class of plant secondary metabolite to validate their prospect as GLP-1R agonists. The use of medicinal plants to treat diabetes mellitus incorporates modes of action such as β-cell regeneration and improved pancreatic insulin secretion, which are accomplished through GLP-1 regulating activities [71]. Although a number of studies have proposed these terpenoids' mode of action in the treatment of diabetes mellitus, none of them have so far been identified as GLP-1 receptor agonists.

Due to its capacity to regulate blood glucose levels, enhance insulin sensitivity, and lower the risk of complications from diabetes, sitosterol may have potential as an antidiabetic medication, according to research. These studies include blood glucose management in diabetic animal models, where alpha-sitosterol administration reduced blood glucose and enhanced insulin sensitivity [72]. Research indicates that by influencing particular signaling pathways in pancreatic beta-cells, it may promote insulin release by improving beta-cell activity [73]. Numerous pathways are involved in the mechanisms behind the antidiabetic effects of sitosterol, including the activation of important signaling molecules like AMP-activated protein kinase (AMPK) and Peroxisome Proliferator-Activated Receptor Gamma (PPAR-γ), which are crucial for insulin sensitivity and cellular energy regulation [74]. However, its agonism of the GLP-1 receptor, however, has not yet been studied. According to Table 5, the BFE of α-sitosterol bound to GLP-1R (-55.61 kcal/mol) indicates a negligibly lower energy of binding than that of danuglipron (-56.47 kcal/mol), a known agonist bound to GLP-1R. For A1236 to bind strongly in the receptor's active site, it must have hydrogen bonds and hydrophobic interactions. A pi-alkyl contact with TRP168 where the Van der Waal energy is higher than -3.5 kcal/mol; a stable hydrogen bond with SER2 and THR172; an alkyl hydrophobic (pi-cation) interaction with PHE342; and a characteristic hydrophobic interaction with TRP4 which is also found in danuglipron. Friedelan-3-one has been shown to significantly lower blood glucose levels in previous studies using experimental animal models of diabetes, including those induced by streptozotocin (STZ) [75]. This suggests that people with type 2 diabetes may benefit from its ability to lower blood sugar levels. Friedelan-3-one also improves insulin sensitivity in diabetic rats, according to research, most likely via changing crucial molecular pathways including the PI3K/Akt signaling pathway. Friedelan-3-one also improves insulin sensitivity in diabetic rats, according to research, most likely via changing crucial molecular pathways including the PI3K/Akt signaling pathway. This pathway is crucial for insulin signaling, and enhancing it may reduce insulin resistance, a hallmark of type 2 diabetes, according to [76]. Friedelan-3-one may help reduce the systemic inflammation associated with diabetes because of its anti-inflammatory qualities. AMP-activated protein kinase (AMPK), a crucial regulator of lipid and glucose metabolism, may be triggered by Friedelan-3-one, one of the main ways it may have antidiabetic benefits. AMPK activation enhances glucose uptake and utilization, which aids in blood glucose regulation [77]. Friedelan-3-one scavenges free radicals and boosts the activity of natural antioxidant enzymes to lessen oxidative damage to tissues, especially pancreatic cells [78]. Alpha amyrin, a naturally occurring triterpenoid found in a variety of plants, has garnered attention for its potential antidiabetic properties, which are attributed to its ability to regulate blood glucose levels, shield the pancreatic beta-cells from harm, and improve insulin sensitivity [79,80]. By inhibiting the α-glucosidase enzyme, which is involved in the intestinal breakdown of complex carbs into simple sugars, alpha amyrin may also have antidiabetic effects [81]. Alpha amyrin has been shown to control lipid metabolism, which is often broken down in diabetic conditions. By increasing HDL cholesterol and decreasing LDL, triglycerides, and total cholesterol, it alters lipid profiles. De Melo et al., [73] claim that this effect concurrently improves insulin sensitivity and reduces the risk of cardiovascular issues in diabetics. Beta-amyrin activates AMP-activated protein kinase (AMPK), a crucial regulator of the energy balance within cells. Because AMPK activation enhances insulin sensitivity and muscle cell uptake of glucose, it lowers blood sugar. In diabetic rats, beta-amyrin treatment has been shown to significantly reduce fasting blood glucose levels [79]. By altering significant signaling pathways, such as the PI3K/Akt pathway, which is essential for insulin signaling, beta-amyrin has also been reported to increase insulin sensitivity. Enhancing the uptake and utilization of glucose reduces insulin resistance in organs such as the liver and muscles, as well as hyperglycemia [82]. Beta-amyrin inhibits the enzyme alpha-glucosidase, which converts carbs into glucose in the intestines. Thengyai et al., [81] claim that by delaying the absorption of glucose, this inhibition lowers the rise in blood sugar that occurs after meals. Furthermore, by raising the amounts of natural antioxidant enzymes like catalase and superoxide dismutase (SOD), beta-amyrin decreases blood sugar levels. These enzymes help reduce oxidative stress, which exacerbates diabetes problems like neuropathy and retinopathy [83]. Lupeol's antidiabetic potential is demonstrated by its capacity to control insulin secretion, improve insulin sensitivity, and lessen inflammation and oxidative stress. It is a promising option for diabetes management because of its anti-inflammatory properties and impacts on important signaling pathways like AMPK and PPAR-γ. In diabetic rats, lupeol therapy raised insulin levels via enhancing the activity of pancreatic islet cells, based on a study by [84]. According to this action, lupeol may increase insulin secretion by boosting the activation of particular signaling pathways related to insulin production and glucose metabolism. Additionally, research indicates that lupeol may have an antidiabetic impact via activating the AMPK pathway, which improves insulin sensitivity and glucose uptake—two critical factors in the management of diabetes [85]. Activation of peroxisome proliferator-activated receptor gamma (PPAR-γ) is another significant way via which lupeol may signal insulinotropic effects. According to research by Giacoman-Martínez et al., [80] lupeol improves insulin sensitivity and glucose metabolism via activating PPAR-γ. Reactive oxygen species (ROS) and chronic inflammation are recognized to play a role in the development of type 2 diabetes and insulin resistance [84,86]. Along with pi-hydrophobic bonds with LEU166 and LYS 162, as well as alkyl-hydrophobic interactions with LEU3, LEU182, LEU107, LEU106, TYR110, and TRP179 amino residues, a strong van der Waal bond was identified with TRP4 for the lupeol (B1245). The conventional H-bonds with GLN 186 and CYS 191 are among the other noteworthy interactions. The positive control of this experiment, Danuglipron, was shown to interact with ASP 163 and GLU103 through hydrogen bonding; it created a pi-pi staking hydrophobic link with the residue of TRP4 (33) and a pi-alkyl hydrophobic bond with ALA74 and PHE342 (381). The binding energy of MM/GBSA is ΔG = -56.47 kcal/mol. The main structural components of danuglipron are the piperidine ring, benzyl ether, fluoro-pyrimidine, and benzimidazole. An ideal small molecule agonist is a fused aromatic heterocyclic core with pi-pi stacking interactions with GLP-1R residues; carboxylic acid groups provide improved activity and PK properties; aromatic heterocycles can form π - π stacking interactions with TRP4, which stabilizes the ECD binding conformation; and linkers facilitate conformations of small molecules like aliphatic rings and alkanes [10].

Considering that large peptides like GLP-1 interact with the extracellular and transmembrane domains of the GLP-1R [58], and obviously seeing that our investigated compounds showed molecular interactions with many amino acid residues located in these domains, it is plausible to suggest these terpenoids as agonists of GLP-1R. Essentially, TRP4 as a critical amino acid promotes cAMP production in GLP-1R-expressing cells, thereby aiding the pharmacological action of small-molecule GLP-1R agonists [13]. Based on the strong binding of lupeol (B1245) in the active site of GLP-1R and interaction with TRP4 similar to danuglipron, its molecular pharmacology may require TRP4 in the GLP-1R for signaling, hence prompting cAMP accumulation in GLP-1R-expressing cells leading to insulinotropic effect. Our findings establish a strong theoretical basis and can directly guide experimental efforts; essentially a focus on triterpenes possibility of agonism potential against GLP-1R, which has a paucity of information in literature. Specifically, the predicted binding poses and energetics highlight key residues and interactions that can be prioritized for future in vitro and in vivo validation.

The TRP4, ECD, ECL1, TM1 and TM6 are reported to be functionally important to the activation of GLP-1R [11]. The interaction of TRP4 of the ECD with ECL1, TM1 and TM6 buried Lupeol and stabilized it in the binding pocket (Figure 10). Lupeol (B1245) could be a promising agonist of GLP-1R as it interacts with TRP4(33) similarly with Danuglipron via van der Waal interaction of –3.4 kcal/mol and a MM/GBSA shows binding energy of ΔG = -57.23kcal/mol lower than that of PC recorded to be ΔG = -56.47 kcal/mol. Though compound B1245 is not structurally like Danuglipron, it interacts with significant residues at the binding pocket on GLP-1R and can activate it with the resultant effect of cyclic AMP production, and consequent GLP-1 secretion. Remarkably, a non-peptide agonist in a prior study demonstrated an aromatic interaction with TYR110, which is similar to the way B1245 interacts with GLP1R when it binds agonistically [12]. B1245 interacts with TM1 residues like LEU106 and ECL1 residues like CYS191 and GLN186 that form strong hydrogen bonds; a similar interaction was previously documented for LY3502970, another non-peptide agonist [12]. This may provide further insight into the process by which B1245 activates GLP-1R. Lupeol or its analogues may be able to be developed as novel GLP-1R agonists due to all of these interactions with GLP-1R and conformations that are ready to stabilize its binding pocket.

4. Conclusions

Our study involved the screening of different classes of phytochemicals against glucagon-like peptide-1 receptor (GLP-1R). Among these compounds, lupeol possesses the most promising effect with binding free energy higher than the standard compound, danuglipron. The lupeol showed conformational stability in the active site of GLP-1R, with no significant difference from the RMSD of danuglipron. GLP-1 receptor agonism potential of lupeol is demonstrated in this study by its interaction with amino acid residue TRP4, extracellular and transmembrane domains. These are characteristic interactions observable in many other known agonists of GLP-1R. This current work has shown the potential of lupeol and other terpenoids as GLP-1R agonist. Therefore, this study has opened up the opportunity for further experimental works to validate the GLP-1R agonistic properties of these and more terpenoids for possible antidiabetic drug design.