In silico approach: Prediction of ADMET, molecular docking, and QSPR of secondary metabolites in Mangroves

Mangrove plants are known to produce various secondary metabolites (SMs) such as polyisoprenoids (dolichol and polyprenol), emodin, and luteolin to show anticancer, anti-inflammatory, antiviral, and antibacterial activities. This study aimed to predict the multiple activities of the SMs of mangroves based on several enzymes that utilize the in silico method. The properties of absorption, distribution, metabolism, excretion, and toxicity for emodin, luteolin, polyprenol C80, dolichol-17 (C85), and dolichol-20 (C100) displayed variation. Prediction for Lipinski’s rule showed that emodin and luteolin have lower molecular weight than polyprenol C80, dolichol C85, and dolichol C100. Emodin and luteolin were further docked with cyclooxygenase-2, beta-lactamase, CYP450-dependent 14-alpha demethylase, 3C-like protease, and P-glycoprotein as protein targets using the Molegro Virtual Docker Ver.5.5 approach. In comparison to celecoxib, ketoconazole, clavulanic acid, remdesivir, and verapamil, emodin and luteolin had higher rerank scores. The quantitative structure–property relationships depicted that the electronic parameter, highest occupied molecular orbital ( E HOMO ), was the physical chemistry parameter that influenced the total clearance. The present findings emphasized that emodin and luteolin from the SMs of mangroves have multiple activities as potent inhibitors of cancer cells and bacterial infections.


INTRODUCTION
Currently, computational designs based on the structure of the target receptors responsible for a compound's toxicity and activity in the body can aid in the discovery of novel drugs. This method of computation is known as "structure-based drug design" or "computer-aided drug design" (CADD) (Song et al., 2009). The use of CAAD in finding new drug candidates provides information on the molecular properties of the active compound (Dara et al., 2021). Active compounds as ideal new drug candidates have bioavailability (absorption, distribution, metabolism, and elimination) and toxicological data. The value of oral drug absorption depends on the dose, solubility, and permeability according to Lipinski's rule of five. Lipinski's rule of five aims to determine the amount of a compound that is absorbed from the gastrointestinal tract through passive diffusion. The fulfillment of Lipinski's five rules for the active compound shows its solubility in water. Solubility in water is related to faster metabolic processes, elimination time, and side effects of drug accumulation. The appropriate activity test determines if the active compound is sufficiently soluble (Duchowicz et al., 2009).
Each active compound has physicochemical properties related to its structure. The predictive use of the chemical properties of the active compound structure has lipophilic, electronic, and steric parameters. In this study, the quantitative relationship between the structure and chemical properties was determined using the quantitative structure-property relationships (QSPR) model. In addition, the QSPR model can be used to explain chemical properties that affect pharmacokinetic values, especially elimination, which is related to the determination of drug dose values. The elimination parameter used in this study was the total clearance (CL TOT ) of the active compound (Praditapuspa et al., 2021). The QSPR model has the advantage of being more economical and reducing the need for time-consuming and expensive experiments (Duchowicz et al., 2009;Sandfort et al., 2020).
The active site of the protein that binds to the drug compound is identified using information from the target protein's structure. Compounds that are expected to bind to the target protein and have biological activity can be designed based on the prediction of the active site. The target protein's structure can be modeled using crystal structure data obtained from databases such as the Protein Data Bank (PDB) or by performing comparative modeling to predict the 3D structure of the unknown target protein based on its similarity (>40% similarity) with the known target protein if the protein structure is not found in the available database. The availability of data is critical to the success of drug discovery and development campaigns (Song et al., 2009). Drug development relies heavily on natural products. Natural products have become increasingly popular in pharmacological research during the last 20 years, including in antibacterial medication development. Terpenes, alkaloids, flavonoids, and phenols are only a few of the secondary metabolites (SMs) produced by typical plant metabolic processes that have antibacterial properties. By enhancing membrane permeability, inhibiting enzyme synthesis, or blocking biochemical reactions, natural compounds can address multiple bacterial targets at the same time (Guo et al., 2020).
Indonesia is known as the world's most important mangrove area with a diverse range of species . Several studies have reported, as presented in Table 1, that mangrove plants contain SMs, i.e., polyisoprenoids (for example, dolichol and polyprenol), anthraquinone (for example, emodin), and flavonoids (for example, luteolin). In addition, they provide properties that are anticancer, anti-inflammatory, antioxidant, antiviral, and antibacterial (Ho et al., 2007;Illian et al., 2018;2019May Zin et al., 2017;Momtazi-Borojeni et al., 2013;Odontuya et al., 2005;Qurrohman et al., 2020;Sumardi et al., 2018;Wu et al., 2008). The chemical structures of these compounds are shown in Figure 1. Due to the wide-ranging biological activities of polyisoprenoids, anthraquinone, and flavonoids in mangrove plants, these SMs are regarded as potential natural resources for pharmacological and other properties (Table 1). Thus, these observations in concert prompted us to implement the in silico method to predict the activities and value of the bond energy between emodin, luteolin, dolichol, and polyprenol on the enzymes of cyclooxygenase-2 (COX-2), betalactamase, CYP450-dependent 14-alpha demethylase, 3C-like protease (3CL-pro), and P-glycoprotein (P-GP).

Prediction for Lipinski's rule of five
According to "the rule of five" in the discovery setting, when there are more than 5 H-bond donors and 10 H-bond acceptors, the molecular weight surpasses 500, and when the computed LogP surpasses 5, it is more likely that there will be poor absorption or penetration (Lipinski, 2004). The physicochemical parameters were predicted by the ChemBio 3D Ultra ver.13 software.

Determination of QSPR
The determination of QSPR aimed to determine the quantitative relationship between the structure and physicochemical parameters (Shoombuatong et al., 2017). The physicochemical parameters included the lipophilic, electronic, and steric parameters. The lipophilic parameters included LogP and LogS. The electronic parameters included E TOT , E HOMO , and E LUMO . The steric parameters included molecular weight and molecular volume. The prediction of CL TOT value as a parameter of drug elimination was made by pkCSM. SPSS ver. 26 was used for the determination of the QSPR equation.

Preparation of ligands
ChemBioDraw Ultra ver. 13.0 (Cambridge Soft) was used to draw the structures of emodin and luteolin. ChemBio 3D 13.0 was used to convert the compound's 2D structures to 3D structures. The Merck molecular force field (MMFF94) method was used to optimize the molecule and minimize the geometry of the ligand, and the results were saved in SYBYL.mol2 format, which the Molegro Virtual Docker (MVD) software could read.
The pharmacokinetic data in Table 2 revealed that emodin and dolichol-20 (C100) do not inhibit CYP3A4 enzymes, suggesting that it is safe to use emodin and dolichol-20 (C100) with a wide range of CYP3A4 substrates (including alprazolam, atorvastatin, vincristine, halothane, hydrocortisone, zidovudine, carbamazepine, codeine, cortisol, caffeine, lidocaine, lovastatin, midazolam, nifedipine, paracetamol, tacrolimus, tamoxifen, testosterone, phenytoin, cyclosporine, cyclophosphamide, erythromycin, R-warfarin, and S-warfarin). This is due to the fact that the isoenzyme type CYP3A4 is the primary metabolizer of approximately 75% of all medicines metabolized by CYP (Zhang et al., 2021). According to a prior study examining the 200 most widely prescribed medications in the United States, the CYP isoenzyme type CYP3A4 metabolizes around 46% of pharmaceuticals (Sychev et al., 2018). Another pharmacokinetic predictor is toxicity. Based on the Ames test toxicity data, it is known that the five compounds are not mutagenic, skin-sensitizing, or hepatotoxic. Emodin and luteolin compounds are not inhibitors of hERG I or hERG II associated with potassium channels. Meanwhile, polyprenol (C80) and dolichol (C85 and C100) compounds are seen from the inhibition parameters of hERG (the human Eth-er-ago-go-related). These three compounds are known to inhibit hERG II. hERG is a gene that encodes a protein potassium channel that contributes to heart rate activity. Other parameters, such as the maximum tolerated dose, LD 50 , and LOAEL values, can be classified as good and having a low maximum tolerance dose. However, these five compounds are known to be toxic to Tetrahymena pyriformis (Fadlan et al., 2021).

Determination of QSPR
One of the seven regression equations with the best correlation coefficient and significance was chosen (see Table 5). The following was considered the best regression equation for the review (Widiyana et al., 2016): 1/CL TOT = 0.255 E HOMO + 3.110.
The equation was chosen since it has the best correlation coefficient (R = 0.985), the smallest significance of 0.002 (<0.05), the smallest standard error (SE = 0.071), and the largest F-value (95.555). The physical chemistry parameter that affected the CL TOT , according to the equation, was the electronic parameter E HOMO (as shown in Tables 5 and 6).

Molecular docking studies
Emodin and luteolin were further docked with betalactamase, CYP450-dependent 14-alpha demethylase, COX-2, 3CL-pro, and P-GP as protein targets using the MVD Ver. 5.5 docking simulation approach (Widiyana et al., 2021). Based on two criteria, such as ligand-binding position and fitness function score comparison, this algorithm determined the best docking score. The RMSD was used to determine the most optimal ligandbinding location. Table 7 presents the interactions of emodin and luteolin with COX-2, beta-lactamase, CYP450-dependent 14-alpha demethylase, 3CL-pro, and P-GP.
Molecular docking is a computer simulation of a ligand binding to a receptor that aids in predicting the protein target (in which the molecule will form a bond) and assessment of the binding affinity. Furthermore, the pose organizer can rotate the hydrogen molecules in both the receptor and the ligand to their ideal position. It is also possible to utilize it to rerank the ligands (using a rank score). The middle panel enables the recalculation of scoring functions, including the MolDock score, binding affinity score, and reranking score, while analyzing docking data. If the poses are imported from an MVD results file, these scoring function values have already been computed. Recalculating the scores for each of the three metrics is done by pressing the "recalculate energies" button (using the coefficients specified in the files for the binding affinity and reranking scores). Using a more computationally sophisticated scoring algorithm than the one used during docking, it is feasible to rerank poses. Numerous components (such as van der Waals forces, electrostatic interactions, and solvent terms) can be manually adjusted in this reranking score function. In MVD, the MolDock, MolDock GRID, and rerank scores are not expressed in chemically meaningful units. However, the MVD can provide an approximation estimate of the binding affinity (Piramanayagam and Lisina, 2014). In numerous types of activity, the proteins mentioned above are shown to be overexpressed. Since the target molecules' X-ray crystallographic structures are known, structure-based rational drug design methodologies can be used to build novel lead compounds based on them. The active site of the target compounds was discovered and molecular docking was performed using ChemBio 3D 13.0 containing 3D structures of chemicals obtained from marine systems. MVD handles every step of the docking procedure, from preparing the compounds to determining the target protein's probable binding sites and predicting the ligands' binding modes. MVD provides high-quality docking by combining a revolutionary optimization technique with a user interface that prioritizes usability and efficiency. The computational experiment yielded the identification of a 5-second metabolite that docked perfectly into the target's active site. For each ligand, a visual evaluation of the top-ranking molecules (hits) generated a list of roughly five small molecules. The docked postures were listed along with the rerank scores that correlate with them. The  MolDock score calculated for "hits" was significantly higher than the MolDock score calculated for the original ligand-protein docking. According to the energy principle, the lowest score was used since a molecule is most stable at the lowest energy level. As a result, it was worthwhile to test the actual binding affinities of these small molecules to the target protein to see if the computer predictions corresponded to the biological scenario (Widiyana, 2021;Zhang et al., 2021).
The receptor structure used in this study was obtained from the Research Collaboratory for Structural Bioinformatics (RCSB) PDB. Each protein provides information about the active site and the composition of the secondary metabolite-protein binding space as well as the opportunity to use the enzyme or protein in a functional conformation. In the current investigation, we tested five activities, including antibacterial activity with a betalactamase inhibitor mechanism, antifungal activity with a CYP450dependent 14-alpha demethylase inhibitor mechanism, antiinflammation activity with a COX-2 inhibitor mechanism, antiviral activity with a 3CL-pro inhibitor mechanism, and anticancer activity with a P-GP inhibitor mechanism (Puspaningtyas, 2013). The validation stage was carried out to calibrate the docking method in silico by repeating the prediction of the native ligand's binding ability to proteins and then comparing it with the position of the original native ligand, which is innate to the receptor on the ligand-binding site pocket. The assessment parameter used to provide an assessment of validity was the RMSD value. The RMSD value was close to zero, indicating that the position of the copy ligand after the superimposition was approaching that of the native ligand, implying that the method utilized would be more precise. Validation was carried out at the pocket ligand-binding site with 10 replications. To ensure that the orientation of the ligand matches the behavior in the ligand-receptor binding model and the position obtained from the docking study using the MVD program, the parameters of this docking method first needed to be validated by the structure of the receptor used, namely PDB, which was replicated 10 times. The rerank score is a value that represents the binding energy required to form a bond between the ligand and the receptor, and it predicts compound activity. It also improves the ligand-receptor interaction. Emodin and luteolin were docked as ligands with COX-2, beta-lactamase, CYP450-dependent 14-alpha demethylase, 3CL-pro, and P-GP as protein targets using the MVD Ver. 5.5 docking simulation approach. This algorithm determined the best docking score based on two criteria: ligandbinding position and fitness function score comparison. The lower value of the rerank indicated the strength of the ligand-receptor bond. The most optimal ligand-binding location was determined using the RMSD.
Docking analysis was performed for antibacterial activity using the mechanism of inhibition toward beta-lactamase using PDB code 1LLB, which contains the inhibitor crystal structure of AmpC beta-lactamase from Escherichia coli in a complex with ATMO-penicillin (https://www.rcsb.org/structure/1llb), and PDB code 6W33, which contains the crystal structure of Class A beta-lactamase from Bacillus cereus in a complex with the betalactamase inhibitor of clavulanate acid (https://www.rcsb.org/ structure/6W33). At this point, the test compounds' interactions with the 1LLB and 6W33 receptors were compared to interaction with clavulanic acid as the standard. Interactions with the 1LLB receptor had a higher rerank score of −74.04 than interaction with clavulanic acid (as presented in Table 7), and since emodin and luteolin had the lowest rerank scores (−85.78 and −88.47 kcal/ mol, respectively), this indicated that emodin and luteolin have higher binding energies than the others. Furthermore, interaction with the 6W33 receptor had a higher rerank score of −79.53 than interaction with clavulanic acid (see Table 7), and since luteolin had the lowest rerank score (−82.83 kcal/mol), luteolin revealed higher binding energies than the others. These results predicted that emodin and luteolin might have antibacterial activity and are related to overcoming beta-lactamase enzyme-mediated multidrug resistance through the mechanism of beta-lactamase inhibition.
Docking analysis was carried out for antifungal activity using the mechanism of inhibition toward CYP450-dependent 14-alpha demethylase (PDB code 1JIP), which contains ketoconazole-induced conformational changes in the active site of cytochrome mammalian P450 enzymes (https://www.rcsb.org/structure/1JIP). At this point, the test compounds' interactions with the 1JIP receptor were compared to interaction with ketoconazole as the standard. Interaction with the 1JIP receptor in comparison to interaction with ketoconazole had a higher rerank score of −135.08 (see Table 7), and since emodin and luteolin had the lowest rerank scores (−71.62 and −73.81 kcal/mol, respectively), emodin and luteolin exhibited lower binding energies than the others. These results predicted that emodin and luteolin might lack antifungal activity through the mechanism of CYP450 inhibition. Further investigation into the antifungal activity of emodin and luteolin using the other mechanisms is needed.
Docking analysis was carried out for anti-inflammation activity using the mechanism of inhibition toward COX-2 (PDB code 3LN1), which contains the structure of celecoxib bound at the COX-2 inhibitor active site (https://www.rcsb.org/ structure/3LN1). At this point, the test compounds' interactions with the 3LN1 receptor were compared to interaction with celecoxib as the standard. Interaction with the 3LN1 receptor had a higher rerank score of −137.80 than interaction with celecoxib (see Table 7), and since emodin and luteolin had the lowest rerank scores (−85.02 and −106.84 kcal/mol, respectively), emodin and luteolin exhibited lower binding energies than the others. These results predicted that emodin and luteolin might lack anti-inflammatory activity through the mechanism of COX-2 inhibition (Park et al., 2016). Further investigation into the anti-inflammatory activity of emodin and luteolin using the other mechanisms is needed.
Docking analysis was carried out for antiviral activity using the mechanism of inhibition toward the 3CL-pro inhibitor (PDB code 6LU7), which contains the crystal structure of COVID-19 main protease in a complex with an inhibitor N3 (https://www.rcsb.org/structure/6LU7). At this point, the test compounds' interactions with the 6LU7 receptor were compared to interaction with remdesivir as the standard. Interaction with the 6LU7 receptor in comparison to interaction with remdesivir had a higher rerank score of −128.35 (see Table 7), and since emodin and luteolin had the lowest rerank scores (−83.99 and −96.22 kcal/ mol, respectively), emodin and luteolin exhibited lower binding energies than the others. These results predicted that emodin and luteolin might lack antiviral activity against SARS-CoV-2 through the mechanism of 3CL-pro inhibition (Mandal et al., 2021). Further investigation into the antiviral activity of emodin and luteolin using the other mechanisms is needed.
Docking analysis was carried out for anticancer activity using the mechanism of inhibition toward P-GP using PDB code 1MV5, which contains the crystal structure of the LmrA ATPbinding domain (https://www.rcsb.org/structure/1mv5), and PDB code 4DGN, which contains the crystal structure of maize CK2 in a complex with the inhibitor luteolin (https://www.rcsb.org/ structure/4dgn). At this point, the test compounds' interactions with the 1MV5 and 4DGN receptors were compared to interaction with verapamil as the standard. Interaction with the 1MV5 receptor had a higher rerank score of −74.28 (see Table 7) than interaction with verapamil, and since emodin had the lowest rerank score (−77.94), emodin exhibited higher binding energies than the others. However, luteolin, which had a rerank score of −73.94, can show anticancer activity since it has a value similar to the rerank score of verapamil. Furthermore, interaction with the 4DGN receptor had a higher rerank score of −80.97 (see Table 7) than verapamil, and since emodin and luteolin had the lowest rerank scores (−92.56 and −101.57 kcal/mol, respectively), emodin and luteolin exhibited higher binding energies than the others. These results predicted that emodin and luteolin might show anticancer activity and are related to overcoming P-GP-mediated multidrug resistance through the mechanism of P-GP inhibition.

CONCLUSION
The present study suggested that emodin and luteolin have lower toxicity than the standard, while the standard has a higher affinity for protein enzymes than emodin and luteolin. The findings emphasized that these compounds may show multiple activities as potent inhibitors of cancer cells as well as bacterial infections. Therefore, for further investigation to ascertain the effectiveness against cancer cells and bacterial infections as well as fungal and viral infections, synthesis and in vitro evaluation are required.

AUTHORS' CONTRIBUTIONS
DNI, APW, and MB contributed to conceptualization; APW, ARH, DNI, HM, and MB contributed to methodology; ARH and APW used the software; APW, ARH, SSAM, and DNI contributed to validation; MB collected the resources; DNI contributed to original draft preparation; DNI, APW, ARH, HM, SSAM, and MB contributed to review and editing; MB supervised the study; MB contributed to funding acquisition. All authors have read and agreed to the published version of the manuscript.