Natural products have continued to be used as sources of potential anticancer agents over the years. Almost, half of the approved anticancer agents are from natural products or derivatives exhibiting broad pharmacological activities  . Emodin is an anthraquinone, is naturally distributed and well known to possess broad pharmacological activities  . Emodin (D1) has been reported to inhibit the proliferation of many cancer cells which include; Colon cancer   , breast cancer   , gall bladder cancer   , pancreatic cancer  , lung cancer   and human cervical cancer  . Poor bioavailability and toxicity in vivo have been documented to limit emodin as cancer chemotherapy  . Natural emodin glycoside derivatives have been reported to possess high antitumor activities than emodin  . It has been shown that modification of the glycoside by addition of sugar chain at carbon-3 hydroxide (C3-OH) increases solubility and antitumor activity of emodin  .
Recently, Xing and co-workers  , reported a novel derivatives of emodin (Table 1) which strongly inhibited anti-proliferative activities on different human cancer cell lines  . Thus, possession of broad anti-proliferative activities makes them become potential anticancer agents  . Despite of the fact that the anticancer activity of these derivatives on various cancer cells has been reported, their biological targets, molecular mechanisms of action and pharmacokinetic profile remain uncovered. Recently, the use of in silico approach has been found to help addressing the biological targets and molecular mechanism of action of small molecules with protein or enzymes. In this regard, molecular docking, modeling and pharmacokinetic (ADMET) studies were carried out on emodin derivatives specifically to three key enzymes involved in cancer formation namely; N-Acetyltransferase 2 (NAT2), Cyclooxygenase 2 (COX2) and Topoisomerase 1 (TOP1).
N-Acetyltransferase is an important enzyme known to catalyse the transfer of acetyl groups from acetyl CoA to Arylamines  . N-Acetyltransferase exist in two isozymes; NAT1 and NAT2. The two forms of enzymes are polymorphic and catalyse both O-acetylation (activation) and N-acetylation (deactivation) of aromatic and heterocyclic amine carcinogen  . The later isoform, Arylamine N-Acetyltransferase 2 (NAT2) is known to metabolize arylamine and hydrazine moieties present in several chemicals, carcinogens and therapeutic drugs  . The polymorphic gene encoding of NAT2 often results to slow or rapid acetylator phenotypes  , which in-turn, causes drug-induced toxicities as well as risks of developing cancers of colon, bladder and lung  . Studies have shown a positive correlation of NAT2and rapid acetylator phenotype with colon cancer  , while slow NAT1 has been positively correlated with urinary bladder cancer  . Heterocyclic aromatic amines (HAAs) are family of mutagenic compounds produced in meat subject to high temperature cooking. Consumption of HAAs is highly associated with the risk of colon, breast, lung, skin, liver cancers   . Dietary of HAAs are transformed by the polymorphic NAT2 enzymes to carcinogen   . Dietary NAT2 enzyme is known to express in genotypes dependent in colon epithelium. Increased NAT2 enzyme activity activates HAAs in colon, and thus, increasing the risk of colon cancer. Inhibition of the activity of NAT2 enzyme play an important role in colon cancer prevention and treatment.
Cyclooxygenase (COX) enzymes exist in two isoforms, cyclooxygenase 1
Table 1. Binding energies (kcal/mol).
(COX1) and cyclooxygenase 2 (COX2). COX is known to catalyze the conversion of arachidonic acid to prostaglandins which play a vital role in the proliferation of cancer cells  . COX2 overexpress in human breast tumor cells and is positively correlated to the development of breast cancer  . Consequently, its inhibition is essential towards treatment of breast cancer.
Human topoisomerase 1 (TOP1) enzyme is an important drug target in cancer treatment  . The enzyme is responsible for catalysing the breaking and re-joining of phosphodiester of the DNA strand throughout cell cycle. It relaxes super coiled DNA during DNA replication and transcription   . TOP1 inhibitors work by blocking the ligation step of the cell cycle by generating single and double strands breaks which in turn harm the integrity of the genome; the break results to apoptosis cell death  . Apparently, TOP1 inhibitors play a great role towards cancer treatment.
Having showed broad anticancer activities on various cancer cells  , emodin derivatives (Table 1) was therefore investigated in silico to reveal whether NAT2, COX2 and TOP1 enzymes are the molecular targets they exert their action. Furthermore, their pharmacokinetic profiles are predicated to unveil their ADMET properties.
2. Materials and Methods
2.1. Proteins and Ligands Accession
The three dimension (3D) crystal structure of NAT2 complexed with Coenzyme A, COX2 and TOP1 enzymes that were used for evaluation were obtained from RCSB protein data bank (PDB ID: 2PFR, PDB ID: 3NTG, PDB ID: 1A35) respectively,    . Emodin derivatives were searched from the literature  .
2.2. Ligands and Protein Preparation
Ligand and Protein used in this study were prepared as previously reported   . Ligand was generated using Chem3D Pro 12, geometry were optimized and energy minimized by using Vega ZZ and ArgusLab using AM1 force field   to obtain coordinates with minimum energy and stable conformation. The optimized ligand was served in .pdb file format. The binding sites of NAT2, COX2 and TOP1 were analyzed by using 3D Ligand Site predication server  . The three proteins were separately prepared and energy minimized using Vega ZZ.
2.3. Docking Procedures
Docking experiments were done as reported by   . Firstly, for NAT2 enzyme, docking experiment was validated by redocking the acetyl CoA. The acetyl CoA was extracted from the binding site and redocked again; the acetyl showed to bind in a similar pocket interacting with amino acid residue as before it was docked which validated the docking method. Docking was done by using PyRx- virtual screening tool, with AutoDockVina docking option based on scoring functions  . The energy interaction of NAT2, COX2 and TOP1 with ligands was assigned as grind point. A blind docking with flexible conformation was allowed to allow ligands to search the binding site. The parameters were set as default, except for energy of interaction between the derivatives and NAT2, COX2 and TOP1 which were evaluated using atomic affinity potentials computed on a grid.
2.4. ADMET Analysis
The pharmacokinetics (ADMET) were predicted using ADMET predictor ver. 8.0  , The predicated ADMET properties are; Blood-Brain Barrier (BBB), Human Intestinal Absorption (HIA), Caco-2 cell permeability, Pgp inhibition, CYP 450 inhibitor and substrate (CYP1A2, 2C9, 2D6, and 3A4) and hERG inhibition.
2.5. Physicochemical Descriptors and SAR Correlation Analysis
Various physicochemical descriptors (Table 3) were correlated to their structural activities as NAT2 inhibitors using ADMET Modeler ver. 8.0  and SPSS ver.16. The observed and predicated activities (Table 4) were further used to study correlated analysis of the activity. The two sets; training set, a set of data used to discover potentially predictive relationship, and a test set which was used to assess the strength and utility of a predictive relationship were used to compute the observed and predicated activities using ADMET modeler ver. 8.0. The correlation analysis was computed using SPSS ver. 16.
3. Results and Discussion
3.1. Molecular Docking Studies
Drug design and development are costly and time consuming. Computer Aided Drug Design (CADD) provides a valuable alternative to time and cost in designing and developing drugs  . Emodin derivatives were docked to NAT2, COX2 and TOP1 enzymes to evaluate their inhibition. Docking results of emodin derivatives and standard drugs against the three enzymes are presented in Table 1. These derivatives were selected for docking studies due to their broad anticancer activities they have shown  . Possession of broad activities motivated us to further investigate their interaction with selected enzymes (NAT2, COX2 and TOP1) as target for cancer chemotherapeutics.
A docking study was done by using PyRx-virtual screening tool. The lowest docking pose energy with lower root mean square deviation was selected as the docking score with best protein affinity (Table 1). Docking study of emodin derivatives to target enzymes was done along with standard drugs of the respective enzyme inhibitor. Docking of emodin derivatives to NAT2 enzyme showed interesting results. Derivatives D2, D5, D6 and D9 strongly inhibited NAT2 enzyme with binding affinity of −10.6, −10.2, −10.2 and −9.7 kcal/mol, respectively. Derivative D9 was in agreement to reported experimental results  which strongly inhibited various cancer cells, NAT2 could be the target enzyme of D9 which could work in a similar manner as described by  . The other derivatives D2, D5 and D6 probably inhibit NAT2 in a different way as D9 work. When compared to emodin, derivatives D2, D5, D6 and D9 inhibited strongly NAT2 than emodin which had −9.5 kcal/mol. Furthermore, when compared to known NAT2 inhibitor (capecitabine), derivatives D2, D5, D6 and D9 strongly inhibited NAT2 enzymes than capecitabine which had −8.3 kcal/mol (Table 2). Docking of emodin derivatives to COX2 enzymes indicated that, emodin had high inhibition than all derivatives. Derivatives D2 (−8.3 kcal/mol), D4 (−8.2 kcal/mo), D5 (−8.2 kcal/mol) and D6 (−8.1 kcal/mol) indicated lower inhibition to COX2 enzyme compared to emodin (−8.8 kcal/mol). Further, all derivatives had lower inhibition when compared to celecoxib a known COX2 potent inhibitor (Table 1). Docking results for these derivatives suggests that COX2 is not their target enzyme. Emodin derivative did not strongly inhibit TOP1 enzyme when compared to camptothecin a known TOP1 inhibitor. Docking results further suggested that, derivative D5, D6 and D8 can inhibit TOP1 enzyme. Experimental results  have shown that D8 and D9 strongly inhibit various human cancer cell lines in vitro. In the present study, molecular docking results are agreement to experimental finding  . It can further be suggested that D8 is an
Table 2. Pharmacokinetics (ADMET) profile of selected emodin derivatives.
*The site of metabolism is shown in Figure 4, **Properties were predicated using Pre-ADMET server. aLipinski rule of five used to evaluate lipophilicity, it requires molecular a to have MlogP < 5, bpredicated Madin-Darby canine kidney (MDCK) cell permeability recommended range (<25 poor and >500 great), cBlood brain barrier recommended range (−3 to 1). dPredicated Human intestinal absorption recommended range. ePlasma protein binding, fPermeability to carcinoma cell recommended range (<5 low and >100 high), gPredicated inhibition concentration to human ether a-go-go-related gene, the potential risk for inhibitors ranges 5.5 - 6.
inhibitor of TOP1 enzyme, while D9 is an inhibitor of NAT2 enzyme. Docking analysis of all emodin derivatives to three key cancer chemotherapeutic target suggest that, NAT2 is the target enzyme for these derivatives with D2, D5, D6 and D9 having better binding energy (Table 1). Thus, derivatives D2, D5, D6, D8 and D9 were further investigated for their pharmacokinetic properties.
Emodin derivatives D2, D5, D6, D8 and D9 were further analysed for their interaction with NAT2 and TOP1 enzymes, as they strongly inhibited these enzymes compared to COX2 enzyme. Interaction analyses are presented in Figures 1(a)-(i). Interaction of D2 with NAT2 involved five hydrogen bonds (Figure 1(a), Figure 1(b)), the bond formed involved Ser287---OH (3.28 Å), the hydroxyl group (OH) which formed the hydrogen bond with amino acid Ser287 was from the ring B (Figure 2(a)). The other four hydrogen bond were contributed from the rhamnoside moiety which includes; Ser287---OH (3.16 Å), Thr289---OH (3.12 Å), Ser129---OH (3.02 Å) and Gly126---OH (3.05 Å). Interaction of D5 with NAT2 formed seven hydrogen bond (Figure 2(c), Figure 2(d)), all hydroxyl groups in the rings formed hydrogen with NAT2 (Figure 2(b)). The residues involved in forming hydrogen bonds includes; Ser216---OH (3.19 Å), Phe93---OH (3.34 Å), Ser125---OH (3.05 Å) of ring near sugar moiety (ring C), Thr289---OH (3.36 Å) and Gly125---OH (3.38 Å) both from ring B (Figure 2(b)), Gly125---OH (2.98 Å) and Ser129---OH (2.86 Å) both from ring A (Figure 2(b)). Derivative D6 formed five hydrogen bonds with NAT2, the amino acid residues involved in forming hydrogen bonds were; Lys188---OH (3.08 Å) from sugar moiety, Tyr94---O (2.87 Å) from ring B (Figure 2(c)), Gln163---O (3.10 Å), Ile---OH (3.07 Å) and Arg165---OMe (3.09 Å). Derivative D9 formed six hydrogen bonds with NAT2 enzyme, one hydrogen bond was from ring C (Figure 2(d)). Other hydrogen bonds were from the sugar moiety which contributed to its increase in activity. The interaction of D9 with NAT2 by forming hydrogen bond with amino acid residue Try94 could further explain the increased activity of D9. The hydrogen bond interaction with amino acid residue was Tyr94---OH (3.36 Å) from ring C, Lys188---MeOMe (3.43 Å), Ile290---OH (2.99 Å), Thr289---OH (2.60 Å) (Figure 1(g), Figure 2(d)). The interaction of D8 with TOP1 involved several hydrogen bond, the sugar moiety were responsible in forming many hydrogen bonds (Figure 1(h) and Figure 1(i)). It was interesting to note that, like known TOP1 inhibitors such as camptothecin which interacts and forms hydrogen bonds with Arg364, Lys532 and Asn722  , D8 also interacted and formed hydrogen bonds with such amino acid residues. Other amino acid residue formed hydrogen bond with D8 are Asp440---OH (3.14 Å) from ring A, Lys443---O (2.86 Å) from ring B and Lys443---OH (3.04 Å) ring C. Derivative D8 could inhibit TOP1 enzymes via a similar mechanism of action as camptothecin.
Docking analysis further showed that D2, D5 D6 and D9 with different conformation bind to similar binding pocket at different conformation of NAT2 (Figure 3). Binding of ligands at the same pocket with different conformation may be due to changes of receptor specific conformation as explained by  . Binding to same pocket may further suggest similar mechanism of action for these derivatives.
3.2. Pharmacokinetic (ADMET) Analysis
Many drugs fail to enter into clinical markets due to poor pharmacokinetics. In the present study, we explored the pharmacokinetics of the selected compounds
(a) (b)(c) (d)(e) (f)(g)(h) (i)
Figure 1. (a)-(g) molecular interaction of D2, D5, D6 and D9 with NAT2 enzyme, (h)-(i) interaction of D8 with TOP1 as viewed by PyMol. Interaction was stabilized by presence of hydrogen bonds indicated by yellow dotted line.
Figure 2. 2D representations of the interaction of D2, D5, D6 and D9 with NAT2 shown in figure (a), (b), (c) and (d), respectively.
using ADMET PredictorTM. The pharmacokinetics parameters (Table 2) analyzed includes: MlogP, S + logP, S + logD, BBB, Caco-2, HIA, S + MDCK, PPB, Pgp inhibition, CYP 450 inhibitor and substrate (2D6 and 3A4), hERG inhibition, S + Sw, skin permeability and rat acute toxicity, CYP 3A4 intrinsic clearance, CYP 2D6 intrinsic clearance and kinetic parameters for CYP 2D6 and 3A4 enzymes. Lipophilicity defined as the ability of a chemical compound to dissolve in fats, oils, lipids and non-polar solvents, is a physicochemical property which affects drug transport through lipid structure and drug interaction and with the target protein  . Lipophilicity and water solubility was predicated using different models: Moriguchi model of octanol-water partition coefficient (MlogP), octanol-water partition coefficient (S + logP) and octanol-water distribution
Figure 3. Positions of derivatives in NAT2 active site after docking, (a) Superposing of D2, D6 and D9 showing similar position (b) D5 showing different binding position in the same pocket. Binding in different position could be explained by the acetyl chain on the sugar ring causing it to bind in different position, unlikely for D2, D6 and D9.
coefficient (S + logD) calculated from S + pKa and S + logP. Water solubility was predicated using native water solubility model (S + Sw, (mg/mL)).
All compounds showed lipophicility and aqueous water solubility to be in an acceptable range (Table 2). The Lipinski rule of five was further used to evaluate lipophilicity and water solubility, the rule requires that, for a compound to have good lipophilicity it should have no more than 5 MlogP value (MlogP ≤ 5). The predicated BBB filter indicated low for all derivatives suggesting that they may not permeate to the brain and thus not causing damage to the central nervous system, the computed log BB fell with the recommended range (−3 to 1) for log BB  . Oral bioavailability was predicated by using Madin-Darby canine kidney (MDCK) cell permeability, the recommended range for MDCK is <25 poor and >500 great. Results indicated that, except derivative D9 all other derivatives were indicated to have poor oral bioavailability (Table 2). Permeability to carcinoma cell (Caco-2) was studied (Table 2), all the derivatives indicated moderately permeability with values falling within the recommended range (<5 low and >100 high). The predicated plasma protein binding was below 90 (<90%) suggesting that even though they bind to protein, there is available fraction of molecules required for exerting therapeutic effects (Table 2). It is well acknowledged that, metabolism plays an important role in drug-drug interaction and bioavailability. Cytochrome CYP 450 is known to metabolize drugs during phase I metabolism  . Metabolism of the selected compounds was evaluated using the developed model for cytochrome P450 (CYP) site of metabolism, kinetic parameters (estimated Michael-Menten Km constant for predicated sites of metabolism (µM), estimated Michael-MentenVmax constant for predicated sites of metabolism (nmol/min/nmol enzyme) and estimated intrinsic clearance (Clint) for predicated sites of metabolism (μL/min/mg HLM protein) and inhibition. CYP 2D6 and CYP 3A4 are important isoform enzymes involved in drug metabolism, CYP 3A4 is known to metabolize half of the drugs. Thus, prediction of CYP inhibitor is of great importance in drug development. In this study, the kinetic parameters and inhibition of CYP 2D6 and CYP 3A4 was predicated using the model developed in ADMET PredictorTM (Table 2).
Results showed that, all compounds were non inhibitors of CYP 2D6 and non- substrate for all kinetic parameter (Table 2). Results suggest that, these compounds could be well metabolized by CYP 2D6, easy cleared from the body and thus causing little toxicity. However, it was further noted that, D8 and D9 inhibited CYP 3A4 enzyme (Table 2), while D2, D5, and D6 were nonsubstrate to CYP 3A4 Km, Vmax and Clint, suggesting to be cleared without causing significant toxicity. Further analysis showed that, derivatives D8 and D9 were predicted to be substrate to CYP 3A4 (Figure 4). Human toxicity was evaluated us- ing the model human ether-a-go-go related gene inhibition (hERG_pIC50). hERG blockers are known to prolong the QT interval and results to fatal cardiac arrhythmia known as Torsades de pointes. Thus, potential blockers of hERG potassium channel need to be investigated early during drug development and has recently been a major concern in pharmaceutical industries. The model filter provides Yes or NO to question “is the compound potential hERG inhibitor?”, the model further predicts hERG_pIC50 (Table 2). The hERG risk value for hERG_pIC50 starts at 5.5 to 6. The pIC50 for all compounds was below 5.5 (<5.5) indicating compounds to be non-blocker of the hERG, thus non-toxic. Generally, the predicated ADMET properties for derivatives suggest that these derivatives can be possible inhibitors of NAT2 and TOP1 enzymes with desirable pharmacokinetic properties.
3.3. SAR Correlation Analysis
An attempt to correlate the activity of compounds with their structure or property descriptors as NAT2 inhibitors was done using computational approaches. Computational approaches are nowadays used to classify compounds and correlate their activities using physicochemical descriptors. In the present study, different physicochemical descriptors (Table 3) were calculated using computational tools and online software and correlated to their activities pIC50 (Table 3). Predicated and observed values (Table 4) were also correlated which gave r2 = 0.56 upon removal of the outlier (D4 and D8) from the training set. The correlation results of the physicochemical descriptors with the pIC50 are presented in Table 5. The index of reflectance (IoR) correlated positively and significantly contributed to enhancement of activities (P < 0.05) however, the lipophilicity (logP) as well as molar reflectance (MR) negatively correlated to activities of the compounds.
In the present study, emodin derivatives were in silico valuated for their inhibitory activities and pharmacokinetics on NAT2, COX2 and TOP1 enzymes as agents for colon and other forms of cancer. Docking studies suggested that D8 to be a target inhibitor of TOP1 while D5, D6 and D9 targets inhibitors of NAT2 enzymes. Pharmacokinetics suggested that these compounds can be potential anticancer agents. Physicochemical parameter correlated to the compounds activities.
Figure 4. CYP 3A4-Site of metabolism for D2, D6, D8 and D9. Display of atomic properties for the first step of metabolic oxidation by a CYP 450. Results are generated for likely substrate (sites in red) and likely non-substrate (sites in grey).
Table 3. Physicochemical descriptors.
aIC50 in µg/mL of A549 cells  , data for this cell type were chosen because emodin derivatives showed better inhibition to this cell unlikely for other cancer cells. bDensity it related the size and bulk of the substituent, it is calculated by ACD/Lab, cSurface tension calculated by ACD/Lab, dIndex of reflectance calculated by ACD/Lab, eThelogP was calculated using ADMET PredictorTM. fMolar reflectance was calculated using ACD/Lab, gBalaban distance connectivity index of the hydrogen suppressed molecular graph, hBalaban mean square distance index of the hydrogen-suppressed molecular graph.
Table 4. Comparison of predicted and observed activities for the training and test set.
jResidual = Observed – Predicated.
Table 5. Correlation analysis showing correlation of physicochemical descriptors and inhibitory activities.
*Correlation is significant at the 0.05 level, **predicted IC50 values from the IC50  .
D.M.S. thanks Dr. Michael Lawless (Simulation Plus, UK) for giving an access to use an evaluation copy of ADMET PredictorTM for this study.