Int J Pharm Pharm Sci, Vol 7, Issue 6, 77-91Original Article



Principal K. M. Kundnani College of Pharmacy, Cuffe Parade, Mumbai 400005, India

Received: 09 Mar 2015 Revised and Accepted: 05 Apr 2015


Objective: Identifying new inhibitors of Epidermal Growth Factor Receptor (EGFR) by virtual screening using a pharmacophore model followed by docking.

Methods: A pharmacophore model was developed using a dataset of 77 chemically diverse EGFR inhibitors using PHASE. Statistically valid Three Dimensional Quantitative Structure Activity Relationship (3D-QSAR) equations were generated for the pharmacophore model. This was followed by database screening to obtain probable hits. Docking of the probable hits into the crystal structure of EGFR was used as a second filter. Docking studies were carried out using GLIDE. Calculation of ADME properties of the probable hits arising out of docking further reduced the number of hits.

Results: A five-point pharmacophore was generated for EGFR inhibitors reported in literature. The pharmacophore indicated that the presence of two aromatic ring features (R), one acceptor feature (A), one donor feature (D) and one hydrophobic feature (H) is necessary for potent inhibitory activity. The generated pharmacophore yielded statistically significant 3D-QSAR model, with a correlation coefficient r2 of 0.9905 and q2 of 0.8764. Virtual screening using the best pharmacophore model resulted in 372 hits. Docking studies as a second filter reduced the hits to 8. Application of drug-likeness as a third filter gave 6 leads.

Conclusion: 6 leads with satisfactory pharmacokinetics properties were identified as potential EGFR inhibitors. This study may facilitate development of some new potential EGFR inhibitors.

Keywords: EGFR inhibitors, 3D-QSAR, Pharmacophore, Docking, Virtual screening.


A majority of the therapeutically used anticancer agents are associated with a major drawback of being nonspecific and non-targeted in their action. However, advances in biochemistry and molecular biology have lead to identification of an array of pathways which are disregulated in cancer [1]. Targeting disregulated pathways in an attempt to control cancerous growth is an approach which is proposed to be less toxic to the normal cells. The epidermal growth factor receptor (EGFR) is a member of the receptor family of tyrosine kinases [2]. Receptor tyrosine kinases are involved in a range of cancer cell behaviors, including excess growth, invasion and blood vessel formation [3]. EGFR is known to be over-expressed in several human solid tumors [4]. Therefore, EGFR has received much attention as a target for anticancer drugs [5, 6]. Erlotinib, Gefitinib and Lapatinib are some of the important drugs belonging to the class of EGFR inhibitors in clinical practice [7, 8].

Small molecules that inhibit the kinase activity of EGFR are of considerable interest as new therapeutic antitumor agents for the treatment of EGFR mediated cancers. These molecules act by binding, either reversibly or irreversibly, to the C-terminal tyrosine kinase domain of EGFR. Anilinoquinazolines are among the most explored compounds reported as EGFR inhibitors [9, 10]. These compounds, however, may not be very specific for inhibiting EGFR mutants because of their size and pattern of H-bonding interaction. This calls for exploration of other classes of compounds.

Pharmacophore modeling is considered as one of the important approaches used in the ligand-based drug discovery [11-13]. A pharmacophore model is a very useful tool for the screening of the database, and if combined with docking at a later stage, it can combine an advantages of both ligand-based and structure-based approaches for drug discovery. In view of these facts, the present study involved development of pharmacophore model, QSAR equations and use of these data for the screening of commercially available databases in search of novel chemical scaffolds.


Selection of dataset

A total of 77 EGFR inhibitors (table 1) reported in the literature [14] belonging to different chemical classes were used as a dataset for present pharmacophore modelling studies. The activity values were expressed as pIC50 ranged from 4 to 9 spanning five orders of magnitude. The ligands were considered as actives if pIC50 value was greater than ≥8.5 and considered as inactives if pIC50 value was ≤ 5. This resulted in a set of 23 highly active molecules and 17 inactive molecules. The dataset was divided into a training set (54 compounds) and a test set (23 compounds; table 1). Since the usefulness of the QSAR model depends on the diversity of structures and on effective spanning of the activity range, training set was selected after considering these parameters. The inhibitor, Erlotinib (Compound No.3), was also included in the training set.

Table 1: Structures and biological activities of training and test set compounds

Compound No. Structure Actual pIC50 Predicted pIC50 Residual
1 8 7.89 0.11
2 (T) 7.658 8.64 -0.982
3 8.699 8.98 -0.281
4 8.523 8.75 -0.227
5 8.699 8.58 0.119
6 8.699 8.5 0.199
7 8.699 8.92 -0.221
8 8.699 8.69 0.009
9 8.523 8.64 -0.117
10 8.523 8.53 -0.007
11 8.699 8.62 0.079
12 8.699 8.72 -0.021
13 (T) 8.046 8.63 -0.584
14 9 9.11 -0.11
15 8.187 8.32 -0.133
16 8.921 8.52 0.401
17 8.495 8.58 -0.085
18 8.337 8.67 -0.333
19 (T) 8.201 8.65 -0.449
20 9.097 8.91 0.187
21 8.585 8.74 -0.155
22 (T) 8.194 8.62 -0.426
23 (T) 8.18 7.72 0.46
24 (T) 7.997 8.73 -0.733
25 8.066 7.78 0.286
26 8.78 8.86 -0.08
27 8.569 8.79 -0.221
28 (T) 8.167 7.98 0.187
29 8.602 8.71 -0.108
30 8.347 8.09 0.257
31 8.602 8.8 -0.198
32 8.018 8.1 -0.082
33 (T) 8.347 7.78 0.567
34 8.745 8.75 -0.005
35 (T) 8.347 8.95 -0.603
36 (T) 8.367 7.95 0.417
37 8.721 8.45 0.271
38 8.569 8.38 0.189
39 (T) 8.222 8.74 -0.518
40 (T) 8.097 7.92 0.177
41 (T) 8.046 8.26 -0.214
42 8.046 7.75 0.296
43 8.097 8.34 -0.243
44 9 8.87 0.13
45 (T) 8.046 7.58 0.466
46 8.222 8.15 0.072
47 (T) 8.222 7.88 0.342
48 8.699 8.24 0.459
49 (T) 5.131 4.66 0.471
50 (T) 5.06 4.79 0.27
51 5.149 4.77 0.379
52 4 4.31 -0.31
53 4 4.06 -0.06
54 (T) 5.222 4.67 0.552
55 5.051 5.21 -0.159
56 5 5.06 -0.06
57 (T) 5.155 5.04 0.115
58 (T) 5.125 5.05 0.075
59 5 5.11 -0.11
60 (T) 5.237 5.12 0.117
61 5 4.82 0.18
62 5 5.12 -0.12
63 5 4.92 0.08
64 5 5.05 -0.05
65 5.244 5.29 -0.046
66 (T) 5.201 5.15 0.051
67 4 3.84 0.16
68 4 4.05 -0.05
69 4 3.9 0.1
70 4 4.05 -0.05
71 4.301 4.36 -0.059
72 4.301 4.17 0.131
73 4.986 5.2 -0.214
74 (T) 5.268 6.32 -1.052
75 4.862 4.77 0.092
76 5.282 5.3 -0.018
77 5.125 5.39 -0.265

(T) indicates compounds in test set

Generation of common pharmacophore hypotheses

PHASE, version 3.3, Schrodinger, LLC, New York, NY, 2011 implemented in the Maestro 9.3 software package was used to generate pharmacophore and three-dimensional quantitative structure–activity relationship (3D-QSAR) models for EGFR inhibitors. The two-dimensional (2D) chemical structures of the compounds were sketched in Maestro graphical interface and converted into corresponding standard three-dimensional (3D) structures using the best conformation model generation method with OPLS-2005 force field and mixed Monte Carlo multiple minimum/low mode (MCMM/lMOD) search algorithm to reduce redundant conformers. Conformations with energy higher than the global minimum by 10kcal/mol were rejected. Whenever the structure contained a chiral carbon, all possible isomers were generated for the study.

A set of six built-in pharmacophoric features including hydrogen-bond (H-bond) acceptor (A), H-bond donor (D), hydrophobic group (H), hydrophobic aromatic rings (R), positively ionizable group (P) and negatively ionisable group (N) are provided in Phase. This set of features was used to generate pharmacophoric sites for all the compounds. Various pharmacophore hypotheses were generated by taking into consideration active compounds using a tree-based partitioning technique that groups together similar pharmacophores based on their intersite distances. A tree depth of four was used in the present work. The terminal box size was 1 A °. The generated hypotheses were subjected to a stringent scoring process. The hypotheses which survived the scoring process were taken up for building QSAR models. Scoring with respect to actives was conducted using default parameters for site, vector, and volume terms.

Building QSAR models

Phase provides the means to build 3D-QSAR models for a set of ligands which are aligned to a selection of hypotheses and to visualize these models along with the ligand structures and the hypotheses. The QSAR models were developed for ligands belonging to the training set with a range of activities. The molecules of the training set were placed in a grid box of cubes, with each cube dimension corresponding to 1A °. Each cube was allotted a binary number based on the presence or absence of the pharmacophoric features present in it. This representation was then used to generate a 3D-QSAR model using partial least squares (PLS) method. All hypotheses successfully generated and scored were then used to build pharmacophore-based 3D-QSAR models with grid spacing 1A °, random seed value of one, and 5 PLS factors. The quality of the generated QSAR models was assessed by examining the associated statistical parameters such as standard deviation (SD), root mean square error (RMSE), variance ratio (F) and correlation coefficients: r2 and q2.

Validation of pharmacophore and 3D-QSAR models

The generated 3D-QSAR models were validated by leave-one-out method and prediction of activities of the test set compounds using the QSAR model. Result of the QSAR validation was also included as criterion for selection of best pharmacophore hypothesis.

Database screening using pharmacophore hypothesis

The best pharmacophore hypothesis was used as a 3D query for screening ZINC database. The purpose of this screening was to retrieve potential hits suitable for further development. Conformers were generated for each molecule in the database using Confgen. The database screening was carried out using phase virtual screening protocol implemented in Schrodinger with Best/Flexible Search option. The molecules retrieved from the database mapped all the five features of the pharmacophore.

Molecular docking

These studies were done using GLIDE (Maestro, version 8.5, Schrodinger, LLC, 2008) software and the crystal structure of EGFR (PDB code: 1M17), complexed with 4-anilinoquinazoline inhibitor, Erlotinib. The crystal structure was cleaned by deleting the ligand and the cofactors. This was followed by adding hydrogen atoms in their standard geometry, adjusting the bond orders and formal charges. The crystal structure was then refined, and the geometries were optimized with the OPLS-2005 force field using standard protocol and parameters. The inhibitor was extracted from the complex and re-docked (fig. 1).

The final docked conformation of the inhibitor was aligned to the original conformation and root mean square deviation (RMSD) calculated. RMSD value of 1.12 confirmed the accuracy of the docking program. The docking studies were carried out using extra precision mode of Glide using default parameters. The active site was defined by the generation of a grid box such that the co-crystallized ligand occupied the center of the box.

Fig. 1: Validation of docking protocol by re-docking 1M17 associated ligand, Erlotinib into the active site of the receptor

All the hits obtained from pharmacophore screening were subjected to molecular docking using standard precision mode of Glide which gave G-score value corresponding to each compound. The G-score value was calculated by taking into consideration the favorable van der Waals, coulombic, lipophilic and hydrogen-bonding interactions and penalizing for steric and buried polar clashes. Compounds with lower G-score were eliminated and the remaining compounds were subjected to extra precision docking. This was followed by interaction-based selection of 8 hits.

Qik prop descriptors

These hits obtained from the database were also subjected to further filter via Lipinski’s rule of five to identify compounds with potent EGFR inhibitory activity and favorable absorption, distribution, metabolism and excretion (ADME) properties. The ADME properties were calculated using QikProp. In the present study, QikProp was run in normal processing mode with default options.


Pharmacophore model

The active ligands were used to identify common pharmacophore hypotheses by following a tree-based partitioning technique that groups together similar pharmacophores according to their inter site distances. Six top-ranked five featured hypotheses were chosen for building QSAR equations (table 2).

Table 2: Summary of phase 3D-QSAR statistical results of six top-ranked hypotheses

Hypotheses SD r2 F RMSE q2 Pearson-R
ADHRR.442 0.2491 0.9851 635.1 0.6786 0.7718 0.9206
ADHRR.256 0.2131 0.9891 871.2 0.7893 0.6913 0.9008
ADHRR.409 0.1641 0.9935 1475.9 0.6066 0.8177 0.9251
ADHRR.187 0.1988 0.9905 1002.2 0.4994 0.8764 0.9507
AADHR.359 0.2043 0.99 948.5 0.5982 0.8227 0.9307
AAHRR.313 0.1952 0.9909 1040.1 0.6427 0.7954 0.9059

The very first pharmacophore model for the ATP-competitive inhibitors of EGFR was proposed [15] way back in 1999. This model highlights four important structural features, a hydrogen bond acceptor, a donor, a hydrophobic site and a ring. The crystal structure of EGFR tyrosine kinase domain, along with ATP-competitive inhibitor, 1M17 published later focuses on the acceptor forming a hydrogen bond with Met769 along with hydrophobic interactions. Other pharmacophores reported in the literature for EGFR inhibitors include AADHR reported for the class of phenylureas [16], AAADRR for pyrrolopyrimidines [17] and AHHHR for a diverse class containing anilinoquinazolines, thiazolo and thienopyrimidines [18]. All these pharmacophores highlight the importance of an acceptor and a ring structure apart from other hydrophobic interactions.

3D-QSAR results

With known experimental activities, a 3D-QSAR model was created for each hypothesis, using ligands that are aligned to the associated pharmacophore on five points. Since the dataset consists of significant diversity, a pharmacophore-based 3D-QSAR model was developed which considered the sites on a molecule that is matched to the hypothesis. Statistical validity and robustness of the QSAR models can be evaluated using various parameters including SD, F, Pearson-R, r2 values of the training set and q2 values of the test sets. We started the evaluation of QSAR equations with a look at SD, r2 and F values. These values are associated with training set results. Table 2 shows six pharmacophore hypotheses which indicated satisfactory results of the training sets. Pearson-R and q2values, which are associated with the test set, were taken into consideration in the next step. A stringent criterion of the difference between r2 and q2 values of less than or equal to 0.2 was applied to filter the hypotheses further. This resulted in an elimination of hypothesis ADHRR.256. The final filter of RMSE value less than 0.5 gave the best pharmacophore model, ADHRR.187. Fig. 2a and 2b shows a plot of predicted activity versus actual activity for training and test set, respectively, based on the best hypothesis ADHRR.187.

Fig. 2: Scatter plot for the predicted activity and actual activity values for EGFR inhibitors in (a) Training set and (b) Test set

The best predictive five featured hypothesis (ADHRR.187) developed in the present study consisted of one hydrogen-bond acceptor, one hydrogen-bond donor, one hydrophobic feature (H) and two aromatic ring features. The pharmacophore hypothesis and inter-site distances between the pharmacophoric sites is depicted in fig. 3 and shown in table 3.

Fig. 3: Intersite distance (Å) in the geometry of the pharmacophore hypothesis (ADHRR.187). Red spheres with vectors represent acceptor feature, blue spheres with vectors represent donor feature, green spheres represent hydrophobic feature and orange rings represent aromatic feature

When the highest active compound was mapped on the generated hypothesis (fig. 4), the quinazoline nitrogen mapped on the acceptor feature, secondary amine nitrogen mapped on the donor whereas Fluorine mapped on the hydrophobic feature. A benzene ring and a quinazoline ring was mapped on the ring residues.

Table 3: Distances of all pharmacophore hypothesis sites of 3D-QSAR model ADHRR.187

Site 1 Site 2 Distance (A °)
A2 D5 3.236
A2 H6 6.807
A2 R10 3.689
A2 R11 4.035
D5 H6 5.117
D5 R10 3.68
D5 R11 3.356
H6 R10 8.763
H6 R11 3.133
R10 R11 6.561

Fig. 4: Top ranked pharmacophore model mapped on the highest active compound (Compound 91)

Database screening and docking

The generated pharmacophore hypothesis was used as a 3D querry for searching the ZINC database. This resulted in identification of 372 hits from the ZINC database. Docking was used as a second filter to arrive at probable hits.

To validate the docking protocol, Docking began with redocking of the inhibitor from the co-crystallized ligand. The RMSD between the docked pose and the crystallized pose of the same ligand was found to be 1.12, thus validating the docking protocol. This was followed by docking of the hits obtained using pharmacophore screening into the crystal structure of EGFR (PDB code: 1M17). The criteria for selection of best 10 compounds out of the 372 hits were based on obtaining a good dock score (>7.0) along with the docked structure retaining most of the important interactions which are seen in the crystal structure. Hydrogen bonding interaction with Met769 was defined as the most important interaction [19] and the hits in which this interaction was absent were rejected. The procedure gave eight compounds which fitted in these criteria. Fig. 5 shows structures of the identified compounds from ZINC database and fig. 6 shows the docked images of the same. Table 4 shows interaction of these eight compounds with the active site of the tyrosine kinase domain of EGFR.

Fig. 5: Structures of ZINC database leads

Fig. 6: Two dimensional representations of lead molecules in their binding orientation in the active site of EGFR

Table 4: Summary of docking results for 8 leads

Compounda Glide score Interacting residues
ZINC44012931 -7.86 Met769, Leu768, Gly772, Leu694, Leu820, Pro770, Ala719, Ile720, Ile765
ZINC32650290 -7.20 Met769, Leu768, Gly772, Leu694, Leu820, Pro770, Ala719, Ile720, Ile765, Ile753
ZINC32544573 -8.54 Met769, Asp831, Leu768, Gly772, Leu694, Leu820, Pro770, Ala719, Ile720, Ile765, Leu753
ZINC017025205 -9.22 Met769, Asp831, Lys721, Thr830, Leu768, Gly772, Leu694, Leu820, Pro770, Ala719
ZINC04822304 -9.48 Met769, Glu738, Thr766, Leu768, Gly772, Leu694, Leu820, Pro770, Ala719
ZINC01680588 -9.94 Met769, Asp831, Thr766, Leu768, Gly772, Leu694, Leu820, Pro770, Ala719
ZINC00611711 -7.16 Met769, Leu768, Gly772, Leu694, Leu820, Leu753, Pro770, Ala719, Ile720, Ile765
ZINC12804309 -7.03 Met769, Leu768, Gly772, Leu694, Pro770, Ala719

aLigand Ids for ZINC database

Qik Prop descriptors

To assess drug likeness of the identified compounds, these compounds were subjected to calculation of a variety of ADME properties. These included partition coefficient, water solubility, percent human oral absorption, permeation through the blood brain barrier, and cell permeability. The results as listed in table 5 indicate that compounds ZINC04822304 and ZINC01680588 show a significantly lower oral absorption. These compounds belong to the class of pyrrolopyrimidine and imidazopyrimidine respectively. Both these compounds contain a free amino group in their structures. The other identified hits belong to chemical classes of chalcone, pyrazolopyrimidine, anilinoquinazoline, imidazopyrimidine-2,4-dione and benzofuranopyrimidine.

Table 5: ADME properties of the 8 lead molecules by using QikProp

Compounda QPlogPo/wb QPlogSc QPlogBBd QPPCacoe QPPMDCKf % human oral absorptiong
ZINC44012931 2.935 -3.181 -0.157 1995.746 2639.456 100
ZINC32650290 3.273 -4.391 -0.123 2456.422 1306.837 100
ZINC32544573 2.664 -4.842 -0.694 968.785 478.039 95.992
ZINC017025205 -0.241 -1.949 -0.982 305.972 236.556 70.025
ZINC04822304 -1.082 -2.127 -1.901 43.63 21.247 49.956
ZINC01680588 -0.469 -2.402 -1.433 94.18 94.17 59.53
ZINC00611711 0.62 -3.137 -1.328 173.587 74.535 70.658
ZINC12804309 4.484 -5.399 0.306 5457.062 3096.815 100

aZINC database ligand IDs, bPredicted octanol/water partition co-efficient log p (acceptable range: −2.0 to 6.5), cPredicted aqueous solubility; S in mol/l (acceptable range: −6.5 to 0.5), dPredicted BBB permeability (acceptable range: −3 to 1.2), ePredicted Caco-2 cell permeability in nm/s (acceptable range: <25 is poor and>500 is great), fPredicted apparent MDCK cell permeability in nm/s (acceptable range: <25 is poor and>500 is great), gPercentage of human oral absorption (<25% is poor and>80% is high).


In conclusion, we have developed a pharmacophore model for ATP-competitive inhibitors of EGFR, using a diverse set of ligands. A QSAR model was developed using this pharmacophore hypothesis. The development of the pharmacophore model and QSAR model involved use of stringent statistical criteria coupled with knowledge of the structure of the receptor. This was then used for screening the ZINC database to obtain probable hits. The hits were further filtered by docking and calculation of ADME properties to arrive at six molecules with probable EGFR inhibitory activity.


The authors sincerely acknowledge DBT (BT/PR14373/MED/ 30/530/2010) for funding received.


The authors declare no conflict of interest


  1. Sawyers C. Targeted cancer therapy. Nature 2004;432:294-7.
  2. Arteaga C. Targeting HER1/EGFR: a molecular approach to cancer therapy. Paper presented at the Seminars in oncology; 2003.
  3. Yarden Y. The EGFR family and its ligands in human cancer: signalling mechanisms and therapeutic opportunities. Eur J Cancer 2001;37:3-8.
  4. Ritter CA, Arteaga CL. The epidermal growth factor receptor-tyrosine kinase: A promising therapeutic target in solid tumors. Paper presented at the Seminars in oncology; 2003.
  5. Mass RD. The HER receptor family: a rich target for therapeutic development. Int J Radiation Oncol Biol Physics 2004;58:932-40.
  6. Lurje G, Lenz HJ. EGFR signaling and drug discovery. Oncol 2010;77:400-10.
  7. Eberhard DA, Johnson BE, Amler LC, Goddard AD, Heldens SL, Herbst RS, et al. Mutations in the epidermal growth factor receptor and in KRAS are predictive and prognostic indicators in patients with non-small-cell lung cancer treated with chemotherapy alone and in combination with erlotinib. J Clin Oncol 2005;23:5900-9.
  8. Burris HA, Hurwitz HI, Dees EC, Dowlati A, Blackwell KL, O'Neil B, et al. Phase I safety, pharmacokinetics, and clinical activity study of lapatinib (GW572016), a reversible dual inhibitor of epidermal growth factor receptor tyrosine kinases, in heavily pretreated patients with metastatic carcinomas. J Clin Oncol 2005;23:5305-13.
  9. Aparna V, Rambabu G, Panigrahi SK, Sarma JARP, Desiraju GR. Virtual screening of 4-anilinoquinazoline analogues as EGFR kinase inhibitors: importance of hydrogen bonds in the evaluation of poses and scoring functions. J Chem Inf Model 2005;45:725-38.
  10. Liu LT, Yuan T, Liu H, Chena S, Wu Y. Synthesis and biological evaluation of substituted 6-alkynyl-4-anilinoquinazoline derivatives as potent EGFR inhibitors. Bioorg Med Chem Lett 2007;17:6373–7.
  11. Kubinyi H. ed. 3D QSAR in drug design: volume 1: Theory methods and applications. Vol. 1. Springer Science & Business Media; 1993.
  12. Dixon SL, Smondyrev AM, Rao SN. PHASE: a novel approach to pharmacophore modeling and 3d database searching. Chem Biol Drug Des 2006;67(5):370-2.
  13. Kurogi Y, Guner OF. Pharmacophore modeling and three-dimensional database searching for drug design using catalyst. Curr Med Chem 2001;8(9):1035-55.
  14. La Motta C, Sartini S, Tuccinardi T, Nerini E, Da Settimo F, Martinelli A. Computational studies of epidermal growth factor receptor: docking reliability, three-dimensional quantitative structure-activity relationship analysis, and virtual screening studies. J Med Chem 2009;52:964-75.
  15. Traxler P, Green J, Mett H, Urs Se´quin, Furet P. Use of a pharmacophore model for the design of EGFR tyrosine kinase inhibitors: isoflavones and 3-phenyl-4(1H)-quinolones. J Med Chem 1999;42:1018-26.
  16. Rakesh JG, Shaheera B, Aarumugam P, Vaiyshnavi R, Thulasibabu R. Emerging drug discovery paradigm in non small cell lung cancer: pharmacophore modeling, atom-based 3D-QSAR and virtual screening of novel EGFR inhibitors. J Drug Discovery Ther 2014;2(23):9-17.
  17. Sudha A, Srinivasan P, Rameshthangam P. Exploration of potential EGFR inhibitors: a combination of pharmacophore-based virtual screening, atom-based 3D-QSAR and molecular docking analysis. J Recept Signal Transduction Res 2014;35(2):137-48.
  18. Gupta AK, Bhunia SS, Balaramnavar VM, Saxena AK. Pharmacophore modelling, molecular docking and virtual screening for EGFR (HER 1) tyrosine kinase inhibitors. SAR QSAR Environ Res 2011;22(3-4):239-63.
  19. Hou T, Zhu L, Chen L, Xu X. Mapping the binding site of a large set of quinazoline type EGF-R inhibitors using molecular field analyses and molecular docking studies. J Chem Inf Comp Sci 2003;43(1):273-87.