Wnt signaling should be blocked to authorized natural heart regeneration; hence, finding potential small molecules inhibitor for one of these pathway components like TNKS could activate heart regeneration after injury. In this regard, virtual screening was performed over approved drugs, regulated chemicals, and herbal isolates from the SWEETLEAD database over TNKS1 and TNKS2. The structures were sorted and ranked based on the docking binding energy (
Tables 1 and
2).
Identification of the Active site
Tankyrase proteins stabilized the β-catenin during the degradation process with Axin protein (
29). Tankyrase belongs to the poly(ADP)-ribose polymerases (PARPs) family, which utilizes NAD
+ to catalyze ADP-ribose (ADPr) on target proteins. Tankyrase has two isoforms: tankyrase 1 (TNKS1) and tankyrase 2 (TNKS2), also known as ADP-ribosyltransferases ARTD5 and ARTD6, respectively. Tankyrase contains three particular domain, five ankyrin repeat (ANK) domain with regulating role in protein-protein interaction, the sterile alpha motif domain (SAM) for oligomerization located N terminal and the catalytic ADP-ribosyltransferase (ARTD) domain at the C-terminus (
30-
32). TNKS1 and TNKS2 indicate 85% similarity in amino acid sequences and have 1327 and 1166 residues, respectively. The catalytic domain shares 89% overall sequence identity and contains three conserved central amino acids: His1184, Tyr1213, and Glu1291 in TNKS1, and the residue His1031, Tyr1060, and Glu1138 in TNKS2. Except for the common domain, as noted before, TNKS1 has a special His, Pro, Ser rich (HPS) region in the N-terminus with unknown function
Figure 2 (
33).
CAVER 3.0.1. software was used to predict the possible tunnels of the TNKS1 and TNKS2. The results indicate that for each protein, one tunnel leads into the key residues
Figure 3. The tunnel identified by CAVER for TNKS1 is much longer than the predicted tunnel of the TNKS2. Despite the sequence similarity between TNKS1 and TNKS2, the binding groove shapes make a difference in virtual screening results.
Docking based virtual screening
Virtual screening was performed over the 9127 ligands of the SWEETLEAD database over TNKS1 and TNKS2, and the obtained structures were ranked based on the docking binding energy (
34). The top ten screened compounds were listed in
Tables 1 and
2. Docking binding energy has shown the score ranging from -12.84 to -11.5 kcal/ mol and -11.70 to -11.08 kcal/ mol for TNK1 and TNKS2.
The result indicates that in TNKS1, the ligands with longer carbon side chains and higher volume were placed well in the binging groove, but in TNKS2, the ligands with shorter side chains were placed following the binding track volume were able to fit in the binding groove. Interestingly, on average, the binding groove of TNKS1 indicates more hydrophobic properties than TNKS2, and the compounds with the higher value of miLogP could fit well in the binding groove.
Three top ligands of each screening were investigated for their interactions. In general, lipophilic interactions and van der Waals interaction play an essential role in binding. The amino acids involved in binding phenothrin and TNKS1 are critical residue Tyr1213 with van der Waals and π-alkyl interactions, Tyr1224 with π-sigma, and π-alkyl interaction. Lots of π-alkyl or alkyl interactions with His1184, Ile1192, Tyr1203, Ile1204, Ala1215, Ile1212, and critical residue Glu1291 is also forming van der Waals interactions. Regarding bornyl trans-cinnamate, hydrogen bond with Ile1204 and π-π stacking with Tyr1224 and several π-alkyl interactions with His1184, Phe1188, Ile1192, Phe1197, His 1201, Tyr1203, and Ile1212 were observed. Eletriptan indicates four hydrogen bonds with crucial residue His1184, Gly1185, His1201, and Gly1211. Ile1212 form π-σ and π-π stacking interactions were seen with Tyr1224 and π-alkyl interactions with Ala1215 and Lys1220.
TNK2 interact with ethyl rosinate (Abietic acid ethyl ester) with van der Waals interaction of Pro1034 and several and π-alkyl interactions with key residue His1031, Phe1035, Tyr1050, essential residue Tyr1060, Ala1062, Lys1067, and Tyr1071. Bornyl trans-cinnamate forms hydrogen bond with essential residue Tyr1060, π-π stacking with Tyr1071, and lots of π-alkyl interactions with key residue His1031, Pro1034, Phe1035, Ala1062, Lys1067, and Ile1075. Allylestrenol indicates hydrogen bond with Ala 1049 and several π-alkyl interactions with critical residue His1031, Pro1034, Tyr1050, essential residue Tyr1060, Ala1062, Lys1067, and Tyr1071
Figure 4.
Molecular dynamics simulation
The top compounds of each list, phenothrin, and ethyl rosinate in complexes with the TNKS1 and TNKS2, were further evaluated for stability by molecular dynamics simulations for 100 ns, and the structural deviations were calculated. The backbone atom deviation was applied to determine the stability of both complexes. RMSD plot shows that both complexes stabilized after 20 ns, but the TNKS1 in complex with phenothrin indicated slight deviations and reached a steady station faster
Fig 5A. The same trend is seen for the RMSD plot of ligand, and the phenothrin reached a stable state very soon and maintained stability in the entire simulation period
Figure 5B. The average RMSD value of phenothrin in complex with TNKS1 was 0.2 nm, close to the obtained average RMSD value of 0.136 nm for the ligand ZINC28852318 (N-(4-morpholinophenyl)-3-(3-oxopiperazin-1-yl)propanamide) reported recently (
35).
RMSF profile indicates that the active site’s vital amino acids do not have flexibility in both complexes, and the occupation of the ligands reduced the fluctuation of the active site. However, just as the RMSD plot showed that the complexes of TNKS1 with phenothrin had stabilized earlier, there is less fluctuation in the RMSF plot indicating the ligand’s ability to deprive the protein of flexibility
Figure 6The calculated number of the hydrogen bonds shows that both ligands had at least one hydrogen bond during the simulation period except the 30 ns of simulation of the complexes of phenothrin. However, the number of hydrogen bonds has reached three at the end of simulations, but the ethyl rosinate established only one hydrogen bond during the simulation (
Figure 7).
MM-PBSA binding free energy calculations
G_mmpbsa package was used to calculate the binding free energy, which is an important indicator that accounts for the potential affinity of the ligand to the receptor. In general, complexes with lower binding free energy can be considered more stable, and their ligands are expected to have a higher inhibitory effect and potency. Each conformation of the MD simulation trajectory was used to calculate the binding free energies by the MM-PBSA method.
Accordingly, the binding free energy for the two complexes of phenothrin in complexes with the TNKS1 and the ethyl rosinate in complexes with the TNKS2 were employed over the MmPbSaStat.py python script (
28). The total free binding energy of each complex component, i.e., the energy of the complex, protein, and ligand, was calculated by this script.
Furthermore, the cumulative sum of the molecular mechanic’s potential energy in a vacuum and the free energy of solvation includes the polar (electrostatic) and nonpolar (non-electrostatic) solvation energy are the components of free binding energy calculations. The nonpolar solvation energy was usually calculated by the model of solvent-accessible surface area (SASA). All types of energies and the value of standard deviation were calculated by the g_mmpbsa package. Then to obtain the average of the free energy of each component, they were summed together. Finally, the total binding free energy is obtained by subtracting the receptor and ligand’s total free energy from the complex’s total free energy.
Table 3 summarizes the interaction energies and the binding free energy for the two complexes.
The results of MM-PBSA calculation of the Free Gibbs energy of the phenothrin in complexes with the TNKS1 indicate the slightly higher binding affinities to the TNKS1 than the ethyl rosinate in complexes with the TNKS2. Phenothrin in complexes with the TNKS1 was better than ethyl rosinate in complexes with the TNKS2 in all the calculated energy formats except (electrostatic energy). Phenothrin average binding free energy reached –115.18 kJ/mol, while ethyl rosinate average binding free energy reached –111.87 kJ/mol, which is similar to the trend of docking results in binding energy. The overall results of the dynamic simulations supported our design concept and validated the entire virtual screening approach; they also emphasized the potential inhibitory effect of phenothrin on TNKS1.
ADME properties
SwissADME server was used to determine the drug-likeness of selected compounds phenothrin and ethyl rosinate (
36). As
Table 4 shows, both compounds fulfill the criteria of Lipinski’s Rule of five (
37) in terms of hydrogen bond donors (HBD) and hydrogen bond acceptors (HBA) logP and molecular weight (MW). The solubility class of phenothrin is poorly soluble, but ethyl rosinate indicates moderate solubility. Both ligands show topological polar surface area (TPSA) values lower than 140 Å
2, which means good cell permeability. Therefore, both compounds indicate the potential of drug-likeness properties.