WaterSwap Analysis, a Computation-based Method for the Discovery of Effective and Stable Binding Compounds for Mutant EGFR Inhibition

doi:10.21203/rs.3.rs-981656/v1

Download PDF

Research Article

WaterSwap Analysis, a Computation-based Method for the Discovery of Effective and Stable Binding Compounds for Mutant EGFR Inhibition

https://doi.org/10.21203/rs.3.rs-981656/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this latest preprint version

The epidermal growth factor receptor of the tyrosine kinase family has been largely targeted in mutations associated with non-small cell lung cancer. EGFR inhibitors have been produced that bind allosterically to the C797S mutant EGFR enzyme. Here, the Waterswap tool has been used for the interpretation and visualization of trajectories of mutant EGFR-ligand complexes. Virtual screening of the generated compounds has been carried out along with its molecular docking and ADMET analysis. Out of the generated library of compounds, the top 15 have been selected. Waterswap calculated the binding free energies of the compounds and thermodynamic properties of the enumerated compounds were compared with that of standard EAI045. It was observed that styrylquinoline stabilized better than EAI045. The results show that Waterswap analysis offers a promising new path in the hunt for improved tools for analyzing and visualizing molecular driving forces in protein-ligand complex simulations.

Molecular Genetics

EGFR

Waterswap

Virtual screening

Molecular docking

ADMET

Styrylquinoline

Cancer is the second biggest cause of death in the world, and it is one of the primary causes of morbidity and mortality. In 2015, cancer claimed the lives of 8.8 million people. Cancer is responsible for about one-sixth of all fatalities worldwide [1]. Despite decades of attempts to find and develop small molecule anticancer medications, the demand for novel antitumor medicines with enhanced tumour efficiency, selectivity, and safety remains critical [2, 3].

One of the most significant enzymes in the tyrosine kinase family is the Epidermal Growth Factor Receptor. EGFR is a transmembrane receptor that is found on the cell surface. The EGFR enzyme is classified as (ErbB1, HER1), ErbB2 (HER2, neu in rodents), ErbB3 (HER3), and ErbB4 (HER4) based on structural similarities [4]. EGFR and its members have become a popular target for anticancer therapy, particularly in the treatment of non-small cell lung cancer [5]. Activating mutations in the EGFR gene had a significant impact on non-small cell lung cancer therapy approaches (NSCLC). The recent identification of novel EGFR inhibitors was critical in overcoming such mutations [6].

Drugs' affinity for their targets is described by thermodynamic characteristics. As a result, accurate thermodynamic property calculations have long been a goal of computational chemistry. Hence, accurate thermodynamic property calculations have long been a computational chemistry goal. In the industry, study of solvation patterns surrounding pharmaceuticals and binding sites has become a major goal for molecular modelling [7–9]. According to an estimate, hydrophobic surfaces cover 75% of surface drug sites [10], implying that the hydrophobic effect is one of the most energetic contributions to ligand binding and recognition. Individual water molecules that serve a crucial stabilising role, for example, by building bridges between the ligand and the receptor, have been shown by high-resolution co-crystallized structures of protein-ligand complexes [11, 12]. Improved binding affinity has also been linked to a network of H-bond-enriched water molecules around the ligand [13–15].

The water-swap reaction coordinate was described by Woods et al. as a new reaction coordinate (WSRC) [16]. The ligand and water are exchanged from being bound to the protein to being solvated by a periodic box of bulk water in this approach. Because no molecules are moved into a vacuum, the enormous decoupling free energy associated with double decoupling procedures is avoided. Furthermore, a new approach for restricting the identity of solvent clusters within bulk has been established, allowing the swap to be accomplished without the requirement for any restrictions on the ligand or swapped water, and this strategy applies to both buried and solvent-exposed binding sites [16]. The concept behind the water-swap reaction coordinate is that it permits the identity of water clusters to be limited during simulation. During a simulation, this identity constraint [17] works by modifying how water molecules are labeled. We used computation-based methods such as virtual screening, molecular docking, and WaterSwap analysis to find an effective compound with stable binding and mutant EGFR inhibition effect in this study.

Insilico methodology

ChemT tool was used for library generation and enumeration. The molecular docking tool, AutoDock Vina, was used for docking studies of all compounds on the EGFR enzyme. SwissADME tool was used to determine ADMET studies. Waterswap study was performed on trial version of Cresset’s Flare software.

ChemT and Library Enumeration

ChemT, a Biochemcore-based software, was used to enumerate the library. In Marvin Sketch, the R-group position on the template styrylquinoline. The styrylquinoline was then assigned to the identified R-groups through the enumeration model of the software. EAI045 input determined the standard template. Additionally, few available pharmacokinetic properties such as molecular weight, H-bond donors, and Lipinski rule of five was set for designing drug-like compounds [18].

Molecular Docking

To investigate the binding pocket, residue interactions, and predicted binding affinity of the enumerated compounds, molecular docking protocol was carried out using AutoDock tools. Through Chimera modeler, the protein was fetched from the open-source PDB libraries. Explicit hydrogen bonds were added followed by the removal of water molecules. Using Gastegier charges, ions were removed and finally, the protein was minimized using the Amber94 force field. The total number of loops was set to 1000 to ensure maximum minimization of the protein. The minimization was terminated when the energy converged or the root mean square deviation (RMSD) reached a maximum cut-off of 0.30 Å. Ramachandran plot was set to check the total number of disallowed residues. Through AutoDock 4.0 GUI, the minimized PDB was converted to a PDBQT file followed by generation of the grid. To perform flexible docking, torsions were chosen in the selected residues and saved as flex.PDBQT. The configuration file, which was used to outline the grid area as well as the protein and ligand input, was defined and saved as a text document. All the ligands were prepared using Gastegier charges and saved as a PDBQT file. The files saved in the fold were called through the command line under the vina environment to start the docking process. The results for each input ligand were saved as a text file consisting of docking scores and RMSD values of each ligand conformation which was used for the analysis of data. The output file was then converted into Maestro readable extension for image saving [19–22].

WaterSwap methodology

The PDB: 5D41 was downloaded initially and full preparation on proteins and ligands was carried out. Under flare protein preparation wizard, at pH 7.0 small side chains movement were allowed, atoms were removed from residues with incomplete backbone, along with filling of 1 and 2 residue gaps at 6.00 A° active site size. Then the protein was minimized under the normal calculation method with a gradient cutoff of 0.200 kcal/mol/A with 1300 maximum iterations.

The ligands were loaded in the sdf format for docking. Docking was performed under the normal calculation method under Autodock flexible docking protocol. The ligands to dock were already selected, the template ligand was defined, the co-crystallized ligand was selected for docking grid preparation, the minimized protein was selected and docking calculation was started [16, 23].

For WaterSwap calculation, the protein along with docked ligand was prepared again and minimized as per the above method. The calculation method was kept custom and no changes were made. The small molecule force field was kept as AMBER with GAFF version with AM1-BCC as the charge method. The solvate box buffer was kept at 10.0 A. Open MM platform was GPU, minimized energy tolerance as 0.25 (kcal/mol), the system was equilibrated for 100 ps. For Waterswap, the iterations were kept as 1000 -1300 (1300 for docked ligand), maximum threads were kept on Auto, with Water monitor distance as 7.00 A, the protein was selected, the ligand was selected along with the proper chain, and the calculation was started.

Computational Methods

Database building and Library Screening

The most potent substituted styrylquinoline with the highest antitumour activity was taken from the literature [24], and by modification, at R₁ position a library of around 2000 compounds was generated virtually. Based on the obtained fitness score using ChemT, the top 1000 hits were screened for further pre-ADMET analysis. All 1000 hits were screened through Lipinski’s Rule of 5 and around 100 hits were chosen for further study fitting with the criteria. In the virtual screening workflow protocol, all 100 hits were subjected to a flexible docking approach using AutoDock. A score cutoff of -9.20 kcal/mol for docking discarded lower-scoring compounds selecting 15 highest scoring compounds. Top 15 highest scoring, best pocket fitting, and possessing similar binding interactions with that of standard EAI045, novel and potent compounds were obtained through virtual screening and were selected for further studies.

Molecular Docking

Molecular docking was carried out on mutant EGFR enzyme (PDB ID: 5D41) to check the binding pattern and interactions. It was seen that all 15 molecules docked correctly into the allosteric binding pocket of the protein with common interacting amino acid residues LYS745, LEU788, ASP855 and PHE856. From Figures 1a and 1b, it was observed that all compounds fit into the hydrophobic pocket of the protein interacting with the common amino acid residues. The combined docking scores along with structures are shown in Table 1.

Figure 1b Docked images of top 15 compounds over T790M/C797S mutant EGFR (PDB: 5D41). Blue lines in the figures indicate π-π stacking interaction and pink lines indicate hydrogen bonding interaction.

ADMET studies

SwissADME website (http://www.swissadme.ch/index.php) was used to quantify the physicochemical and pharmacokinetic properties of the synthesized derivatives with respect to EAI045. Table 2 shows the analysis of the pharmacokinetic parameters needed for ADMET study of compounds KSK (1-15).

Table 2

ADMET study of top 15 compounds with EAI045.
Compound	#HB acceptor	#HB donors	TPSA^a	MLOGP^b	GI absorption	BBB permeant	Violation of Lipinski’s Rule	Bio -availability Score	Synthetic Accessibility
EAI045	5	2	110.77	1.93	High	No	0	0.55	3.28
KSK-1	2	2	70.83	3.65	High	No	0	0.55	3.37
KSK-2	5	2	88.68	0.81	High	No	0	0.55	3.99
KSK-3	1	1	42.15	3.75	High	No	0	0.55	3.33
KSK-4	2	1	56.4	3.15	High	No	0	0.55	4.16
KSK-5	2	1	36.36	3.85	High	Yes	0	0.55	3.05
KSK-6	2	0	58.5	3.78	High	No	0	0.55	3.94
KSK-7	3	2	75.27	3.38	High	No	0	0.55	3.34
KSK-8	3	1	100.88	3.91	Low	No	0	0.56	4.39
KSK-9	3	1	70.89	3.41	High	No	0	0.55	4.47
KSK-10	2	1	41.05	3.54	High	Yes	0	0.55	3.54
KSK-11	2	0	29.02	3.18	High	No	0	0.55	3.25
KSK-12	3	1	45.59	3.41	High	Yes	0	0.55	4.19
KSK-13	2	0	67.82	3.98	High	No	0	0.55	3.28
KSK-14	3	1	36.36	3.75	Low	No	0	0.55	3.24
KSK-15	2	1	55.04	3.33	High	Yes	0	0.55	3.41
^aTPSA 40Å - 160Å
^bMLOGP must be less than 4.15

The drug must pass Lipinski’s rule to show good oral absorption in humans [25, 26]. All 15 compounds do not show Lipinski’s rule violation, passing the criteria. The MLOGP values were also within the acceptable range, showing the compounds' correct drug absorption and distribution within the body. Almost all the compounds show good GI absorption and some of them were blood-brain barrier permeant. Through bioavailability scores, it was indicated that a potent compound is available at the site of action. Overall, the compounds showed satisfactory pharmacokinetic parameters range defined for human use along with good docking scores [27–29].

Waterswap

When dealing with the cavitation and large-value problem of double decoupling, we used Waterswap, an absolute binding free energy method. Continuum solvent methods do not include the molecular detail of protein–water, ligand–water, and protein–water–ligand interaction interactions because they use an explicit model of water [23].

It uses a single simulation to replace the ligand-bound protein with bulk water molecules of the same shape and volume as shown in Figure 2. An increase in protein-ligand interactions occurs when the protein's resulting cavity is filled with a cluster of water molecules. With only one reaction coordinate, it's possible to get a more accurate free energy value because it's not the product of two huge numbers [30].

Moreover, the beneficial waters around the protein-ligand complex are also identified to obtain a stable pocket. These beneficial waters also become a potential bridge in creating interactions between likely residues and ligand. Both proteins favourable and ligand favourable waters are important to identify to assess the pocket interactions and their stability [31]. Two different types of colours are given for the protein depending on its preference i.e., the protein present in green colour states that it prefers ligand, and protein present in red colour states that it prefers water as indicated in Figure 3. The water box identified by the system further provides the required binding free energy for the protein and ligand [16].

To analyze Waterswap studies, binding free energy (∆G), shown in Figure 4 comes into play describing the stability of the waters around the pocket. If ∆Gprotein in the complete simulation decreases, it does not stabilize with the set number of iterations, the possibility of the existence of unfavourable water distorting the pocket intensity is high. In such cases, the ligand would exit the simulation box resulting in an inactive compound [32]. The similar identity of ∆Gwater is identified through drastic increase or instability of the waters resulting in distortions in the protein. This could also result in terminating the protein-ligand interactions.

EGFR protein becomes the perfect example of identifying and comparing the potential waters for the binding affinity and checking the potency of the complexed ligand. A mutant EGFR protein with known standard EAI045 is initially run for 1300 iterations using Waterswap to identify all the favourable waters and their binding free energy. Then the ligand is swapped to check the potential binding affinity and stability of the designed compounds concerning the identified favourable waters. If the binding free energy results obtained during the complete iterative simulations are similar or better than the results for standard drug EAI045, the compounds designed are the potential to withstand the waters and show potency.

Figure 5 describes the relative binding free energy graph during the 1300 iterations performed for EAI045 compared to the styrylquinoline compound identified from the literature. The blue line in the graph shows the average binding free energy change for the number of iterations. The thermodynamic properties compared in the Figure showed that the styrylquinoline stabilized better than the standard compound EAI045. Towards the end of the iterations, the graph showed stability after a drop in relative binding energy in styrylquinoline. Overall, all the attributes taken into account concludes the potential of styrylquinoline compound as an EGFR inhibitor.

The results show that Waterswap analysis offers a promising new path in the hunt for improved tools for analyzing and visualizing molecular driving forces in protein-ligand complex simulations.

The binding free energy calculations using the Waterswap analysis for interpreting and visualizing trajectories of EGFR-ligand complexes have been demonstrated. These visualizations provide precise molecular-level information that should help medicinal chemists construct better-binding ligands. Although we tested these decompositions on mutant EGFR, we believe they should apply to other protein-ligand systems as well. The thermodynamic properties of the enumerated compounds were compared with that of standard EAI045, and it was observed that styrylquinoline stabilized better than EAI045. The results show that Waterswap analysis offers a promising new path in the hunt for improved tools for analyzing and visualizing molecular driving forces in protein-ligand complex simulations.

Acknowledgement

The authors K. S. K and P. S. W. are thankful to Indian Council of Medical Research (ICMR) for providing funding in support of Senior Research Fellowship with Ref No: 3/2/2/33/2018/Online OncoFship/NCD-III. The authors are also thankful to The Head, Department of Chemical Technology, Dr. Babasaheb Ambedkar Marathwada University, Aurangabad 431 004 (MS), India, for providing the laboratory facility and Cresset for providing us with the free trial version of the Flare^TM WaterSwap software.

Declaration of Interest

Funding: The project was funded by Indian Council of Medical Research (ICMR) in support of Senior Research Fellowship with Ref No: 3/2/2/33/2018/Online OncoFship/NCD-III.

Conflicts of interest/Competing interests: The authors declare no competing ﬁnancial interest.

Availability of data and material: The data that support the findings of this study are available from the corresponding author, P.S.W, upon reasonable request.

Code availability: Not applicable.

Authors' contributions: All the mentioned authors contributed equally in the preparation of manuscript.

Hoelder S, Clarke PA, Workman P (2012) Discovery of small molecule cancer drugs: successes, challenges and opportunities. Mol Oncol 6:155–176
Dayyani F, Gallick GE, Christopher CL, Corn PG (2011) Novel Therapies for Metastatic Castrate-Resistant Prostate Cancer. J Natl Cancer Inst Monogr 103:1665–1675
You WK, Sennino B, Williamson CW, Falcón B, Hashizume H, Yao LC, Aftab DT, McDonald DM (2011) VEGF and c-Met blockade amplify angiogenesis inhibition in pancreatic islet cancer. Cancer Res 71:4758–4768
Wieduwilt MJ, Moasser MM (2008) The epidermal growth factor receptor family: Biology driving targeted therapeutics. Cell Mol Life Sci 65(10):1566–1584
Sharma SV, Bell DW, Settleman J, Haber DA (2007) Epidermal growth factor receptor mutations in lung cancer. Nat Rev Cancer 7(3):169–181
Santhi MP, Bupesh G, SenthilKumar V, Meenakumari K, Prabhu K, Sugunthan S, Manikandan E, Saravanan K (2016) Anticancer activity and drug likeliness of quinoline through insilico docking against cervical and liver cancer receptors. Indian J Med Res Pharm 3:83–93
Mason JS, Bortolato A, Weiss DR, Deflorian F, Tehan B, Marshall FH (2013) High End GPCR Design: Crafted Ligand Design and Druggability Analysis Using Protein Structure, Lipophilic Hotspots and Explicit Water Networks. Silico Pharmacol 1:23–46
Bodnarchuk MS (2016) Water, Water, Everywhere··· It’s Time to Stop and Think. Drug Discov 21(7):1139–1146
Graves AP, Wall ID, Edge CM, Woolven JM, Cui G, Le Gall A, Hong X, Raha K, Manas ES (2017) A Perspective on Water Site Prediction Methods for Structure Based Drug Design. Curr Top Med Chem 17(23):2599–2616
Snyder PW, Lockett MR, Moustakas DT, Whitesides GM (2014) Is It the Shape of the Cavity, or the Shape of the Water in the Cavity? Eur Phys J Spec Top 223(5):853–891
Lu Y, Wang R, Yang CY, Wang S (2007) Analysis of Ligand-Bound Water Molecules in High-Resolution Crystal Structures of Protein– Ligand Complexes. J Chem Inf Model 47(2):668–675
García-Sosa AT (2013) Hydration Properties of Ligands and Drugs in Protein Binding Sites: Tightly-Bound, Bridging Water Molecules and Their Effects and Consequences on Molecular Design Strategies. J Chem Inf Model 53(6):1388–1405
Michel J, Tirado-Rives J, Jorgensen WL (2009) Energetics of Displacing Water Molecules from Protein Binding Sites: Consequences for Ligand Optimization. J Am Chem Soc 131:15403–15411
Krimmer SG, Cramer J, Betz M, Fridh V, Karlsson R, Heine A, Klebe G (2016) Rational Design of Thermodynamic and Kinetic Binding Profiles by Optimizing Surface Water Networks Coating Protein-Bound Ligands. J Med Chem 59(23):10530–10548
Snyder PW, Mecinovic J, Moustakas DT, Thomas SW, Harder M, Mack ET, Lockett MR, Heroux A, Sherman W, Whitesides GM (2011) Mechanism of the Hydrophobic Effect in the Biomolecular Recognition of Arylsulfonamides by Carbonic Anhydrase. Proc Natl Acad Sci U S A 108(44):17889–17894
Woods CJ, Malaisree M, Hannongbua S, Mulholland AJ (2011) A water-swap reaction coordinate for the calculation of absolute protein–ligand binding free energies. J Chem Phys 134:1–13
Tyka MD, Sessions RB, Clarke AR (2007) Absolute Free-Energy Calculations of Liquids Using a Harmonic Reference State. J Chem Phys 111:9571–9580
Abreu RMV, Froufe HJC, Daniel POM, Queiroz MJRP, Ferreira ICFR (2011) ChemT, an open-source software for building template-based chemical libraries. SAR QSAR Environ Res 22(5-6):603–610
Karnik KS, Sarkate AP, Tiwari SV, Azad R, Wakte PS (2021) Free energy perturbation guided Synthesis with Biological Evaluation of Substituted Quinoline derivatives as small molecule L858R/T790M/C797S mutant EGFR inhibitors targeting resistance in Non-Small Cell Lung Cancer (NSCLC). Bioorg Chem 115:105226–105238
Trott O, Olson AJ (2010) AutoDock Vina: improving the speed and accuracy of docking with a new scoring function, efficient optimization and multithreading. J Comput Chem 31:455–461
Forli S, Huey R, Pique ME, Sanner M, Goodsell DS, Olson AJ (2016) Computational protein-ligand docking and virtual drug screening with the AutoDock suite. Nat Protoc 11:905–919
Seeliger D, Groot BL (2010) Ligand docking and binding site analysis with PyMOL and Autodock/Vina. J Comput Aided Mol Des 24:417–422
Woods CJ, Malaisree M, Michel J, Long B, McIntosh-Smith S, Mulholland AJ (2014) Rapid decomposition and visualisation of protein–ligand binding free energies by residue and by water. Faraday Discuss 169:477–499
Magda AA, El-Sayed, Walaa M, El-Husseiny, Abdel-Aziz NI, ElAzab AS, Abuelizz HA, Abdel-Aziz ALM (2018) Synthesis and biological evaluation of 2- styrylquinolines as antitumour agents and EGFR kinase inhibitors: molecular docking study. J Enzyme Inhib Med Chem 33(1):199–209
Angelis ID, Turco L (2011) Caco- 2 cells as a model for intestinal absorption. Curr ProtocToxicol 47:20–26
Lipinski CA, Lombardo F, Dominy BW, Feeney PJ (2012) Experimental and computational approaches to estimate solubility and permeability in drug discovery and development settings. Adv Drug Deliv Rev 64:4–17
Lokwani DK, Azad R, Sarkate AP, Reddanna P, Shinde DB (2015) Structure Based Library Design (SBLD) for new 1,4-dihydropyrimidine scaffold as simultaneous COX-1/COX-2 and 5-LOX inhibitors. Bioorg Med Chem 23:4533–4543
Bhosle MR, Khillare LD, Mali JR, Sarkate AP, Lokwani DK, Tiwari SV (2018) DIPEAc promoted one-pot synthesis of dihydropyrido[2,3-d:6,5-d′]dipyrimidinetetraone and pyrimido[4,5-d]pyrimidine derivatives as potent tyrosinase inhibitors and anticancer agents: in vitro screening, molecular docking and ADMET predictions. New J Chem 42:18621–18632
Tiwari SV, Seijas JA, Vazquez-Tato MP, Sarkate AP, Karnik KS, Nikalje APG (2018) Ionic Liquid-Promoted Synthesis of Novel Chromone-Pyrimidine Coupled Derivatives, Antimicrobial Analysis, Enzyme Assay, Docking Study and Toxicity Study. Molecules 23:1–23
Woods CJ, Malaisree M, Pattarapongdilok N, Sompornpisut P, Hannongbua S, Mulholland AJ (2012) Long Time Scale GPU Dynamics Reveal the Mechanism of Drug Resistance of the Dual Mutant I223R/H275Y Neuraminidase H1N1-2009 Influenza Virus. Biochemistry 51:4364–4375
Kurdekar V, Giridharan S, Subbarao J, Nijaguna MB, Periasamy J, Boggaram S, Bharatham K, Potluri V, Shivange AV, Sadasivam G, Padigaru M, Venkitaraman AR Structural features underlying the activity of benzimidazole derivatives that target phosphopeptide recognition by the tandem BRCT domain of the BRCA1 protein. bioRxiv preprint. https://doi.org/10.1101/555623
Liu J, Liu L, Tian Z, Li Y, Shi C, Shi J, Wei S, Zhao Y, Zhang C, Bai B, Chen Z, Zhang, H In silico discovery of a small molecule suppressing lung carcinoma A549 cells proliferation and inducing autophagy via mTOR pathway inhibition. Mol Pharm. https://doi.org/10.1021/acs.molpharmaceut.8b00996

Due to technical limitations, table 1 is only available as a download in the Supplemental Files section.

Download PDF

Version 1

posted

You are reading this latest preprint version

WaterSwap Analysis, a Computation-based Method for the Discovery of Effective and Stable Binding Compounds for Mutant EGFR Inhibition

Status:

Version 1

Abstract

Figures

Introduction

Materials And Methods

ChemT and Library Enumeration

Molecular Docking

WaterSwap methodology

Result And Discussion

Computational Methods

Database building and Library Screening

Molecular Docking

ADMET studies

Waterswap

Conclusion

Declarations

References

Tables

Supplementary Files

Status:

Version 1