Insights from molecular modeling, docking and simulation of imidazole nucleus containing chalcones with EGFR kinase domain for improved binding function.

EGFR Kinase domain is a crucial role player cell surface receptor protein activated by specific binding of its ligand like EGFR. Importance of this protein as a therapeutically important drug target towards treating various cancer types has been proven elsewhere in previous literature. In this present study, we have designed a novel series of five compounds and computationally evaluated their potential to act as inhibitors of EGFR kinase domain towards anti-cancer activity. Our docking study shows compounds have the potential to dock into the active site of the EGFR Kinase domain with a binding energy in a range of -5.46 to - 7.32 Kcal/mol, Among all the compounds, compound 2 was found to be the lead like molecule with the binding energy of -7.32 kcal/mol with predicted IC50 value of 4.33 micro molar level. Molecular dynamic simulation studies for this compound 2 in complex with EGFR kinase domain has revealed several interesting molecular interactions with some of the important residues present at the active binding site of EGFR Kinase domain. Conclusively, novel designed compound 2 of the present study have shown promising anti-cancer potential worth considering for further evaluations.

Background

EGFR is a cell surface protein that binds to epidermal growth factor. Binding of the protein to a ligand induces receptor dimerization and tyrosine autophosphorylation and leads to cell proliferation [1]. Cancer is one of the major diseases associated with EGFR. Receptor tyrosine kinase binding ligands of the EGF family and activates several signaling cascades to convert extracellular cues into appropriate cellular responses. Known ligands include EGF, epigen/EPGN, TGFA/TGF-alpha, BTC/betacellulin, amphiregulin, HBEGF/heparin-binding EGF and epiregulin/EREG. Ligand binding triggers receptor homoand/ or heterodimerization and autophosphorylation on key cytoplasmic residues [2]. Rao et all has concluded that EGFR is the major targeted domain for cytotoxic activity [3]. Ke-Ming Qiu et al has proposed that a nitrogenous ring pyrazole showed high binding affinity towards EGFR domain [4]. Er Bin Yang et al. reported that these are the major inhibitors of egfr kinase domain due to their high hydrogen bonds and hydrophobic interactions [5]. Yun-Yun Xu et all have stated that curcumin analogues are also a potential egfr inhibitors while showing their antiproliferative activity against tumor cells [6]. Based on the above study, the rationalized chalcone molecules might be show better anticancer activity against the targeted domain EGFR.

Computational methods:

Software and program

Schrodinger’s maestro visualization program v9.5 [7] is utilized to visualize the receptors, ligand structures, hydrogen bonding network, to calculate length of the bonds and to render images. Chemsktech was used to draw the ligand compounds. Autodock is the primary docking program used in this work for the semi-flexible docking studies. Preparation of the ligands and protein receptors in pdbqt file and determination of the grid box size were carried out using Auto-Dock Tools version 1.5.6 [8]. Schrodinger’s Desmond module [9] was used for molecular dynamic simulation studies. Protocol used for the docking and simulation is described elsewhere in detail [10,11, 12,13].

Results & Discussion

Docking of the compounds with EGFR Kinase domain active site

Firstly, we have performed the docking studies for the five compounds (Figure 1 in supplementary document) with the EGFR Kinase domain protein targeting its active binding site in order to know the binding energy involved in this complex formation and to know the molecular interactions responsible for this target specific inhibition. Docking results are tabulated in table 1. All the five compounds studied in this present work have shown to be successfully docking inside the active site of EGFR kinase domain with a binding energy in a range of -5.46 to -7.32 Kcal/mol. As per the molecular docking results, it was revealed that Compound 2 has the best estimated -7.32 Kcal/mol of binding energy (Table 1) for the EGFR kinase domain inhibited complex formation by forming two hydrogen bond with LYS721 and CYS773 residues. Apart from hydrogen bonds, this compound was found to be forming hydrophobic interactions with MET769; LEU694; PRO770; PHE771; LEU768; ALA719; ILE720; LEU764; LEU753; ILE765 and MET742 along with hydrophilic interactions with THR 830; THR766 and LYS704 residues (Figure 1). Apart from above interactions picationic, pi-pi interactions and other non-covalent bonds have thought to play key role in stabilizing the docked compounds in the active site of EGFR Kinase domain towards forming this stable complexes. Interaction snapshots for rest of compounds have been present in Figure S1 in the supplementary document.

Figure 1

A) Docking snapshot of compound-2 (3-[4-(benzyloxy)phenyl]-1-[4-(1H-imidazol-1-yl)phenyl]prop-2-en-1-one) at the active binding site of EGFR kinase domain; b) 3D docking snapshot showing compound-2 (3-[4-(benzyloxy)phenyl]-1-[4-(1Himidazol- 1-yl)phenyl]prop-2-en-1-one) forming a hydrogen bond with LYS704 residue and c) represents the 2D interactions between compound-2 (3-[4-(benzyloxy)phenyl]-1-[4-(1H-imidazol-1-yl)phenyl]prop-2-en-1-one) and EGFR Kinase domain.

Table 1

Docking energies of compounds 1-5 with EGFR Kinase domain

Compound Name	Compound	Docking binding energy in Kcal/mol	Predicted IC50 value in micro molar	No
(1-[4-(1H-imidazol-1-yl)phenyl]-3-[2-(trifluoromethyl)phenyl]prop-2-en-1-one)	#1	-5.46	99.96	1
(3-[4-(benzyloxy)phenyl]-1-[4-(1H-imidazol-1-yl)phenyl]prop-2-en-1-one)	#2	-7.32	4.33	2
(1-[4-(1H-imidazol-1-yl)phenyl]-3-[2-(trifluoromethyl)phenyl]prop-2-en-1-one)	#3	-5.6	78.91	3
(3-(2-chloro-6-fluorophenyl)-1-[4-(1H-imidazol-1-yl)phenyl]prop-2-en-1-one)	#4	-6.89	8.96	4
(3-[4-(diethylamino)phenyl]-1-[4-(1H-imidazol-1-yl)phenyl]prop-2-en-1-one)	#5	-6.53	16.27	5

IC50 prediction

In order to understand the plausible experimental anti-cancer activity of the present studied compounds, we have carried out the half maximal inhibitory concentration (IC50) value predictions. IC50 value is a useful parameter to quantitatively measure the effectiveness of compound to inhibit a given biological process by half and is universally used to symbolize the inhibitory effect of compounds [14]. Table 2 shows the predicted IC50 value for the compounds 1-5 were in a range of 4.33 to 99.96 micro molar. Among which the compound 2 has shown the best possible inhibitory potential with 4.33 micro molar. IC50 values obtained clearly demonstrated plausible high inhibitory potential of present studied compounds with kinase domain of EGFR.

Table 2

Molecular dynamic simulations trajectory statistical data analysis of EGFR Kinase domain protein in its apo state; in complex with compound-2 (3-[4-(benzyloxy) phenyl]-1-[4-(1H-imidazol-1-yl)phenyl]prop-2-en-1-one)

Parameter	EGFR apo		EGFR+ compound 2 complex EGFR apo
	Mean	Range	Mean	Range
Total Energy	-9003.562	-9780.786, -8243.045	-9359.584	-10044.489, -479.591
VDW Energy	-323.745	-437.542, -205.494	-320.960	-433.675, -186.172
Coulomb’s energy	-5423.000	-5759.860, -5082.624	-5554.231	-5963.839, -5114.732
RMSD	2.363	0.000, 3.291	1.667	0.000, 2.269
Intra H-Bonds	210.609	181.000, 232.000	219.812	197.000, 243.000
Inter H-Bonds	N/A	N/A	0.027	0.000, 1.000
Radius of Gyration	20.191	19.843, 20.579	20.003	19.723, 20.272

MD simulations of EGFR Kinase domain protein in its apo state (without any ligand) in comparison with in complex with compound 2

In order to analyze the stability and conformational changes along with the aim of revealing their underlying molecular interactions at atomic level of EGFR Kinase domain in its native state without any ligand; in comparison to in complex with our present investigated compound 2. Before starting the analysis part of the MD simulations trajectories, firstly we have analyzed the simulation quality parameters, which is an important issue to deal with to conform that the molecular dynamic simulations were carried out under our given temperature, pressure and simulation box volume conditions throughout the simulated time. For all the molecular dynamic simulations, input parameters given such that the temperature is at 300 k; pressure is at 0 atmospheric and volume of the simulations water box is kept at 365000 Å3 throughout the simulated time. As evident with Figure S2 in the supplementary material, the given input simulation parameters were maintained thoroughly throughout the simulated time of 10 ns each, thus conforming the quality of the simulations carried out.

Simulation event analysis

After confirming the simulation quality parameters, we have studied the Root mean square deviation (RMSD), Root mean square fluctuations (RMSF), Radius of gyration (ROG), intra molecular hydrogen bonds; total energy and its total secondary structure elements (SSE) contributions with time dependant function of MD simulations in order to understand the stability and conformational changes of the EGFR Kinase domain protein in its apo state and in complex with compound 2 at each frame of the MD simulations trajectory. The statistical data for the analysis of the MD trajectories has been presented in Table 2 and the relevant graphs were represented in Figure 2.

Figure 2

Event analysis graphs of EGFR in complex with compound-2 (3-[4-(benzyloxy) phenyl]-1-[4-(1H-imidazol-1- yl)phenyl]prop-2-en-1-one). Red color lines represent EGFR kinase domain in its native state without any ligand and blue color lines represent EGFR kinase domain in complex with compound -2 (3-[4-(benzyloxy)phenyl]-1-[4-(1H-imidazol-1-yl)phenyl]prop-2-en-1-one).

Molecular interactions of EGFR Kinase domain-compound 2 complex during MD simulations

We have used Simulation interactions diagram program integrated within Desmond module of Schrödinger for studying the detailed inter-molecular interactions between EGFR Kinase domain and compound 2. There were about 28 contacts found in between EGFR Kinase domain and compound 2 in total among which one hydrogen bond was observed with HIS781 residue and most of other contacts were found to be hydrophobic contacts followed by water bridging and ionic bonds. Table 3 shows the molecular interactions profile of compound 2 with EGFR kinase domain. Hydrophobic contact were found to forming with residues LEU694; SER696; PHE699; VAL702; ALA719; LYS721; MET742; LEU753; LEU764; GLN767; LEU768; MET769; PHE771; HIS781 and LEU820. Whereas, LYS704 was found to be forming an Ionic bond with compound 02. Several hydrogen bonds were observed facilitated by water brides with the residues LYS692; LEU694; LYS704; CYS751; THR766; MET769; PRO770; PHE771; ASP776; TYR777; ARG779; GLU780; HIS781; THR830 and ASP83. Apart from above molecular interactions LYS721; PHE771 and HIS781 residues were found to be forming pi-pi stackings with compound 2 towards stabilizing the protein-ligand complex in its active site (Figure 3).

Table 3

Molecular Interaction profile of Compound -2 (3-[4-(benzyloxy) phenyl]-1-[4-(1H-imidazol-1-yl) phenyl] prop-2-en-1-one) with EGFR Kinase domain during MD simulation

No	Drug Target	Compound	Molecular interactions
			H-Bonds	Hydrophobic interactions	Ionic bonds	Water bridging residues	Pi-Pi stacking
1	EGFR Kinase domain	#2	HIS781	LEU694; SER696; PHE699; VAL702; ALA719; LYS721; MET742; LEU753; LEU764; GLN767; LEU768; MET769; PHE771; HIS781; LEU820.	LYS704	LYS692; LEU694; LYS704; CYS751; THR766;MET769; PRO770; PHE771; ASP776; TYR777; ARG779;GLU780; HIS781; THR830; ASP831.	LYS721; PHE771; HIS781.

Figure 3

A) A schematic of detailed compound -2 (3-[4-(benzyloxy)phenyl]-1-[4-(1H-imidazol-1-yl)phenyl]prop-2-en-1-one) atom interactions with the EGFR Kinase domain protein residues. Interactions that occur more than 0.0% of the simulation time in the selected trajectory (0.00 through 10.00 nsec), are shown; B) Protein-ligand interactions (or 'contacts') of EGFR Kinase domain- Compound -2 (3-[4-(benzyloxy)phenyl]-1-[4-(1H-imidazol-1-yl)phenyl]prop-2-en-1-one) complex. The stacked bar charts are normalized over the course of the trajectory: for example, a value of 0.7 suggests that 70% of the simulation time the specific interaction is maintained.

Conclusion

The present study provides a rationalization to the ability of present studied novel compounds as a valuable small ligand molecule with strong binding affinity towards EGFR Kinase domain for plausible anti-cancer activity. Our computational analysis evidenced that the large value of binding energy is involved in binding of present investigated five compounds with the EGFR Kinase domain consolidating their complex’s thermodynamic stability; moreover predicted IC50 values further substantiated our hypothesis that these compounds have the potential to inhibit EGFR Kinase domain. Further, de novo simulations for 10 ns suggest that ligand interaction with the residues of EGFR, all or some of which fall under catalytic active site important residues for its structurally stability and/or functionality, could be critical for its inhibitory activity.