The Discovery of a Novel IκB Kinase β Inhibitor Based on Pharmacophore Modeling, Virtual Screening and Biological Evaluation

Abstract

Background:

IκB kinase β (IKKβ) plays a pivotal role in the NF-κB signaling pathway and is considered a promising therapeutic target for various diseases.

Materials & methods:

The authors developed and validated a 3D pharmacophore model of IKKβ inhibitors via the HypoGen algorithm in Discovery Studio 2019, then performed virtual screening, molecular docking and kinase assays to identify hit compounds from the ChemDiv database. The compound with the highest inhibitory activity was further evaluated in adjuvant-induced arthritis rat models.

Results:

Among the four hit compounds, Hit 4 had the highest IKKβ inhibitory activity (IC₅₀ = 30.4 ± 3.8), and it could significantly ameliorate joint inflammation and damage in vivo.

Conclusion:

The identified compound, Hit 4, can be optimized as a therapeutic agent for inflammatory diseases.

Plain language summary

This research paper focuses on the development and validation of IκB kinase β (IKKβ) inhibitors. IKKβ is a crucial enzyme that plays an important role in the NF-κB signaling pathway, which is involved in many diseases such as inflammatory diseases and cancers. The researchers used computer-aided drug design strategies to identify potential IKKβ inhibitors. First, they used a model to screen a large database of chemical compounds. Then, they conducted further tests to pinpoint the ones that could effectively inhibit IKKβ. Out of all the tested compounds, one referred to as ‘Hit 4’ showed the highest inhibitory activity. It was even able to significantly reduce joint inflammation and damage in rat models. Although many drugs targeting IKKβ have been developed, none are commercially available yet due to issues with efficacy or safety. Therefore, the findings of this study are significant and could lead to the development of new effective therapeutic agents for inflammatory diseases.

Graphical abstract

Keywords::

• IKKβ inhibitor
• molecular docking
• pharmacophore modeling
• virtual screening

IKKβ is a versatile serine kinase that plays a pivotal role in the well-characterized NF-κB signaling cascade, as well as being engaged in multiple NF-κB-independent biological processes. IKKβ, together with another catalytic subunit, IKKα, and a regulatory subunit, NEMO, constitute the IKK complex, which is responsible for NF-κB activation [1]. The canonical NF-κB pathway is rapidly initiated upon various stimuli, leading to IKKβ phosphorylation and activation in a NEMO-dependent manner. This results in K48-linked ubiquitination and proteasomal degradation of IκBα [2], allowing NF-κB to translocate into the nucleus and regulate target gene transcription [3]. The target genes of the canonical pathway are generally involved in inflammation, pro-proliferation, anti-apoptosis and angiogenesis [1,4–6]. Thus, dysregulation of NF-κB signaling has been an important driver of various inflammatory diseases, most solid cancers [7–9] and hematopoietic malignancies [10]. Moreover, IKKβ also phosphorylates substrates such as p53 [11,12], FOXO3a [13] and RIPK1 [14,15] independent of NF-κB, facilitating cancer cell survival and metabolism.

The crystal structures of IKKβ have revealed key features, including the kinase domain (KD), ubiquitin-like domain (ULD) and scaffold dimerization domain [16–18]. IKKβ in mammalian cell extracts exist predominantly as a dimer, with the scaffold dimerization domain mediating the dimerization process. However, higher-order assemblies require conformational changes, where KD–KD interactions can be realized through a V-shaped interface, and the activation loop becomes accessible for transphosphorylation. It is well recognized that the critical step for IKKβ’s final activation lies in the specific phosphorylation on two serine residues (Ser181 and Ser177), while the process of how IKKβ turns from inactive to active still lacks comprehensive interpretation. Recent studies have highlighted the indispensable role of NEMO in relaying upstream ubiquitination signals and directly inducing IKKβ activation; thus, blocking NEMO–IKKβ interaction using competitive peptides inhibits IKKβ activation in vitro [19].

Many selective IKKβ inhibitors have been developed, mostly ATP-competitive ones such as MLN120B [20]. Some ATP noncompetitive inhibitors such as BMS-345541 [21,22] have also been designed to target unique allosteric sites or key cysteine residues. However, they have failed clinical trials due to efficacy or safety issues. Representative IKKβ inhibitors that have entered clinical trials include CSP-1103 (showed good safety in phase I but phase II was terminated), IMD-0354 (phase I showed poor efficacy), SAR-113945 (phase I demonstrated suboptimal efficacy), BMS-345541 (showed toxicity in phase I). Hence, currently no IKKβ inhibitor is commercially available.

Computer-aided drug design provides an efficient strategy for discovering new inhibitors. In this study, the authors conducted ligand-based virtual screening against a commercial compound database (ChemDiv’s 300k Representative Compounds Library) and subsequently performed molecular docking to identify novel IKKβ inhibitors. Four commercially available hit compounds were identified, which then underwent in vitro bioactivity assays. Further validation of the effects of Hit 4 was carried out using complete Freund’s adjuvant-treated arthritic rats. In conclusion, the results of the study may facilitate subsequent efforts to optimize and develop novel IKKβ inhibitors.

Materials & methods

Pharmacophore model generation

Molecule collection

A reliable pharmacophore model requires appropriate IKKβ inhibitors as input. The following criteria should be met [23]: wide activity range, inclusion of the most potent inhibitor in the training set and structural diversity.

56 chemically diverse IKKβ inhibitors with experimental IC₅₀ values ranging from 3 to 19, 100 nM were collected from available literature resources [20,24–40]. 20 representative inhibitors constituted the training set, while the rest formed the test set. The training set was used to generate the pharmacophore model, which was evaluated by the test set. Inhibitory activities (IC₅₀) of both sets spanned over four orders of magnitude, according to which the selected molecules were classified as highly active (IC₅₀ <10 nM, ++++), active (10 nM ≤ IC₅₀ <100 nM, +++), moderately active (100 nM ≤ IC₅₀ <1000 nM, ++) and inactive (1000 nM ≤ IC₅₀, +). Their 2D structures were drawn in ChemBioDraw Ultra 14.0 and then converted to 3D conformations in Discovery Studio 2019, from Accelrys (CA, USA). Subsequently, the ‘Calculate Principal Components’ module was applied to convert these compounds to discrete points in the 3D coordinate system. Discrete points taking up more than 70% of the descriptor space indicates that the collected set is representative.

Molecule preparation

Prior to virtual screening and molecular docking, compounds were prepared using the ‘Ligand Preparation and Full Minimization’ modules. ‘Ligand Preparation’ generates low-energy conformations by processing ligands via the removal of irrelevant atoms, addition of hydrogen atoms and generation of ionization states. In this study, all input atoms were retained, with consideration of all potential protonation states, generating up to 32 conformers per compound. ‘Full Minimization’ further optimizes the ligand conformations from ‘Ligand Preparation’ using the ‘CHARMm forcefield’, with 1000 steps of minimization applying a maximum step size of 0.001. Default parameters were used in this step. With the above two-step preparation, optimized low-energy conformations were obtained for subsequent virtual screening and docking procedures. In this paper, the authors adopted the 300k Representative Compounds library from the ChemDiv database (www.chemdiv.com) for virtual screening, and the collected 207,525 compounds were prepared in the same fashion as above.

3D quantitative structure–activity relationship pharmacophore modeling

The authors utilized the HypoGen algorithm in Discovery Studio’s 3D quantitative structure–activity relationship (QSAR) pharmacophore modeling protocol to generate a ligand-based pharmacophore model. First, known inhibitory activities (experimental IC₅₀ values against IKKβ) and uncertainty values were assigned to all compounds by adding new attributes. The uncertainty value was set to 2.0 Å, since the IC₅₀ values of both the training and test sets spanned the required four orders of magnitude [41]. Next, up to 255 conformers were generated for each compound, with the energy threshold at 10 kcal/mol and minimum interfeature distance at 1.5 Å, using the best conformational search option. All other parameters were kept default. The ‘Feature Mapping’ protocol was then applied to the training set to identify distinct pharmacophoric features, including hydrogen bond acceptor, hydrogen bond donor, positive ionizable, aromatic ring and hydrophobic groups. The pharmacophore feature maps of the training compounds were utilized to construct 3D-QSAR pharmacophore models. In this study, the top ten hypotheses from HypoGen were selected for further analysis based on their cost values (null, fixed and total costs), root mean square (RMS), correlation coefficient (R²) and fit values. The null cost represents the highest cost of a hypothetical pharmacophore without any selected features. The fixed cost equals the cost of the simplest model that perfectly fits all the input data. A model’s reliability depends on whether its total cost value is distant from the null cost and close to the fixed cost. If the cost difference between null and total costs is 40–60 units, the model has a 75–90% predicted correlation probability [42]. Additionally, the RMS indicates the difference between estimated and experimental IC₅₀ values of the training compounds, reflecting the model’s predictive capability. The correlation coefficient is calculated by linear regression of the geometric fit index and is expected to approach one.

Pharmacophore model evaluation

An optimal pharmacophore model should satisfy two criteria: effectively predict the activity of compounds within and outside the training set and achieve statistical significance. The authors utilized two methods to validate the generated hypotheses: external validation and the Fisher randomization test.

External validation

To evaluate the predictive capability, the authors tested the hypothesis against 36 compounds from the test set. A reliable model should accurately classify these molecules based on activity. Low energy conformations were generated using ‘Conformation Generation’ with the FAST algorithm. ‘Ligand Pharmacophore Mapping’ was applied to the conformers using the ‘BEST’ search option, with ‘Maximum omitted features’ set to -1. Other parameters were kept default.

Fisher randomization test

The Fisher randomization test assesses the statistical significance and structure–activity correlation of a hypothesis. In this validation, random hypotheses were constructed using randomized activity data of the training set compounds. The random hypotheses possessed the same features and parameters as the original hypothesis. Statistical significance was calculated as $S=\left[1-\left(\frac{1+X}{Y}\right)\right]\times 100\%$, where X is the total number of hypotheses having a total cost lower than the original Hypo X and Y is the total number of initial HypoGen runs plus random runs. With 95% confidence level, 19 random spreadsheets were generated automatically. If the random hypotheses showed similar or better cost, RMS deviation and correlation, it implies the original hypothesis was generated by chance [43].

Drug-likeness evaluation & virtual screening

To filter out compounds with poor drug-likeness, the authors first applied Lipinski’s rule of five, which defines criteria for compounds with good bioavailability: molecular weight ≤500 Da, hydrogen bond donors ≤5, hydrogen bond acceptors ≤10, partition coefficient (logP) ≤5 and rotatable bonds ≤10.

After initial screening, absorption, distribution, metabolism, excretion, toxicity (ADMET) properties of the remaining compounds were predicted using the ‘ADMET Descriptor’ protocol in Discovery Studio. Key pharmacokinetic parameters were analyzed for hit compounds, including human intestinal absorption, aqueous solubility, blood–brain barrier permeability, plasma protein binding and CYP2D6 enzyme inhibition [44]. This step was crucial for drug discovery.

Compounds passing the above steps were mapped to Hypo 1 to obtain those matching all pharmacophoric features, which were then used in molecular docking studies.

Molecular docking

Molecular docking involves sampling ligand poses within the receptor binding site, and is a powerful technique to predict optimal binding modes of hit compounds and estimate their affinities using scoring functions. Compound 1 from the test set, with high activity (IC₅₀ = 4 nM), was chosen as the reference to compare other compounds. The docking procedure contains three main steps: prepare the receptor–ligand system, perform binding calculations and analyze the interaction results.

The crystal structure of human IKKβ (hIKKβ, Protein Data Bank [PDB] ID: 4KIK) was downloaded from the PDB as the receptor. The ‘Protein Preparation’ module was utilized to remove crystallographic water and optimize the structure by adding missing hydrogen atoms with CHARMm forcefield. Binding site definition is critical for accurate docking and scoring. The potential IKKβ binding sites were identified using the ‘Define and Edit Binding Site’ module.

Hits from the database were refined by LibDock and CDOCKER (CHARMm-based DOCKER). LibDock rigidly docks ligands into hotspots on the receptor based on structural features. Polar and apolar hotspots match ligand polar and apolar atoms, respectively. With user-specified settings and max hits 10, LibDock rapidly screened compounds and prioritized those with LibDockScore higher than the reference. These were passed to CDOCKER, which enables full ligand flexibility and the CHARMm forcefields [45]. Cluster radius was set to 0.5 Å. CDOCKER ranks outputs using nonbonded interaction energy between protein and ligand, with more negative values indicating more favorable binding.

In vitro kinase activity assay

Measurement of the selected compound’s inhibitory activity on IKKβ was performed in 96-well plates (Sigma Aldrich, Castle Hill, Australia) of the ADP-Glo™ analysis kit (Promega Corporation, WI, USA). The kinase reaction buffer contained 3 μl standard buffer (167 mM 4-[2-hydroxyethyl]-1-piperazine ethanesulfonic acid [HEPES]-NaOH pH 7.5, 10 mM MgCl₂, 10 mM MnCl₂), 1.5 μl kinase dilution buffer (50 mM HEPES-NaOH pH 7.5, 0.25 mg/ml PEG 20,000 and 1 mM dithiothreitol), 80 ng IKKβ protein kinase, 2.5 μl HEPES (50 mM) and 0.24 μg RBER IRStide. Finally, 2 μl ATP (19.5 μM) and the test compound were added for a total volume of 10 μl. After a 40-min incubation at 37° C, the unconsumed ATP was depleted using ADP-Glo reagent at room temperature for 40 min. Subsequently, the kinase detection reagent was added to convert ADP to ATP, while luciferase and luciferin were also introduced to detect ATP. They were incubated at room temperature in a dark room for 60 min, during which the luminescence was recorded by PerkinElmer EnVsion^® multimode reader (MA, USA). After the reaction was terminated, the change in ADP content was detected to reflect the IKKβ enzyme activity level. Finally, inhibition curves with different compound concentrations were generated to determine the IC₅₀ values against IKKβ through curve fitting.

In vivo pharmacological effects

Studies have demonstrated that IKKβ activation participates in the pathogenesis of osteoarthritis. Inhibiting IKKβ activation alleviates cartilage damage in osteoarthritis mouse models [46,47]. Moreover, complete Freund adjuvant (CFA) activates IKKβ and the NF-κB pathway to induce inflammatory arthritis models. Therefore, the authors utilized a complete Freund’s adjuvant-induced arthritis (AIA) rat model to observe the IKKβ inhibitory effects of selected compounds as follows.

30, 160–180 g, female Sprague–Dawley rats were housed at 23–25°C, with 40–60% humidity and a 12-h light–dark cycle under protocols. The rats were randomly divided into five groups (N = 6 per group): normal control group, AIA group, indomethacin (10 mg/kg) + AIA group, Hit 4 (10 mg/kg) + AIA group and Hit 4 (30 mg/kg) + AIA group. Rats in the experimental groups received intracutaneous injection of 1 ml CFA into the left hind foot. The sham control group received an equal volume of paraffin oil. 10 days postimmunization, all groups showed first signs of arthritis, including visible redness and swelling. Meanwhile, 10 days after CFA administration, the Hit 4 (10 mg/kg) + AIA group and Hit 4 (30 mg/kg) + AIA group received 10 and 30 mg/kg Hit 4, respectively, for 14 days via gastric gavage, while the positive control group received 10 mg/kg indomethacin for 14 days. 24 days after CFA injection, all animals were sacrificed. Then, the inflamed contralateral ankle joints were collected, fixed in 4% paraformaldehyde, decalcified with 5% formic acid and paraffin-embedded. Hematoxylin and eosin staining allowed light microscopic examination of histopathological changes.

For macroscopic observation, the authors also estimated the volume of the paw and body weight. Additionally, arthritis severity was graded by arthritis index (AI) based on the following parameters: no swelling (grade 0), erythema in a finger or mild swelling (grade 1), swelling in fingers (grade 2), swelling of the ankle or wrist (grade 3) and severe arthritic swelling in wrist and fingers (grade 4). The maximum value of AI was 16 for each mouse [48].

Results & discussion

3D quantitive structure–activity relationship pharmacophore modeling

Supplementary Figure 1 presents the structures and biological activities of the collected compounds. A training set containing 20 molecules with IC₅₀ values ranging from 3 to 19,100 nM and a test set containing 36 molecules with IC₅₀ values ranging from 4 to 10,000 nM were prepared for pharmacophore modeling. Principal component analysis was carried out to demonstrate the structural diversity of the selected compounds. As shown in Supplementary Figure 2, the training set compounds (colored yellow) and the test set compounds (colored blue) covered over 70% of the descriptor space, demonstrating their representativeness.

Using the pharmacophoric features of hydrogen bond acceptor, hydrogen bond donor, positive ionizable, aromatic ring and hydrophobic, pharmacophore modeling was performed via 3D-QSAR pharmacophore generation as described previously. Supplementary Table 1 summarizes the top ten quantitative models ranked by decreasing cost difference. Hypo 1 has the largest cost difference (72.48), lowest RMS (1.653) and highest correlation coefficient (0.880). It shows total and null costs of 106.203 and 178.683 bits, respectively. Compared with other hypotheses, the large cost difference of 72.48 bits suggests over 90% probability that Hypo 1 can effectively correlate experimental and estimated data. Its correlation coefficient and RMS indicate small deviation of predicted from experimental activities. Thus, Hypo 1 is the best ligand-based pharmacophore model for further analysis. As depicted in Figure 1, Hypo 1 comprises four targeted pharmacophoric features: one hydrogen bond acceptor, one hydrogen bond donor, one hydrophobic feature and one hydrophobic-aromatic feature. Figure 1 also shows its spatial arrangement and geometric parameters.

Figure 1. The best pharmacophore model generated by HypoGen.

(A) Pharmacophoric features (hydrogen bond acceptor, green; hydrogen bond donor, magenta; hydrophobic, cyan; hydrophobic-aromatic ring, blue). (B) 3D spatial arrangement and geometric parameters of Hypo 1.

HBA: Hydrogen bond acceptor; HBD: Hydrogen bond donor; HY-AR: Hydrophobic-aromatic ring; HY: Hydrophobic.

Supplementary Table 2 lists the estimated inhibitory activities of 20 training set compounds using Hypo 1, demonstrating its ability to estimate most compounds’ biological activities and correctly classify them by activity scale. Although predictions diverged for moderately active and inactive compounds (compound 9 was overestimated as active, compound 15 underestimated as inactive), almost all highly active and active training compounds were correctly classified by activity scale with predictions very close to actual activities. Regression analysis within the training set further verified Hypo 1’s predictive accuracy, with the correlation coefficient between estimated and experimental IC₅₀ values being 0.93 (Supplementary Figure 3). Furthermore, as shown in Figure 2, where the most and least active training set compounds were aligned to Hypo 1, training set compound 1 could be matched to Hypo 1 with a fit value of 8.51. In conclusion, this preliminary analysis indicates Hypo 1’s validity.

Figure 2. Alignment of Hypo 1 to the compound of the training set.

(A) Mapping of training set compound 1 (IC₅₀ = 3 nM) in Hypo 1. (B) Mapping of training set compound 20 (IC₅₀ = 19,100 nM) in Hypo 1.

Pharmacophore validation

External validation

Prediction of test sets was the first step to further validate Hypo 1. The test set comprised 36 structurally diverse molecules with inhibitory activity ranging from 4 to 10,000 nM. These compounds were also categorized into four groups based on their IC₅₀ values. The ‘Ligand Pharmacophore Mapping’ protocol with the ‘Best Flexible Search’ option was then applied to map each compound, and the results are presented in Supplementary Table 3 & Figure 3. Most compounds were anticipated to fall within the same activity range as their experimental values, with the exception of two moderately active compounds whose activities were overestimated. Furthermore, Hypo1 exhibited a significant correlation between the predicted and actual biological activities of the test set compounds (R² = 0.87; Supplementary Figure 3). This result validated that Hypo1 was capable of predicting the inhibitory activities of both the internal training set and the external test set compounds.

Figure 3. Chemical structures of the final four hit compounds.

Fischer randomization test

Fischer randomization testing further evaluated Hypo 1’s statistical significance. Comparing total costs of these random models with the original Hypo 1 showed that none of the top 19 random hypotheses scored lower total costs than the original. Supplementary Figure 4 shows the cost difference between HypoGen and Fischer randomizations. These cross-validation results strongly evidence that Hypo 1 was not generated by chance but represents a true correlation in the training set.

Virtual screening

Virtual screening is a convenient approach widely used for the discovery of innovative lead compounds from large databases. In this study, virtual screening was performed according to the designed strategy as presented in the following.

Drug-likeness filtration

Database ChemDiv 300k Representative Compounds Library was initially screened using Lipinski’s rule of five and ADMET property calculations in Discovery Studio 2019. As a result, 6914 out of 207,525 selective natural compounds fulfilled the above criteria and were subjected to the subsequent screening. This step successfully narrowed the scope of screening by filtering out unfavorable compounds that were poorly absorbed in human bodies and lacking in drug-likeness.

Pharmacophore mapping

Subsequently, Hypo 1 was employed as a 3D query to search the database using Ligand Pharmacophore Mapping protocol with flexible search option; 117 hit compounds with fit values higher than 7.72 (the highest fit value of the training set compounds, test 24) were quickly obtained and were chosen for further docking studies.

Before the docking studies, we verified the reliability of the docking programs by redocking the reference ligand, test 2 (an inhibitor from the test set with IC₅₀ = 8 nM) [25], to the active binding sites of human IKKβ protein (PDB ID: 4KIK). The RMS deviation was used to compare the docking postures with the initial postures as shown in Supplementary Table 4, and it was found that the docking positions were almost identical in both programs (Supplementary Figure 6).

Molecular docking

Molecular docking was applied to further refine the hit compounds and, in the best possible way, mitigate the risk of false-positive and false-negative errors in the previous pharmacophore mapping; 117 hits from the L6020 database together with the compounds in the test set docked at the active sites of 4KIK using the LibDock module. Within the test set, test 2, with the phenyl-(4-phenyl-pyrimidin-2-yl)-amine moiety as the molecular core [34], had the highest LibDockscore (117.9 kcal/mol). Therefore, it was selected as the reference inhibitor for further screening. A total of 19 compounds with higher LibDockscore than that of test 2 were obtained, and all of them were subjected to the CDOCKER protocol, docking at the binding sites again with higher resolution. Based on the scoring function of -CDOCKER_INTERACTION_ENERGY, finally the authors obtained four ideal hit compounds, as displayed in Figure 3, each with a greater value of CDOCKER interaction energy than that of test 2 (51.99 kcal/mol). The four hit compounds were named Hit 1, Hit 2, Hit 3 and Hit 4 according to the decreasing order of calculated -CDOCKER_INTERACTION_ENERGY values (60.74, 54.39, 52.80 and 52.56 kcal/mol, respectively). Their interactions with 4KIK are illustrated in Figure 4 and discussed as follows.

Figure 4. 2D diagrams of the docking results.

Hit 1 (compound13_439), Hit 2 (compound10_288), Hit 3 (compound14_8), Hit 4 (compound10_19), reference inhibitor (compound 2 of the test set). Bound compounds are shown in ball-and-stick models, with hydrogen, nitrogen, oxygen, sulfur and fluorine colored in white, blue, red, yellow and cyan, respectively. Nonbonded interactions are represented as dotted lines, with salt bridge colored in orange, carbon–hydrogen bond colored in light green, conventional hydrogen bond colored in green, alkyl or π–alkyl colored in pink, π–σ colored in purple, halogen (fluorine) interactions colored in cyan, π–sulfur colored in yellow. The solvent-accessible surface of an interacting residue is surrounded by a blue halo, and the diameter of the circle is proportional to the solvent-accessible surface.

Noncovalent molecular interactions underlie molecular recognition and subsequent biological processes between receptors and ligands. Comprehensive analysis of these interactions is therefore crucial for drug discovery and design. Coulomb’s law explains the basis of intermolecular interactions, elucidating key factors such as atomic distances, charge quantities and the dielectric constant of the medium. Notably, these weak atom–atom interactions can be grouped into five main categories: hydrogen bonds and electrostatic, hydrophobic, halogen and other miscellaneous interactions. Each category can further be divided into subtypes – for example, hydrogen bonds consist of traditional and nontraditional ones; hydrophobic interactions include alkyl, π–alkyl, π–σ and π–sulfur interactions; and salt bridges, with their relatively stronger noncovalent nature, belong to both hydrogen bonds and electrostatic interactions. The CDOCKER program results are visualized in the 2D graphs in Figure 4, with the binding details enumerated in Table 1.

Table 1. Molecular interaction details of the hit compounds with IκB kinase β.

Compound	Hydrogen bond	Hydrophobic	Electrostatic	Other	Score (kcal/mol)^†
	H-donor:H-acceptor (length/Å)
Hit 1	Lys106:HZ1-O21(2.48) H51-O22(2.10, intermolecular) Gly102:HA1-O21(2.67) Gly102:HA2-O22(2.34)	Lys44 3 Ile165 Val74 3 Val29 2 Ala42 Leu21 Met96	Lys106	2 Met96:SD	60.74
Hit 2	H33-Glu97: O(2.70) Gly102:HA2-N26(2.32) Lys106:HZ2-O29(2.27) Lys147HZ3:F20(2.48) Lys147:HE1-F20(2.89) H42-O2(2.274, intermolecular)	Lys44 2 Ile165 2 Val29 Leu21		2 Met96:SD Asn150:OD1-F20	54.39
Hit 3	Lys44:HZ2-O20(1.81) Lys106:HZ2– O11(1.81) H23-Glu97: O(1.95) H44-Asp166: OD1(2.05)	Leu21 Val29 Lys44 Met96 Ile165			52.80
Hit 4	Cys99:HN-O16 Asp166:HN-O7 H36-Glu61:OE2 H41-Glu97:O Tyr98: HA-O16	Ile165 Val29 Val152 Leu21:HB1			52.56
Reference inhibitor (test 2)	Asp103:HN-O24(2.78) H39-Asn150:OD1(2.77) H40-Asn150:OD1(2.38) H41-Asp103:OD1(2.06) H42-Leu21:O(2.28) Gly102:HA1-O24(2.54) Gly102:HA2-O24(2.87) Gly24:HA1-O27(2.63) Gly24:HA2-O26(2.79) Gly24:HA2-O27(2.71) H29-Glu97:O(2.31)	2Leu21 2Val29 Cys99 Ala42 2Val152 2Ile165		Met96:SD	51.99

† Score is calculated according to the -CDOCKER_INTERACTION_ENERGY scoring function.

All small molecules that bind to the kinase domain of the unphosphorylated Chain A interact with both the N-terminal (1–109) and C-terminal lobe (110–307). The simplified 2D graph and table provide deeper insights into their docking results. All hit compounds formed hydrogen bonds and hydrophobic interactions with 4KIK’s active site, while only Hit 1 formed a salt bridge and Hit 2 uniquely showed halogen (fluorine) interactions. Hit 1 scored highest, though differences with other hits were small. Specifically, Hit 1’s hydroxyl oxygen (O21) and carboxylic acid oxygen (O22) serve as hydrogen acceptors, with 4KIK’s Gly102 polarized carbons as donors, forming two 2.67 and 2.34Å carbon hydrogen bonds. Stability also stems from hydrophobic interactions like alkyl and π–alkyl with Val29, Ile165 and Lys44. A π–sulfur interaction with Met96 surrounds its benzopyrrole ring. Further, the carboxylate anion–Lys106 contact creates a relatively strong electrostatic interaction absent in other receptor–ligand simulations, potentially explaining Hit 1’s top ranking. However, unlike other hits, Hit 1 lacked conventional hydrogen bonds with 4KIK. Lys106, located at one salt bridge end, is a solvent-accessible interacting residue denoted by a blue halo in the 2D diagram. Thus, Hit 1’s inhibitory ability in solution remains experimentally unvalidated. Hit 2 also formed various 4KIK contacts through short-range interactions like halogens (fluorine) and conventional hydrogen bonds. Hit 3 only showed 4KIK hydrogen bonds and π–alkyl hydrophobic interactions but had more conventional hydrogen bonds. Notably, Hit 4 formed the most hydrogen bonds with 4KIK, with one Hit 4 forming four conventional hydrogen bonds, one carbon–hydrogen bond and π–σ/other hydrophobic interactions. The lack of salt bridges, halogens or other interactions may explain its lowest score.

IKKβ inhibitory activity in vitro

The Promega ADP-Glo analysis kit (Promega Corporation) measured the inhibitory potency on IKKβ of the four hit compounds and the known inhibitor ertiprotafib (IC₅₀ = 400 nM) [40], and the results are summarized in Table 2. The reported IC₅₀ value of ertiprotafib against 4KIK was 412.8 ± 20.5 nM, validating the test conditions. All hits displayed considerable inhibitory activity. Contrary to docking results, the top-scoring Hit 1 exhibited the lowest inhibition in solution, even below ertiprotafib. Hit 4 was strongest, while the ranks of Hits 2 and 3 were unchanged. Discrepancies could stem from nonbonded interactions and residue roles in IKKβ activation.

Table 2. IC₅₀ values of hit compounds and reference inhibitors against IκB kinase β.

Compound	IC50 ± SD^† (nM)
Hit 1	420.4 ± 24.6
Hit 2	83.2 ± 7.3
Hit 3	101.6 ± 11.4
Hit 4	30.4 ± 3.8
Ertiprotafib	412.8 ± 20.5

† Concentration of compound that inhibits 50% kinase inhibition activity, presented as mean ± SD, performed each test three-times.

First, docking showed Hit 1 formed favorable salt bridges with human IKKβ (hIKKβ). However, in solution, the positively charged Lys106 can be interfered by highly electronegative water oxygen atoms, diminishing the salt bridge’s binding stability contribution. Although remaining hydrophobic/hydrogen bonds guide Hit 1 to the active site, these interactions are relatively weak. In contrast, Hit 4 uniquely contacted several key hIKKβ residues – Glu61, Glu97, Tyr98, Cys99 and Asp166. Activated hIKKβ has the C-helix’s Glu61 interact with the β3 strand’s catalytic Lys44 to stabilize ATP binding, while Hit 4 hydrogen bonding with Glu61 weakens hIKKβ’s catalysis. Glu97, Tyr98 and Cys99 in the hinge connect the kinase lobes to aid ATP/substrate binding. Hit 4’s hydrogen bonds with these residues indicate ATP site selectivity. Hit 4 also recognizes the conserved DFG motif in the activation loop by hydrogen bonding Asp166, potentially preventing ASP from contacting all three ATP phosphates.

Pharmacological effects of Hit 4 in vivo

The pharmacological effects of Hit 4 in vivo were further confirmed by evaluating the rat AIA model through both macroscopic and microscopic observations. Macroscopic evaluation of AIA rats showed paw edema and inflammatory polyarthritis induced in all experimental animals (Figure 5). Secondary inflammation appeared around day 10, so therapeutic administration of Hit 4 (10 and 30 mg/kg) began then. Compared with the model group, Hit 4 (30 mg/kg)-treated groups significantly reduced rat paw swelling in a concentration-dependent manner, with efficacy similar to 10 mg/kg of indomethacin.

Figure 5. Macroscopic observation of the paws.

Representative pictures of swelling joints are presented (A). Paw volume (B). Arthritis index (C) was evaluated following the published standard. Body weight (D) was measured.

Subsequently, histological sections of ankle joints helped investigate Hit 4’s effects on pathological changes in AIA rats. As shown in Figure 6, in normal rats, the ankle joint structure was intact, with no inflammatory cell infiltration. The synovial cells were neatly arranged, the cartilage surface was smooth and no effusion was observed in the joint cavity. In contrast, the model group displayed obvious synovial cell proliferation, inflammatory cell infiltration, pannus formation and cartilage destruction. These symptoms were alleviated in AIA rats treated with Hit 4 in a dose-dependent manner. Although Hit 4 appeared to relieve paw swelling to a similar extent as indomethacin on gross examination, histological results showed that compared with indomethacin, high-dose Hit 4 (30 mg/l) could significantly reduce inflammatory cell infiltration and hyperemia of the synovium.

Figure 6. Histological changes of the ankle joints.

Representative sections of ankle joints stained with hematoxylin and eosin. Cartilage damage, ×2, and synovial inflammation, ×200, were examined. In the Hit 4-treated samples, significant recovery of chondrocytes and cartilage was observed, indicating that Hit 4 promoted cartilage regeneration.

Conclusion

To summarize, utilizing several molecular modeling approaches, a new scaffold for IKKβ inhibitors was developed. A valid 3D-QSAR pharmacophore model was constructed, which was then utilized to screen a ChemDiv database based on pharmacophore matching. The virtually screened compounds were filtered by docking them to the IKKβ structure with PDB code 4KIK. After the molecular modeling procedures, the four hit compounds were evaluated biologically against IKKβ. The approximated IC₅₀ values of Hits 1–4 were 420.4, 83.2, 101.6 and 30.4 nM, respectively. Moreover, an in vivo study verified that Hit 4 mitigates inflammation in arthritic rats treated with complete Freund adjuvant. Therefore, Hit 4 represents a novel molecular framework which could be further enhanced and employed to create promising potent IKKβ inhibitors.

In conclusion, our study provides valuable insights into the identification and development of potential IKKβ inhibitors. However, we acknowledge a few limitations in our work. First, the selectivity of these kinase inhibitors was not assessed in this study. Considering the ATP-competitive nature of the inhibitors, it is crucial to evaluate their selectivity among different kinases, which we plan to address in future work. Second, many of the compounds investigated contain sulfamide, primary sulfonamide, carboxylate or phenolic moieties, which are known to interact with a broad range of metalloenzymes. This aspect was not fully considered in our study, and we recognize the importance of investigating these potential interactions in subsequent studies. Also, it should be noted that the reasons for the clinical failure of previous IKKβ inhibitors are related to the complexity of the IKKβ pathway. Systemic inhibition of IKKβ may affect normal immune responses.

Despite these limitations, we believe our findings are a significant step forward in the quest for potent IKKβ inhibitors. We look forward to addressing these limitations in future research, which will involve rigorous selectivity testing and thorough investigation of potential interactions with metalloenzymes. Future drug development may consider inhibitors with higher specificity and lower side effects or may explore local administration. A deeper understanding of the precise mechanisms of the IKKβ pathway in diseases is still needed.

Summary points

• Hypo 1, comprising one hydrogen bond acceptor, one hydrogen bond donor, one hydrophobic feature and one hydrophobic-aromatic feature, was the best pharmacophore model generated by HypoGen based on the training set compounds.

• Hit compounds (n = 117) were filtered out from ChemDiv 300k Representative Compounds Library (n = 207,525), and four of them were finally selected for their molecular docking results against human IKKβ (Protein Data Bank ID: 4KIK).

• Hit 4 showed excellent inhibitory activity against IKKβ (IC₅₀ = 30.4 ± 3.8), probably because of its hydrogen bonding with the C-helix’s Glu61, as well as the hinge area’s Glu97, Tyr98 and Cys99 of IKKβ.

• Hit 4 (30 mg/kg) was also able to reduce paw swelling and ameliorate joint inflammation in adjuvant-induced arthritis rats, with efficacy slightly higher than that of indomethacin (10 mg/kg).

Author contributions

L Li: data collection, manuscript preparation, kinase assays, animal experiments, data analysis; S Gong: research conception, manuscript revision. Authors contributed to the article and approved the submitted version.

Ethical conduct of research

All animal experiments were carried out following the principles outlined in the Declaration of Helsinki and approved by the Institutional Animal Care and Use Committee of Xi’an Jiaotong University.

Data sharing statement

The results, data and figures in this manuscript have not been published elsewhere, nor are they under consideration by another publisher. Data will be made available on request.

Acknowledgments

The authors would like to acknowledge scholars’ contribution to the ChemDiv for the availability of the data. The authors are grateful to the editors and reviewers for their helpful comments on this paper. Meanwhile, the authors acknowledge the efforts of all researchers in developing Discovery Studio as well as ChemDraw software. Financial support via Xi’an Science and Technology Plan and the Key Research and Development Project of Shaanxi Province is acknowledged.

Supplementary data

To view the supplementary data that accompany this paper please visit the journal website at: www.tandfonline.com/doi/suppl/10.4155/fmc-2023-0261

Financial disclosure

This work was supported by the National Natural Science Foundation of China (grant no. 81903268), the Shaanxi Province Key R&D Program (grant no. S2023-YF-YBSF-0261), Wu Jieping Medical Fund and the Shaanxi Province Key R&D Program (grant no. 2022SF-110). The authors have no other relevant affiliations or financial involvement with any organization or entity with a financial interest in or financial conflict with the subject matter or materials discussed in the manuscript apart from those disclosed.

Competing interests disclosure

The authors have no competing interests or relevant affiliations with any organization or entity with the subject matter or materials discussed in the manuscript. This includes employment, consultancies, honoraria, stock ownership or options, expert testimony, grants or patents received or pending or royalties.

Writing disclosure

No writing assistance was utilized in the production of this manuscript.