A New Theobromine-Based EGFR WT and EGFR T790M Inhibitor and Apoptosis Inducer: Design, Semi-Synthesis, Docking, DFT, MD Simulations, and In Vitro Studies

: The essential pharmacophoric structural properties were applied to design a new derivative of theobromine as an antiangiogenic EGFR inhibitor. The designed candidate is a (para-nitrophenyl)acetamide derivative of the natural alkaloid, theobromine ( T-2-PNPA ). The potentialities of T-2-PNPA to inhibit the EGFR protein were studied computationally in an extensive way. Firstly, the molecular docking against EGFR WT and EGFR T790M demonstrated T-2-PNPA ’s capabilities of binding with the targeted receptors. Then, the MD experiments (for 100 ns) illustrated through six different studies the changes that occurred in the energy as well as in the structure of EGFR– T-2-PNPA complex. Additionally, an MM-GBSA analysis determined the exact energy of binding and the essential residues. Furthermore, DFT calculations investigated the stability, reactivity, and electrostatic potential of T-2-PNPA . Finally, ADMET and toxicity studies conﬁrmed both the safety as well as the general likeness of T-2-PNPA . Consequently, T-2-PNPA was prepared for the in vitro biological studies. T-2-PNPA inhibited EGFR WT and EGFR T790M with IC 50 values of 7.05 and 126.20 nM, respectively, which is comparable with erlotinib activities (5.91 and 202.40, respectively). Interestingly, T-2-PNPA expressed cytotoxic potentialities against A549 and HCT-116 cells with IC 50 values of 11.09 and 21.01 µ M, respectively, which is again comparable with erlotinib activities (6.73 and 16.35, respectively). T-2-PNPA was much safer against WI-38 (IC 50 = 48.06 µ M) than erlotinib (IC 50 = 31.17 µ M). The calculated selectivity indices of T-2-PNPA against A549 and HCT-116 cells were 4.3 and 2.3, respectively. This manuscript presents a new lead anticancer compound ( T-2-PNPA ) that has been synthesized for the ﬁrst time and exhibited promising in silico and in vitro anticancer potentialities.


Introduction
According to estimates of the WHO, cancer will be the second leading cause of death in the next few years [1]. Cancer therapy is a big challenge for medicinal chemists who aim to develop drugs that inhibit cancer growth by interacting with specific molecular targets and causing damage to the cancer cells [2]. There are various studies linking tumor growth and survival to the increased vascularity (angiogenesis) in cancer cells [3]. Thus, a potential strategy to combat cancer is to counteract this phenomenon [4]. Cancer cells grow and undergo angiogenesis through epidermal growth factor receptors (EGFRs) [5]. A cascade of downstream signaling pathways is triggered by EGFR's overexpression of, which stimulates proliferation and differentiation, as well as survival of cells [6]. The link between the promotion of solid tumor growth and the elevated levels of EGFR was discovered in multiple types of cancer [7]. Additionally, a reduction in survival rate has been associated with increased EGFR expression in several types of cancer, and EGFR expression in these cancers has served as a strong prognostic indicator [8]. Fortunately, researchers were able to target cancer cells with EGFR receptors due to their overexpression in cancer cells as compared to normal cells [9].
Natural products were historically the most essential base of drug development [10,11]. Presently, as many as one-third of the new FDA-approved drugs in the last two decades are derived from natural sources [12]. In the field of natural anticancer drug discovery, xanthine derivatives showed a variety of antimutagenic activities against ovarian [13], prostate [14], leukemia [15], and breast [15] cancers. Theobromine is a xanthine alkaloid that was first discovered in 1841 [16]. In addition to being found in chocolate, it is an ingredient of other foods such as leaves from the tea plant [17]. Through in vitro and in vivo testing, theobromine showed promising anticancer activities by inhibiting DNA synthesis and inhibiting the growth of glioblastoma multiforme (U87MG cell line) [18]. Similarly, theobromine inhibited the angiogenesis in lung cancer [19]. It also exhibited promising antiangiogenic activity via inhibition of VEGF in both in vivo and in vitro testing of ovarian cancer [20].
The study of the structure-activity relationship leads to the discovery and repurposing of more potent drugs and enhances drug-likeness, pharmacokinetics, and pharmacodynamics of new drugs [21]. Computational chemistry is used for researching chemistry problems using software employing many techniques. Computational chemistry was widely used in the pharmaceutical industry for studying how potential drugs interact with biomolecules [22]. The computational chemistry applications starts from the molecular and drug design [23] and extends to the structure-based approaches as docking and MD simulations [24] as well as the ligand-based approaches such ADMET [25,26], DFT [27], structure similarity [28], and pharmacophore assessment [29].
There are three generations of the FDA-approved anticancer drugs targeting EGFR. The first generation (such as erlotinib I [30]) was very active against the non-small-cell lung cancer. The second generation (such as afatinib II [31]) overcomed EGFR resistance due to mutation. The third generation (such as olmutinib III [32] and avitinib IV [33]) overcome the second generation's induced toxicity and showed improved activities towards EGFR T790M over the EGFR WT . In the literature, some 1H-pyrazolo [3,4-d]pyrimidine derivatives such as compounds V [34] and VI [35] were reported to have good EGFR inhibitory activity. Moreover, our team discovered some 1H-pyrazolo [3,4-d]pyrimidine derivatives (such as compound VII [36]) and thieno [2,3-d]pyrimidine derivatives (such as compound VIII [37]). These compounds were proved to have promising EGFR kinase inhibitory activity with good cytotoxic activity.

Rationale
As shown in Figure 1, four essential features should be present in the chemical structure of EGFR inhibitors in order to bond correctly with the active site [38]. These features are a flat hetero aromatic system, an NH spacer, a terminal hydrophobic head, and a hydrophobic tail. These moieties have the potentialities to occupy the adenine-binding pocket, the linker region, the hydrophobic region I, and the hydrophobic region II, respectively of the EGFR active site [39][40][41].
aromatic systems are composed of pyrimidine moiety fused with a five-membered ring. Depending on the ligand-based drug design approach, we designed a theobromine derivative that kept the main features of the lead compounds (III-VIII).
As shown in Figure 2, we have designed a new compound (2-(3,7-Dimethyl-2,6-dioxo-2,3,6,7-tetrahydro-1H-purin-1-yl)-N-(4-nitrophenyl)acetamide, T-2-PNPA) that has a theobromine moiety as a hetero aromatic system to occupy the adenine-binding pocket. The theobromine moiety is composed of a pyrimidine moiety fused with the five-membered 1H-imidazole ring, meaning that it is considered a good ring equivalent of the thieno [3,2-d]pyrimidine, 7H-pyrrolo [2,3-d]pyrimidine, and 1H-pyrazolo [3,4-d]pyrimidine moieties. Additionally, the acetamide moiety was designed as a linker which may increase the flexibility of the designed compound and consequently increase its affinity against the target protein. The p-nitrobenzene moiety was designed with a hydrophobic head to occupy the hydrophobic region I. Furthermore, the methyl group at the 7-position of theobromine moiety was utilized as a hydrophobic tail to occupy the hydrophobic region II (Figure 2).  Observing the structures of Olmutinib III, avitinib IV, compounds V, VI, VII, and VIII, there are diverse hetero aromatic systems such as thieno [3,2-d]pyrimidine, 7H-pyrrolo [2,3d]pyrimidine, and 1H-pyrazolo [3,4-d]pyrimidine. These hetero aromatic systems can occupy the adenine-binding pocket of the EGFR active site. Interestingly, all these aromatic systems are composed of pyrimidine moiety fused with a five-membered ring. Depending on the ligand-based drug design approach, we designed a theobromine derivative that kept the main features of the lead compounds (III-VIII).
As shown in Figure 2, we have designed a new compound (2-(3,7-Dimethyl-2,6dioxo-2,3,6,7-tetrahydro-1H-purin-1-yl)-N-(4-nitrophenyl)acetamide, T-2-PNPA) that has a theobromine moiety as a hetero aromatic system to occupy the adenine-binding pocket. The theobromine moiety is composed of a pyrimidine moiety fused with the five-membered 1H-imidazole ring, meaning that it is considered a good ring equivalent of the thieno [3,2d]pyrimidine, 7H-pyrrolo [2,3-d]pyrimidine, and 1H-pyrazolo [3,4-d]pyrimidine moieties. Additionally, the acetamide moiety was designed as a linker which may increase the flexibility of the designed compound and consequently increase its affinity against the target protein. The p-nitrobenzene moiety was designed with a hydrophobic head to occupy the hydrophobic region I. Furthermore, the methyl group at the 7-position of theobromine moiety was utilized as a hydrophobic tail to occupy the hydrophobic region II (Figure 2).  We decided to proceed with a molecular docking study of T-2-PNPA aga ATP-binding site of wild EGFR (EGFR WT ; PDB: 4HJO) using erlotinib as a standa validation step generates an acceptable RMSD value of 1.40 (Supplementary Ma Furthermore, the binding mode of erlotinib has shown the essential hydrogen bo interaction with the key amino acid residue, Met 769, and this was consistent with described in the literature [42] as presented in Figure 3.  We decided to proceed with a molecular docking study of T-2-PNPA against the ATP-binding site of wild EGFR (EGFR WT ; PDB: 4HJO) using erlotinib as a standard. The validation step generates an acceptable RMSD value of 1.40 (Supplementary Materials). Furthermore, the binding mode of erlotinib has shown the essential hydrogen bond (HB) interaction with the key amino acid residue, Met 769, and this was consistent with the one described in the literature [42] as presented in Figure 3.  Similar to how erlotinib binds, T-2-PNPA was able to employ its 6-oxo group of the 1H-purine moiety to anchor the adenine pocket via a key HBing formation with Met793. Additionally, in the same pocket, the whole nucleus (3,7-dimethyl-2,6-dioxo-2,3,6,7-tetrahydro-1H-purine moiety) was involved in three hydrophobic bond (HI)s with Leu694 and Cys773. Regarding the hydrophobic pocket, the 4-nitrophenyl moiety was engaged in two pi-pi and two electrostatic interactions with Ala719, Val702, Lys721, and Asp831, respectively. As previously stated, our assumption was that the nitro group would interact with the molecular target via a charge-reinforced HB or electrostatic charge ( Figure 4). Similar to how erlotinib binds, T-2-PNPA was able to employ its 6-oxo group of the 1H-purine moiety to anchor the adenine pocket via a key HBing formation with Met793. Additionally, in the same pocket, the whole nucleus (3,7-dimethyl-2,6-dioxo-2,3,6,7tetrahydro-1H-purine moiety) was involved in three hydrophobic bond (HI)s with Leu694 and Cys773. Regarding the hydrophobic pocket, the 4-nitrophenyl moiety was engaged in two pi-pi and two electrostatic interactions with Ala719, Val702, Lys721, and Asp831, respectively. As previously stated, our assumption was that the nitro group would interact with the molecular target via a charge-reinforced HB or electrostatic charge ( Figure 4).

Molecular Docking Studies against Mutant EGFR Protein (EGFR T790M )
To reinforce the above findings, the docking results of T-2-PNPA against the mutant type of EGFR (EGFR T790M , PDB: 3W2O) were discussed. Initially, aligning the co-crystalized ligand (TAK-285) in the mutant EGFR T790M generated an RMSD value of 1.11. This value validated the docking procedure (Supplementary Materials). The binding profile of TAK-285 to the EGFR T790M active site is shown in Figure 5. The obtained results and the reported data harmonized completely. The designed compound, T-2-PNPA, was stacked onto the EGFR T790M catalytic sit similarly to the native ligand, as shown in Figure 6. In the adenine pocket, the 3,7-dime thyl-2,6-dioxo-2,3,6,7-tetrahydro-1H-purine moiety achieved eight HIs with Leu844 Leu718, Ala743, and Met793 as well as two HBs with Thr845 and Met793. Moreover, two electrostatic interactions were observed in the hydrophobic pocket via the interaction o p-nitrophenyl moiety with Lys745 and Glu762. The designed compound, T-2-PNPA, was stacked onto the EGFR T790M catalytic site similarly to the native ligand, as shown in Figure 6. In the adenine pocket, the 3,7-dimethyl-2,6-dioxo-2,3,6,7-tetrahydro-1H-purine moiety achieved eight HIs with Leu844, Leu718, Ala743, and Met793 as well as two HBs with Thr845 and Met793. Moreover, two electrostatic interactions were observed in the hydrophobic pocket via the interaction of p-nitrophenyl moiety with Lys745 and Glu762.

MD Simulations
The MD analyses performed on the production run of 100 ns show that, overall, the system is stable. The RMSD plot ( Figure 7A) shows equilibrated values for EGFR and EGFR-T-2-PNPA complex (blue and green curve, respectively) through the whole trajectory with an average value of 2.9 Å . Interestingly, the RMSD of the T-2-PNPA (red curve) shows very large fluctuations during the trajectory which indicates its large free motion inside the binding pocket, however, its distance from the center of mass of EGFR remains nearly constant with an average of 11.62 Å ( Figure 7F). The RoG ( Figure 7B), and SASA ( Figure 7C) show a stable protein fluctuation with an average of 19.79 Å , 15755 Å 2 , respectively. H-bonds formed between the T-2-PNPA and EGFR ( Figure 7D) shows that at least one H-bond is formed in 214 frames out of 1000 frames (21.4%). The amino acids' fluctuation was depicted in the RMSF plot ( Figure 7E) showing low fluctuation value (less than 2 Å ) except for the G695:G700 (reaching 3.6 Å ), and E842:V852 loop regions (reaching 3.6 Å ) and the free C-terminal (reaching 7 Å ).

MD Simulations
The MD analyses performed on the production run of 100 ns show that, overall, the system is stable. The RMSD plot ( Figure 7A) shows equilibrated values for EGFR and EGFR-T-2-PNPA complex (blue and green curve, respectively) through the whole trajectory with an average value of 2.9 Å. Interestingly, the RMSD of the T-2-PNPA (red curve) shows very large fluctuations during the trajectory which indicates its large free motion inside the binding pocket, however, its distance from the center of mass of EGFR remains nearly constant with an average of 11.62 Å ( Figure 7F). The RoG ( Figure 7B), and SASA ( Figure 7C) show a stable protein fluctuation with an average of 19.79 Å, 15755 Å 2 , respectively. H-bonds formed between the T-2-PNPA and EGFR ( Figure 7D) shows that at least one H-bond is formed in 214 frames out of 1000 frames (21.4%). The amino acids' fluctuation was depicted in the RMSF plot ( Figure 7E) showing low fluctuation value (less than 2 Å) except for the G695:G700 (reaching 3.6 Å), and E842:V852 loop regions (reaching 3.6 Å) and the free C-terminal (reaching 7 Å).   38 kcal/mol, respectively. Additionally, we performed a decomposition study (Figure 9) to identify the amino acids that contribute within 1 nm of the EGFR-T-2-PNPA complex, which have the most significant contribution to the binding. L694, V702, K721, and L820 are the amino acid residues that have contributed with an energy value less than −1 kcal/mol (−2.69, −1.97, −1.22, and −1.13 kcal/mol, respectively). On the other hand, D776 (+0.46 kcal/mol), and D831 (+0.7 kcal/mol) show a positive free energy (non-favorable) contribution to the binding to the EGFR protein.

PLIP Study
As already known, MD simulation can produce a very large number of frames that can have a minimal difference. Many methods to reduce the number of frames by selecting the ones that are conformationally different have been developed such as clustering and

PLIP Study
As already known, MD simulation can produce a very large number of frames that can have a minimal difference. Many methods to reduce the number of frames by selecting the ones that are conformationally different have been developed such as clustering and

PLIP Study
As already known, MD simulation can produce a very large number of frames that can have a minimal difference. Many methods to reduce the number of frames by selecting the ones that are conformationally different have been developed such as clustering and wavelet analysis [43]. In this work, we have utilized the clustering technique to reduce the data dimensionality and to obtain a representative frame for every produced cluster. As clarified in the experimental section (Supplementary Materials), cluster numbers were automatically calculated by the elbow method to produce four successive clusters. For each one, the PLIP webserver was employed to identify both number and types of interactions between T-2-PNPA and EGFR (Table 1). As can be seen, HIs are the common type of interactions between T-2-PNPA and EGFR with a total of 11 HIs compared to only three H-bonds. Additionally, a .pse file was generated to understand the 3D conformation of the T-2-PNPA interactions inside EGFR (Figure 10).
Processes 2022, 10, x FOR PEER REVIEW 12 of 23 wavelet analysis [43]. In this work, we have utilized the clustering technique to reduce the data dimensionality and to obtain a representative frame for every produced cluster. As clarified in the experimental section (Supplementary Materials), cluster numbers were automatically calculated by the elbow method to produce four successive clusters. For each one, the PLIP webserver was employed to identify both number and types of interactions between T-2-PNPA and EGFR (Table 1). As can be seen, HIs are the common type of interactions between T-2-PNPA and EGFR with a total of 11 HIs compared to only three Hbonds. Additionally, a .pse file was generated to understand the 3D conformation of the T-2-PNPA interactions inside EGFR ( Figure 10).

DFT Studies
Density functional theory (DFT) is one of the best computational techniques for quantitatively predicting and rationalizing the reactivity of biomolecules. The predicted properties of reasonable accuracy can be obtained by DFT. Many quantum chemical parameters, such as the lowest unoccupied molecular orbital energy, ELUMO, and the highest occupied molecular orbitals energy, EHOMO, gap's energy (ΔE) EHOMO-ELUMO, global electrophilicity (ω), dipole moment (D), electronegativity (χ), softness (σ), and chemical hardness (η), all affect how surface atoms interact electronically with the target.

DFT Studies
Density functional theory (DFT) is one of the best computational techniques for quantitatively predicting and rationalizing the reactivity of biomolecules. The predicted properties of reasonable accuracy can be obtained by DFT. Many quantum chemical parameters, such as the lowest unoccupied molecular orbital energy, ELUMO, and the highest occupied molecular orbitals energy, EHOMO, gap's energy (∆E) EHOMO-ELUMO, global electrophilicity (ω), dipole moment (D), electronegativity (χ), softness (σ), and chemical hardness (η), all affect how surface atoms interact electronically with the target.

Geometry Optimization
The level of theory used for all calculations was DFT/B3LYP/6-311G++(d,p). Without any restrictions, T-2-PNPA optimization was carried out and depicted in Figure 11. The intermediate and reactant were bound via the C14-N2 bond. The C14-N2 bond length and the two angles on both sides were recorded on the optimized structure, Figure 11. The dipole moment (Dm), which indicates the interaction within the compound, also determines the molecule's polarity. According to Table 2, the Dm value is 11.321 Debye and the measured ground total energy (TE) is −34,966.4 eV.

Geometry Optimization
The level of theory used for all calculations was DFT/B3LYP/6-311G++(d,p). Without any restrictions, T-2-PNPA optimization was carried out and depicted in Figure 11. The intermediate and reactant were bound via the C14-N2 bond. The C14-N2 bond length and the two angles on both sides were recorded on the optimized structure, Figure 11. The dipole moment (Dm), which indicates the interaction within the compound, also determines the molecule's polarity. According to Table 2, the Dm value is 11.321 Debye and the measured ground total energy (TE) is −34,966.4 eV.

Frontier Molecular Orbital Analysis (FMO) Analysis
All chemical reactions greatly benefit from orbital analysis. In reality, molecular orbitals (MOs) are crucial for a deeper comprehension of chemical processes and electrical electronic characteristics. Around the year 1952, Fukui proposed the border orbital theory, which connects reactivity to the characteristics of the HOMO and LUMO molecular orbitals.
The computed energy gap (Egap) between T-2-PNPA's HOMO and LUMO is −4.260 eV according to the FMO study, and the sketch diagram is shown in Figure 12. The calculated Egap value is low which indicates a small HOMO-LUMO gap in the border orbitals. The molecule under investigation, T-2-PNPA, showed high polarizability and chemical reactivity according to the calculated Dm and Egap values [44]. Table 2 shows the computed energies for HOMO (EHOMO), LUMO (ELUMO), and Egap.

Frontier Molecular Orbital Analysis (FMO) Analysis
All chemical reactions greatly benefit from orbital analysis. In reality, molecular orbitals (MOs) are crucial for a deeper comprehension of chemical processes and electrical electronic characteristics. Around the year 1952, Fukui proposed the border orbital theory, which connects reactivity to the characteristics of the HOMO and LUMO molecular orbitals.
The computed energy gap (E gap ) between T-2-PNPA's HOMO and LUMO is −4.260 eV according to the FMO study, and the sketch diagram is shown in Figure 12. The calculated E gap value is low which indicates a small HOMO-LUMO gap in the border orbitals. The molecule under investigation, T-2-PNPA, showed high polarizability and chemical reactivity according to the calculated Dm and E gap values [44]. Table 2   The global reactivity indices provides information about the reactivity or inhibitory ability of a compound from chemical point of view. It can be efficiently obtained via quantum calculations [45]. Theoretical static calculations are used to derive the reactivity indices involving molecular orbitals [46]. Global reactive indices that characterize a molecule's capability to be stable or reactive include global chemical hardness (η) and softness (σ). When a molecule's softness value (σ) is high and its hardness value (η) is low, it is found to be a strong inhibitor. The ionization potential (IP), energy change (∆E), maximal charge acceptance (Nmax), electrophilicity (ω), chemical potential (µ ), electronegativity (χ), global chemical softness (σ), global chemical hardness (η), and electron affinity (EA) of the molecule were all calculated using Koopmans' theory:

Global Reactive Indices and Total Density of State (TDOS)
The global reactivity indices provides information about the reactivity or inhibitory ability of a compound from chemical point of view. It can be efficiently obtained via quantum calculations [45]. Theoretical static calculations are used to derive the reactivity indices involving molecular orbitals [46]. Global reactive indices that characterize a molecule's capability to be stable or reactive include global chemical hardness (η) and softness (σ). When a molecule's softness value (σ) is high and its hardness value (η) is low, it is found to be a strong inhibitor. The ionization potential (IP), energy change (∆E), maximal charge acceptance (N max ), electrophilicity (ω), chemical potential (µ), electronegativity (χ), global chemical softness (σ), global chemical hardness (η), and electron affinity (EA) of the molecule were all calculated using Koopmans' theory: Due to the possibility of quasi-degenerate levels in the border region, LUMO and HOMO could not give an accurate description of FMO. Although FMO analysis is important, the HOMO-LUMO theory has some limitations and new efforts are necessary to understand when and how the HOMO energy works and does not work [47]. The occupied states' number per unit volume is determined at equilibrium by a chemical system's probability distribution function and the density of its states. This quantity is used to investigate the physical characteristics of chemical systems. Figure 13 shows that T-2-PNPA has a considerable narrow energy gap, as determined by the TDOS analysis, and the orbitals under the HOMO orbital had the highest reported intensities. Such findings highlighted the potency of the prospective inhibitor.

Egap = ELUMO -EHOMO
Due to the possibility of quasi-degenerate levels in the border region, LUMO and HOMO could not give an accurate description of FMO. Although FMO analysis is important, the HOMO-LUMO theory has some limitations and new efforts are necessary to understand when and how the HOMO energy works and does not work [47]. The occupied states' number per unit volume is determined at equilibrium by a chemical system's probability distribution function and the density of its states. This quantity is used to investigate the physical characteristics of chemical systems. Figure 13 shows that T-2-PNPA has a considerable narrow energy gap, as determined by the TDOS analysis, and the orbitals under the HOMO orbital had the highest reported intensities. Such findings highlighted the potency of the prospective inhibitor. Figure 13. Molecular frontier orbital spectrum, their energies and TDOS were calculated at B3LYP/6-311++G(d,p).

TED and ESP Surface Potential Maps
The total electron density surface (TED), also known as the potential (electron + nuclei) mapping, shows the distribution, size, shape, and dipole moments of the compound's ESP. Figure 14 illustrates the TED and ESP of T-2-PNPA graphically. Green, blue, yellow and red sections indicate, respectively, the neutral electron, electron deficient, slightly rich, and electron-rich zones. The most negative potential zones are found around the carbonyl and nitrate groups, while the slightly electron-rich, yellow-colored regions are found over the terminal phenyl and N9 atom of the purine moiety. The purine substituted methyl (C13) and C8H27 of purine are blue patches that indicate an electron deficient zone. Figure 13. Molecular frontier orbital spectrum, their energies and TDOS were calculated at B3LYP/6-311++G(d,p).

TED and ESP Surface Potential Maps
The total electron density surface (TED), also known as the potential (electron + nuclei) mapping, shows the distribution, size, shape, and dipole moments of the compound's ESP. Figure 14 illustrates the TED and ESP of T-2-PNPA graphically. Green, blue, yellow and red sections indicate, respectively, the neutral electron, electron deficient, slightly rich, and electron-rich zones. The most negative potential zones are found around the carbonyl and nitrate groups, while the slightly electron-rich, yellow-colored regions are found over the terminal phenyl and N9 atom of the purine moiety. The purine substituted methyl (C13) and C8H27 of purine are blue patches that indicate an electron deficient zone.

ADMET Profiling Study
A compound becomes a drug based on its pharmacokinetic assessment and biological activity, so to prevent the late withdrawal of a new compound it should be studied early in the discovery process for ADMET properties [48]. Though ADMET identifies absorption, distribution, metabolism, excretion, and toxicity, there are a number of in vitro

ADMET Profiling Study
A compound becomes a drug based on its pharmacokinetic assessment and biological activity, so to prevent the late withdrawal of a new compound it should be studied early in the discovery process for ADMET properties [48]. Though ADMET identifies absorption, distribution, metabolism, excretion, and toxicity, there are a number of in vitro and in vivo studies that can explore the ADMET descriptors, but in silico studies are much more advantageous due to limitations of time, effort, cost, and to avoid the firm regulations governing animal testing [49]. Using Discovery studio ® , these ADMET parameters were predicted computationally for T-2-PNPA based on its reference molecule, erlotinib. The ADMET results of T-2-PNPA comparing erlotinib ( Figure 15) showed an acceptable degree of druglikeness as it was predicted to have the very low ability to pass the BBB and is expected to be non-hepatotoxic and non-inhibitor against the cytochrome P-450, CYP2D6. Additionally, as shown in Table 3, T-2-PNPA exhibited good aqueous solubility and medium intestinal absorption levels.

ADMET Profiling Study
A compound becomes a drug based on its pharmacokinetic assessment and biological activity, so to prevent the late withdrawal of a new compound it should be studied early in the discovery process for ADMET properties [48]. Though ADMET identifies absorption, distribution, metabolism, excretion, and toxicity, there are a number of in vitro and in vivo studies that can explore the ADMET descriptors, but in silico studies are much more advantageous due to limitations of time, effort, cost, and to avoid the firm regulations governing animal testing [49]. Using Discovery studio ® , these ADMET parameters were predicted computationally for T-2-PNPA based on its reference molecule, erlotinib. The ADMET results of T-2-PNPA comparing erlotinib ( Figure 15) showed an acceptable degree of drug-likeness as it was predicted to have the very low ability to pass the BBB and is expected to be non-hepatotoxic and non-inhibitor against the cytochrome P-450, CYP2D6. Additionally, as shown in Table 3, T-2-PNPA exhibited good aqueous solubility and medium intestinal absorption levels.

In Silico Toxicity Studies
For the development of drugs, the early estimation of toxicity is crucial to minimize failures in later development and clinical stages [50]. The use of in silico methods in toxicity prediction played a crucial role in reducing time wasting and ethical issues in the traditional in vitro and in vivo testing [51]. As the important principle of in silico prediction of toxicity is the structure-activity relationship (SAR)-predictive toxicity, the software compares the chemical structural characteristics of the molecule with those of thousands of molecules, which have reported safety or toxicity [52]. Computing eight parameters of toxicity using Discovery studio (Table 4) indicated the general safety of T-2-PNPA comparing erlotinib.

Chemistry
Scheme 1 described the synthetic pathway of T-2-PNPA. Initially, theobromine 1 was treated with alcoholic KOH while being continuously stirred to produce potassium salt 2 [53]. On the other hand, the 4-nitroaniline was treated with chloroacetylchloride in DMF to afford the target intermediate 4, 2-chloro-N-(4-nitrophenyl)acetamide. In order to finalize T-2-PNPA, the formed potassium salt 2 was then combined with 2-chloro-N-(4nitrophenyl)acetamide 4 in DMF using a catalytic amount of KI.

In Vitro EGFR Inhibition
To verify the conducted computational results, the potentialities of T-2-PNPA against both wild subtype, EGFR WT , and the mutant subtype EGFR T790M of EGFR protein were investigated in vitro. T-2-PNPA firmly inhibited both proteins with IC50 values of 7.05 and 126.20 nM, respectively. Those results ( Figure 16) were bordering erlotinib's value (5.91 nM) against EGFR WT and much stronger than that against EGFR T790M (202.40 nM). T-2-PNPA's IR spectrum was emphasized by the appearance of absorption bands at 1709 and 1661 cm −1 corresponding to carbonyl groups. In addition, the 1 H NMR detected an amidic proton at 11.16 ppm. Additionally, a distinctive singlet signal for the CH 2 group could be seen at 4.76 ppm. The validity was further supported by the 13 C NMR spectrum, which showed distinct peaks at around 43.4, 42.5, 33.7, and 29.9 ppm for the two CH 2 and two CH 3 moieties, respectively.

In Vitro EGFR Inhibition
To verify the conducted computational results, the potentialities of T-2-PNPA against both wild subtype, EGFR WT , and the mutant subtype EGFR T790M of EGFR protein were investigated in vitro. T-2-PNPA firmly inhibited both proteins with IC 50 values of 7.05 and 126.20 nM, respectively. Those results ( Figure 16) were bordering erlotinib's value (5.91 nM) against EGFR WT and much stronger than that against EGFR T790M (202.40 nM).

In Vitro EGFR Inhibition
To verify the conducted computational results, the potentialities of T-2-PNPA against both wild subtype, EGFR WT , and the mutant subtype EGFR T790M of EGFR protein were investigated in vitro. T-2-PNPA firmly inhibited both proteins with IC50 values of 7.05 and 126.20 nM, respectively. Those results ( Figure 16) were bordering erlotinib's value (5.91 nM) against EGFR WT and much stronger than that against EGFR T790M (202.40 nM).  The preceded computational studies (design, molecular docking, MD simulations, MM-GPSA, and DFT calculations) predicted the potentialities of T-2-PNPA to bind correctly and inhibit both types of EGFR protein (EGFR WT and EGFR T790M ). These in silico results were highly in consonance with the obtained results of the in vitro assay of T-2-PNPA against EGFR WT and EGFR T790M .

Cytotoxicity
Because T-2-PNPA showed good inhibitory (in silico and in vitro) potentialities against EGFR WT and EGFR T790M , it is predicted to express promising anticancer effects. Accordingly, the cytotoxicity of T-2-PNPA was assessed in vitro against A549 and HCT-116 malignant besides WI-38 normal cell lines using erlotinib as a reference drug. Table 5 shows that T-2-PNPA expressed anticancer properties almost the same as erlotinib, but with much higher safety. In detail, T-2-PNPA inhibited the growth of the malignant cells with IC 50 values of 11.09 and 21.01 µM, respectively. Additionally, it was determined that T-2-PNPA demonstrated an excellent level of safety against WI-38 cell lines (Figure 17), exhibiting a high IC 50 value of 48.06 µM and very high selectivity indices (SI) against cancer cell lines of 4.3 and 2.3, respectively.