Quinoxalinones as A Novel Inhibitor Scaffold for EGFR (L858R/T790M/C797S) Tyrosine Kinase: Molecular Docking, Biological Evaluations, and Computational Insights

Combating acquired drug resistance of EGFR tyrosine kinase (TK) is a great challenge and an urgent necessity in the management of non-small cell lung cancers. The advanced EGFR (L858R/T790M/C797S) triple mutation has been recently reported, and there have been no specific drugs approved for this strain. Therefore, our research aimed to search for effective agents that could impede the function of EGFR (L858R/T790M/C797S) TK by the integration of in silico and in vitro approaches. Our in-house quinoxalinone-containing compounds were screened through molecular docking and their biological activity was then verified by enzyme- and cell-based assay. We found that the four quinoxalinone-containing compounds including CPD4, CPD15, CPD16, and CPD21 were promising to be novel EGFR (L858R/T790M/C797S) TK inhibitors. The IC50 values measured by the enzyme-based assay were 3.04 ± 1.24 nM; 6.50 ± 3.02 nM,10.50 ± 1.10 nM; and 3.81 ± 1.80 nM, respectively, which are at a similar level to a reference drug; osimertinib (8.93 ± 3.01 nM). Besides that, they displayed cytotoxic effects on a lung cancer cell line (H1975) with IC50 values in the range of 3.47 to 79.43 μM. In this proposed study, we found that all screened compounds could interact with M793 at the hinge regions and two mutated residues including M790 and S797; which may be the main reason supporting the inhibitory activity in vitro. The structural dynamics revealed that the screened compounds have sufficient non-native contacts with surrounding amino acids and could be well-buried in the binding site’s cleft. In addition, all predicted physicochemical parameters were favorable to be drug-like based on Lipinski’s rule of five, and no extreme violation of toxicity features was found. Altogether, this study proposes a novel EGFR (L858R/T790M/C797S) TK inhibitor scaffold and provides a detailed understanding of compounds’ recognition and susceptibility at the molecular level.


Introduction
Cancers pose a public health burden and challenge to all nations and remain a major cause of mortality and morbidity. Lung cancer is the most common type of cancer worldwide. It has been reported that more than 2.2 million were new cases; among these non-small cell lung cancer (NSCLC) accounts for 84% of all newly diagnosed cases by phenyliodine (III) diacetate (PIDA)-induced oxidative C-N bond coupling. Compounds exhibiting a comparable level of binding affinity to osimertinib were then selected to be tested experimentally. The potential compounds were predicted for their drug-likeness properties and toxicity features. Then, the MD simulations were employed to gain insights into binding recognition and susceptibility, guiding the mechanism of inhibitory action at the atomic level. Altogether, our proposed quinoxalinone derivatives could become viable candidates to be further developed as the effective fourth generation of EGFR inhibitors.

Docking-Based Virtual Screening
To rapidly identify the potential compounds from 30 in-house quinoxalinones capable of binding to the TK domain of EGFR (L858R/T790M/C797S), molecular docking was performed by employing the Autodock Vina XB software package [34,35]. As shown in Figure 2, the predicted binding energy of quinoxalinone-containing compounds was in a range −7.8 to −5.4 kcal/mol whilst the osimertinib drug was predicted to be −7.4 kcal/mol. Selection of potential compounds was based on the predicted binding energy, which is theoretically proportional to the dissociation constant (kd) and is widely used in several computer-aided drug discovery campaigns. For this purpose, compounds showing binding energy lower than −7.0 kcal/mol were selected which resulted in four quinoxalinonecontaining compounds including CPD4, CPD15, CPD16, and CPD21. Note that their chemical structures are shown in Figure 3.  [4,33]. (adapted from Refs. [4,33] and created with https://BioRender.com, accessed on 30 October 2022).

Docking-Based Virtual Screening
To rapidly identify the potential compounds from 30 in-house quinoxalinones capable of binding to the TK domain of EGFR (L858R/T790M/C797S), molecular docking was performed by employing the Autodock Vina XB software package [34,35]. As shown in Figure 2, the predicted binding energy of quinoxalinone-containing compounds was in a range −7.8 to −5.4 kcal/mol whilst the osimertinib drug was predicted to be −7.4 kcal/mol. Selection of potential compounds was based on the predicted binding energy, which is theoretically proportional to the dissociation constant (k d ) and is widely used in several computer-aided drug discovery campaigns. For this purpose, compounds showing binding energy lower than −7.0 kcal/mol were selected which resulted in four quinoxalinonecontaining compounds including CPD4, CPD15, CPD16, and CPD21. Note that their chemical structures are shown in Figure 3.

Drug-Likeness
The potent quinoxalinone-containing compounds against EGFR (L858R/T790M/C797S) were predicted the drug-likeness character based on Lipinski's rule of five by inspecting their physicochemical properties (molecular weight (MW), the numbers of hydrogen bond donors (HBD) and acceptors (HBA), rotatable bond (RB), polar surface area (PSA) and Log P). According to Table 1, the predictive results revealed that all compounds could confer the drug-likeness property and obey the acceptable value within the criteria of the rules as follows: (i) molecular weight ≤ 500 Da, (ii) hydrogen bond donors ≤5, and hydrogen acceptors ≤10, (iii) rotatable bond ≤10, (iv) polar surface area ≤140 Å 2 and (v) lipophilicity (expressed as Log P) ≤5. Therefore, we believed that these four compounds were likely to be developed as promising novel mutated EGFR inhibitors.

Drug-Likeness
The potent quinoxalinone-containing compounds against EGFR (L858R/T790M/C797S) were predicted the drug-likeness character based on Lipinski's rule of five by inspecting their physicochemical properties (molecular weight (MW), the numbers of hydrogen bond donors (HBD) and acceptors (HBA), rotatable bond (RB), polar surface area (PSA) and Log P). According to Table 1, the predictive results revealed that all compounds could confer the drug-likeness property and obey the acceptable value within the criteria of the rules as follows: (i) molecular weight ≤ 500 Da, (ii) hydrogen bond donors ≤5, and hydrogen acceptors ≤10, (iii) rotatable bond ≤10, (iv) polar surface area ≤140 Å 2 and (v) lipophilicity (expressed as Log P) ≤5. Therefore, we believed that these four compounds were likely to be developed as promising novel mutated EGFR inhibitors.

SASA, Number of H-Bonds, and Contact Atoms
To investigate the buried capacity of each compound within the cleft of the binding site, the solvent-accessible surface area (SASA) was computed by inspecting the amounts of water within the 5.0 Å from the ligand. As shown in Figure 5, the averaged SASA in the last 20 ns (considered reaching an equilibrated state, Figure S4) of each complex was in the order of CPD4 (674 Å 2 ) < CPD15 (730 Å 2 ) < osimertinib (736Å 2 ) < CPD21 (768 Å 2 ) < CPD16 (866 Å 2 ). The smaller area of occupied solvents could imply the better-fitting ligand conformation within the binding pocket and the binding events might not be interfered with by the solvents. Therefore, CPD4 was the most well-buried to binding pocket; however, the intermolecular interactions of CPD16 and surrounding amino acids might be probably reduced. Apart from that, the numbers of H-bonds were analyzed during 80-100 ns. It was found that quinoxalines could form fewer H-bonds when compared to osimertinib (1-2 bonds), suggesting that this type of interaction was not dominantly responsible for the quinoxalinone compounds' recognition since they comprise gradually lower numbers of hydrogen bond donors and acceptors when compared to osimertinib (total numbers of HBD and HBA are 3-4 while there are 7 for osimertinib, Table 1).

Hot-Spot Residues
To investigate the key binding residues within the focused site of EGFR tyrosine kinase, the per-residue decomposition energy (ΔG ) was elucidated via the MM/PBSA approach. The positive and negative ΔG values represent the ligand destabilization and stabilization, respectively. We noted that only residues exhibiting ΔG value We also identified the number of contact atoms within the 5.0 Å from the ligand. As shown in Figure 5, the native contacts of quinoxalinones were less than osimertinib. Nevertheless, the non-native contacts demonstrated a similar range falling onto 50-80 atoms (last 20 ns), indicating that the screened compounds could interact with the surrounding amino acids when the dynamics was taken into account. In the case of osimertinib, due to its structure being larger than quinoxalinones (1.5 to 2-fold larger molecular weight, Table 1), it could plausibly expose to a larger number of proximate amino acids when compared to quinoxalinones.

Hot-Spot Residues
To investigate the key binding residues within the focused site of EGFR tyrosine kinase, the per-residue decomposition energy (∆G bind residue ) was elucidated via the MM/PBSA approach. The positive and negative ∆G bind residue values represent the ligand destabilization and stabilization, respectively. We noted that only residues exhibiting ∆G bind residue value lower than −0.10 kcal/mol were elucidated and compared with other residues. As seen in Figure 6, quinoxalinones showed somewhat a similar pattern of binding amino acids compared to osimertinib. The common amino acids participating in the ligand binding include L718, F723 (except CPD15), V726 (except CPD21), A743, M790, L792, G796, S797, L844, and T854 (except CPD21), which agreed well with the previous report [36]. The amino acids largely contributing to ligand binding consist of both G loop and hinge region. The L718 and V726 (except CPD21) were found to noticeably influence to the binding in almost all screened compounds and osimertinib, which agreed well with the previous findings [36]. In addition, the M793 lining at the hinge region was detected in all screened compounds' binding, which has been highlighted as a satisfied condition of successful EGFR TK inhibitors by controlling the entry of inhibitors to the ATP-binding pocket [36,37]. Interestingly, all screened compounds could interact with two mutated residues including M790 and S797, which are lining in the ligand binding site. It is worth noting that CPD16 showed higher binding energy to M793 than osimertinib, revealing a satisfied condition of being successful EGFR TK inhibitors. Altogether, this finding could provide an explanation in terms of structural basis that support the great inhibitory activity in vitro of all screened compounds.
In addition, the energy contributions to each amino acid were observed. All screened compounds shared somewhat a similar pattern to osimertinib in which the most amino acids were recognized through van der Waals interaction energies (∆E vdW + ∆G nonpolar solv ) while the electrostatic term (∆E elec + ∆G polar solv ) was not preferably dominant toward the binding (Figure 7). This finding agreed well with the previous study that suggested a deep hydrophobic pocket within the ATP-binding site of EGFR TK that is controlled by the residue at 790 position (known as a gatekeeper) [2]. However, CPD16 showed a higher electrostatic contribution to M793 than other quinoxalinones and even osimertinib, suggesting this compound could be well recognized during the binding event. Interestingly, the residue T854 demonstrated a higher contribution of van der Waals interaction energies to all screened compounds when compared to osimertinib, implying the main distinct character of inhibitory actions of quinoxalinones. In addition, the energy contributions to each amino acid were observed. All screened compounds shared somewhat a similar pattern to osimertinib in which the most amino acids were recognized through van der Waals interaction energies (∆EvdW ∆G ) while the electrostatic term (∆Eelec ∆G ) was not preferably dominant toward the binding (Figure 7). This finding agreed well with the previous study that suggested a deep hydrophobic pocket within the ATP-binding site of EGFR TK that is controlled by the residue at 790 position (known as a gatekeeper) [2]. However, CPD16 showed a higher electrostatic contribution to M793 than other quinoxalinones and even osimertinib, suggesting this compound could be well recognized during the binding event. Interestingly, the residue T854 demonstrated a higher contribution of van der Waals interaction energies to all screened compounds when compared to osimertinib, implying the main distinct character of inhibitory actions of quinoxalinones. Figure 6. Binding free energy decomposed to each residue (only some key influential amino acids were shown). (A) The level of residue contribution was shaded in a red-blue-white color scale, which was ranging from higher to lower contribution, respectively. Note that the mutated amino acids were shown in red word. (B) A close-up view of key influential amino acids contributing to ligand binding in complex with screened compounds and osimertinib.

Toxicity Prediction
Toxicity is one of the crucial factors that has been taken into consideration prior to continuing the next phase of drug discovery campaigns. Herein, we employed the avail able databases to predict the main crucial features necessary for inspecting compounds toxicity. As listed in Table 2, all screened compounds showed no predicted mutagenicity tumorigenicity, irritant, and reproductive effects (except CPD4; low level) as well as the

Toxicity Prediction
Toxicity is one of the crucial factors that has been taken into consideration prior to continuing the next phase of drug discovery campaigns. Herein, we employed the available databases to predict the main crucial features necessary for inspecting compounds' toxicity. As listed in Table 2, all screened compounds showed no predicted mutagenicity, tumorigenicity, irritant, and reproductive effects (except CPD4; low level) as well as the predicted LD 50 values were in the range of slightly toxic (5000 mg/kg > LD50 > 500 mg/kg). For osimertinib, it was predicted to be highly prone to mutagenicity and was moderately toxic (500 mg/kg > LD50 > 50 mg/kg). Hence, all screened quinoxalinones were predicted to be safe and might reduce the rate of failure caused by toxicity.

Effect of CPD Derivatives on the Cell Viability
The cytotoxic effect of CPD derivatives on H1975 lung cancer cells expressing L858R/ T790M EGFR as well as Vero normal kidney cells was evaluated using MTT assay. Osimertinib, a third-generation EGFR TKi used to treat adenocarcinoma patients carrying EGFR mutation [22], was used as the positive control. As shown in Table 3, all CPD derivatives as well as osimertinib inhibited the growth of H1975 cells in a concentration-dependent manner. Of all CPD derivatives, CPD4 exhibited the most potent cytotoxicity toward H1975 cells (IC 50 = 3.47 ± 2.20 µM), followed by osimertinib (IC 50  Interestingly, all CPD derivatives, especially CPD15, CPD16, and CPD21 exhibited low cytotoxic effects on Vero normal kidney cells (IC 50 > 100 µM). These results suggested that CPD derivatives are promising anticancer agents with high safety profiles. We noted that the cell viability of H1975 and Vero was shown in Figure S5 and Figure S6, respectively.

Docking-Based Virtual Screening
In this study, protein-ligand docking was performed using the crystal structure of EGFR (L858R/T790M/C797S) in a complex with the known drug, osimertinib (PDB ID: 6LUD [38]). The missing residues (No. 748-751 and 991-999) were homologically built using the SWISS-MODEL server [39]. The newly constructed protein structure was then validated by plotting a Ramachandran diagram using PROCHECK [40]. In addition to the protein structure, the quinoxalinones synthesized by Wimonsong et.al., [41] were manually constructed into 3D structures by using the Gaussian 09 program [42].
For docking validation, the crystalized ligand was defined as a center in the binding site, and redocked into the same binding site, which was in X, Y, and Z coordinates as of −48.818, −2.761, and −18.500, respectively, while the box size was set to 40 Å for every dimension. The crystallized and redocked poses were then aligned to observe the verification of the docking protocol used as shown in Figure S3. Accordingly, the validated docking protocol was employed for all ligands during virtual screening. The docking-based screening was executed by the Autodock Vina XB software package (Sirimulla Researcg Group@University of Texas at EI Paso, EI Paso, TX, USA).

EGFR Tyrosine Kinase Inhibition
The inhibitory activity of each candidate compound towards EGFR triple mutations was performed using the luminescent ADP-Glo™ kinase assay (Promega Corporation, Madison, WI, USA) as previously described [43,44]. According to the reaction recipe, 8 µL of buffer (40 mM Tris-HCl pH 7.5, 20 mM MgCl 2 and 0.1 mg mL −1 bovine serum albumin) was added to a 384-well plate, followed by 5 µL of enzyme (1.25 ng µL −1 ). Then, 2 µL of inhibitors and 10 µL of a mixture containing 25 µM ATP and 2.5 µM poly(glu-tyr) were added and the reaction was incubated for 1 h at room temperature, followed by adding 5 µL of the ADP-Glo reagent and the reaction was further incubated for 40 min. Subsequently, 10 µL of kinase detection reagent was added and incubated at room temperature for 30 min to convert ADP into ATP. The generated luminescence by ATP in a luciferase reaction was then measured using a microplate reader (Infinite M200 microplate reader, Tecan, Männedorf, Switzerland), which corresponds to kinase activity. All assays were performed in triplicate and the obtained data were represented as the relative inhibition (%) of inhibitors compared to the control with no inhibitor as shown in Equation (1) [45];

Molecular Dynamic (MD) Simulations
MD simulations were applied to gain insights into the key binding residues guiding the mechanism of inhibitory action in the dynamical and physiological system. The initial structure of each EGFR (L858R/T790M/C797S) TK and focused ligand was retrieved from the best docking pose (lowest E binding ). Then, the constructed complex was subjected to run under the periodic boundary condition with the isothermal-isobaric (NPT) scheme (310 K, 1 atm), according to the previous studies [35,36,46,47]. The AMBER ff14SB force field and the generalized AMBER force field version 2 (GAFF2) were applied to treat bonded and non-bonded interaction parameters of the simulated complex [48] whilst the water model was used to solvate the system [49]. To neutralize the overall charge of the system, either sodium or chloride ion was then randomly introduced according to the calculated numbers of charges. Minimization of the hydrogen atoms and water molecules was carried out using 500 steps of steepest descent (SD) followed by 1500 steps of conjugated gradient (CG) methods while the rest of the molecules were fixed. Next, the minimization of the proteinligand complex (constrained solvents) and the whole complex system was minimized again by using the same procedure, according to the previous publications [30,50,51].
The particle mesh Ewald summation approach [52] was used to handle electrostatic interactions, while all hydrogen atoms were constrained by the SHAKE algorithm [53]. The temperature was heated from 10 to 310 K by using a Langevin thermostat [54] with a collision frequency of 2 ps −1 whilst the pressure was created by the Berendsen barostat [55]. The MD production lasted to 100 ns by the 2 fs increment of a time step. The MD outputs were elucidated through the cpptraj module, while the per-residue decomposition energy (∆G residue binding ) was computed by MM/PBSA.py implemented in AMBER16.

Pharmacokinetic Properties Prediction
The focused pharmacokinetic properties including physicochemical characteristics and drug-likeness were predicted by using the web-based program SwissADME (SIB Swiss Institute of Bioinformatics, Écublens, Switzerland, www.swissadme.ch/, accessed on 30 October 2022) [56]. The screened compounds' features were calculated compared to the known drug, osimertinib.
The studied cell lines were seeded into 96-well plates (cell density of 6000 cells/well) Scheme 1. Synthesis of CPD21.
The studied cell lines were seeded into 96-well plates (cell density of 6000 cells/well) and were incubated overnight. Then, cells were treated with screened compounds and the known drug with various concentrations (125, 62.5, 31.25, 15.63, 7.81, 3.91, 1.95, and 0.98 µM) for 72 h. After incubation, the effect of compounds on cell viability was measured using the MTT assay as follows: 10 µL of MTT solution (5 mg/mL) was added and then incubated for 3 h. The medium was removed, followed by adding 50 µL of DMSO to dissolve the crystal formazan product. The purple solution of crystal formazan was then measured at 570 nm using a microplate spectrophotometer (Infinite M200 microplate reader, Tecan, Männedorf, Switzerland). The complete medium with 1.0% DMSO (no inhibitor) was used as a control in accordance with the reaction containing 1.0% DMSO of all wells.

Statistical Analysis
The data was represented as mean ± standard error of the mean (SEM). Differences between groups were compared using one-way ANOVA (Tukey's test) for multiple comparisons. The differences in means were determined at the confidence level p ≤ 0.05.

Conclusions
Finding pre-clinical drug candidates against acquired drug resistance of EGFR tyrosine kinase is a great challenge and necessary. Herein, synthesized quinoxalinone-containing compounds were explored through in silico and in vitro methods. From screening to biological assay, CPD4, CPD15, CPD16, and CPD21 showed great inhibitory activity in enzymatic assay with IC 50 values down to a nanomolar scale as being a comparable level when compared to the known drug, osimertinib. They could also inhibit the NSCLC cell growth in a concentration-dependent manner with safety profiles. Based on computational modeling, we proposed some common amino-acid residues which are responsible for the binding of all screened compounds and the drug include L718 and V726 at the G loop as well as M790, M793, and S797 at the hinge region. Moreover, all screened compounds could be intermolecularly recognized by two previously mutated amino acids including M790 and S797, which is the main reason to advocate the inhibitory capability against EGFR (L858R/T790M/C797S) TK. In addition, all predicted physicochemical parameters based on Lipinski's rule of five had met the criteria of drug-likeness and were predicted to show no highly toxic features. Altogether, we proposed that CPD4, CPD15, CPD16, and CPD121 could become viable candidates for development as the fourth generation of EGFR (L858R/T790M/C797S) TK. Animal testing is encouraged to investigate further.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/molecules27248901/s1, Table S1: Chemical structures of all synthesized quinoxalinonecontaining compounds and their original code in the previously published paper, Figure S1: 1H spectrum of CPD21 (CDCl3, 400 MHz), Figure S2: 13C NMR spectra of CPD21 (CDCl3, 100 MHz), Figure S3: Alignment of the re-docked pose (green) and available crystallized osimertinib (grey) showing the validated docking protocol used to perform docking-based virtual screening, Figure S4: Calculated backbone root-mean-square displacement (RMSD) within 5 Å of all screened compounds and osimertinib, Figure S5: H1975 cell viability against treated screened compounds and osimertinib at different doses. Data were represented as mean ± SEM with three replicates, Figure S6: Vero cell viability against treated screened compounds and osimertinib at different doses. Data were represented as mean ± SEM with three replicates.