Computational study on Cross-talking Cancer Signalling Mechanism of Ring Finger Protein 146, Axin and Tankyrase Protein Complex

Lakshmanan Loganathan, Keerthana Natarajan & Karthikeyan Muthusamy

To cite this article: Lakshmanan Loganathan, Keerthana Natarajan & Karthikeyan Muthusamy (2019): Computational study on Cross-talking Cancer Signalling Mechanism of Ring Finger Protein 146, Axin and Tankyrase Protein Complex, Journal of Biomolecular Structure and Dynamics, DOI: 10.1080/07391102.2019.1696707
To link to this article: https://doi.org/10.1080/07391102.2019.1696707
Full Terms & Conditions of access and use can be found at https://www.tandfonline.com/action/journalInformation?journalCode=tbsd20

Computational study on Cross-talking Cancer Signalling Mechanism of Ring Finger Protein 146, Axin and Tankyrase Protein Complex
Lakshmanan Loganathan, Keerthana Natarajan, Karthikeyan Muthusamy* Department of Bioinformatics, Alagappa University, Karaikudi, Tamilnadu, India.
*Corresponding Author: Dr. Karthikeyan Muthusamy, Assistant Professor, Department of Bioinformatics, Alagappa University, Karaikudi, Tamil Nadu, India, 630 003. Email: [email protected], Phone No: +91-4565-223344, FaxNo:+91-4565-225202

Abstract
Cancer is distinguished by uncontrolled cell growth and it is regulated by several environmental and genetic factors. The Wnt β-Catenin signaling pathway has been considered as the most significant colon cancer-targeted pathway. AXIN plays a major regulatory role in the Wnt signaling mechanism. The AXIN after PARsylated by TNKS is ubiqutinated by RNF146 through its WWE domain that leads to degradation of AXIN protein. Several studies have been proposed highlighting the inhibition of the PARsylation mechanism that mediates the degradation of AXIN and improves β-catenin stability. The present study focused on the identification of potential inhibitors for the inhibition of RNF146-TNKS complex through identifying potential RNF146 inhibitors to prevent ubiquitination of AXIN, further to confirm the regulatory role and inhibition mechanism of RNF146-AXIN and RNF146-TNKS. The docked complex was then evaluated using various computational analysis. Molecular interactions analysis was performed to observe the interacting residues between the protein complex. The compounds from various databases were docked with the RNF146 and complex proteins. Both the protein complex and ligand were analyzed for the confirmation of structural stability using molecular dynamics simulations. Selected compounds’ atomic configuration and electron profile were analyzed through DFT calculations and ADME/T (Physico-chemical) properties. As a result, we found several common lead compounds for RNF146, TNKS protein inhibition. Therefore, the docked compounds may act as a better antagonist molecule for RNF146, TNKS and associated signaling molecules. Further, experimental validations are required to prove the potency of the identified compounds.

Keywords: Cancer, Wnt β-Catenin signaling pathway, Molecular Docking, RNF146-AXIN, RNF146-TNKS, Molecular Dynamics.
List of Abbreviations: DFT – Dendity Functional Theory; HOMO – Highly Occupied Molecular Orbitals; LUMO – Lowest Unoccupied Molecular Orbitals; MD – Molecular Dynamics; MESP – Molecular Electrostatic Potentials; NCI – National Cancer Institute; PAR – poly-ADP-ribosylation; PDB – Protein Data Bank; RMSD – Root Mean Square Deviation; RMSF – Root Mean Square Fluctuation; RNF146 – Ring Finger Protein 146; SASA – Solvent Accessible Surface Area; TNKS – Tankyrase; TCM – Traditional Chinese Medicine; VSW – Virtual Screening Workflow; XP – Extra Precision;

Introduction
Worldwide, Cancer is a major public health issue, nearly 1 in 6 premature death was caused by cancer(Plummer et al., 2016). There are over 100 different types of cancer, and each is classified by the type of cell that is initially affected. According to the report, USA cancer mortality and morbitity rate is increased, about 1,685,210 and 5,95,690 respectively(Siegel, Miller, & Jemal, 2019). By 2030, the estimated burden of cancer incidence rate likely to be increased with 23.6 million new cases. Increasing incidence rate is not a regional specific issue, invariably it affects all type of racial populations(Bray et al., 2018). The Indian Council of Medical Research projected new cancer cases and death cases in India are 37,250 and 13,820 and also estimated that about 1.45 million patients will develop cancer with the number expected to rise to 1.73 million by 2020(Imran, Waseem, & Kishwar, 2011). The choice of treatment for colorectal cancer can depend on several factors, including the patient’s health, the size of the tumor, and its location. Surgery is the most common treatment option, and the type of surgery used, again, depends on variables such as the location of the cancer and the existence and extent of metastasis.

Wnt signaling mechanism is one of the major regulatory signal of colon cancer. It is also plays a major role in embryonic development and also controls homeostatic self-renewal in various tissues(Kirubakaran, Arunkumar, Premkumar, & Muthusamy, 2014; Steinhart & Angers, 2018; Zhan, Rindtorff, & Boutros, 2017). It causes several types of hereditary diseases and cancer through germline and somatic cell mutations. The Wnt canonical pathway regulates the production and stability of β-catenin protein that has crucial roles in cancer development, whereas the non- canonical pathways have no role in β-catenin stabilization(Mariotti, Pollock, & Guettler, 2017). Inhibition of the key proteins in Wnt β-catenin pathway has a chance of greater therapeutic value in treating cancer patients. Ring finger protein 146 (RNF146), a novel PARsylation-directed ring finger E3 ubiquitin ligase, is located at 6q22.1-q22.33 of human chromosome(Kang et al., 2011). RNF146 is a really interesting new gene (RING) finger protein that contains a RING and WWE domains. RNF146 was renamed Iduna and was recently shown to possess a PAR binding motif (PBM) and protects against parthanatos via binding to PAR. Ring finger proteins contain ring fingers, which are considered to be the functional module for E3 ubiquitin ligase activity(Zhou, Chan, Xiao, & Tan, 2011).

RNF146 was confirmed to be the novel PARsylation-directed E3 ubiqitin ligase to induce ubiquitination of PAR modified axin and positively regulate Wnt signaling. WWE domain-containing proteins are tightly linked with PARsylation and ubiquitination signalling pathways. RNF146 acts in cooperation with tankyrase proteins (TNKS), which mediate PARsylation of target proteins. Parsylated substrate is likely held in the RNF146-TNKS complex via its interaction with TNKS, as the RNF146-TNKS interaction is required for ubiquitination of substrate in cells. RNF146 strongly induced both K48- and K63-linked ubiquitylation of TNKS protein, whereas AXIN showed predominantly K48-linked ubiquitin. Literatures concluded that RNF146 directly interacts with the PAR polymerase tankyrase (TNKS). TNKS parsylates RNF146 and in reciprocal RNF146 ubiquitinates TNKS(Zhang et al., 2011).

The literature study suggested that RNF146 ubiquitinates PARsylation dependent AXIN that leads to the degradation of β-catenin(DaRosa, Klevit, & Xu, 2018). In this study, targeting the RNF146, TNKS and Axin simultaneously to distract and inhibit the co-expression and complex formation (Figure 1). The figure 1, elucidate the concept in a diagrammatic representation which explains the identified inhibitors will be directed to inhibit the RNF146 mediated Ubiquitination and TNKS mediated PARsylation of Axin. Thereby, the axin level increase to form a destructive complex that performs proteosomal degradation the β-catenin. The computational power was utilized to model the RNF146-Axin, RNF146-TNKS complex through moldecular docking approach and thereby finding the potent inhibitors that are highly specific to the binding cavity of the complex. Finding a common inhibitor that has specificity to all three molecules will improve the inhibition of PARsylation and Ubiquitination in Wnt-β-catenin signaling cascade.
Material and Methods Protein structure datasets

The structure of the RNF146 (WWE domain) (PDB ID: 3V3L), AXIN (RGS domain) (PDB ID: 1DK8) and TNKS (Catalytic domain) (PDB ID: 4UFY) was retrieved from the RCSB Protein Data Bank. Usually, the PDB structures may contain co-crystallized ligands, water molecules, metal ions, cofactors and also it may have some missing side-chains, along with the improper charges and bond orders(Madhavi Sastry, Adzhigirey, Day, Annabhimoju, & Sherman, 2013). In order to get a proper structure, the obtained proteins were prepared using protein preparation wizard module of Maestro, Schrodinger 2018-4.

Ligand datasets and Preparation
Chemical compounds from Enamine, May Bridge, National Cancer Institute (NCI), and Traditional Chinese Medicine (TCM) database were included in this study. Approximately, 15,34,000 compounds were taken for three-dimensional (3D) conversion and minimization of chemical compounds using Ligprep, Maestro, Schrodinger 2018-4. Conformers were generated using a rapid torsion angle search approach followed by minimization of each generated structure using the OPLS-3e force field, with 30 implicit GB/SA solvent models(Roos et al., 2019). A maximum of 32 conformers were generated per structure using preprocess minimization of 100 steps and post process minimization of 50 steps. Each minimized conformer was filtered through a relative energy than window of 10 kcal/mol and a minimum atom deviation of 1.00 Å. This value (10 kcal/mol) sets an energy threshold relative to the lowest- energy conformer(Pérez-Benito, Keränen, Van Vlijmen, & Tresadern, 2018). Finally, 25,20,148 compounds including the generated conformers from the four database were taken for molecular docking studies.

Active Site Prediction
The active site of the target proteins (RNF146, TNKS and AXIN) contains the catalytic and binding sites are essential for molecular docking analysis. The target structures were taken for the active site prediction using SiteMap module of Schrodinger2018-4(Halgren, 2009). In oder to predict or identify the selective amino acids that participates in the protein complex interaction. The amino acids, which are located in hydrophobic, hydrogen bonds acceptor and donor rich sites will be preferable to locate the active site of the protein. The sitemap provides five top ranked sites in the protein, among them best site were selected based on the sites score and the area coverage. The active site centre obtained by active site determination programs was used as a grid center in the grid parameter file, rest all the parameters were set to default values.

Protein-protein docking
Protein-protein docking methods is a significant approach to understand the protein networks and their non-covalent inter-protein interactions. For each new configuration, the energy of interaction is calculated based on terms such as surface complementarities, electrostatic interactions, van der Waals interaction and additional terms depending on the method developed. In this study, RNF146 with TNKS and AXIN were docked independently using HADDOCK server version 2.2 to study the complex protein interactions and conformational changes. We have utilized the easy interface of HADDOCK server for defining the active and passive residues of both receptor and ligand protein. The interaction interface obtained from the Sitemap program was used to define the active residues. The active and passive residues were also verified with the reported structural studies of the proteins. The results of the Haddock run provide several clusters of structure to validate the minimized protein complex. The best protein complex was selected based on the Haddock score with less RMSD score and good z-score(Van Zundert et al., 2016). A docking method that incorporates experimental data to dock the two protein structures was developed by Domingues et al., 2003(Dominguez, Boelens, & Bonvin, 2003).

Structure Based Virtual Screening
The structure based virtual screening workflow is used to screen the large databases of chemical compounds to filter out the possible lead molecules. The prepared compounds from different databases such as Enamine, Maybridge, NCI and TCM (25,20,148) were screened in the active site of the RNF146 (PDB ID: 3V3L), RNF146-Axin and RNF146-TNKS complex using Glide-VSW module. The virtual screening workflow is a well defined protocol to filter the chemical compounds by screening and docking, initially Lipinski’s and Qikprop filter was applied, followed by HTVS (High Throughput Virtual Screening), SP (Standard Precision) and finally XP (Extra Precision) docking programs was used to screen the potential ligands. Glide HTVS and SP docking uses a series of hierarchical filters to find the best possible ligand binding locations in a previously built receptor grid space. Further, the poses generated from Glide SP were again refined in the final process called Glide XP(Friesner et al., 2006). In each stage of HTVS, SP and XP the top 10% compounds were only allowed to next stage. During pulling the screened compounds to the next stage, the HTVS retains all the good state and best
scoring state of ligand poses, whereas the SP retain only good scoring state and XP retain hardly the best scoring state respectively. For this study, RNF146, RNF146-TNKS and RNF146-Axin complex was used as the receptors. Identifying the common inhibitor for all the complex molecules would be the potential one to propose a promising antagonist molecule to inhibit the complex signaling pathways. Therefore, the database compounds were screened individually with all the three receptors and the results were scrutinized to find the common inhibitors. The top-hit compounds were ranked according to the Glide XP score, Glide energy and interacting residues of protein with the compound. The Glide XP score is the total GScore or final score to identify the best ligand. The identification of the best pose is calculated through the nonbonded interactions, whereas the Glide energy is the score of modified Coulomb-van der Waals interaction energy. The Glide energy is useful to compare the binding affinities of different ligands. The top hit ten ligands from all the four database were further taken for various analysis(Loganathan, Muthusamy, Jayaraj, Kajamaideen, & Balthasar, 2018).

Binding Energy calculation
The binding energy calculations was done for the top ten protein-ligand complex taken from the top-hit list in structure based virtual screening using Prime/MM-GBSA approach(Cappel et al., 2016). MM-GBSA, a post scoring approach helps in the evaluation of molecular docking and in the validation of the accuracy of different charges. This method is helpful to predict the binding energy for a set of ligands poses to the receptor. Following equations were used for the calculation of the binding energy: ΔGbind = ΔE + ΔGsolv + ΔGSA (1) ΔE = Ecomplex –Eprotein -Eligand (2) Where Ecomplex, Eprotein, and Eligand are the minimized energies of the protein-inhibitor complex, protein, and inhibitor, respectively; ΔGsolv is generalized born electrostatic solvation energy of the complex; and ΔGSA is a non-polar contribution to the solvation energy due to the surface area(Nagamani, Kesavan, & Muthusamy, 2012).

Molecular dynamics of protein-ligand complex and protein-protein complex
Molecular dynamics simulation was performed on a cluster of 8 processors running on CentOS platform. To obtain the structural stability and conformational features on protein-ligand and protein-protein complexes, the molecular dynamics simulation was carried out usin GROMACS (GROingen Machine for Chemical Simulation) version 5.0.4 with GROMACS96 43al force field(Loganathan et al., 2019). Topologies of small molecules were generated by using PRODRG server with full charges. The complexes were relaxed to eliminate unwanted atomic contacts and then embedded in dodecahedron water box of simple point charge216 (spc216) water model(Loganathan & Muthusamy, 2018). The system was neutralized by adding appropriate Na+ and Cl- ions. Steepest descent algorithms have been used to minimize the system to remove unwanted van der Waals contacts and to reduce the net force through motion of each atom before taking the same to the molecular dynamic simulation studies. The particle-mesh Ewald (PME) method was used to enlighten the system in a normal cutoff of 0.9 nm. Then the solvated system is subjected to the constant temperature and pressure bath. The constant temperature of 300k and a pressure of 1 bar; the Berendesen thermostat were maintained using a coupling time of 0.1 and 1.0 ps (Pron et al., 2013). Under equilibrium period, a free molecular dynamics simulation of 10 ps was simulated with a force constant of 1000 kj mol-1 nm-2 were applied. All covalent bond lengths involving hydrogen atoms were constrained by the LINCS algorithm. Finally, a free molecular dynamics simulation of 20 ns and 50 ns run has been performed for the protein-ligand and protein-protein complexes respectively. The statistical analysis, such as RMSD, RMSF, hydrogen bond analysis and their plots of each complex was analyzed from the output file generated from GROMACS using Origin pro version 9.8 (Essman et al., 1995).
ADME/T prediction

ADMET property is an essential physicochemical properties helpful in the drug development process. The quality of experimental data varies enormously and is often limited in quantity and chemical diversity. As the capacity for biological screening and chemical synthesis have dramatically increased, so we have demands for large quantities of early information on absorption, distribution, metabolism, excretion (ADME)(Chatterjee, Cutler, Doerksen, Khan, & Williamson, 2014). The identified best lead compounds were studied for their absorption, distribution, metabolism and excretion (ADME) properties using Qikprop module, Schrodinger2018-4. The required principle and physiochemical properties of drug compounds are predicted from the defined dataset. It also evaluates the acceptability of the entire known and screened compounds on the basis of Lipinski’s rule of 5(Nagamani, Kesavan, & Muthusamy, 2018). The toxicity of the compounds was predicted using ProTox-II server. The compounds were individually submitted for the prediction of carcinogenicity, mutagenecity, cytotoxicity and hepatotoxicity. The compounds are classified in to six groups based on the toxicity scale(Banerjee, Eckert, Schrey, & Preissner, 2018; Drwal, Banerjee, Dunkel, Wettig, & Preissner, 2014). Each classified group presume the pharmacokinetic nature of the identified lead compound.

Molecular electrostatic potential (MESP) calculation
The top-hit ligand molecules which have been examined as the attractive compound as a common inhibitor for all targets identified through molecular docking studies were taken as input to the density functional theory (DFT) calculations. More specifically, we have taken the TCM_18652, TCM_12889917, NCI_109568, NCI_31491 and Maybridge_22818 was selected as the attractive ligand molecules. All DFT computations were carried out using Jaguar module, Schrodinger. Complete geometry optimization of the compound was carried out using hybrid density functional theory with Becke’s three-parameter exchange potential and the Lee-Yang- Parr correlation functional (B3LYP), using the basis set 6-31G**++ level. To simulate physiological conditions, energy calculations were performed in an aqueous environment using the adaptive Poisson-Boltzmann solver. The quantum chemical descriptors, including MESP, HOMO, LUMO and aqueous solvation energy were computed using Jaguar(Bochevarov et al., 2013). The electrostatic potentials were sampled over the entire accessible surface to provide a measure of the electrostatic potential at roughly the van der Waals surface of the molecules and over the space extending beyond the molecular surface, providing a measure of the charge distribution from the point of view of an approaching reagent. The regions of the positive electrostatic potential indicate the excess positive charged test probe. Negative potential indicate areas of excess negative charge, which is attracted of the positively charged probe(Chinnasamy, Chinnasamy, Nagamani, & Muthusamy, 2014).

The 3D iso surfaces of the MESPs at the van der Waals contact surface represent electrostatic potentials superimposed onto a surface of constant electron density (0.01 e/au3). These color coded iso surface values indicate the overall molecular size and the location of negative or positive electrostatic potentials. The most negative and positive electrostatic potentials are colored deepest red and blue respectively.

Results and Discussion
Virtual Screening of RNF146
The chemical compounds from four different databases (25,20,148) were screened for inhibitors against RNF146 protein using VSW in Glide, Schrodinger2018-4. The compounds for screening were initially filtered out with a Lipinski’s rule of five and other important physico- chemical properties using Qikprop module. The active site was generated based on the co- crystallized ligand that was present in the RNF146 protein crystal structure(PDB ID: 3V3L). The residues occupy the active site are Tyr-107, Tyr-110, Trp-114, Tyr-144, Gln-153, Arg-161, Arg- 163, Lys-175, and Arg-181 (Table 1). These residues were considered for the screening analysis as an important factor for inhibition of RNF146 protein. The 2-D interaction pattern of the docked compounds is displayed in Figure 2.

The best hit compounds listed in the Table 1, are ranked based on their good glide XP score, glide energy, binding free energy and interacting residues. From the results obtained it is shown that Traditional Chinese Medicine (TCM) database screened compounds have a greater binding affinity against RNF146 protein with glide XP score ranges from -14.248 to -13.925 kcal/mol along with this NCI compounds were shown glide XP score ranges from -13.892 to -13.257 kcal/mol. Both the database covers the natural compounds and anti-cancer compounds and it clearly emphasizes the drug selectivity for the target. The interacting residues of RNF146 protein with chemical compounds are Tyr107, Arg110, Trp114, Gln153, Arg163, Tyr144, Lys175, and Arg181. Among these, Tyr107 and Tyr144 are forming pi- pi interaction with the compounds which significantly increase the stability and binding affinity. The glide energy ranges from -61.874 to -44.719 kcal/mol for the top hit compounds, it clearly reveals that the listed compounds were high binding affinity towards the receptor. Among them, 18652_TCM, 18551_TCM, 12889917_TCM and 12889905_TCM were similar in their molecular weight (420.21 g/mol). Although it share common pharmacophore features, the binding conformation and reactive groups differ slightly. Therefore, the binding affinity and Glide score of the hit compounds vary very little.

Binding free energy was calculated for the top hit compounds using MM-GBSA of Prime, Schrodinger. Here, we observed the identified compounds have more binding efficacy and free energy with the RNF146 protein. It is clearly shown that these hit compounds could be the promising RNF146 inhibitors. Since RNF146 is responsible for ubiquitination of AXIN and inhibiting RNF146 can avoid degradation of AXIN, which plays an important role in colon cancer pathogenesis. The top hit compounds are very similar in their chemical structure and composition. The chemical features of the top hit compounds are responsible for their activity.

AXIN active site analysis
The active site of the AXIN protein was predicted using Sitemap, Schrodinger 2018-4. The site 4 has the highest site score 0.891 was taken as active site and the residues involved are Pro-166, Glu-171, Arg-174, Leu-175, Thr-201, Phe-204, Ile-205, Gly-207, Lys-211, Gln-212,
Leu-213, Ile-214, Asp-215, Ala-217, Met-218, Phe-219, Gln-221 and Ala-222 respectively. The colored region displayed in the figure 3, represent the active site of the protein obtained from Sitemap. The active site region consist of yellow mesh region indicates hydrophobic region in the protein that is responsible for the flexible and stable binding of the ligand (Figure 3). The regions apart from hydrophobic represents as a donor (blue mesh) or acceptor (red mesh) hydrogen bonds.

RNF146-AXIN protein docking
RNF146 was docked with AXIN using the HADDOCK web server 2.2. The similar active site residues were denoted as the active residues required for docking and the residues around the active site was considered as passive residues for the protein-protein docking. The best HADDDOCK score is -146.4 +/- 13.6 and was taken as an excellent model. The score details are mentioned above in the (Table 2). The interaction between RNF146 and AXIN was calculated using protein interaction analysis module in Schrodinger. Based on the interaction between the proteins, haddock algorithm provides the best cluster with higher haddock score and it was selected as the best model for further studies. The 3-D view of ligand interaction with the receptor is displayed in Figure 4.

The interactions between two proteins are Trp-114 have a pi-pi stacking with Phe-204 and Asp-215 responsible for the stability of the structure. Asp-215 has salt bridge interaction with Lys-175. Lys-163 has salt bridge interaction with Glu-158. Glu-165 has salt bridge and hydrogen bond interaction with His-159, Arg-161. Asp-220 has both hydrogen bond and salt bridge interaction with Arg-163 (Table 3).
From the earlier studies, it has been reported that lysine residue has significant role in ubiquitination and from the above results it is proved that lysine has been involved and acts a interacting residue when RNF146 binds with AXIN(Zhang et al., 2011). The complex structure was taken for the further studies to identify the inhibitors for the inhibition of RNF146-AXIN complex which is responsible for ubiquitination process.
Ligand docking with RNF146-AXINcomplex

The best RNF146-AXIN complex structure was docked with chemical compounds from different databases using the Ligand docking module, Schrodinger. This study was designed to confirm whether the compounds found as a best hit for RNF146, does have a significant inhibiting role when the RNF146 was complexed with AXIN. To perform this analysis the best hit compounds of four databases were taken from the RNF146 screened results. Therefore, totally forty compounds were proceeded for docking analysis to find the binding mechanism of each compound in the complex structure. Top hit compounds were selected based on the ranking of the Glide XP score (Table 4). But the result has only two database compounds which cover the top 10 places. Therefore, top ten best hits from the screened list of each database was considered as the promising inhibitors. The ligand interaction pattern of TCM-12889905 compound in 2-D is represented in Figure 5.

In RNF146-AXIN complex docking with selected ligands Traditional Chinese Medicine has a greater binding efficiency with glide score ranges from -7.302 kcal/mol to -6.675 kcal/mol and NCI compounds glide score ranges from -6.549 kcal/mol to -6.262 kcal/mol. The glide energy ranges from -35.346 to -25.565 kcal/mol, which confirms the binding affinity of the docked molecules. The interaction residues between the protein complexes with ligands are Arg- 119, Leu-126, Gly-141, Glu-231, and Gln-227. In order to validate the inhibitors it is very essential to analysis the dynamic nature of the molecules. The complex structure with top hit compound is further taken for molecular dynamic simulation studies to check the stability of the complex.

Molecular Dynamics Simulations
The attractive compounds identified from RNF146, RNF146-AXIN complex and RNF146-AXIN complex docking studies were taken for the Molecular Dynamics Simulation using GROMACS for the time period of 20ns and 50ns respectively.

Structural Stability Analysis
The RMSD of backbone atoms was calculated with respect to the initial state of the structure as a function of time in order to evaluate the degree of conformational drift in the protein. RMSD plots of the RNF146 with top hit compounds, RNF146-AXIN complex and RNF146-AXIN complex with ligands are shown in the Figure 6. Each simulation was characterized by an initial rise in RMSD.

The TCM_18652 is an attractive among the top hits, it has relatively very good stabilization within 0.2 nm. Moreover, all the five best hit compounds are not much deviated and stabilized within 0.4 nm. The complex dynamics also observed to have better RMSD value in an average of
0.5 nm. The results strongly represent the screened compounds and docked complex are more reliable and further it can be proceed for experimental studies.

Hydrogen Bond Analysis
In order to understand the binding of the molecules in the active site of RNF146 with ligands, RNF146 binding with AXIN protein and RNF146-AXIN protein complex with ligands, the time dependence of hydrogen bonds between receptor and ligands were monitored which shows that the calculated intermolecular hydrogen bonds are correlated with the biological activities of the lead molecules (Figure 7). All the molecules showed constant H-bond interaction throughout the simulation period and showed that top hit compounds have more stability in RNF146 protein.

RMSF Analysis
Fluctuations in the residues during the dynamics was monitored to understand the residue fluctuation during the motion of the system. The residues Arg119, Gly141, Leu126, and Glu231 of RNF146 protein are the majorly observed interaction in the docking studies. All of these residues were observed carefully in the plot for the fluctuation, it was observed all these interactions were less fluctuated with only less than 0.3 nm (Figure 8). Apart from the selective residues, most of the residues having interaction with ligand molecule are well stable and the regions fluctuated over 0.3 nm are identified as the loop and glycine rich regions of the protein. In the complex protein dynamics plot, both the protein residues was plotted individually and observed the significant fluctuation in the interacting residues. The RMSF plot of complexed protein clearly explains the stability of the proteins and ligands that provides an extensive results to move forward.

Protein-Protein docking of RNF146 and TNKS
RNF146 was docked with TNKS using Commercial tool Schrodinger and non- commercial tool HADDOCK web server 2.2. To dock the RNF146 with TNKS protein, prepared structure of the proteins was imported. The residues involved in the poly(ADP-ribose) (PAR) binding site at the tail of WWE domain of RNF146 protein is Phe-144 to Arg-167 and residues that are responsible for PARsylation in TNKS are Tyr1050, Tyr1060, Tyr1071, His1040, Asp1045, Glu1138, Ser1068, Gly1032, His1031, Ser1033, Phe1035, Gly1043, and His1048 residues. Using these important residues protein-protein docking analysis was done through protein-protein docking module, Schrodinger. In the RNF146-TNKS docking the interactions between the proteins was calculated by protein interaction analysis, Schrodinger (Figure 9).

The interactions between two proteins are His-90 have a pi-pi stacking with Trp-24 responsible for the stability of the structure. Arg-20 has hydrogen bond, salt bridge interaction with Glu- 88 and Glu-187. Gly-95 has hydrogen bond, salt bridge interaction with Glu-68 and Lys-156, Ala-91 has a hydrogen bond with Arg-73 (Table 5). The residues Histidine, Glutamine, Tryptophan is very much important in PARsylation and from the obtained results proves that all three residues has been involved and acts an interacting residue when RNF146 binds with TNKS. Comparing both docking, Schrodinger has the best score and interaction present. The complex structure was taken for the further studies to identify the inhibitors for the inhibition of RNF146-TNKS complex which undergoes an allosteric interaction.

Protein complex-Ligand complex
The best RNF146-TNKS complex structure was docked with chemical compounds from different databases using Protein-Ligand docking module, Schrodinger. This study was designed to check for the inhibitors to inhibit RNF146-TNKS complex. All the four database compounds full fills the top 10 places (Table 6). Hence, top ten best hits from the screened list of each database were taken. Therefore, totally forty compounds were proceeded for docking analysis to find the binding mechanism of each compound in the complex structure. Initially the compounds were screened, filtered out with a Lipinski’s rule of five and other important physico-chemical properties using Qikprop module. Top hit compounds were selected based on the ranking of the Glide XP score. The best hit molecule NCI-109568 is represented as 3-D interaction pose in Figure 10.

In RNF146-TNKS complex docking with selected ligands NCI has greater binding efficiency docking score ranges from -8.854 kcal/mol to -8.560 kcal/mol followed by May bridge docking score from -8.044 kcal/mol. The interaction residues between the protein complexes with ligands are Arg-29, Asp-87, Gly-85, Phe-77, His-90, Ser-75, Tyr-102. The TCM compound TCM_18652 is already found to be listed in top hit in docking studies of RNF146 and its complex with AXIN. This compound has good number of interactions and favorable glide XP score against RNF146-TNKS complex. In order to validate the inhibitors it is very essential to analysis the dynamics of the molecules. The complex structure with top hit compound is further taken for molecular dynamic simulation studies to check the stability of the complex.

Molecular Dynamics Simulations of RNF146-TNKS complex
The top hit compounds from RNF146-TNKS complex and RNF146-TNKS complex with ligands was taken for the Molecular Dynamics Simulation in GROMACS for the time period of 50ns and 10ns respectively. The results were analysed using Origin Pro 8.0.

Structure Stability Analysis
The RMSD of the backbone atoms was calculated with respect to the initial state of the structure as a function of time in order to evaluate the degree of conformational drift in the protein. RMSD plots of the RNF146 with top hit compounds, RNF146-TNKS complex are shown in the Figure 11. Each simulation was characterized by an initial rise in RMSD. The complex dynamics observed to have better RMSD value in an average of 0.4 nm for last 15ns. The results strongly represent the screened compounds and docked complex are more reliable and further it can be proceed for experimental studies.

RMSF Analysis
Fluctuations in the residues during the dynamics were monitored to understand the residue fluctuation during the motion of the system. In the complex protein dynamics plot, both the protein residues was plotted individually and observed the significant fluctuation in the interacting residues (Figure 12). The RMSF plot of complex protein clearly explains the stability of the proteins and provides the extensive results to move forward. The plot describes the RNF146 protein residues (Red color) has considerable fluctuations in key residues which indicate the protein conformational changes induced by TNKS protein. A similar observation was found in TNKS (Green color) has a maximum fluctuation of 0.5 nm for residues ~160-180. Overall, the RNF146-TNKS complex residues have significantly stable and important residue fluctuations are under acceptable range.

SASA Analysis
Surface area was monitored during dynamics by SASA program of Gromacs. The solvent accessible surface area of the protein plays a major role in the function and conformational changes of the protein. The SASA analysis affords a significant results describes the variation in the area and surface of the proteins (Figure 13). Protein-protein complex has considerable changes during the simulation period. The changes observed to retain similar throughout the simulation has offered that the system is well stable and has good binding affinity

Density Functional Theory Analysis
The top hit compounds among the attractive compounds which share specific properties responsible for the good binding affinity was selected for DFT analysis. Among them the TCM compounds have very favorable results in aforesaid analysis and were taken for the electrostatic potential surface analysis to find the electro negativity and positive potential regions. The results show that all the three compounds share specific properties which are highly noticed to have a good binding affinity. The attractive compounds found in docking studies have taken for molecular electrostatic potential surface analysis to find the electronegativity and positive potential regions on the selected compounds. The MESP maps of the two compounds were displayed in Figure 13, indicates that the electropositive potential regions (blue color) near the aryl chain linked with an oxygen atom to the benzene ring. The similar pattern was identified for the other two compounds of NCI. The most electronegative potential regions were found near the carboxyl groups present in the compounds. However, among the selected compounds electro positive potential regions were more occurred in the specific sites which are having more interactions with active site residues. On the other hand, the orbital energies associated with frontier orbitals which provide the useful means of characterizing electron donor and acceptor of the selected molecules. Most of the chemical reactions were governed by these orbitals which actually plays a central role. The energy oof the HOMO is related directly to the ionization potential, whereas the LUMO is related with electron affinity and both are the indicators of possible electrophilic and nucleophilic attack sites on the molecules respectively. The HOMO and LUMO sites are plotted on the molecular surfaces for the compounds were displayed in Figure 13. The HOMO and LUMO energies of the compounds under study are small, ranging from -0.20455 to -0.21333 eV and -0.02589 to -0.02457 eV respectively (Table 7), which clearly indicates the fragile nature of the bound electrons. The small energy gap of HOMO and LUMO energies varies between -0.178 and -0.188 eV. It quietly signifies less stability and makes both rapid electron transfer and exchange equally possible by making these compounds very reactive.

ADME/T Property Analysis
Qikprop module of Schrodinger was used to analyze the drug likeness (Lipinski’s rule of five) and ADME (absorption, distribution, metabolism, excretion) evaluated the identified 10 RNF146 inhibitors that demonstrated physically significant descriptors and pharmaceutically relevant properties. The drug likeness descriptor such as molecular weight, QPlogS, QPlogHERG, QPCaco, QPPMDCK, Human Oral Absorption, QplogPo/w, Rule of five (RoF), QplogKp were analyzed and found that all compounds are in the range of Lipinski;s rule of five. The drug likeness descriptors such as molecular weight, QPPCaco, QPLog HERG, LogP Po/w (octanol/water), LogP MDCK, and percentage of human oral absorption based on Lipinski’s rule of 5 are given in Table 8, which indicate that all the compounds were within an acceptable range of Lipinski’s rule of 5.

Moreover, we analyzed the ADME properties of the attractive and top listed compounds; QPlogPo/w is the partition coefficient ranging from 0.21 to 1.95, which is crucial for the absorption and distribution estimation of drugs within the body, while QPPCaco is a cell permeability factor that determines the metabolism of a drug. The range of QPPMDCK is from 0.21 to 3.41. Human oral absorption of 6 compounds ranged from 21 to 63 %. Overall, the pharmacokinetic properties of identified compounds were of an acceptable range likely having drug-like nature for human use. Therefore, these compounds might act as druglike molecules.

Toxicity Prediction
The toxicity of the selected compound was predicted using ProTox server. The results obtained helps to determine the toxicity in different levels specifically hepatotoxicity, Immunotoxicity, Cytotoxicity, Carcinogenecity and Mutagenecity. All the compounds are classified into Class II to IV category which is acceptable and may require few alteration if adverse events found (Table 9). The TCM_18652, TCM_1288917, NCI_31491 and Maybridge_22818 are found to have a better toxicity profile among the six compounds. Most of the compound has predicted to have mutagenecity profile, except TCM_1288917 and Maybridge_22818. The prediction has the probability ranges from 0.58 to 1.0. Overall, the toxicity results conclude that the selected compounds has the average level of toxic response. Further, in vitro testing are required for the confirmation of bearable toxic records.

CONCLUSION
Molecular interactions of the proteins RNF146, TNKS and Axin regulates the major role in the Wnt signalling mechanism for the degradation of the beta-catenin. The Axin protein plays a critical role in regulating ON/OFF mechanism of the Wnt pathway. This Axin protein is activated by the ubiquitination by RNF146 protein in a specific manner. Inhibition of this step would help to prevent the degradation of axin and inhibits the beta-catenin entry into the nuclear receptor (Figure 1). In order to identify the potent small molecule inhibitors to stop the ubiquitination process of Axin and RNF146 proteins, the structure based virtual screening was carried out. The top hit compounds of TCM and NCI databases shows excellent binding affinity with RNF16. The binding free energy analysis was also analyzed for the above identified compounds. Arg119, Gly141, Leu126, and Glu231 are the key interacting residues observed in the active site which interacting with the above said ligands. The RNF-Axin complex was made through protein-protein docking to study the complex mechanism of the proteins. RNF146 strongly induce both K48- and K63-linked ubiquitylation of TNKS protein, whereas AXIN showed predominantly K48-linked ubiquitin. Further, to support the cross talking of the RNF146, TNKS and AXIN protein, RNF146-TNKS complex were prepared using protein- protein docking complex. Subsequently, the screened inhibitors of four databases against RNF146 were used against the TNKS to the confirm the dual inhibition nature.

Totally, forty compound was docked in the RNF146-Axin interacting sites and the compounds TCM_18652; TCM_18551 shows good binding affinity and interaction. Protein-ligand docking was performed and found that compounds from NCI and Maybridge database showed good binding affinity with the RNF146-TNKS complex protein. The complex structure was proceeded for molecular dynamics simulation studies to analyze the structural stability of the protein and ligand. The results clearly show that both the ligand and protein complex were very well stable and less fluctuated in the throughout simulation process. The compounds TCM-18652, NCI-109658, NCI-12889905 and Maybridge-22818 shows good binding efficiency with RNF146 protein, RNF146-AXIN and RNF146_TNKS complex. These two small molecules are specifically inhibits all the three types of complex, which helps in understanding of the molecular mechanism between the cross talking proteins as well as for deregulation of the interactive complex. The top listed compounds from TCM and NCI database was analyzed for ADME/T properties and observed all the compounds has an acceptable range of properties and can be considered as drug like molecules. As a result of this study, above said two compounds can be further proceeding to the experimental analysis will confirm its potency against RNF146 and Axin proteins.

CONFLICT OF INTEREST
There is no conflict of interest.

ACKNOWLEDGEMENTS
The authors gratefully acknowledge the Department of Biotechnology, Govt. of India for the Major Research Project (Ref. No. BT/PR8491/BID/7/459/2013). The authors are thankful to the MHRD-RUSA 2.0 – F.24/51/2014-U, Policy (TNMulti-Gen), Dept. of Edn. Govt. of India for the support.

References
Banerjee, P., Eckert, A. O., Schrey, A. K., & Preissner, R. (2018). ProTox-II: A webserver for the prediction of toxicity of chemicals. Nucleic Acids Research, 46(W1), W257–W263. https://doi.org/10.1093/nar/gky318
Bochevarov, A. D., Harder, E., Hughes, T. F., Greenwood, J. R., Braden, D. A., Philipp, D. M.,
… Friesner, R. A. (2013). Jaguar: A high-performance quantum chemistry software program with strengths in life and materials sciences. International Journal of Quantum Chemistry, 113(18), 2110–2142. https://doi.org/10.1002/qua.24481
Bray, F., Ferlay, J., Soerjomataram, I., Siegel, R. L., Torre, L. A., & Jemal, A. (2018). Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA: A Cancer Journal for Clinicians, 68(6), 394–424. https://doi.org/10.3322/caac.21492
Cappel, D., Hall, M. L., Lenselink, E. B., Beuming, T., Qi, J., Bradner, J., & Sherman, W. (2016). Relative Binding Free Energy Calculations Applied to Protein Homology Models. Journal of Chemical Information and Modeling, 56(12), 2388–2400. https://doi.org/10.1021/acs.jcim.6b00362
Chatterjee, A., Cutler, S. J., Doerksen, R. J., Khan, I. A., & Williamson, J. S. (2014). Discovery of thienoquinolone derivatives as selective and ATP non-competitive CDK5/p25 inhibitors by structure-based virtual screening. Bioorganic and Medicinal Chemistry, 22(22), 6409– 6421. https://doi.org/10.1016/j.bmc.2014.09.043
Chinnasamy, S., Chinnasamy, S., Nagamani, S., & Muthusamy, K. (2014). Identification of potent inhibitors against snake venom metalloproteinase (SVMP) using molecular docking and molecular dynamics studies. Journal of Biomolecular Structure and Dynamics, 33(September 2014), 1–12. https://doi.org/10.1080/07391102.2014.963146
DaRosa, P. A., Klevit, R. E., & Xu, W. (2018). Structural basis for tankyrase-RNF146 interaction reveals noncanonical tankyrase-binding motifs. Protein Science : A Publication of the Protein Society, 27(6), 1057–1067. https://doi.org/10.1002/pro.3413
Dominguez, C., Boelens, R., & Bonvin, A. M. J. J. (2003). HADDOCK: A protein-protein docking approach based on biochemical or biophysical information. Journal of the American Chemical Society, 125(7), 1731–1737. https://doi.org/10.1021/ja026939x
Drwal, M. N., Banerjee, P., Dunkel, M., Wettig, M. R., & Preissner, R. (2014). ProTox: A web server for the in silico prediction of rodent oral toxicity. Nucleic Acids Research, 42(W1), W53-8. https://doi.org/10.1093/nar/gku401
Friesner, R. A., Murphy, R. B., Repasky, M. P., Frye, L. L., Greenwood, J. R., Halgren, T. A., … Mainz, D. T. (2006). Extra Precision Glide: Docking and Scoring Incorporating a Model of Hydrophobic Enclosure for Protein−Ligand Complexes. Journal of Medicinal Chemistry, 49(21), 6177–6196. https://doi.org/10.1021/jm051256o
Halgren, T. A. (2009). Identifying and characterizing binding sites and assessing druggability.
Journal of Chemical Information and Modeling, 49(2), 377–389. https://doi.org/10.1021/ci800324m
Imran, A., Waseem, W., & Kishwar, S. (2011). Cancer Scenario with Future Perspectives.
Cancer Therapy, 8, 56–70. Retrieved from https://pdfs.semanticscholar.org/2129/3a0c4f983f9e2ba5920acf5edbb96efbb019.pdf
Kang, H. C., Lee, Y.-I., Shin, J.-H., Andrabi, S. A., Chi, Z., Gagné, J.-P., … Dawson, T. M. (2011). Iduna is a poly(ADP-ribose) (PAR)-dependent E3 ubiquitin ligase that regulates DNA damage. Proceedings of the National Academy of Sciences of the United States of America, 108(34), 14103–14108. https://doi.org/10.1073/pnas.1108799108
Kirubakaran, P., Arunkumar, P., Premkumar, K., & Muthusamy, K. (2014). Sighting of tankyrase inhibitors by structure- and ligand-based screening and in vitro approach. Molecular BioSystems, 10(10), 2699–2712. https://doi.org/10.1039/c4mb00309h
Loganathan, L., Gopinath, K., Sankaranarayanan, V. M., Kukreti, R., Rajendran, K., Lee, J.-K., & Muthusamy, K. (2019). Computational and Pharmacogenomics Insights on Hypertension Treatment: Rational Drug Design and Optimization Strategies. Current Drug Targets, 20. https://doi.org/10.2174/1389450120666190808101356
Loganathan, L., & Muthusamy, K. (2018). Investigation of drug interaction potentials and binding modes on direct renin inhibitors. A computational modeling studies. Letters in Drug Design & Discovery, 15. https://doi.org/10.2174/1570180815666180827113622
Loganathan, L., Muthusamy, K., Jayaraj, J. M., Kajamaideen, A., & Balthasar, J. J. (2018). In silico insights on tankyrase protein: A potential target for colorectal cancer. Journal of Biomolecular Structure and Dynamics, 11(0), 1–29.

https://doi.org/10.1080/07391102.2018.1521748

Madhavi Sastry, G., Adzhigirey, M., Day, T., Annabhimoju, R., & Sherman, W. (2013). Protein and ligand preparation: Parameters, protocols, and influence on virtual screening enrichments. Journal of Computer-Aided Molecular Design, 27(3), 221–234. https://doi.org/10.1007/s10822-013-9644-8
Mariotti, L., Pollock, K., & Guettler, S. (2017, December). Regulation of Wnt/β-catenin signalling by tankyrase-dependent poly(ADP-ribosyl)ation and scaffolding. British Journal of Pharmacology. Wiley-Blackwell. https://doi.org/10.1111/bph.14038
Nagamani, S., Kesavan, C., & Muthusamy, K. (2012). E-Pharmacophore mapping and docking studies on Vitamin D receptor (VDR). Bioinformation, 8(15), 705–710. https://doi.org/10.6026/97320630008705
Nagamani, S., Kesavan, C., & Muthusamy, K. (2018). Atom-based and Pharmacophore-based 3D – QSAR Studies on Vitamin D Receptor (VDR). Combinatorial Chemistry & High Throughput Screening, 21(5), 329–343. https://doi.org/10.2174/1386207321666180607101720
Origin(Pro), Version 2019 9.8; OriginLab Corporation, Northampton, MA, USA.
Pérez-Benito, L., Keränen, H., Van Vlijmen, H., & Tresadern, G. (2018). Predicting Binding Free Energies of PDE2 Inhibitors. The Difficulties of Protein Conformation. Scientific Reports, 8(1), 4883. https://doi.org/10.1038/s41598-018-23039-5
Plummer, M., de Martel, C., Vignat, J., Ferlay, J., Bray, F., & Franceschi, S. (2016). Global burden of cancers attributable to infections in 2012: a synthetic analysis. The Lancet Global Health, 4(9), e609–e616. https://doi.org/10.1016/S2214-109X(16)30143-7
Roos, K., Wu, C., Damm, W., Reboul, M., Stevenson, J. M., Lu, C., … Harder, E. D. (2019). OPLS3e: Extending Force Field Coverage for Drug-Like Small Molecules. Journal of Chemical Theory and Computation, 15(3), 1863–1874. https://doi.org/10.1021/acs.jctc.8b01026
Siegel, R. L., Miller, K. D., & Jemal, A. (2019). Cancer statistics, 2019. CA: A Cancer Journal for Tegatrabetan Clinicians, 69(1), 7–34. https://doi.org/10.3322/caac.21551
Steinhart, Z., & Angers, S. (2018). Wnt signaling in development and tissue homeostasis.
Development, 145(11), dev146589. https://doi.org/10.1242/dev.146589
Van Zundert, G. C. P., Rodrigues, J. P. G. L. M., Trellet, M., Schmitz, C., Kastritis, P. L.,
Karaca, E., … Bonvin, A. M. J. J. (2016). The HADDOCK2.2 Web Server: User-Friendly Integrative Modeling of Biomolecular Complexes. Journal of Molecular Biology, 428(4), 720–725. https://doi.org/10.1016/j.jmb.2015.09.014
Zhan, T., Rindtorff, N., & Boutros, M. (2017). Wnt signaling in cancer. Oncogene, 36(11), 1461–1473. https://doi.org/10.1038/onc.2016.304
Zhang, Y., Liu, S., Mickanin, C., Feng, Y., Charlat, O., Michaud, G. A., … Cong, F. (2011).
RNF146 is a poly(ADP-ribose)-directed E3 ligase that regulates axin degradation and Wnt signalling. Nature Cell Biology, 13(5), 623–629. https://doi.org/10.1038/ncb2222
Zhou, Z., Chan, C. H., Xiao, Z., & Tan, E. (2011). Ring finger protein 146/Iduna is a poly(ADP- ribose) polymer binding and PARsylation dependent E3 ubiquitin ligase. Cell Adhesion & Migration, 5(6), 463–471. https://doi.org/10.4161/cam.5.6.18356