For this study, a set of 14 tetrahydroquinoline derivatives obtained from Fathy et al. (2020) as potential anti-EGFR agents with inhibitory activities (IC₅₀) in μmol L^–1 (Table 1) was investigated. All the studied compounds were optimized using Spartan 14 software (Waziri et al., 2023) to achieve a stable conformation with the least amount of energy. Molecular mechanics force field (Vanommeslaeghe et al., 2014) was used to remove the strain energy from the molecules, and density functional theory was used to execute the optimization using the standard 6-31+G* (d, p) basis set, which includes Becke’s (1993) gradient, exchange correlation and the Yang et al. (2005), Parr et al. (1999) correlation functional (i.e., B3LYP).

Table 1.
Chemical structure, IUPAC name and IC₅₀ of studied compounds.

5-FU: 5-Flurouracil.

Fathy et al. (2020)

The obtained descriptors were extracted from the optimized compounds together and were used in developing valid and reliable QSAR model. The extracted descriptors were used as independent variable and the experimental inhibition concentration (IC₅₀) was used as dependent variable. The QSAR model development was achieved using genetic function approximation (GFA) via material studio software (Puzyn et al., 2010) and the MLR-GFA equation for the model is shown in Eq. 1 below:

(1)

where α is the regression constant, X., X..... X. are the descriptors and β_1, β.…β.are the coefficient of the corresponding descriptors. The developed QSAR model was validated by considering series of statistic factors such cross validation R² (CV.R² ) and Adjusted R²(R²_a) (Eqs. 2 and 3) and the developed model will be considered valid when the calculated value for R²_a and CV.R² were ≥ 0.6 and ≥ 0.5 respectively.

(2)

where 𝑌𝑜𝑏𝑠 = experimentally observed IC₅₀, 𝑌𝑐𝑎𝑙 = calculated IC₅₀ and 𝑌̅𝑜𝑏𝑠 = average of the experimentally observed IC₅₀

(3)

where N = no of compounds observed, P = no of molecular descriptors used in the QSAR model and R. = R- squared value obtained from the QSAR model.

Crystal structure of the EGFR (PDB ID: 4LRM) (Yasuda et al., 2013) was retrieved from the Protein Data Bank (https://www.rcsb.org/) and employed as target receptors in this study. Also, the tetrahydroquinoline derivatives (Table 1) synthesized by Fathy et al. (2020) were used as ligands. The 2D structures of the ligands were drawn with the help of Chem Professional 15.0 and saved as structure-data files (SDFs).

Chimera 1.14 was employed to prepare the protein by eliminating water molecules, numerous ligands, nonprotein component, and other extraneous substances downloaded together with the proteins (Pettersen et al., 2004). The ligands were saved as SDF files and treated protein was converted to PDBQT format using Autodock tool 4.2. Grid box size x = 18.49 Å, y = 22.59 Å, size z = 17.09 Å and grid center dimensions x = 41.77, y = 361.60 and z = 17.09 were set for 4LRM.

The prepared ligands were docked into the binding site of the receptor which was determined by using CastP online server) and the docking calculation was performed using Autodock Vina from PyRX workspace (Trott and Olson, 2010). The inhibitors were treated as flexible throughout the docking simulations. A force-field-based energy scoring function was used to score the ligand orientations, and the highest-scoring binding structure was chosen. Finally, 3D views of the protein-ligand complex were analyzed by using UCSF Chimera 1.14 and Discovery Studio 2020 was used to create the 2D images of the molecular interactions (Capra et al., 2009; Oyeneyin et al., 2022).

The biological action of medications and their metabolic fate in an organism are intimately related to their absorption, distribution, metabolism, and excretion (ADME) and toxicity (ADMET) qualities. In silico predictive models were used to determine the ADMET properties of the test substances. The ADME properties of the substances were determined using the SwissADME online server (http://www.swissadme.ch/index.php/) (Daina et al., 2017). The acute toxicity class, LD50, hepatotoxicity, carcinogenicity, mutagenicity, cytotoxicity, and immunotoxicity of the substances were all predicted using the ProTox-II web server (Banerjee et al., 2018).

The Swiss ADME and ProTox-II online servers take one or more query molecules in canonical SMILES format as input and use a large database to accurately predict the physicochemical properties, pharmacokinetics, solubility, lipophilicity, drug-likeness, bioavailability score, therapeutic properties, and toxicity of compounds.

The calculated descriptors from optimized compounds were screened and for further processing. The calculated descriptors were the energy of the highest occupied molecular orbital (E_HOMO), lowest unoccupied molecular orbital energy (E_LUMO), dipole moment (DM), volume, area, lipophilicity (log P), polar surface area (PSA), number of hydrogen bond donor (HBD), hydrogen bond acceptor (HBA). The calculated descriptors also helped in QSAR study. As shown in Table 2, C3 possess higher E_HOMO value and according to report by Semire et al. (2017), compound with high E_HOMO value possess greater tendency to donate electron to the neighboring compounds and it is expected to interact with the target; therefore, C3 with –4.29 eV is expected to interact well with neighboring compounds. As reported by Oyebamiji et al. (2018), the lower the E_LUMO value the better the tendency of the compound to interact with the target; therefore, C8 proved to have a tendency to interact than other studied compounds (Table 2).

Table 2.
Selected descriptors for the QSAR model.

DM: dipole moment; MW: molecular weight; Vol: volume; PSA: polar surface area; POL: polarizability; HBD: hydrogen bond donor and HBA: hydrogen bond acceptor.

Five distinct models were developed and the best model was chosen and reported among them, since it met the minimum conditions for the evaluation of a valid QSAR model, as reported by Veerasamy et al. (2011). Tables 2–5 as well as Figs. 1 and 2 present the findings of the QSAR study and the developed QSAR model was presented in Eq. 4. The selected and reported model is given by the Eq. 4 with the following validation terms: R. = 0.9525, adj.R. = 0.9314 and CV.R. = 0.9719.

(4)

Figure 1.
Scatter plot of predicted IC₅₀ versus the experimental IC₅₀ for the reported model.

Figure 2.
Scatter plot of residual versus the experimental IC₅₀ for the reported model.

The negative coefficients of the Log P and HBA descriptors clearly suggested that they have a negative contribution to the compounds inhibitory actions. This suggests that reduction in the amount of these independent descriptors improves the inhibition concentration of the studied compounds, and vice versa. The positive coefficients of dipole moment and molecular weight in the model, on the other hand, indicated that these independent descriptors contributed positively to the inhibitory effects of the ligands under study. It suggests that by increasing the amount of these descriptors in the structures with potential anti-EGFR therapeutic agents, the inhibition concentration these drug-like compounds increases and vice versa (Ibrahim et al., 2020; Oke et al., 2022).

Also, as reported by Oyebamiji et al. (2021), inhibition concentration of any developed QSAR model is not enough to judge the potency of such model and this therefore called for validation of the developed QSAR model. The calculated factors considered for validation were cross validation (CV.R2) and adjusted squared correlation coefficient (adj.R2).

Table 3.
Statistical parameters and limit needed for the QSAR model assessment.

The validation factors revealed in this study matched the validation parameters in Table 3, indicating that the model created is both predictive and resilient.

For the compounds under examination, Table 4 displays the experimental IC₅₀, predicted IC₅₀, and residual values. The minimal residual values observed between the experimental and predicted IC₅₀ in the table validated the model’s high anticipated power (meaning that the reported model was reliable with high predicted power). Figure 1 also shows a plot of the predicted IC₅₀versus experimental IC₅₀ for the compounds, and the distribution of the predicted IC₅₀ and experimental IC₅₀ of the compounds along the line reaffirmed the model’s dependability. The R2 values of the internal validation (0.9525) and the plot (0.9944) agreed with one another, confirming the reported model’s stability and reliability (Grisoni et al., 2018). Furthermore, Fig. 2 shows a scatter plot of residuals against the experimental IC₅₀ which the remarkable appearance of both sets of residuals on the plot's top and lower sides confirms that the reported model was free of methodological mistake (systematic deviations).

Table 4.
The experimental IC₅₀, predicted IC₅₀ and the residual values for the studied compounds.

Table 6 shows the correlation statistical analysis of the independent descriptors in the presented model, which revealed that there is no relationship between the descriptors in the model. This demonstrated the good performance of the descriptors utilized in creating the reported model.

The variation inflation factor (VIF) data were also examined to see if the descriptors employed in the model had any multicollinearity issues. In general, a VIF value of 1 or a value between 1 and 5 indicates that there is no intercorrelation between the descriptors. Nevertheless, if the VIF value is larger than 10, the produced model is unstable, and the model should be rechecked if desired (Beheshti et al., 2016). Table 5 shows that the VIF values for each descriptor were less than 5, indicating that the descriptors were highly orthogonal to each other and that there was no intercorrelation between them. This demonstrates that the descriptors used to construct the provided model do not have a multicollinearity problem.

A one-way analysis of variance was used to determine the association between the descriptors and each compound’s biological activity (ANOVA). Table 5 shows the probability value of each descriptor at the 95% confidence level (p = 0.05). As a result, the alternative hypothesis stating that there is a direct relationship between the biological activity of each compound and the descriptor swaying the built model is accepted; however, the null hypothesis stating that there is no direct relationship between the biological activity of each compound and the descriptor swaying the created model is rejected.

Table 5.
Statistical parameters that affect the model.

Table 6.
Coefficient of Pearson’s correlation for descriptor in QSAR model.

3.3 Molecular Docking Studies

Molecular docking technique was used to recognize the compounds which bind well to the EGFR pocket. All tested compounds had strong binding affinity with the receptor. They perfectly fit into the active site of the receptor by interacting with the amino acid residue in the active site of the receptor. The protein residues involved in the interaction, types of interaction involved, distance and binding affinity obtained for each complex formed by the studied compounds and standard were displayed in Tables 7 and 8. The protein residues involved in the interaction and type of interactions observed for C14-EFGR complex was also displayed in Fig. 3.

The complexes formed by 4LRM receptor showed good interactions (Fig. 3). The tested compounds showed binding affinity ranging from –4.6 to –10.1 kcal mol^–1. According to Adepoju et al. (2022), the lower the binding affinity value of any compound, the better the inhibiting ability of such compound; thus, Compound C14 with –10.1 kcal mol^–1 was observed to have highest tendency to inhibit the studied receptor than other studied compounds. Pi-anion interaction was observed in all the complexes except in compounds C5, C7, C8, and C9. A salt bridge was formed between Asp858, and the nitro group of the ligand in complex formed by compound C12. Complexes formed by compound C4, C6, C8, C11, and C14 showed a pi-sulfur interaction between Cys800 and the respective ligands. Pi-sigma bond was observed for all complexes except for complexes formed by C2, C5, C9, C12 and C13. only the complex formed by C14 had unfavorable acceptor-acceptor interaction. Complexes formed by compound C2, C5, C9, C10, C11, C12, and C14 showed a conventional hydrogen bond while others do not. Other interactions observed in the complexes are carbon hydrogen bond, pi-pi stacked, alkyl and pi-alkyl (Table 7).

Table 7.
Binding Interactions and Binding Affinity between Ligands and Receptor 4LRM.

Figure 3.
(a) 3D and (b) 2D views of the molecular interactions of amino acid residues of receptor 4LRM with compound 14.

Regardless of the compounds inhibitory capability against target proteins, the ADMET properties of the compounds are a critical element in determining their pharmacological activities such as drug-likeness, lipophilicity, pharmacokinetics, and toxicity profile when used as drug. In silico ADMET prediction is a quick and low-cost method of determining whether a molecule will be well distributed to its intended site of action, metabolized efficiently, and quickly removed from the body without causing hazardous side effects (Ntie-Kang, 2013). The prediction output for drug likeness, lipophilicity pharmacokinetics, and toxicity profiles of test molecules are displayed on Tables 9–11 respectively.

Drug-likeness analysis is a qualitative evaluation of oral bioavailability based on structural or physicochemical examination of molecules. In this study the rule base filters which include Lipinski’s et al. (2013), Veber et al. (2002), and Muegge and Ghose (2001) rules were employed to determine if a molecule is impermeable or badly absorbed. If a molecule does not break more than two rules, it is deemed orally bioavailable. As indicated on Table 8, the molecular weights of the compounds range from 130 (5-FU) to 590 g mol^–1 (C14). The HBA ranged from 3 to 8 and HBDs ranged from 0 to 2; number of rotatable bonds are between 0 (5-FU) to 9 (C11). The TPSA are between 55.14 (5-FU) and 134.93 (C13 and C14). MR ranged from 27.64 (5-FU) to 150.11 (C14). The less the drug candidate can infiltrate, the higher the TPSA and molecular weight values are, and vice versa.

A target molecule’s synthetic accessibility (SA) score is used to evaluate it before it is synthesized. It is a parameter used to rate the ease or difficulty of synthesizing a compound. The SA score runs from 1 (very easy) to 10 (very tough) (Swierczewska et al., 2015). The studied compounds possess good SA, with SA scores ranging from 1.52 to 4.13 as shown in Table 8.

Table 8.
Drug-likeness prediction output of test compounds.

5-FLU: 5-flurouracil; MW: molecular weight; MR: molar refractivity; HBA: hydrogen bond acceptor; HBD: hydrogen bond donor; TPSA: topology polar surface area; MR: molar refractivity; SA: synthetic accessibility.

The partition coefficient between n-octanol and water (log Po/w) is the fundamental descriptor for lipophilicity (Log P). The tested compounds’ log P values, water solubility, and bioavailability are reported in Table 9. Lipophilicity and water solubility are two important physicochemical qualities that influence a drug’s behavior. An orally delivered medicine should be lipophilic enough to pass through the intestinal lining, permeate target cell membranes, and be hydrophilic enough to move through the aqueous blood. The higher a compound’s log P value, the more lipophilic it is and the less water soluble it is (Johnson et al., 2021). The log P values for the studied compounds ranged from –0.13 (5-FU) to 5.77 (C4), which revealed that C4 is not as soluble as the referenced drug (Table 9).

The Silicos-IT prediction model was used to calculate water solubility as the logarithm of molar solubility in water (log S). The log S value of a drug affects its capacity to dissolve, and the lower the value, the better (Darvas et al., 2002). The compounds had log S values ranging from –0.76 (poorly soluble, C4) to –3.5 (soluble, C10) indicating. This trend is in accordance with the values of log P earlier discussed.

The bioavailability score employs the total charge, the TPSA, and the Lipinski filter to provide a semiquantitative indication of the compounds’ likelihood of being effective oral medicines (Testa and Krämer, 2009). The less the drug candidate can penetrate, the higher the TPSA and molecular weight values are, and vice versa. The bioavailability score of the test compounds is 0.55, which means that these compounds have a 55% chance of a minimum of 10% oral absorption in rat or human colon carcinoma (Caco-2) permeability (Testa and Krämer, 2009) indicating that they could be useful as oral drug.

Table 9.
Predicted lipophilicity (Log P) values, water solubility and bioavailability of the studied compounds.

5-FLU: 5-flurouracil; BS: bioavailability score.

Table 10 showed the result of pharmacokinetics prediction of the test molecules. P-glycoprotein (P-gp), a well-studied plasma membrane ATP-binding cassette transporter, is responsible for the active efflux of xenobiotics through biological membranes to protect the body against foreign substances. This efflux pump contributes to drug resistance by preventing some medications from entering sensitive regions. Except for C13 and C14, all the chemicals had a high likelihood of being absorbed in the gastrointestinal tract. This shows that when these chemicals are taken orally, they have the potential to be absorbed in the gastrointestinal tract (Lynch and Price, 2007).

The enzyme cytochrome P450 monooxygenase plays an important role in the metabolism and removal of drugs. It is a group of isoenzymes that catalyzes several processes in the first phase of drug metabolism (Hollenberg, 2002). The suppression of five main isoforms (CYP1A2, CYP2C19, CYP2C9, CYP2D6, and CYP3A4), which are substrates of 50–90% of medications and primary sources of pharmacokinetics-related drug-drug interactions. Because these enzymes are not inhibited, the chemicals have a high chance of being converted and, as a result, of becoming bioavailable when taken orally. Inhibition of the CYP isomers by the substances, on the other hand, can result in poor bioavailability due to failure to be metabolized and severe side effects due to their retention (Sahin and Benet, 2008). Compounds C1, C2, and C9 are potential inhibitor of CYP1A2 among the tested compounds. All the compounds are potential inhibitor of CYP2C19 except compound C2, C10, and standard 5-FU. As shown in table 10, only compound C2 and 5-FU (referenced drug) were non-inhibitors of CYP2C9 while other studied compounds proved to possess the ability to inhibit CYP2C19. Compounds C1, C2, C3, C6, C9, and C12 are potential inhibitor of CYP2D6 while others are not. All compounds are potential inhibitors of CYP3A4, except compounds C4, C7, C8, 5-FU. It was observed that the standard 5-FU was is not an inhibitor of any cytochrome P450, which suggests that the standard is more likely to be converted and bioavailable after oral delivery than the compounds investigated.

The skin is a selective barrier that allows diverse compounds to permeate through at different rates based on their physicochemical properties. Hence, the skin permeability (LogKp) is a key parameter for the evaluation of molecules that might require transdermal administration (Dixon et al., 2006). The more negative the log Kp (with Kp in cm s^–1), the less skin permeant is the molecule. As highlighted in Table 10 log Kp (cm s^–1) of the test compounds ranged from –7.73 (least permeant) to –4.62 (most permeant).

The skin acts as a selective barrier, allowing various molecules to infiltrate at varied rates depending on their physicochemical qualities (Cheng and Dixon, 2003). As a result, skin permeability (LogKp) is an important characteristic to consider when evaluating compounds that may need to be administered transdermally. The lower the log Kp (in cm s^–1), the less permeant the molecule is to the skin. The test compounds’ log Kp (cm s^–1) ranged from –7.73 (least permeant) to –4.62 (most permeant), as shown in Table 9. The values indicate that the compounds show low skin permeability.

Table 10.
Pharmacokinetics prediction output of test compounds.

5-FLU: 5-flurouracil; GI: gastrointestinal; and P-gp: P-glycoprotein.

The toxicity prediction output of the compounds is shown in Table 11. Compounds C8, C12, C13, and C14 fall under the oral toxicity class 5, while other tested compounds fall under oral toxicity class 4. Human hepatotoxicity (H-HT) is a term that describes and illustrates numerous types of liver damage that can lead to the organ failing or death (Cheng and Dixon, 2003). Compounds C10, C11, C12, C13, and C14 were hepatotoxic. Carcinogenicity compounds can damage the genome, the disrupt cellular metabolic processes and cause cancer. Compounds C2, C3, C4, C5, C6, C11, C12 and the standard are carcinogenic. Immunotoxicity is defined as the negative impact of toxic compounds on the functioning of both the local and systemic immune systems. Only compounds C1, C9, C10 and the standard are nonimmunotoxic. The mutagenicity test is used to identify substances that can cause mutations or malignant growth in humans (Loving et al., 2009). Compounds C2, C3, C4, C5, C6, and C12 are mutagenic, and only compound C2 is cytotoxic among the test compounds. The study clearly revealed that most of the compounds possess one form of toxicity except for compounds C1 and C9, which did not show any tendency for hepatotoxicity, carcinogenicity, immunotoxicity, mutagenicity and cytotoxicity, indicating that they are safe as possible therapeutic agents. These two compounds did not violate any of the base rule filter, they are moderately soluble and not P-gp substate.

Table 11.
Toxicity profiles of test compounds.

5-FLU: 5-flurouracil.

abeloyebamiji@gmail.com