Virtual screening
In the present study, the main objective is the identification of potential new drugs as E6 inhibitors. A total of 6340 molecules that already have anticancer activity are processed by setting up drug-like filters to evaluate the similarity of these molecules with drugs, a high-throughput screening was reread through PyRx software (Dallakyan & Olson, 2015). The chemical database was subjected to a simple filter to remove non-drug molecules to generate targeted databases for HPV16 E6. Accordingly, using Lipinski based filter, we have removed 962 non-suitable molecules chemicals allowing to narrow the database to 5378 drug-like compounds.
The prepared chemical database was subjected to high-throughput screening to identify new candidate molecules potentially suitable for E6/E6AP/P53 inhibition, using AutoDoc vina involved in the PyRx tool which generated 9 distinct conformations for each ligand, ranked by binding affinity (kcal/mol).
Molecular Docking
Autodoc vina and MGL tools were used to dock the 28 candidate molecules for HPV16 E6 treatment. This resulted in a final list 3 compounds; F0679-0355, F3345-0326 and F3374-0275; showed the highest affinities: -10.4Kcal/mol, -10Kcal/mol and -9.9Kcal/mol, respectively. Total energy, hydrogen bonds (HBond) and other interactions of the 3 selected compounds are reported in Table I. From the molecular docking results of the three candidate compounds, we found that the majority of these compounds share interactions with the same amino acids responsible for the formation of the E6/E6AP/p53 complex, and are involved in Van Der Waals interaction and pi-alkyl and pi-sigma interactions with residues Val31, Tyr32, Leu50, Cys51, Val53, Val62, Leu67, Tyr70, Ile73, Leu96, Cys97, Asp98, Leu99 and Leu100 of E6 contribute to the LxxLL binding pocket. Whereas, Gln6, Glu7, Arg8, Arg10, Gln14, Glu18, Arg40, Glu41, Val42, Tyr43, Asp44, Phe45, Ala46, Phe47, Asp49, Leu50, Leu100, Cys106, Gln107, Lys108, Pro109, Leu110, Cys111, Pro112, Glu113, Gln114 and Lys115 of E6 contribute to the p53 binding pocket which is consistent with recent studies (table 1).
Table 1. Molecular Docking result analysis of top 3 compounds
Compound name
|
Molecular
Formula
|
PubChem CID
|
Hydrogen bonds
|
Other interactions
|
Binding
affinity
|
F3345-0326
|
C24H19N5O2
|
4903840
|
Ser(71)
|
Val(53), Val(62), Leu(50), Cys(51), Phe(45), Val(31), Tyr(32), Ser(74), Ile(73), Arg(129), His(78), Tyr(70), Arg(131), Leu(67), Gln(107), Arg(102)
|
-10.4
|
F3374-0275
|
C22H13N3O3
|
17016408
|
Cys(51),
Arg(102)
|
Gln(107), Leu(67), Tyr(32), Val(31), Val(31), Val(62), Phe(45), Val(53), Leu(50), Gln(130), Trp(132), Leu(100), Arg(131)
|
-10
|
F0679-0355
|
C19H20FN5O2
|
3155624
|
Arg(131),Gln(107),Tyr(32)
|
His(78),Ser(74),Ser(71),Leu(67),Arg(77)Ile(73),Tyr(70),Leu(50),Val(62),Val(31),Cys(51),Ala(61),Val(53),Phe(45)
|
-9.9
|
Molecular docking of the selected molecule F0679-0355 on HR E6 oncoprotein is reported in Fig 1 and showed that this ligand forms three strong hydrogen bonds with the amino acids Arg(131), Gln(107) and Tyr(32) of the HR E6 protein via 8 van der Waals bonds with the amino acids His(78), Ser(74), Ser(71), Arg(77), Tyr(70), Val(53), Val(31) and Phe(45). The docking study showed also the formation of 2 Halogen bonds between the ligand and the amino acids Cys(51) and Ala(61), and interacting via a π -Alkyl with the amino acids Ile(73) and Leu(50).
Results of molecular docking of the selected molecule F3345-0326 on HR E6 oncoprotein is reported in Fig 2 and clearly showed the formation of one strong hydrogen bond with the amino acids Ser(71) of the HR E6 protein and 9 van der Waals bonds with the amino acids Arg(129), His(78), Ser(74), Arg(131), Tyr(70), Leu(67),Val(53),Val(31),Gln(107), and 3 interacting via a π -Alkyl with the amino acids Val(62), Leu(50), Ile(73). The docking study also showed the formation of 1 π-Sulfur bonds between the ligand and the amino acids Cys(51)
Molecular docking of the selected molecule F3374-0275on HR E6 oncoprotein is reported in Fig 3 and showed the formation of two strong hydrogen bonds with the amino acids Cys(51) and Arg(102) of the HPV16 E6 protein and via 6 van der Waals bonds with the amino acids Gln(107), Leu(67), Tyr(32), Val(31), Phe(45) and Gly(130), and interacting via a π -Alkyl with the amino acids Val(53), Val(62), Leu(50) and Arg(131), and via π –sigma with Leu(100). The docking study showed also the formation of one carbon-hydrogen bond between the ligand and the amino acid Trp(132).
In silico pharmacokinetics (PK)/pharmacodynamics (PD) evaluation
The purpose of the ADMET preclinical study is to remove poorly performing molecules and focus on the best drug candidates. The present study investigated the suitability of the five proposed compounds for use as anti-cancer drugs based on properties, absorption, distribution, metabolism, excretion, and toxicity, which are critical elements in drug development.
These pharmacokinetic properties were obtained using the admet SAR approach evaluating blood-brain barrier (BBB) permeability, human intestinal absorption (HIA), Caco-2 cell permeability and the AMES assay that are the main elements used to determine drug capacities and results are reported in Table 2. Crossing the blood-brain barrier (BBB) is a crucial property widely used to determine the usefulness of the chemical as a drug, indicating whether drugs can cross the blood-brain barrier (Upadhyay, 2014). In this study, the 3 evaluated compounds have shown low BBB index and were therefore considered poorly distributed in the brain.
ADMET analysis showed also that the 3 ligands can be absorbed by the human gut, are not metabolized by most cytochrome P450 monooxygenases and don’t exhibit any acute toxicity, mutagenic effect or carcinogenic potential.
Table 2. Analysis of the pharmacokinetic properties and toxicity of the 3 candidate molecules
Compound name
|
F0679-0355
|
F3345-0326
|
F3374-0275
|
Absorption and Distribution
|
blood-brain barrier
|
0.5177
|
0.6933
|
0.9385
|
human gut absorption
|
0.9953
|
0.8780
|
0.9868
|
Caco-2 permeability
|
0.6748
|
0.7412
|
0.6142
|
Substratglycoprotéine P
|
Yes
|
Yes
|
Yes
|
inhibitor of glycoprotein P
|
No
|
No
|
No
|
Metabolism
|
CYP450 2C9 Substrate
|
No
|
No
|
No
|
CYP450 2D6 Substrate
|
No
|
No
|
No
|
CYP450 3A4 Substrate
|
Yes
|
Yes
|
Yes
|
CYP450 1A2 Inhibitor
|
No
|
No
|
No
|
CYP450 2C9 Inhibitor
|
No
|
Yes
|
No
|
CYP450 2D6 Inhibitor
|
No
|
No
|
No
|
CYP450 2C19 Inhibitor
|
No
|
Yes
|
No
|
CYP3A4 Inhibitors
|
No
|
No
|
No
|
Excretion and Toxicity
|
Hepatotoxicity
|
No
|
No
|
No
|
Carcinogens
|
No
|
No
|
No
|
AMES mutagenicity
|
No
|
No
|
No
|
Molecular dynamics
In this study, the in-depth pharmacokinetic analysis identified three candidate compounds for HPV16 E6 inhibition with the most advantageous binding interactions and best pharmacokinetic profiles. The Desmond simulation trajectories were analyzed, and the root mean square deviation (RMSD) Fig 4, root mean square fluctuation (RMSF) Fig5 and protein-ligand contacts were calculated from the MD trajectory analysis.
The left Y-axis shows the variation of protein RMSD through time. The right Y-axis shows the variation of ligand RMSD through time. A) RMSD of F0679-0355_4xr8 complex. B) RMSD of F3374-0275_4xr8 complex. C) RMSD of F3345-0326_4xr8 complex
As reported in Fig 4, the evolution of RMSD values with time for the C-alpha atoms of three complexes F0679-0355_4xr8, F3374-0275_4xr8 and F3345-0326_4xr8, indicates that the 3 complexes reach stability at 10 ns. From then, changes in RMSD values remain within 1.0 Å for the target (4xr8) during the simulation period, which is quite acceptable. Ligands fit to protein RMSD values fluctuate within 1.5 Å till 50 ns after being stable. These indicate that the ligands remained stably bound to the binding site of the receptor during the simulation period. Of particular interest, F0679-0355_4xr8 complex showed better results comparatively.
Figure 5 manifests the RMSF value per residue of the target bound to the ligands. Residues showing larger peaks belong to the loop regions, as determined from the MD trajectories (Fig 6), or to the N- and C-terminal regions. The target bound to F0679-0355 showed a comparatively better RMSF, as shown in Figure 5A. The low RMSF values of the binding site residues show the stability of the ligand binding to the protein.
The majority of the essential ligand-protein interactions determined by MD are hydrogen bonds and hydrophobic interactions, as illustrated in Fig 4. The stacked bar graphs have been normalized over the pathway: therefore, a value of 1.0 indicates that for 100% of the simulation time the particular interaction was conserved Fig 7. Values above 1.0 are likely to be achieved, as some protein residues could have multiple contacts of the same subtype with the ligand.