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

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-13]. Other important effects shown by thymoquinone were anti-diabetic^[14], hepato-protective^[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.

Materials and Methodology

Selection of Thymoquinone Analogues for Study

Literature survey reveals various analogs of thymoquinones from which, we have shortlisted the molecules which have got IC50 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.

Figure 1: Thymoquinone analogs taken for study Click here to View figure

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-silicomethods 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 perform the docking study of thymoquinone analogs, Autodockversion 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.

Redults and Discussion

Molecules Selection and Physiological Properties

From the literature, twenty-two molecules were selected including thymoquinone which were reported to have IC50 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 THY17 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 IC50 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.

Table 1: Physiological Parameters of Thymoquinones

Compound name	Molecular formula	Mol Weight (g/mol)	HA	HD	TPSA	MR	Log P	Drug-likeness
THY01	C10H12O2	164.20	2	0	34.14	47.52	1.018	-1.1996
THY02	C11H15NO2	193.24	3	1	46.17	55.13	1.502	2.7825
THY03	C12H17NO2	207.27	3	1	46.17	59.94	1.925	2.5661
THY04	C10H13NO2	179.22	3	0	60.16	50.23	1.039	-0.022
THY05	C7H6O2	122.12	2	0	34.14	33.10	-0.047	-0.6921
THY06	C9H11NO2	165.19	3	0	37.38	45.61	0.487	1.4535
THY07	C9H11NO2	165.19	3	0	37.38	45.61	0.487	1.8235
THY08	C13H10FNO2	231.22	3	1	46.17	61.89	1.154	-1.4514
THY09	C10H10Br2O2	321.99	2	0	34.14	63.26	2.608	-7.9858
THY10	C11H14O2	178.23	2	0	34.14	52.07	1.514	-7.1574
THY11	C11H12Cl2O2	247.12	2	0	34.14	61.66	2.668	-5.7822
THY12	C10H12O2	164.20	2	0	34.14	47.26	1.333	-4.2276
THY13	C14H20O2	220.31	2	0	34.14	66.23	2.894	-5.0241
THY14	C8H8O2	136.15	2	0	34.14	37.91	0.134	-1.0759
THY15	C14H20O2	220.31	2	0	34.14	66.23	2.894	-4.5682
THY16	C16H24O2	248.36	2	0	34.14	75.84	4.032	-0.1315
THY17	C20H32O2	304.47	2	0	34.14	94.55	6.26	-9.7535
THY18	C10H14O2	166.22	2	2	40.46	50.03	1.984	-2.3359
THY19	C10H12Br2O2	324.01	2	2	40.46	65.43	3.49	-4.1259
THY20	C11H16O2	180.24	2	2	40.46	54.72	2.396	-5.276
THY21	C8H10O2	138.16	2	2	40.46	40.42	1.016	.2.315
THY22	C7H8O2	124.14	2	2	40.46	35.45	0.835	-2.315

HA = Hydrogen bond acceptor; HB = Hydrogen bond donor; MR = Molar Refractivity; TPSA = Topological Polar Surface Area; Log P =Partition co-efficient (P_o/w)

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 figure 1 and table 2, the drug molecules were found to be lipophilic due to presence of alkyl chains and furthermore, the molecules obey 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.

Table 2: Absorption and Distribution of Thymoquinone Analogs

	THY01	THY02	THY03	THY04	THY05	THY06	THY07	THY08	THY09	THY10	THY11
GI Absorption	High	High	High	High	High	High	High	High	High	High	High
BBB permeability	Yes	Yes	Yes	Yes	Yes	Yes	Yes	Yes	Yes	Yes	Yes
P-GP Substrate	No	No	No	No	No	No	No	No	No	No	No
Log Kp (cms-1)	-5.74	-6.25	-6.07	-6.63	-6.53	-6.55	-6.78	-6.11	-5.89	-5.55	-5.14
	THY12	THY13	THY14	THY15	THY16	THY17	THY18	THY19	THY20	THY21	THY22
GI Absorption	High	High	High	High	High	High	High	High	High	High	High
BBB permeability	Yes	Yes	Yes	Yes	Yes	Yes	Yes	Yes	Yes	Yes	Yes
P-GP Substrate	No	No	No	No	No	No	No	No	No	No	No
Log Kp (cms-1)	-5.86	-4.51	-6.26	-4.51	-4.63	-4.06	-5.23	-5.38	-5.52	-6.43	-6.41

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, CYP2C9, 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.

Table 3: Metabolism of Thymoquinone Analogs

Compound name	CYP450 inhibition
	CYP1A2	CYP2C19	CYP2C9	CYP2D6	CYP3A4
THY01	No	No	No	No	No
THY02	Yes	No	No	No	No
THY03	Yes	No	No	No	No
THY04	No	No	No	No	No
THY05	No	No	No	No	No
THY06	No	No	No	No	No
THY07	No	No	No	No	No
THY08	Yes	Yes	No	No	No
THY09	Yes	Yes	No	No	No
THY10	No	No	No	No	No
THY11	Yes	No	Yes	No	No
THY12	No	No	No	No	No
THY13	Yes	No	Yes	No	No
THY14	No	No	No	No	No
THY15	Yes	No	Yes	No	No
THY16	Yes	Yes	No	No	No
THY17	Yes	No	No	No	No
THY18	Yes	No	No	No	No
THY19	Yes	No	Yes	Yes	No
THY20	Yes	No	No	No	No
THY21	No	No	No	No	Yes
THY22	No	No	No	No	Yes

 

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.

Table 4: Docking Results of Thymoqinone Analogs

Comp name	IC50	Binding energy (ΔG)	Inhibitory constant (Ki)	RMSD	Binding residues with bond length (Å)
	OVCAR-8
THY01	11.6	-5.5	93.32µM	161.39	ASN206-1.916
THY02	12.9	-5.98	41.19µM	161.42	ASN228-1.992 with NH.
THY03	16	-6.33	22.93µM	161.27	ASN206-1.992; ASN228-1.872
THY04	35.1	-5.57	82.11µM	161.84	ASN206-2.141
THY05	13.4	-4.62	409.55µM	162.8	ASN206-1.694, ASN228-1.975
THY06	5.6	-4.93	244.68µM	161.2	ASN206-2.211 with (=O)
THY07	35.6	-4.88	266.88µM	161.14	ASN206-1.974 with (=O)
THY08	32	-6.72	11.81µM	159.54	ASN206-2.186, 1.748
THY09	3.6	-6.29	24.54µM	160.79	ASN206-2.203
THY10	13.5	-5.23	146.16µM	160.86	No interactions
THY11	6	-5.97	41.92µM	160.78	ASN206-2.172
THY12	10	-4.73	338.56µM	161.18	GLM11-2.129
THY13	51.2	-5.48	95.39µM	159.78	No interactions
THY14	5.2	-4.81	298.7µM	160.85	ASN228-1.954
THY15	37.4	-5.39	112.09µM	159.82	No interactions
THY16	54.2	-5.64	74.04µM	159.11	ASN206-1.96
THY17	56.5	-7.21	5.2µM	154.93	ALA12-2.034
THY18	8.9	-5.83	53.19µM	161.45	ASN206-2.074, 2.111; GLN15-2.199
THY19	11.6	-6.78	10.7µM	160.17	ASN206-2.027; GLN15-2.166; GLN11-2.05
THY20	12.2	-5.73	46.4µM	160.44	GLN15-2.124; GLN11-2.194; ASM206- 2.22, 2.027
THY21	8.3	-5.2	154.95µM	159.64	GLY142-1.027; ILE171-1.929; SER178-1.792
THY22	6.2	-4.96	233.28µM	161.54	GLN15-2.02; ASN206-2.091, 2.143
Colchicine		-7.07	6.63µM	156.02	SER140-1.939

 

Figure 2: Binding energies of Thymoquinone analogs Click here to View figure

Among the set of twenty-two molecules, compounds with better pIC50 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^thketo 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).

Figure 3: a) Docking images of Thy06 b) Docking images of Thy09 c)Docking images of Thy11 Click here to View figure

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 figure 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 figure 4c and 4d respectively. However, none of the compounds show the activity as potent as colchicine either in terms of their IC₅₀ or in terms of binding affinity but shows closer value. Since Thy14 is the most active compound, it can be inferred that retaining 4^thketo group and reducing the 1^stketo 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.

Figures 4: a) Docking image of Thy05; b) Docking image of Thy22; c) Docking image of Thy14 and d) Docking image of Thy21 Click here to View figure

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. Figure 5 depicts the inhibitory constant values showing their Ki values.

Figure 5: Inhibitory constant values of Thymoquinone analogs Click here to View figure

Toxicological Parameters Prediction

To obtain toxicology parameters of the compounds, we have used SwissADME tool, wherein, 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.

Table 5: Toxic Parameters of Thymoquinone Analogs

Mol name	Mutagenic	Tumorigenic	Reproductive effect	Mol name	Mutagenic	Tumorigenic	Reproductive effect
THY01	High	None	None	THY12	High	None	None
THY02	None	None	None	THY13	High	None	None
THY03	None	None	None	THY14	High	None	None
THY04	None	None	None	THY15	None	None	None
THY05	High	None	None	THY16	Low	None	None
THY06	Low	None	None	THY17	Low	None	None
THY07	None	None	None	THY18	High	Low	None
THY08	Low	None	None	THY19	None	None	None
THY09	None	None	None	THY20	High	High	High
THY10	Low	None	None	THY21	None	None	High
THY11	High	None	None	THY22	None	High	None

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^thketo group and reducing 1^stketo group can help in developing more potent molecules.

Conclusion

Structure optimization was achieved by performing docking studies of various thymoquinone analogs as well as their pharmacokinetic and toxicological parameters were determined. From the results obtained, all the molecules obeyed Lipinski rule except THY17. It was found that compound THY09 was most stable among all the molecules and was shown to have no toxic effects. Further dynamics studies are required to determine the stability of the binding complex. In future aspect, structural modifications of Thy09 such as retaining of 4^thketo group and reduction of 1^stketo, addition of bulky hetero groups such as pyrazine or phenyl sulfur, etc., on 5^th position of the aromatic ring could lead to the compounds with greater potency.

Acknowledgement

The authors are thankful to Management & Principal of Sultan-ul-Uloom College of Pharmacy for providing research facilities for this work.

Conflict of Interest

The authors declare that they have no conflict of interest.

References

1. Majeed A, Muhammad Z, Ahmad H, Hayat SS, Inayat N, Siyyar S.ActaEcologicaSinica. 2020.
2. Sultana S, Asif HM, Akhtar N, Iqbal A, Nazar H, Rehman RU. J. PharmacognosyPhytother.2015;4(4):103-6.
3. Ghannadi A, Hajhashemi V, Jafarabadi H. J. Med. Food. 2005;8(4):488-93.
   CrossRef
4. Chaieb K, Kouidhi B, Jrah H, Mahdouani K, Bakhrouf A. BMCComplemAltern M. 2011;11(1):29.
   CrossRef
5. Khan AU, Tipu MY, Shafee M, Khan NU, Tariq MM, Kiani MR, Shah SI. Pak Vet J. 2018;38:434-7.
6. Benhaddou-Andaloussi A, Martineau L, Vuong T, Meddah B, Madiraju P, Settaf A, Haddad PS.Evid-Based Compl Alt. 2011;2011.
   CrossRef
7. Adam GO, Rahman MM, Lee SJ, Kim GB, Kang HS, Kim JS, Kim SJ. Asian Pac. J. Trop. Med. 2016;9(3):221-7.
   CrossRef
8. Ali BH, Blunden G. Phytother Res. 2003;17(4):299-305.
   CrossRef
9. Ismail M, Al-Naqeep G, Chan KW. Free Radic. Biol. Med.2010;48(5):664-72.
   CrossRef
10. El Gazzar M, El Mezayen R, Marecki JC, Nicolls MR, Canastar A, Dreskin SC. Int. Immunopharmacol. 2006;6(7):1135-42.
    CrossRef
11. Woo CC, Loo SY, Gee V, Yap CW, Sethi G, Kumar AP, Tan KH. Biochem. Pharmacol. 2011;82(5):464-75.
    CrossRef
12. Attoub S, Sperandio O, Raza H, Arafat K, Al‐Salam S, Al Sultan MA, Al Safi M, Takahashi T, Adem A. Fund ClinPharmacol. 2013;27(5):557-69.
    CrossRef
13. Mahmoud YK, Abdelrazek HM. Biomed Pharmacother. 2019;115:108783.
    CrossRef
14. Fararh KM, Shimizu Y, Shiina T, Nikami H, Ghanem MM, Takewaki T. Res. Vet. Sci. 2005;79(3):219-23.
    CrossRef
15. Daba MH, Abdel-Rahman MS.Toxicol. Lett. 1998;95(1):23-9.
    CrossRef
16. Fröhlich T, Reiter C, Saeed ME, Hutterer C, Hahn F, Leidenberger M, Friedrich O, Kappes B, Marschall M, Efferth T, Tsogoeva SB. ACS Med. Chem. Lett. 2017;9(6):534-9.
    CrossRef
17. Ali SA, Meitei KV. J Pharm Pharmacol. 2011;63(5):741-6.
    CrossRef
18. Johnson-Ajinwo OR, Ullah I, Mbye H, Richardson A, Horrocks P, Li WW. Bioorg Med ChemLett. 2018;28(7):1219-22.
    CrossRef
19. Yusufi M, Banerjee S, Mohammad M, Khatal S, Swamy KV, Khan EM, Aboukameel A, Sarkar FH, Padhye S. Bioorg Med ChemLett. 2013;23(10):3101-4.
    CrossRef
20. Xiang M, Cao Y, Fan W, Chen L, Mo Y. Comb. Chem. High Throughput Screen. 2012;15(4):328-37.
    CrossRef
21. SwissADME: a free web tool to evaluate pharmacokinetics, drug-likeness and medicinal chemistry friendliness of small molecules. Sci. Rep. 2017; 7:42717.
    CrossRef
22. Benedetti MS, Whomsley R, Poggesi I, Cawello W, Mathy FX, Delporte ML, Papeleu P, Watelet JB. Drug Metab. Rev. 2009;41(3):344-90.
    CrossRef

---

This work is licensed under a Creative Commons Attribution 4.0 International License.