Methods Coordinate Preparation
3D coordinates were obtained from the X-ray solved, crystal structures of PARN with RCSB codes: 2A1S and 2A1R. The 2A1S entry is the full length, unbound form of PARN, whereas the 2A1R entry contains the catalytic domain of PARN, in bound form with a 3-mer poly(A) chain. The resolution of both ?X-ray structures is 2.6 A overall. All the important parts of both structures, including the catalytic site and the underlying layer, are very clear in their electron densities. For the purposes of this study, the dimeric form of PARN was used in all calculations.

Energy Minimization
Energy minimizations were used to remove any residual geometrical strain in each molecular system, using the Charmm27 forcefield as it is implemented into the Gromacs suite, version 4.5.5 [45?8]. All Gromacs-related simulations were performed though our previously developed graphical interface [49]. An implicit Generalized Born (GB) solvation was chosen at this stage, in an attempt to speed up the energy minimization process.

Molecular Dynamics Simulations
Molecular systems were subjected to unrestrained Molecular Dynamics simulations (MDs) using the Gromacs suite, version 4.5.5 [45?8]. MDS took place in a SPC water-solvated, periodic environment. Water molecules were added using the ?truncated octahedron box extending 7 A from each atom.Figure 5. DNP-poly(A) is a competitive inhibitor of PARN. (A). Double reciprocal plots 1/v versus 1/[substrate] for PARN activity in the presence or absence ( ) of DNP-poly(A) are shown. The DNP-poly(A) concentrations are 0.1 (&), 0.3 (X) and 0.9 (m) and 1.8 (.) mM. Representative plots of at least three independent experiments. Substrate concentrations range from 0.1?.6 mM poly(A). (B). The slopes (KMapp/Vmax) of the double reciprocal lines are plotted versus the DNP-poly(A) concentration used to calculate the Ki value. The intercept of line on x-axis represents i.Molecular systems were neutralized with counter-ions as required. For the purposes of this study all MDS were performed using the NVT ensemble in a canonical environment, at 300 K, 1 atm and a step size equal to 2 femtoseconds for a total 100 nanoseconds simulation time. An NVT ensemble requires that the Number of atoms, Volume and Temperature remain constant throughout the simulation.

2D Structure Activity Relationships and Statistical Analysis
Structure Activity Relationships (SAR) were calculated based on the coefficients of determination R2 and the Pearson’s contingency coefficients C between the Ki activity and the molecular electronic properties. R2 measures how well a regression line represents the data, whereas C measures the relative strength of association between two variables. The R2 values vary between 1 (strong linear association between the two variables) and 0 (weak linear association). The C values vary between 0 (uncorrelated) and 1 (strong correlation). In order to filter the descriptors with important contribution to the observed biological activity of each inhibitor, descriptors with R2.0.2 and C.0.6 were selected (coloured in red in the full list of descriptors presented in Table S3). Notably, most selected descriptors were quantifying the electronic, steric and hydrophobic properties of the 15 inhibitor compounds. These properties have been previously found to be important characteristics that explain the deadenylase activity of PARN and similar catalytic activities of relative enzymes [16,26]. Data patterns between the different modules were then identified in the filtered, based on the above coefficients, data using hierarchical clustering. Additionally, Principal Components Analysis (PCA) was employed on descriptors with non-zero values. All statistical analysis for the estimation of SAR relationships has been conducted using the R statistical software [63].

Sequence Database Search
A combination of key terms and BLAST searches were employed in order to identify homologous PARN protein sequences. The names and/or accession numbers of the characterized PARNs, including human [9], cattle [17], Xenopus laevis [50] and Arabidopsis thaliana [51] PARN, were used to retrieve their corresponding amino acid sequences from UniProtKB [52]. Subsequently, these sequences were used as probes to search the non-redundant databases UniProtKB [52] and GenBank [53] by applying reciprocal BLASTp and tBLASTn [54]. This process was reiterated until convergence.

Phylogenetic Analysis
The retrieved PARN peptide sequences were searched against the InterPro database [55] to identify the boundaries of the catalytic nuclease domain. In order to optimize the sequence alignment, the predicted core nuclease domain was excised from the full-length protein and was used in our phylogenetic analysis. Subsequently, these trimmed sequences were aligned using CLUSTALW [56]. The resulting multiple sequence alignment was then submitted to ProtTest [57] in order to determine the optimal model for protein evolution. Then, a phylogenetic tree employing a maximum-likelihood method implemented in PhyML [58] was reconstructed using the LG amino acid substitution model [59] with four substitution rate categories; the gamma shape parameter (a) and the proportion of invariable sites were estimated from the data. Bootstrap analysis (500 pseudo-replicates) was performed to test the robustness of the inferred tree. The phylogenetic tree was visualized with Dendroscope [60].

Hierarchical Clustering
Hierarchical clustering with resampling was applied to the filtered data to estimate clusters of compounds based on their correlations structures. The pvclust hierarchical clustering algorithm was employed as implemented within the R package [64]. For each cluster the algorithm calculates p-values via multiscale bootstrap resampling to test the robustness of the inferred clustering and report how strongly the cluster is supported by the data. By default pvclust performs hierarchical clustering K6B times, where K = 10 different data sizes and B = 1,000 denotes the number of bootstrap sample [64]. The algorithm provides two types of p-values, the Approximately Unbiased (AU) which are computed by multiscale bootstrap resampling and the Bootstrap Probability (BP) values which are computed by normal bootstrap resampling. Clusters with AU$95% were selected, which are strongly supported by the data.

Motif Construction
Peptide sequences of the PARN family were aligned and edited by employing Utopia suite’s CINEMA alignment editor [61]. Sequence motifs were excised from this alignment and were submitted to Weblogo [62] in order to generate consensus sequences for these motifs.