In Silico Docking Analysis and Admet Prediction of Thymoquinone Derivatives Against Ovarian Cancer

Thymoquinone, the active constituent of Nigella sativa has been reported to have various biological activities. Due to its significance, various analogues of it have been synthesized and reported for anti-cancer activity. In the present research, we have taken the analogs of thymoquinone and performed docking study with an objective to find the binding pattern of all the molecules. Apart from this, pharmacokinetic parameters were predicted along with their toxicological parameters. From the results, the molecule Thy09 was found to have the optimized structure and further modification on this could lead to more potent compounds.


INTRODUCTION
Nature has produced various medicinal plants for treatment of various ailments and among those is the plant Nigella sativa 1 . It is a plant belonging to ranunculaceae family and is also known with various names (black caraway, black cumin, kalojeera, kalonji) 2 . It is known to be a miraculous plant as it possesses extraordinary activities such as analgesic 3 , anti-bacterial 4 , anti-viral 5 , anti-diabetic 6 , hepatoprotective 7 , etc. These activities are basically due to the presence of an active component called as "thymoquinone" 8 . This active component alone has been reported to produce excellent anti-oxidant as well as anti-inflammatory property 9,10 . Other than these, reports concerning its anti-cancer activity were also reported [11][12][13] . Other important effects shown by thymoquinone were anti-diabetic 14 , hepatoprotective 15 , anti-viral 16 as well as its use in treating autoimmune disorders such as vitiligo 17 . Considering the significance of thymoquinone and its structure, various derivatives of it were synthesized earlier and evaluated for their biological activity 18,19 . It was shown to produce promising results against ovarian cancer and also as anti-malarial activity 18 . To understand the structural activity relationship and to identify the pharmacophore, one has to adhere to computational techniques which shortens the time to determine the pharmacophore 20 . Hence, in this research work, we have taken the thymoquinone analogs having action against ovarian cancer and performed docking studies to understand the binding pattern of all the derivatives. Furthermore, we have also predicted various pharmacokinetic parameters and also its toxicology profile to screen the compounds with high risk which can give insight of the optimized molecule with greater activity.

Selection of thymoquinone analogues for study
Literature survey reveals various analogs of thymoquinones from which, we have shortlisted the molecules which have got IC 50 values against cancer (particularly ovarian cancer) 18 . The list of molecules taken for study were sketched in chemDraw (v14) can be seen in Figure 1.

Calculation of physiological parameters
The physiological parameters of the molecules include molecular weight, h-bond donor, acceptor, its partition co-efficient (LogP), molar refractivity and total polar surface area (TPSA). In order to calculate these values, the molecules were sketched in chembiodraw (version 14, developed by perkinelmer) and saved in ".Sdf" format. These saved structures were used to calculate the drug likeness from freely available online swissadme tool 21 .

Prediction of pharmacokinetic parameters
Drugs are said to be good when it is easily absorbed, distributed, metabolized and excreted from the human body. In order to screen the molecules, In-silico methods are majorly adopted wherein the pharmacokinetic parameters are predicted based on the molecular structure by comparing with the large number of molecules in the database and by similarity searching. SwissADME tools does that function and gives the details of the drug's pharmacokinetic profile which include GI absorption, BBB permeability, P-GP substrate, cyp450 inhibition and Log K p (skin permeation). Therefore, all the molecules in study were used to predict their pharmacokinetic parameters.

Obtaining of protein crystal structure
The thymoquinone analogs were shown to have anti-cancer activity (ovarian cancer in particular). The protein structure corresponding to the said activity was found to be ISA0 from the literature. Therefore, the protein molecule (PDB ID: 1SA0) which is a complex of tubulin with colchicine was downloaded from protein databank (www. rcsb.org) with the resolution of 3.58å. The ligand was manually removed and the active site was observed from pdbsum website (http://www.ebi. ac.uk) and the co-ordinates for generating grid box was observed from protein viewer (x=123.324, Y=97.286 And z=6.838).

Preparation of ligands
As the molecules are new and cannot be available in the database, these were sketched in ChemBio3d ultra version 14.0 Software and minimized using minimize structure option of calculate tool from the toolbar. The molecules after minimizing were saved in ".Mol2" format.

Docking
To p e r fo r m t h e d o ck i n g s t u d y o f thymoquinone analogs, autodock version 4.2.6 Was used. The ligand as well as protein files were converted to "pdbqt" format and proceeded further for autogrid and autodock respectively. It was then analyzed for the results.

Toxicological parameters prediction
Most of the molecules, despite having a desired pharmacological activity possess certain side effects and also sometimes my act as a toxic substance resulting in various reactions within the body. Few online tools are available which can predict the toxic properties of compounds based on the similarity of the structure from the databases. Hence, to accomplish this, we have used swissADME tool which is available for free to use to obtain the toxicological parameters which are discussed in results and discussion.

Molecules selection and physiological properties
From the literature, twenty-two molecules were selected including thymoquinone which were reported to have IC 50 value against ovcar-8 (ovarian cancer cell line). These molecules were sketched and saved in "sdf" format and checked for physiological parameters and screen out the compounds which does not obey lipinski rule of five (mol. Wt ≤ 500, Log P ≤ 5, hd ≤ 5, ha ≤ 10). All the molecules were shown to have obeyed the lipinski rule except thy 17 as its log p value was observed to be greater than 5. All the physiological parameters are shown in the following Table 1. From the values obtained, it can be said that all these molecules can be administered orally except THY17 as its IC 50 values is also the highest (Table 4). All the molecules were small with lipophilic in nature as they obey the lipinski RO5, hence possess high absorption.

Pharmacokinetic parameters Absorption profile of thymoquinone analogues
The drugs with good pharmacokinetic profile produce better results when administered. Once a drug it administered; it needs to undergo adme process 22 . Hence, we have predicted the pharmacokinetic profile of the molecules in our study. As per the Fig. 1 and Table 2, the drug molecules were found to be lipophilic due to presence of alkyl chains and furthermore, the molecules obey  THY01  THY02  THY03  THY04 THY05  THY06  THY07  THY08  THY09  THY10  THY11 GI Note: Log K p = Skin Permeation P-glucoprotein is a cell membrane protein that binds to drug molecules decreasing the pharmacokinetic profiles of those drugs. If drugs bind to p-gp, it pumps back the drug into the lumen which then undergoes metabolism in the liver followed by excretion through proximal tubule of kidney. All the molecules which were predicted for p-gp substrate, show no binding as can be seen in Table 3, therefore, possess good absorption and greater efficacy.

Metabolism profile of thymoquinone analogs
Most of the drugs undergoes metabolism in liver in presence of cytochrome P450 which contains various groups of enzymes. The main function of this is to metabolize the drugs to easily excretable forms. The majorly contributing groups for metabolism include CYP1A2, CYP2C19, CYO2C9, CYP2D6 and CYP3A4. If the drugs not metabolized via these enzymes, it can be said that the drug has inhibitory effect on the CYP450 enzyme groups. Therefore, from adme prediction of the drugs, we have observed the drugs to have inhibitory effect against few of the groups of CYP450. For an instance, more drugs have shown inhibitory action against CYP1A2, while the others have very a smaller number of analogs showing inhibitory effect. This provides and insight that the molecules can be metabolized via CYP450 there by converting them to easily excretable forms. The complete metabolism profile of the drugs can be seen in Table 3.

Docking analysis of thymoquinone analogs
To understand the binding affinity of the analogues with the ovarian cancer protein, docking was performed with a tubulin protein (PDB ID:1SA0, 3.58å) obtained from protein data bank. Gasteiger charges were applied on the protein and docked by following genetic algorithm. Among all the molecules, three molecules were shown to have no interactions which are THY10, THY13 and THY15. The docking analysis results are shown in the following Table 4 and the graph can be seen in Figure 2. lipinski rule of 5. Therefore, the molecules can be considered to have greater GI absorption and have permeability to blood brain barrier (BBB). The results of pharmacokinetic profile are shown in Table 2.

Fig. 2. Binding energies of Thymoquinone analogs
Among the set of twenty-two molecules, compounds with better pIC 50 values were observed to be THY06, THY09, THY11, THY12, THY14, THY21 and THY22 with THY09 being the most active. From compounds THY06, THY09 and THY11, it can be deduced that the presence of electronegative atoms on the aromatic ring produces more activity. For instance, presence of a tertiary amine at 5 th position of THY06 instead of isopropyl group as in case of THY01 shows more potency which can be as a result of lone pair of electrons present on the nitrogen that is donating its electrons resulting in more electronegativity at the 4 th keto group there by forming a better bonding with hydrogen of ASN206 residue (bond length: 2.203) While presence of same dimethyl group at position 6 th show reduced activity. It can be conferred that, substituted heteroatoms at position 5 if present shows more potency than the parent molecule (THY01). Introduction of halogens to the parent molecule (THY01) is observed to have greater activity than THY06 as can be observed in the docking result of THY09 and THY11. However, bromo substitutions at position 3 and 6 has more potency than chloro-substitution. Therefore, presence of low electronegative halogens shows better action. The 4 th position keto group is essential as it has shown binding with the active site residues of the 1SA0 protein. The docking images of THY06, THY09 and THY11 can be seen in the Figure 3 (a-c).
Similarly, in case of THY21 and THY22, which are the reduced forms of THY14 and THY05 respectively, all the four compounds have good activity except the fact that, THY05 shows binding with ASN206 and ASN228 with the 1 st and 4 th position ketone group while THY22 shows binding with GLN15 and ASN206 with the hydroxyl groups as can be observed in Fig. 4a and 4b. Furthermore, compound THY21 which is the reduced form of THY14 shows binding with GLY142, ILE171 and SER178 while THY14 shows binding with ASN228 as can be seen from Fig. 4c and 4d respectively. However, none of the compounds show the activity as potent as colchicine either in terms of their IC 50 or in terms of binding affinity but shows closer value. Since THY14 is the most active compound, it can be inferred that retaining 4 th keto group and reducing the 1 st keto group may provide an insight in developing drugs which can shows more binding with the residues of the active site and as well possess good activity against cancer.   Inhibitory constant is the term used to describe the potency of the drug. It is the lowest concentration at which the drug shows inhibitory effect. Autodock provides inhibitory constant values based on the its binding to the protein active site. Lower the ki value, potent will be its effect theoretically in case of in-silico studies. From the inhibitory constant values obtained as shown in Table  5, the compound which is most active (THY09) gives 24.54µm followed by THY11 with 41.92µm which is four to six times greater than colchicine while the other molecules show higher ki values.

Toxicological parameters prediction
To obtain toxicology parameters of the compounds, we have used swissadme tool, where in, it predicts the toxic effects of the molecules based on its structure either to be mutagenic and tumorigenic. It also predicts if the molecules possess any effect on the reproductive system of human beings. From the results as can be seen in Table 5, most of the molecules show no tumorigenic property except THY21 and THY22 while THY20 and THY21 has shown to produce some effect on the reproductive system. Hence, the molecular properties of THY21 and THY22 are to be neglected from considering to be a part of scaffold.
On the other hand, THY09 being the most active compound among the set of molecules shows none of the considered toxic properties, while molecules THY06, THY11 and THY12 shown to possess high mutagenic property. Therefore, THY09 can be considered as the optimized thymoquinone analogs based on which various modifications can be performed either by increasing the side chain or introducing hetero substitution on the 5 th position, retaining 4 th keto group and reducing 1 st keto group can help in developing more potent molecules.