Structural Investigations on 2-Amidobenzimidazole Derivatives as New Inhibitors of Protein Kinase CK1 Delta

Protein kinase CK1δ (CK1δ) is a serine-threonine/kinase that modulates different physiological processes, including the cell cycle, DNA repair, and apoptosis. CK1δ overexpression, and the consequent hyperphosphorylation of specific proteins, can lead to sleep disorders, cancer, and neurodegenerative diseases. CK1δ inhibitors showed anticancer properties as well as neuroprotective effects in cellular and animal models of Parkinson’s and Alzheimer’s diseases and amyotrophic lateral sclerosis. To obtain new ATP-competitive CK1δ inhibitors, three sets of benzimidazole-2-amino derivatives were synthesized (1–32), bearing different substituents on the fused benzo ring (R) and diverse pyrazole-containing acyl moieties on the 2-amino group. The best-performing derivatives were those featuring the (1H-pyrazol-3-yl)-acetyl moiety on the benzimidazol-2-amino scaffold (13–32), which showed CK1δ inhibitor activity in the low micromolar range. Among the R substituents, 5-cyano was the most advantageous, leading to a compound endowed with nanomolar potency (23, IC50 = 98.6 nM). Molecular docking and dynamics studies were performed to point out the inhibitor–kinase interactions.


Introduction
Protein kinase CK1δ (CK1δ) is a serine/threonine kinase belonging to the casein kinase 1 family (CK1), and it is expressed in all eukaryotic organisms [1].At least six human isoforms of CK1 (termed α, γ1-3, δ, and ε) have been cloned and characterized as displaying a highly conservative kinase domain [1].CK1δ isoform is distributed both in the cytosol and in the nucleus and can be associated with membranes, receptors, transport vesicles, cytoskeleton components, and centrosomes, depending on the conditions of the cell [1][2][3].CK1δ isoform modulates several physiological processes, such as circadian rhythm [4], DNA damage repair, cellular proliferation, and apoptosis [5].Expression of CK1δ can vary depending on the tissue and cell types and according to physiological and pathological conditions.CK1δ overexpression and the consequent hyperphosphorylation of its substrates can lead to a dysregulation of transduction signaling and cell functions.CK1δ plays a role in oncogenesis due to the control it exerts on Wnt/β-catenin-, p53-, Hedgehog-, and Hippo-mediated pathways, whose dysfunction can influence the development of cancer cells [5][6][7][8].A high expression of CK1δ was found to be a factor associated with a poor prognosis in patients affected by glioblastoma, lung cancer, and colorectal cancer [5,6,9].
Clear evidence points out that CK1 isoforms could affect the development of central nervous system (CNS) diseases.Familial advanced sleep-phase syndrome (FASPS), a disorder of the circadian rhythm [10], has been linked to CK1δ mutations at the level of the binding domain of Period 2 protein.In Alzheimer's disease (AD) patients, a higher expression of CK1δ in the hippocampus was found compared to the control group, and high levels of CK1δ were evidenced in the neuritic plaques [11].CK1δ phosphorylates tau protein in vitro, leading to the destruction of microtubules [12].CK1δ also phosphorylates TDP-43 (transactive response DNA-binding protein, 43 kDa), a protein contributing to the onset of neurodegenerative processes.Aberrant TPD-43 phosphorylation seems to occur in the earliest phases of amyotrophic lateral sclerosis (ALS) and frontotemporal lobe degeneration (FTLD) [13,14].Aggregates of phosphorylated TPD-43 were found in patients suffering from ALS, FTLD, and other neurodegenerative diseases, including Parkinson's disease (PD) and AD [15,16].CK1 showed the ability to phosphorylate α-synuclein in both in vitro and in vivo experiments [17,18].α-Synuclein is recognized as the hallmark of PD, being implicated in the pathogenesis of the illness.Phosphorylated α-synuclein is the major component of Lewy bodies in the brain of sporadic PD patients and is also present in some cases of familial AD and dementia [19].CK1 and α-synuclein aggregates were found to be colocalized in the post-mortem brains of patients suffering from PD and Lewy-body dementia [20].Parkin is another protein whose mutation and dysfunction have been linked to inherited and sporadic PD.Parkin phosphorylation by CK1δ leads to reduced solubility of the protein, thus causing its aggregation and inactivation.Higher levels of phosphorylated parkin were found in the caudate nucleus of sporadic PD patients [21] compared to controls.Altogether, these results highlight the fact that targeting CK1δ with inhibitors can be a therapeutic strategy to fight against cancer and neurodegenerative diseases.
Since the early 2000s, much effort has been carried out to develop CK1δ inhibitors.Different classes of derivatives, either naturally occurring or suitably synthesized, showed the ability to bind the kinase ATP site [22].From a structural point of view, an ATP-competitive CK1δ inhibitor typically consists of a central nitrogen-containing heterocycle decorated with various hydrophobic groups.In Figure 1, the structures of some potent CK1δ inhibitors are reported.Compound I (PF-670462) was one of the first potent CK1δ inhibitors identified (IC 50 = 14 nM), but it was also active against CK1ε (IC 50 = 7.7 nM) [23,24].Derivative I was tested in a triple transgenic mouse model of AD, furnishing a proof of concept for CK1δ/ε inhibition to reduce AD-related cognitive deficits [25].The oxazole derivative II [26] and its analog imidazole compound III showed high potency and improved selectivity versus the CK1ε isoform [27,28].Compound IV (SR 2890) [29] is one of the most representative purine derivatives linked to a substituted benzimidazole-2-yl ring.Benzimidazole is the core nucleus of a deeply investigated series within which several potent inhibitors were identified, such as the series developed by Bischof [30].Compound V, named Bischof-5 (Figure 2), emerged as a promising compound.It showed nanomolar CK1δ inhibitor activity and good selectivity versus CK1ε isoform and also demonstrated the ability to induce apoptosis on different tumor cell lines [30].Bischof-5 analogs with reduced molecular size were developed to improve solubility and intracellular availability [31].Poor solubility and limited cell membrane permeability are the concerns shared by several inhibitors showing nanomolar activity in enzyme assays while possessing poor efficacy in cellular and in vivo tests.

Chemistry
The synthesis of the target compounds is shown in Scheme 1, Scheme 2 and Scheme 3. Derivatives 1-6 (Scheme 1) were obtained starting from the regioselective alkylation of ethyl pyrazole-3-carboxylate 33 with 1-phenyl-and 2-methoxyphenyl-2-bromoethanones, which gave compounds 34 [32] and 35, respectively.The reaction was carried out in the conditions described to obtain 34, i.e., by using acetonitrile as the solvent and in the presence of K 2 CO 3 at r.t., which afforded the ethyl 1-alkyl-3-carboxylate derivatives 34 and 35 as the predominant regioisomers.The 1-substituted-3-carboxylate structure of 35 was confirmed using NOESY experiments, which showed a spatial closeness of the pyrazole hydrogen atom to the methylene protons.
tigations were performed to predict the kinase-ligand interactions and drive the d of the new compounds.

Chemistry
The synthesis of the target compounds is shown in Schemes 1-3.Derivative (Scheme 1) were obtained starting from the regioselective alkylation of ethyl pyrazo carboxylate 33 with 1-phenyl-and 2-methoxyphenyl-2-bromoethanones, which gave pounds 34 [32] and 35, respectively.The reaction was carried out in the condition scribed to obtain 34, i.e., by using acetonitrile as the solvent and in the presence of K at r.t., which afforded the ethyl 1-alkyl-3-carboxylate derivatives 34 and 35 as the pred inant regioisomers.The 1-substituted-3-carboxylate structure of 35 was confirmed u NOESY experiments, which showed a spatial closeness of the pyrazole hydrogen ato the methylene protons.

R R1 R R1
Derivatives 27-29 were synthesized as depicted in Scheme 3. The 5-nitro-and 4-nitrosubstituted benzimidazole-2-amines 43 [34] and 44 [35], respectively, were reacted with di-tert-butyl-dicarbonate to give their respective triply N-Boc protected derivatives, 69 [41] and 70, as a mixture of regioisomers as a result of benzimidazole tautomerism.Derivative 69 was a mixture of 5-nitro-and 6-nitro-substituted isomers, and 70 was a mixture of 4-nitro and 7-nitro isomers.The regioisomers were not separated but used as such in the next step.
Catalytic hydrogenation of compounds 69 and 70 yielded the corresponding amino derivatives 71 [41] and 72, each constituted by only one isomer (see Section 3 for details).Compound 71 was reacted with acetyl chloride or benzoyl chloride to give derivatives 73 and 74, respectively [42].
Compound 72 was transformed into derivative 75 by reaction with benzoyl chloride.The N-Boc derivatives 73-75 were deprotected by treatment with trifluoroacetic acid and then transformed in the hydrochlorides 76, 77 [42], and 78.Finally, these derivatives, after treatment with triethylamine, were coupled with pyrazol-3-yl acetic acid hydrochloride 68 in the experimental conditions reported above for 13-26, to give compounds 27-29.Derivatives 30 and 31, bearing a five-membered heterocyclic substituent at position 5, were prepared as shown in Scheme 4. The synthesis of the 5-(1,2,4-triazol-1-yl) derivative 30 started from the preparation of the suitable 1,2-phenylenediamine 79 [43], which was transformed into the corresponding benzimidazol-2-amino derivative 80 by cyclization with cyanogen bromide.The 5-(tetrazol-5-yl) derivative 81 was prepared as previously described [44].Derivatives 80 and 81 were reacted with pyrazol-3-yl acetic acid hydrochloride 68 in the experimental conditions employed for 13-26, to give, respectively, compounds 30 and 31.The N-methylbenzimidazole derivative 32 was synthesized as reported in Scheme 5, starting from 5-chloro-2-nitro-N-methylaniline 82 [45], which was catalytically reduced to give the corresponding amino derivative 83.Cyclization of the latter with cyanogen bromide afforded the benzimidazol-2-amino derivative 84, which was coupled with pyrazol- The N-methylbenzimidazole derivative 32 was synthesized as reported in Scheme 5, starting from 5-chloro-2-nitro-N-methylaniline 82 [45], which was catalytically reduced to give the corresponding amino derivative 83.Cyclization of the latter with cyanogen bromide afforded the benzimidazol-2-amino derivative 84, which was coupled with pyrazol-3yl acetic acid 68 to yield the desired final compound 32.The N-methylbenzimidazole derivative 32 was synthesized as reported in Scheme 5, starting from 5-chloro-2-nitro-N-methylaniline 82 [45], which was catalytically reduced to give the corresponding amino derivative 83.Cyclization of the latter with cyanogen bromide afforded the benzimidazol-2-amino derivative 84, which was coupled with pyrazol-3-yl acetic acid 68 to yield the desired final compound 32.The structures of the benzimidazole-2-amine derivatives, both intermediates and final compounds, are drawn as single tautomers, i.e., not considering the prototropic exchange of the benzimidazole-NH and pyrazole-NH.Thus, the structures drawn in figures and schemes are those with the substituents numbered as low as possible on both rings.

Structure-Activity Relationship (SAR) Studies
The newly synthesized derivatives 1-32 were tested on truncated CK1δ using the luminescent kinase assay Kinase Glo ® KIT (Promega Italy S.r.l., Milan, Italy), and the results are reported in Tables 1 and 2. Derivatives were tested at a fixed dose of 40 μM; then, those showing a residual kinase activity percentage lower than 50% were tested at 10 μM.For compounds exhibiting a residual enzyme activity percentage lower than 50% at the latter concentration, the IC50 values were determined.
The first set of derivatives (1-6) was designed taking Bischof-5 as the lead compound and replacing its thiazole ring with a pyrazole.Moreover, to confer more flexibility to the The structures of the benzimidazole-2-amine derivatives, both intermediates and final compounds, are drawn as single tautomers, i.e., not considering the prototropic exchange of the benzimidazole-NH and pyrazole-NH.Thus, the structures drawn in figures and schemes are those with the substituents numbered as low as possible on both rings.

Structure-Activity Relationship (SAR) Studies
The newly synthesized derivatives 1-32 were tested on truncated CK1δ using the luminescent kinase assay Kinase Glo ® KIT (Promega Italy S.r.l., Milan, Italy), and the results are reported in Tables 1 and 2. Derivatives were tested at a fixed dose of 40 µM; then, those showing a residual kinase activity percentage lower than 50% were tested at 10 µM.For compounds exhibiting a residual enzyme activity percentage lower than 50% at the latter concentration, the IC 50 values were determined.These results prompted us to substantially modify the structure of the compound, and by applying a molecular simplification approach, compounds 7-12 were designed.In these derivatives, the lipophilic substituent on the pyrazole was removed to evaluate whether structurally simpler molecules than 1-6 retained the ability to bind the kinase.On the fused benzo ring, apart from Cl and Me, other groups with different lipophilic and  tive of this set, tested at 40 μM, was able to inhibit CK1δ with a significant potency (IC50 values > 40 μM).Significantly better results were obtained with the third set of derivatives (13)(14)(15)(16)(17)(18)(19)(20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30)(31)(32), where the pyrazole ring was distanced through a methylene spacer from the carboxamide bond (Table 2).This study has led to the identification of some compounds (2, 3, 6, 13-24, 26, 30, and 31) showing CK1δ inhibitor activity in the micromolar and submicromolar range (IC 50 < 15 µM), the most active being derivative 23 (R = CN), endowed with nanomolar activity (IC 50 = 98.6 nM).
The first set of derivatives (1-6) was designed taking Bischof-5 as the lead compound and replacing its thiazole ring with a pyrazole.Moreover, to confer more flexibility to the molecules, the NHCO function linking the aryl pendant was replaced by a CH 2 CO linker.Simple substituents (R) were introduced at the 5-position of the fused benzo ring.The unsubstituted derivative 1 (R = H) was scarcely active (IC 50 > 40 µM), and the presence of a methyl or a chlorine atom enhanced activity, as it appears from the IC 50 values of 2 (R = Me) and 3 (R = Cl) (IC 50 = 14.6 and 1.80 µM, respectively).The chlorine atom is particularly effective in reinforcing the hydrophobic interaction of the molecule with the kinase, and this was confirmed in the set of compounds 4-6, where the introduction of an ortho-OMe group on the benzoyl pendant was detrimental to the inhibitory activity.Derivatives 4-6 were less active than 1-3.However, also in this set, the 5-chloro-substituted derivative 6 was the best (IC 50 = 13.2 µM), whereas compounds 4 (H) and 5 (R = Me) were scarcely active (IC 50 > 40 µM).
These results prompted us to substantially modify the structure of the compound, and by applying a molecular simplification approach, compounds 7-12 were designed.In these derivatives, the lipophilic substituent on the pyrazole was removed to evaluate whether structurally simpler molecules than 1-6 retained the ability to bind the kinase.On the fused benzo ring, apart from Cl and Me, other groups with different lipophilic and electronic properties were probed (CF 3 , OMe, and NO 2 ).As shown in Table 1, no derivative of this set, tested at 40 µM, was able to inhibit CK1δ with a significant potency (IC 50 values > 40 µM).
This structural modification was suggested by preliminary docking studies (see Section 2.3), which evidenced a better accommodation of these molecules to the kinase ATP site.
The increased flexibility of 13-32 was also envisaged to improve their water solubility, compared to that of their more rigid analogs 7-12.As shown in Table 2, most of the compounds 13-32 inhibited CK1δ activity, with IC 50 values in the low µM range, and one (23, R = CN) displayed a nanomolar potency (IC 50 = 98.6 nM).
The first synthesized derivatives of this set were 13-18, bearing the same R substituents of their inferior homologs 7-12.Interestingly, the new compounds were significantly more active (IC 50 = 0.12-6.27µM) than 7-12 (IC 50 > 40 µM), thus pointing out the key role of the methylene spacer in improving the complementarity between the inhibitor and the kinase.Regarding the effect of the R substituent on 13-18, the groups affording the highest activity are lipophilic, the best being the 5-NO 2 (18, IC 50 = 0.12 µM) and the 5-Cl (15, IC 50 = 0.485 µM).
Instead, the presence of the hydrophilic 5-OMe group reduced the CK1δ inhibitor activity (17, IC 50 = 6.27 µM).Compound 16, bearing the 5-CF 3 substituent (IC 50 = 1.74 µM), present in the lead Bishof-5, was as potent as the 5-Me derivative 14 (IC 50 = 1.64 µM).Due to these encouraging results, further modifications were carried out on the benzimidazole-2acetamide series.Movement of the NO 2 group from the 5-to the 4-position gave derivative 19, which was 10-fold less potent than the regioisomer 18, although maintaining CK1δ inhibitor activity (IC 50 = 1.22 µM).Among the new substituents introduced at the 5-position, the most effective was the CN, which led to the most active compound of the series (23, IC 50 = 98.6 nM).
The 5-terbutyl-substituted derivative 20 proved to be able to inhibit CK1δ (IC 50 = 1.00 µM) with higher potency, compared to the unsubstituted compound 13 (IC 50 = 4.21 µM).This result would confirm that hydrophobic groups reinforce the interaction with the kinase and indicate that the pocket where the benzo ring is accommodated is roomy enough to allow for the insertion of the terbutyl residue.Coherent to this is the inhibitory effect of the 5,6-dichloro derivative 22 (IC 50 = 0.98 µM).In contrast, the 5-OCF 3 lipophilic group gave a compound (21, IC 50 = 5.99 µM) whose CK1δ inhibitory activity was comparable to that of the 5-OMe derivative 17.
Compounds 24-29 feature amide substituents of different sizes and properties on the fused benzo ring.The 5-CONH 2 substituent proved to be the best among those evaluated, as it afforded a quite good CK1δ inhibitor (24, IC 50 = 2.53 µM).The acidic 5-SO 2 NH 2 group was chosen for its potential ability to be engaged in an ionic interaction with Lys 38 residue of the selectivity pocket of the kinase.Instead, this group was disadvantageous for CK1δ inhibitor activity.In fact, derivative 25 showed an IC 50 value > 40 µM.Methylation of the sulphonamide group, to give derivative 26 (R = SO 2 NHMe), made a certain activity appear (IC 50 = 10.7 µM).A detrimental effect was achieved by the 5-NHCOMe ( 27), 5-NHCOPh (28), and 4-NHCOPh (29) residues, which dropped the ability to inhibit the kinase (IC 50 > 40 µM).At the 5-position, derivatives 30 and 31 bear, respectively, the 1,2,4triazol-1-yl and the tetrazol-5-yl residues, both able to form H-bonding interaction, and the latter possessing acidic properties.Both compounds, 30 and 31, were active as CK1δ inhibitors (30, IC 50 = 2.59 µM; 31, IC 50 = 1.54 µM).The rationale for the synthesis of the N-methyl derivative 32 was to evaluate whether the blockade of the tautomerism could reinforce the interaction with the kinase.Instead, the N-methyl-6-chloro derivative 32 showed a significantly reduced activity against CK1δ (IC 50 = 20.1 µM) compared to its NH analog 15.

Physicochemical and Pharmacokinetic Parameters
Some physicochemical and pharmacokinetic parameters, such as predicted logP, logS, logBB, blockage of HERG K+ channels, Caco-2 cell permeability, binding to human serum albumin, and oral absorption, were calculated for the most active compounds, 15, 18, 22, 23, and the reference PF-670462 (see Section 3.3.10).The results are reported in the Supplementary Materials (Table S1).On the whole, the calculated parameters suggested good pharmacokinetics properties of the new inhibitors, such as good to high oral absorption and good ability to cross the blood-brain barrier.

Molecular Modeling Studies
In the molecular modeling studies conducted on the synthesized compounds, a computational structure-based approach was carried out, starting from the available crystallographic complexes in the Protein Data Bank [46].A CK1δ 3D structure was selected on the basis of the closest scaffold similarity between the co-crystallized inhibitor and the series under investigation.The selection was made in favor of the CK1δ structure coded as 5OKT, in complex with a 2-aminobenzothiazole-like inhibitor.The docking reliability was validated by redocking the PDB crystallographic ligand (IWP-2) to the protein ATP-binding site.The docking algorithm was able to recreate the crystallo-graphic pose (Figure S1) with a low heavy-atom RMSD value of 2.1 Å (the same range of the crystallographic resolution), reaching the value of 0.86 Å when focusing on the 2-aminobenzothiazole scaffold.
The synthesized compounds were docked at the CK1δ binding site.Two docking runs were made, either neglecting or considering two key crystallographic water molecules in the ATP-binding site.In fact, it is well established that water molecules could play an active role in mediating hydrogen-bonding interactions with the ATP-competitive inhibitors within the binding site [47].The first run of docking simulations has been performed excluding all water molecules during the ligand-posing step.In parallel, a second docking run has been carried out considering two water molecules (labeled w295 and w296 in the following analysis) in direct contact with the amino acid residues Tyr56, Glu52, and Lys38.As already anticipated, these two water molecules are highly conserved inside the ATP-binding cleft, and they directly contact the inhibitor in several crystallographic structures [47].
The selection of the representative binding poses has been carried out through a combined approach based on the analysis of electrostatic and van der Waals interactions between the ligand and its recognition site, and subsequently by filtering the poses through a pharmacophore model built on the strongest interactions of the co-crystallized inhibitor with CK1δ in structure 5OKT (see Section 3.3 for details).In particular, the pharmacophore filtering step retained only the docking poses capable of preserving a bi-dentate hydrogen bond with Leu85 residue, as well as the positions of the aromatic groups of the benzothiazole ring.An overview of the selected poses is provided in the video included in the Supplementary Materials (Video S1).In the video's background, heat maps of electrostatic and hydrophobic energies are depicted for each ligand in the series concerning the residues primarily involved in the binding site.
The general binding mode of the compounds in this series highlights the presence of the aromatic benzimidazole scaffold, which is capable of forming a bidentate hydrogen bond with Leu85 and of interacting with the numerous hydrophobic residues of the binding pocket, such as Ile15, Ile23, Ala36, Leu135, and Ile147.Regarding substitutions at positions 4, 5, and 6 of the benzimidazole ring, the presence of substituents with a strong electrostatic component is emphasized, enabling the mediation of hydrogen-bonding interactions with residues like Tyr56 and Lys38.The mentioned residues are common in other previously published computational works [48][49][50].As anticipated, in addition to these residues, at least two water molecules can play a crucial role in mediating additional hydrogen bonds.
To inspect the plausible structure-activity relationship of the series of inhibitors, we have focalized the attention to compounds 7, 13, and 23, combining the binding mode analysis (Figure 3) with the comparison of their Interaction Energy Fingerprint (IEF) profiles (Figure 4).Briefly, compound 7 is a putative representative of compounds 1-12, and compounds 13 and 23 can help to decipher the inhibitory behaviors of compounds 13-32.
To inspect the plausible structure-activity relationship of the series of inhibitors, we have focalized the attention to compounds 7, 13, and 23, combining the binding mode analysis (Figure 3) with the comparison of their Interaction Energy Fingerprint (IEF) profiles (Figure 4).Briefly, compound 7 is a putative representative of compounds 1-12, and compounds 13 and 23 can help to decipher the inhibitory behaviors of compounds 13-32.The binding pocket of CK1δ is mainly composed of hydrophobic residues, with few residues capable of mediating hydrogen bonds or ionic interactions, such as Tyr56, Lys38, and Glu52.As anticipated, in addition to these residues, there are at least two water molecules that can play a crucial role in mediating additional hydrogen bonds.The high concentration of internal hydrophobic residues helps explain the enhanced activity values observed for compounds with lipophilic substitutions at position 5, as mentioned in the previous section.This holds true for both the first dataset of compounds (1-12) and the second cluster (13)(14)(15)(16)(17)(18)(19)(20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30)(31)(32).As anticipated, compound 7 (IC50 > 40 μM) is representative of compounds 1-12 belonging to the pyrazole-amide class of analogs, characterized by a planar conformation (Figure 3A).The planarity of these structures leads to notable steric clashes inside the binding cavity, reducing the pose adaptability within the pocket.
The introduction of a methylene spacer between the pyrazole moiety and the amide group, as in compounds 13-32, improved the corresponding inhibitory activities of the series; and this is clearly exemplified by the activity of compound 13 (IC50 = 4.21 μM), which is analogous to compound 7, with the unique addition of a methylene spacer at the side chain.This elongation increases the flexibility of the side chain and, consequently, might enhance the adaptability of the inhibitor inside the binding cleft.
In addition to this structural consideration, the role played by the substituent in position 5 of the benzimidazole system appears evident from the structure-activity analysis.Compound 23, chosen as a reference structure, exhibits the highest inhibitory activity in the entire series (IC50 = 0.0986 μM) and shares chemical similarities with compounds 18, 19, and 24.All these derivatives have a substituent that is able to establish hydrogen bonding with the water molecules inside the binding site, improving the inhibitory activities of the compounds.
However, the correct balance between the substituent's ability to interact with the water molecules and its volume (steric hindrance properties) contributes to determining the inhibitory capacity of these analogs, as demonstrated by the activities measured for compounds 25-29.Interestingly, compounds 30 and 31, with triazole and tetrazole, were designed to establish a direct interaction with Lys38.This strategy led to measurable activity in the low micromolar range but not to the expected extent.The displacement of water molecules to achieve a direct interaction, especially water 295, seems to be disfavored compared to mediating a hydrogen bond interaction, emphasizing the importance of water molecules within the binding site.The binding pocket of CK1δ is mainly composed of hydrophobic residues, with few residues capable of mediating hydrogen bonds or ionic interactions, such as Tyr56, Lys38, and Glu52.As anticipated, in addition to these residues, there are at least two water molecules that can play a crucial role in mediating additional hydrogen bonds.
As anticipated, compound 7 (IC 50 > 40 µM) is representative of compounds 1-12 belonging to the pyrazole-amide class of analogs, characterized by a planar conformation (Figure 3A).The planarity of these structures leads to notable steric clashes inside the binding cavity, reducing the pose adaptability within the pocket.
The introduction of a methylene spacer between the pyrazole moiety and the amide group, as in compounds 13-32, improved the corresponding inhibitory activities of the series; and this is clearly exemplified by the activity of compound 13 (IC 50 = 4.21 µM), which is analogous to compound 7, with the unique addition of a methylene spacer at the side chain.This elongation increases the flexibility of the side chain and, consequently, might enhance the adaptability of the inhibitor inside the binding cleft.
In addition to this structural consideration, the role played by the substituent in position 5 of the benzimidazole system appears evident from the structure-activity analysis.Compound 23, chosen as a reference structure, exhibits the highest inhibitory activity in the entire series (IC 50 = 0.0986 µM) and shares chemical similarities with compounds 18, 19, and 24.All these derivatives have a substituent that is able to establish hydrogen bonding with the water molecules inside the binding site, improving the inhibitory activities of the compounds.
However, the correct balance between the substituent's ability to interact with the water molecules and its volume (steric hindrance properties) contributes to determining the inhibitory capacity of these analogs, as demonstrated by the activities measured for compounds 25-29.Interestingly, compounds 30 and 31, with triazole and tetrazole, were designed to establish a direct interaction with Lys38.This strategy led to measurable activity in the low micromolar range but not to the expected extent.The displacement of water molecules to achieve a direct interaction, especially water 295, seems to be disfavored compared to mediating a hydrogen bond interaction, emphasizing the importance of water molecules within the binding site.
To better interpret the binding process of this class of inhibitors, molecular dynamics simulations were performed.In particular, SuMD (Supervised Molecular Dynamics) and TTMD (Thermal Titration Molecular Dynamics) were used for studying the binding and unbinding events.
SuMD is an enhanced sampling Molecular Dynamics-based method designed for studying ligand-target association processes at the nanosecond timescale.It integrates classical molecular dynamic simulations with a tabu-like algorithm to increase the probability of observing a binding event [51].On the other hand, TTMD is a computational approach that conducts classical molecular dynamics simulations at increasing temperatures.It uses temperature to boost the system's kinetic energy and promote events that would otherwise require much longer timescales [52].These techniques provide valuable insights into the binding and unbinding kinetics of compounds with the studied scaffold.Compounds 13 and 23 have been selected as representative examples of this class of inhibitors, in particular to interpret at a molecular level the role of the nitrile group in position 5 and the cooperative role of water molecules in stabilizing the positioning of the ligand within the binding pocket.
These analyses typically include measurements of RMSD (Root Mean Square Deviation) compared with the reference structure, which, in this case, is the docking pose obtained earlier.Additionally, per-residue interaction analyses provide insights into the interactions between the ligand and various amino acid residues throughout the simulation.These analyses, in Figure 5 and Video S2, help to determine how the ligand interacts with the protein structure during the binding process.SuMD is an enhanced sampling Molecular Dynamics-based method designed studying ligand-target association processes at the nanosecond timescale.It integra classical molecular dynamic simulations with a tabu-like algorithm to increase the pro bility of observing a binding event [51].On the other hand, TTMD is a computational proach that conducts classical molecular dynamics simulations at increasing tempe tures.It uses temperature to boost the system's kinetic energy and promote events t would otherwise require much longer timescales [52].These techniques provide valua insights into the binding and unbinding kinetics of compounds with the studied scaffo Compounds 13 and 23 have been selected as representative examples of this class of hibitors, in particular to interpret at a molecular level the role of the nitrile group in po tion 5 and the cooperative role of water molecules in stabilizing the positioning of ligand within the binding pocket.
These analyses typically include measurements of RMSD (Root Mean Square Dev tion) compared with the reference structure, which, in this case, is the docking pose tained earlier.Additionally, per-residue interaction analyses provide insights into the teractions between the ligand and various amino acid residues throughout the simulati These analyses, in Figure 5 and Video S2, help to determine how the ligand interacts w the protein structure during the binding process.The inspection of the ligand-binding trajectory clearly shows a preliminary bind event in proximity to the external portion of the P-loop of the protein, in which a stabi ing interaction with Ile15 takes place.When the ligand enters the binding site, both The inspection of the ligand-binding trajectory clearly shows a preliminary binding event in proximity to the external portion of the P-loop of the protein, in which a stabilizing interaction with Ile15 takes place.When the ligand enters the binding site, both the ligand pose and its interaction profile are perfectly in agreement with those observed in the docking analyses.In particular, a strong stabilizing interaction with Leu85 of the hinge region is observed, together with a hydrogen bond interaction between a water molecule and the nitrogen nitrile group.It is worth underlining that this portion of the binding cleft is occupied by a water molecule during the whole simulation time, and as in Figure 6, it could represent the water w295 already described in the docking simulations.In parallel with the exploration of the binding trajectories of compounds 13 (IC50 = 4.21 ± 1.92 μM) and 23 (IC50 = 0.0986 ± 0.0394 μM), a molecular dynamics study of the unbinding event was also conducted using TTMD (Thermal Titration Molecular Dynamics) to evaluate the stability of a ligand-protein complex [53].For both compounds, five simulations were conducted, starting at 300 K and gradually increasing the temperature until the binding mode was lost, resulting in the ligand exiting the binding pocket.The starting point for the dynamic simulations was the docking poses obtained previously.The analyses of the trajectories are shown in Figure 7 and Video S3.The "titration profile" plot reports the average IFPCS value for each TTMD step as a function of the step temperature.In parallel with the exploration of the binding trajectories of compounds 13 (IC 50 = 4.21 ± 1.92 µM) and 23 (IC 50 = 0.0986 ± 0.0394 µM), a molecular dynamics study of the unbinding event was also conducted using TTMD (Thermal Titration Molecular Dynamics) to evaluate the stability of a ligand-protein complex [53].For both compounds, five simulations were conducted, starting at 300 K and gradually increasing the temperature until the binding mode was lost, resulting in the ligand exiting the binding pocket.The starting point for the dynamic simulations was the docking poses obtained previously.The analyses of the trajectories are shown in Figure 7 and Video S3.The "titration profile" plot reports the average IFP CS value for each TTMD step as a function of the step temperature.
This metric makes it easier to qualitatively compare the stabilities of various proteinligand complexes.A straight line joining the initial and final states of the simulation is also drawn in the graph, and its slope (termed MS coefficient), quantifying the stability of the protein-ligand complex, is reported in the legend.The second graph depicts the "Interaction Fingerprint Similarity", which provides information on how the interaction pattern changes over time as compared to the initial reference frame (a colorimetric scale reports temperature).The third and final plots report the time-dependent evolution of the ligand and protein backbone RMSD value [53].
For compound 23, it is evident that in the reported simulation, similarly to the other replicates, the ligand pose remains stable throughout the temperature gradient explored.The interaction similarity values compared to the initial reference, starting at a value below −0.8, also remain constant throughout the temperature gradient explored until the ligand exits the binding site.This steady trend is also reflected in the ligand's RMSD values, which remain stable throughout the simulation.Interestingly, compound 23 loses its binding state only at temperatures above 560 K, when a significant destructuring of the tertiary structure of the protein begins to be noticed.On the contrary, compound 13 shows a loss of the interaction profile at a lower temperature than compound 23, indicating that this ligand may have a shorter residence time due to the lack of the nitrile substituent.
4.21 ± 1.92 μM) and 23 (IC50 = 0.0986 ± 0.0394 μM), a molecular dynamics study of the unbinding event was also conducted using TTMD (Thermal Titration Molecular Dynamics) to evaluate the stability of a ligand-protein complex [53].For both compounds, five simulations were conducted, starting at 300 K and gradually increasing the temperature until the binding mode was lost, resulting in the ligand exiting the binding pocket.The starting point for the dynamic simulations was the docking poses obtained previously.The analyses of the trajectories are shown in Figure 7 and Video S3.The "titration profile" plot reports the average IFPCS value for each TTMD step as a function of the step temperature.These findings suggest that the presence of specific hydrogen bond-accepting substituents at position 5, like the cyano group, plays a crucial role in enhancing the binding affinity and inhibitory activity of the ligand by mediating interactions with internal water molecules within the binding site.This highlights the significance of considering both the protein-ligand interactions and solvent-mediated interactions for a comprehensive understanding of ligand binding and activity.

Chemistry General Methods
All the commercially available reagents and solvents were used as purchased from Merck without further purification.Analytical silica gel plates 60 F 254 (Merck Life Science S.r.l., Milan, Italy) and silica gel 60 (Merck, 70-230 mesh) were used for analytical TLC and column chromatography, respectively.All melting points were determined on a Gallenkamp melting point apparatus and were uncorrected.Compounds were named following IUPAC rules as applied by ChemDrawUtra 9.0.Elemental analyses were performed on unknown compounds with a Flash E1112 Thermofinnigan elemental analyzer for C, H, and N, and the results were within ±0.4% of the theoretical values.All final compounds revealed purity not less than 95%.Catalytic reductions were performed in a Parr hydrogenation apparatus (Parr Instrument Company, Moline, IL, USA).NMR spectra were recorded on a Bruker Avance 400 spectrometer (400 MHz, Bruker Italia S.r.l., Milan, Italy).The chemical shifts are reported in δ (ppm) and are relative to the central peak of the solvent, which was CDCl 3 or DMSOd6.The following abbreviations are used: s = singlet, d = doublet, t = triplet, q = quartet, m = multiplet, br = broad, and ar = aromatic protons.Scanned 1 H-NMR spectra of compounds 1-32 are reported in the Supplementary Materials.
After cooling at 0 • C, the reaction mixture was diluted with H 2 O (about 25 mL), and the suspended solid was collected by filtration, washed with H 2 O (2-3 mL), dried, and purified by recrystallization.Derivative 16 was purified by column chromatography.In the case of derivatives 14 (R = Me) and 26 (R = SO 2 NHMe), dilution of the reaction mixture with H 2 O did not afford any solid, thus the compounds were isolated by extraction with CH 2 Cl 2 (30 mL × 3).The organic phase was washed with brine (30 mL), anhydrified (Na 2 SO 4 ), and evaporated at reduced pressure to give a solid that was treated with Et 2 O (5-10 mL), collected by filtration, and purified by recrystallization (14) or by column chromatography (26).To isolate derivative 15, the crude compound (460 mg) was treated with boiling EtOH (50 mL), and the solid obtained after cooling at r.t. was collected by filtration and recrystallized.
A total of 1M NaOH solution (95 mL) was added dropwise at 0 • C to a cooled (T = 0 • C) solution of the ester 34 or 35 (19 mmol) in EtOH (200 mL).The mixture was heated at 50 • C for 2 h and then concentrated at about one-third of its volume by evaporation of the solvent at reduced pressure.Then, the solution was diluted with about the same volume of H 2 O and acidified to pH = 2 with 6N HCl.The yellow solid that precipitated was collected by filtration, washed with water, and dried.The compounds were pure enough to be used for the next step without further purification.5-( 6)Terbutyl-1H-benzimidazol-2-amine (45).BrCN (3.37 mmol) was added to a solution of commercially available 4-tertbutyl-1,2-phenylendiamine 58 (3.04 mmol) in a 1:1 mixture of H 2 O and MeOH (50 mL), and the mixture was heated at 60 • C for 5 h.About half of the solvent was eliminated under vacuum, and the mixture was extracted with EtOAc (50 mL × 5).Evaporation at reduced pressure gave a red oil, which was treated with H 2 O (10 mL).The mixture was basified to pH = 8-9 with a NaHCO 3 -saturated solution to obtain a solid, which was collected by filtration and used for the next step without further purification.Yield 45%.mp = 201-203 • C. 1 H NMR (DMSOd 6 ) 10.50 (br s, 1H, NH), 7.11 (s, 1 H, ar), 7.00 (d, 1H, ar, J = 7.6 Hz); 6.91 (s, 1 H, ar, J = 7.6 Hz), 6.01 (s, 2H, NH 2 ), 1.29 (s, 9 H, 3CH 3 ).Elemental analysis calcd for C 11 H 15 N 3 : C, 69.81; H, 4.65; N, 10.76; found C, 62.23; H, 4.37; N, 10.51.
The title compounds were synthesized by reacting the suitable 3,4-diaminobenzensulphonamide 63 [39] or 64 (5.34 mmol) with BrCN (5.34 mmol) in a 1:1 mixture of H 2 O and MeOH (60 mL) heated at 50 • C for 5 h.After evaporation of about half the solvent, the mixture was basified to pH = 8 with NaHCO 3-saturated solution, and the suspension was extracted with EtOAc (30 mL × 3).The organic phase was anhydrified (Na 2 SO 4 ), and the solvent evaporated at reduced pressure to give a solid (50), which was recrystallized from MeCN, or an oily residue (51), which was used as such for the next step.
The new compound 70 was prepared in the same conditions described for 69 [41].Briefly, di-tert-butyl dicarbonate (16.8 mmol) was added to a suspension of the suitable 2-aminobenzimidazole derivative 43 [34] or 44 [35] (5.6 mmol) in anhydrous THF (55 mL).The mixture was stirred at r.t. for 5 h, then small amounts of impurities were filtered off, and a second crop of di-tert-butyl dicarbonate (16.8 mmol) was added to the filtrate, followed by 4-dimethylaminopyridine (0.6 mmol).The mixture was stirred at r.t.overnight, then the solvent was evaporated at reduced pressure to give a solid (70), which was collected and used for the next step, or a brownish oil (69), which solidified on air and was purified by silica gel column chromatography (eluted with a gradient of n-hexane/EtOAc from 10:1 to 5:2).Compounds 69 and 70 were a mixture of two regioisomers (69, 5-and 6-nitro substituted; 70, 4-and 7-nitro substituted).
The newly synthesized compound 72 was prepared as previously described for 71 [41].Briefly, compounds 69 or 70 (4.1 mmol), both used as a mixture of the two isomers, were suspended with 10% Pd/C (0.2 g) in MeOH (100 mL), and the mixture was hydrogenated in a Parr apparatus at 35 psi for 8 h.The catalyst was filtered off, and the solvent was evaporated to give a residue that was taken up with a few Et 2 O (2-3 mL) and collected by filtration.Each amino derivative, 71 and 72, was constituted by one compound, since the less abundant isomer was lost in the mother liquor.Both derivatives, 71 and 72, were pure enough to be used for the next step without further purification.We did not determine whether the compounds were 5-or 6-substituted, because their subsequent deprotection generated one isomeric product.
The suitable acyl chloride (5.3 mmol) and triethylamine (5.3 mmol) were added in sequence to a solution of the amino derivative 71 or 72 (4.5 mmol) in anhydrous CH 2 Cl 2 (20 mL).The mixture was stirred at r.t.overnight and then quenched with H 2 O (20 mL).The organic phase was washed with H 2 O (15 mL × 3) and then anhydrified (Na 2 SO 4 ).Evaporation of the solvent at reduced pressure gave a solid, which was collected, dried, and used for the next step without further purification.
A solution of the suitable BoC-derivatives 73-75 (4.0 mmol) in trifluoroacetic acid (2 mL) and CH 2 Cl 2 (2 mL) was stirred at r.t.overnight.Evaporation of the solvents at reduced pressure gave an oily residue, which was taken up with a saturated methanol solution of hydrochloric acid (2 mL).The solution was stirred at r.t. for 30 min, then the solvent was eliminated at reduced pressure, yielding the desired benzimidazol-2amine hydrochlorides.

CK1δ Activity Assays
The inhibitory activity of synthesized compounds on CK1δ has been carried out following a procedure for a luminescence-based assay (KinaseGlo ® kit, Promega Italy Srl, Milan, Italy) already reported in the literature [50,54].Luminescent assays were performed in white 96-well plates in a 40 µL final volume.The buffer used is a milliq water solution of 50 mM HEPES (pH 7.5), 1 mM EDTA, 1 mM EGTA, and 15 mM MgCl 2 .In a typical assay, 10 µL of inhibitor solution (desired concentrations were obtained by diluting in buffer a starting 10 mM solution in DMSO) and 10 µL of enzyme solution (truncated CK1δ, recombinant human, GST-tagged, aa 1-294, Merck Millipore, Merck KGaA, Darmstadt, Germany) were added to the well, followed by 20 µL of assay buffer containing casein substrate (casein solution from bovine milk, 5% in water, Sigma-Aldrich, Merck KGaA, Darmstadt, Germany) and ATP, at a final concentration of 0.05% and 2 µM, respectively.The final enzyme concentration was 6.5 nM.Compound PF-670462 (IC 50 = 4 nM) [24] at a concentration of 40 µM was used as a positive control for CK1δ, and a DMSO/buffer solution was used as a negative control.The final DMSO concentration in the reaction mixture did not exceed 1%.After 60 min of incubation at 30 • C, the enzymatic reactions were stopped with the addition of 40 µL of KinaseGlo reagent (Promega Italy Srl, Milan, Italy).The luminescence signal (relative light unit, RLU) was recorded after 10 min at 25 • C using a Tecan Infinite M100.First, residual enzyme activity percentage was determined in duplicate at 40 µM for each inhibitor with respect to DMSO/buffer only, and IC 50 values were determined for those compounds showing enzyme residual activity lower than 50% at that concentration.IC 50 values were evaluated using ten different inhibitor concentrations ranging from 100 to 0.026 µM.Data were analyzed using Excel and GraphPad Prism software (version 8.0).

Hardware and Software Overview
Most of the molecular modeling operations were carried out using the MOE suite [55].Docking simulations were performed using PLANTS (Protein Ligand ANT System) [56].The computation of MMFF94 partial charges was conducted through MOPAC [57], which is integrated within the MOE suite.
Analyses of docking poses, including energy calculations, visual inspections, and pharmacophoric filtering, were conducted using the tools available in the MOE suite.
The molecular modeling studies were performed on an 8 CPU Linux workstation (Intel ® Xeon ® CPU E5-1620 3.70 GHz).
Molecular dynamics simulations were carried out on an Linux GPU cluster composed of 20 NVIDIA devices ranging from a GTX2080Ti to an RTX3090.

Structure Preparation
The three-dimensional structure of the selected complex for analysis was retrieved from the Protein Data Bank [46] (PDB ID: 5OKT) and subsequently prepared using the tools provided by the MOE 2022.2 suite [55].Upon importing the structure into the MOE environment, the "Structure Preparation" tool was utilized.This step is essential for addressing various issues within the chosen protein chain, which encompass the addition of missing atoms, the capping of N/C-termini, and the management of tautomeric forms and multiple conformations of specific residues.Furthermore, any gaps in the protein structure were reconstructed using a homology-based approach.The determination of protonation states was conducted using MOE's "Protonate3D" tool, which added hydrogen atoms and assigned ionization states to each residue based on a fixed pH of 7.4.The most probable tautomer was calculated based on hydrogen-bonding contributions that significantly contributed to the structure's stability.Fractional occupancies of residues were selected based on the highest probability.Finally, partial charges were computed using the AMBER14 force field, and the hydrogen atoms added during the reconstruction phases underwent an energy minimization process using the same force field.

Molecular Docking
The compound database underwent a comprehensive preparation process for computational analysis, utilizing various tools from the OpenEye suite [58].Firstly, the most probable tautomeric state was assigned to each ligand within the database, by the usage of Tautomers [59].Following this, three-dimensional coordinates of the various molecules were generated with Omega [60], with a specific focus on minimizing their energy.The protonation states of ionizable groups were then determined with Fixpka [61] to ensure accuracy.Lastly, to facilitate subsequent analysis, partial charges were assigned to the atoms within the ligand database thanks to Molcharge [62].
Molecular docking was performed using PLANTS, with the empirical CHEM-PLP scoring function.The conformational search radius was set at 10 Å relative to the center of mass of the crystallographic ligand (PDB ID: 5OKT).
Following this, atomic partial charges were computed for both ligand poses, using the MMFF94 method, and receptors, using the Amber14EHT force field.The contributions of electrostatic and van der Waals forces to the binding energy were calculated using MOE.

Pharmacophoric Filter
The pharmacophore is a set of crucial features that ensure optimal intermolecular interactions between a ligand and its reference protein structure [63].Here, a pharmacophoric model was built with MOE using the crystallographic complex 5OKT.This resulted in the features describing the aromatic scaffold of the ligand and a bi-dentate hydrogen bond between the ligand and Leu85.In particular, donor and acceptor features were positioned on the ligand's 2-aminobenzothiazole NH and N; and acceptor and donor features were positioned on Leu85 carbonyl and amino groups.In all cases, the dimension of the feature sphere was 2 Å.

Interaction Energy Fingerprints (IEF)
Electrostatic interactions were assessed based on the non-bonded energy term of the force field, utilizing the Coulombic function, and they are expressed in kcal/mol.On the other hand, the hydrophobic contribution was calculated using a hydrophobic interaction term associated with a dimensionless score.This calculation is based on the analysis of contact surfaces carried out using MOE 2022.02software, employing a specific scoring function.
To visually represent the results on a "per residue" basis, Heat Maps were utilized and generated with GNUPLOT.These graphical representations employ specific color scales to visualize the calculated interactions between each amino acid of the molecular target (x-axis) and each ligand (y-axis).For IEele (Electrostatic Interaction Energy), colors ranging from blue to red represent energy values from negative to positive.As for IEhyd (Hydrophobic Interaction Energy), colors from white to dark green depict scores ranging from 0 to positive values.

Docking Video Maker
An internal tool called MMsDocking video maker has been exploited to automatically create a video that displays the most pertinent docking data, including docking poses, per-residue IEhyd and IEele data, and experimental binding data, in order to make the visualization and analysis of the data from the docking simulations easier.After obtaining images using the following method, videos were mounted using MEncoder [64].The backdrop heat maps were created using GNUPLOT 5.2.8 [65], based on per-residue IEhyd and IEele data that were computed using MOE.Using the open-source RDKit [66] cheminformatics tools, two-dimensional representations of the molecules were produced.CHIMERA [67] was utilized to reconstruct docking pose representations within the binding location.

System Setup for MD Simulations and Equilibration Protocol
All protein-ligand complexes were prepared according to previously described procedures and subsequently processed using Visual Molecular Dynamics (VMD) 1.9.2 [68] and the AmberTools22 package [69,70].The protein atoms were parameterized employing the ff14SB [71] force field, and the General Amber Force Field (GAFF) [72] was used to parameterize the ligands.Partial charges were assigned to the ligands using the AM1-BCC method [73].
Each investigated system was solvated within a cubic box with a 15 Å buffer, utilizing the TIP3P [74] water model.The appropriate number of sodium and chloride ions was added to neutralize the system and achieve a salt concentration of 0.154 M. Before commencing molecular dynamics (MD) simulations, each system underwent an energy minimization consisting of 500 steps using the conjugate gradient algorithm to remove collisions.
Subsequently, each minimized system underwent a two-phase equilibration protocol.In the first phase, a 0.1 ns simulation was performed in the canonical ensemble (NVT), with harmonic positional restraints applied to both protein and ligand atoms.The second phase consisted of a 0.5 ns simulation in the isothermal-isobaric ensemble (NPT), applying the same restraints only to the ligand's position and the basic protein structure.
All MD simulations presented in this work, both in the equilibration and production phases, were executed with a 2 fs integration time step, maintaining a constant temperature of 310 K for the SuMD simulations and 300 K for the TTMD simulations, using a Langevin thermostat [75], restraining hydrogen bond lengths through the M-SHAKE [76] algorithm, and calculating long-range electrostatic interactions via the particle mesh Ewald method [77] with cubic spline interpolation.A 9.0 Å cutoff was employed for real-space electrostatic and van der Waals interactions, with a switching distance of 7.5 Å. Simulations in the NPT ensemble were conducted by maintaining constant pressure at 1 atm using a Monte Carlo barostat [78].
All MD simulations were conducted using the ACEMD16 3. Supervised Molecular Dynamics (SuMD) is a computational method that allows for the exploration of the ligand-target molecular recognition pathway on a reduced time scale.This process involves a series of unbiased Molecular Dynamics (MD) simulations, followed by an assessment of the simulation's progress using an algorithm like tabu search.In each of these MD simulations, referred to as "SuMD-steps", the system is simulated in the canonical ensemble at a constant temperature of 310 K.A 500-picosecond (ps) SuMD-step was employed in this work.
At the end of each SuMD-step, the distance between the center of mass of the ligand and the binding site is calculated along the trajectory, and the distance data are then fitted using a linear function.If the slope of the resulting linear fit is negative, indicating that the ligand is moving closer to the binding site, the "SuMD-step" is deemed productive and retained for generating the final trajectory.The concluding state of the simulation is used as the starting point for the subsequent step.Conversely, if the slope is positive, suggesting that the ligand is not approaching the binding site, the "SuMD-step" is considered unproductive and is thus discarded.In this case, the step is repeated by randomly adjusting particle velocities using the Langevin thermostat, with the final coordinates from the end of the previous "SuMD-step" retained.
The supervision algorithm is deactivated when the distance between the two centers of mass falls below a threshold value (in this instance, 10 Å).Beyond this point, the simulation continues for an additional 10 nanoseconds (ns) of classical molecular dynamics.
Ligand-receptor interaction energy alongside the SuMD trajectory was obtained through the "NAMD Energy" plugin for VMD, which exploits the NAMD [80] package to retrieve an estimate of the interaction energy defined as the sum of the van der Waals and electrostatic contribution calculated according to the force field.
Root Mean Square Deviation (RMSD) of the ligand's heavy-atom coordinates was computed during the trajectory, using the docking pose as a reference.
Electrostatic interactions were computed on a "Per Residue" basis and plotted in a heatmap displaying time on the x-axis and the 15 most contacted residues on the y-axis and using a color scale to represent the intensity of the interaction, ranging from blue to red for negative to positive values.
All the geometric analyses were performed making use of the MDAnalysis [81,82] Python library; and all the plots were generated with the Matplotlib [83] module.

Thermal Titration Molecular Dynamics (TTMD) Simulations
Thermal Titration Molecular Dynamics (TTMD) is a molecular dynamics protocol developed for the qualitative assessment of protein-ligand unbinding kinetics [53].The TTMD methodology involves a series of short classic molecular dynamics simulations, termed "TTMD-steps", conducted at progressively increasing temperatures.The temperature increment serves to increase the kinetic energy of the system, thereby reducing the simulation time required to observe protein-ligand unbinding events compared to conventional molecular dynamics simulations.To track the simulation's progress, the preservation of the native binding mode is assessed using a scoring function based on the Protein-Ligand Interaction Fingerprint (PLIF) [84].
In this case, the initial temperature was set at 300 K, and the simulation was interrupted when the ligand native binding mode and its interaction fingerprint were lost.The temperature increased by 10K between each "TTMD-step", and the duration of each simulation window was set at 10 ns.The choice of the temperature ramp was based on the preservation of the protein's native fold throughout the simulation, which was monitored by assessing the backbone Root Mean Square Deviation (RMSD).

Physicochemical and Pharmacokinetic Parameters
Some physicochemical and pharmacokinetic parameters (such as predicted logP, logS, logBB, blockage of HERG K+ channels, Caco-2 cell permeability, binding to human serum albumin, and oral absorption) were computed for the most active compounds, 15, 18, 22, 23 (Table S1 in the Supplemetary Materials), with QikProp of the Schrödinger package [86].
SAR analysis highlighted the structural requirements for targeting CK1δ.Molecular modeling studies at the CK1δ crystal structure were performed to rationalize the observed activity data.Molecular docking and molecular dynamics simulations (specifically SuMD and TTMD) were used to compare the most potent compound (23, IC 50 = 98.6 nM) with 5-CN and the unsubstituted analogue (13, IC 50 = 4.21 µM), and the analysis ascribed the potency increase to interaction between the cyano group and two crystallographic water molecules in the binding site.These results pointed out that both protein-ligand and solvent-mediated interactions should be taken into consideration for a comprehensive understanding of ligand binding and activity.
Based on these findings, further modifications will be achieved on this benzimidazole series, either at the level of the benzofused ring or on the lateral chain, where the pyrazole residue will be replaced, among others, with unsaturated nitrogen containing rings to improve the compound's solubility.

Supplementary Materials:
The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/ph17040468/s1,Video S1: Selected docking poses and relative color-coded per-residue interaction energy fingerprint (electrostatic and hydrophobic) for all the compounds.;

new benzimidazole derivatives 1 - 6 ,
taking Bischof-5 as the lead compound (Fig replacing its thiazole ring with a pyrazole nucleus, as well as replacing the NH to the terminal aryl group with a more flexible CH2CO residue.

Figure 4 .
Figure 4. Interaction Energy Fingerprints (IEFs) comparison between compounds 7, 13, and 23; on the left side, the electrostatic contribution comparison is shown, and the hydrophobic one is shown on the right.

Figure 4 .
Figure 4. Interaction Energy Fingerprints (IEFs) comparison between compounds 7, 13, and 23; on the left side, the electrostatic contribution comparison is shown, and the hydrophobic one is shown on the right.

Figure 5 .
Figure5.On the left (A), the ligand's RMSD values during the simulation can be seen compared its reference (docking pose).On the right (B), the image of the "Per Residue" analysis was perform on the 15 residues most frequently contacted by the ligand during the simulation.

Figure 5 .
Figure5.On the left (A), the ligand's RMSD values during the simulation can be seen compared to its reference (docking pose).On the right (B), the image of the "Per Residue" analysis was performed on the 15 residues most frequently contacted by the ligand during the simulation.

32 Figure 6 .
Figure 6.The binding mode of the last frame of the molecular dynamics simulation is on the left, and its heat maps of IEele and IEhyd are on the right.

Figure 6 .
Figure 6.The binding mode of the last frame of the molecular dynamics simulation is on the left, and its heat maps of IE ele and IE hyd are on the right.

Figure 7 .
Figure 7. (A) Titration profile of compound 23.(B) Titration profile of compound 13.(C) Interaction Fingerprint Similarity and RMSD values of compound 23.(D) Interaction Fingerprint Similarity and RMSD values of compound 13.

Table 1 . Inhibition activity of compounds 1-12 against CK1δ R R IC50 (μM) a 1
a Data represent the mean ± SE of three independent experiments performed in duplicate.

R R IC 50 (µM) a 1
a Data represent the mean ± SE of three independent experiments performed in duplicate.

Table 2 . Inhibition activity of compounds 13-32 against CK1δ. R IC50 (μM) a 13
a Data represent the mean ± SE of three independent experiments performed in duplicate.