Virtual Screening Strategy and In Vitro Tests to Identify New Inhibitors of the Immunoproteasome

Immunoproteasome inhibition is a promising strategy for the treatment of hematological malignancies, autoimmune diseases, and inflammatory diseases. The design of non-covalent inhibitors of the immunoproteasome β1i/β5i catalytic subunits could be a novel approach to avoid the drawbacks of the known covalent inhibitors, such as toxicity due to off-target binding. In this work, we report the biological evaluation of thirty-four compounds selected from a commercially available collection. These hit compounds are the outcomes of a virtual screening strategy including a dynamic pharmacophore modeling approach onto the β1i subunit and a pharmacophore/docking approach onto the β5i subunit. The computational studies were first followed by in vitro enzymatic assays at 100 μM. Only compounds capable of inhibiting the enzymatic activity by more than 50% were characterized in detail using Tian continuous assays, determining the dissociation constant (Ki) of the non-covalent complex where Ki is also the measure of the binding affinity. Seven out of thirty-four hits showed to inhibit β1i and/or β5i subunit. Compound 3 is the most active on the β1i subunit with Ki = 11.84 ± 1.63 µM, and compound 17 showed Ki = 12.50 ± 0.77 µM on the β5i subunit. Compound 2 showed inhibitory activity on both subunits (Ki = 12.53 ± 0.18 and Ki = 31.95 ± 0.81 on the β1i subunit and β5i subunit, respectively). The induced fit docking analysis revealed interactions with Thr1 and Phe31 of β1i subunit and that represent new key residues as reported in our previous work. Onto β5i subunit, it interacts with the key residues Thr1, Thr21, and Tyr169. This last hit compound identified represents an interesting starting point for further optimization of β1i/β5i dual inhibitors of the immunoproteasome.


Introduction
The ubiquitin-proteasome system (UPS) is the main non-lysosomal proteolytic system involved in the intracellular protein turnover in eukaryotic cells [1]. In normal cells, the proteasome proteolytic activity is responsible for regular cell cycle progression, homeostasis control, and immune surveillance [2,3]. Defects and aberrations in UPS can lead to uncontrolled cell proliferation and tumor development. For these reasons, novel approaches to cancer therapy exist based on proteasome inhibition [1,4,5].
The 20S constitutive proteasome is composed of four stacked rings: the two inner β-rings, which contain the catalytic subunit, and two outer α-rings, whose function is structural [6].
The catalytic subunits β1c, β2c, and β5c are responsible for the caspase-like (C-L), trypsin-like (T-L), and chymotrypsin-like (ChT-L) activities of the proteasome, respectively, with the ChT-L activity being currently considered the primary target for the development of novel anticancer agents [7].
In addition to the constitutive proteasome, vertebrates possess a specialized form of the proteasome, known as immunoproteasome, expressed in monocytes and lymphocytes and responsible for cell-mediated immunity and for the generation of major histocompatibility complex (MHC) class I ligands [8,9]. Under the stimuli of IFN-g and TNF-α, the βc core particles are replaced by β1i, β2i, and β5i subunits. While β2i and β5i perform the same type of activities of the β2c and β5c subunits, differently, the replacement of β1c with β1i brings the caspase-like activity to background levels. Prominent levels of immunoproteasome core particles have been detected in many autoimmune diseases, such as rheumatoid arthritis and systemic lupus erythematosus, and in a panel of inflammatory diseases, such as Crohn's disease, inflammatory bowel disease, and ulcerative diseases. At the same time, overexpression of immunoproteasome core particles has been detected in hematological malignancies, such as multiple myeloma [9].
In this context, several β5i and/or β1i immunoproteasome-selective inhibitors have been identified [10][11][12]. Recently, particular attention has been devoted to the development of noncovalent proteasome inhibitors, which are free from side effects related to irreversible inhibition of a human target [13][14][15][16].
In this research area, our research group has been actively involved in the development of novel 20S proteasome inhibitors [17][18][19][20][21][22]. More recently, considering the variety of therapeutic applications of the development of immunoproteasome-selective inhibitors, we developed novel noncovalent immunoproteasome inhibitors.
Among all the synthesized immunoproteasome inhibitors, we identified a panel of active inhibitors with K i values in the low micromolar or sub-micromolar ranges towards the β5i and/or β1i subunits of immunoproteasome [23]. Compound 1, N-benzyl-2-(2oxopyridin-1(2H)-yl)propanamide (Scheme 1), was shown to be the most potent and selective inhibitor with a K i value of 21 nM against the sole β1i subunit. The catalytic subunits β1c, β2c, and β5c are responsible for the caspase-like (C-L), trypsin-like (T-L), and chymotrypsin-like (ChT-L) activities of the proteasome, respectively, with the ChT-L activity being currently considered the primary target for the development of novel anticancer agents [7].
In addition to the constitutive proteasome, vertebrates possess a specialized form of the proteasome, known as immunoproteasome, expressed in monocytes and lymphocytes and responsible for cell-mediated immunity and for the generation of major histocompatibility complex (MHC) class I ligands [8,9]. Under the stimuli of IFN-g and TNF-α, the βc core particles are replaced by β1i, β2i, and β5i subunits. While β2i and β5i perform the same type of activities of the β2c and β5c subunits, differently, the replacement of β1c with β1i brings the caspase-like activity to background levels. Prominent levels of immunoproteasome core particles have been detected in many autoimmune diseases, such as rheumatoid arthritis and systemic lupus erythematosus, and in a panel of inflammatory diseases, such as Crohn's disease, inflammatory bowel disease, and ulcerative diseases. At the same time, overexpression of immunoproteasome core particles has been detected in hematological malignancies, such as multiple myeloma [9].
In this context, several β5i and/or β1i immunoproteasome-selective inhibitors have been identified [10][11][12]. Recently, particular a ention has been devoted to the development of noncovalent proteasome inhibitors, which are free from side effects related to irreversible inhibition of a human target [13][14][15][16].
In this research area, our research group has been actively involved in the development of novel 20S proteasome inhibitors [17][18][19][20][21][22]. More recently, considering the variety of therapeutic applications of the development of immunoproteasome-selective inhibitors, we developed novel noncovalent immunoproteasome inhibitors.
Among all the synthesized immunoproteasome inhibitors, we identified a panel of active inhibitors with Ki values in the low micromolar or sub-micromolar ranges towards the β5i and/or β1i subunits of immunoproteasome [23]. Compound 1, N-benzyl-2-(2oxopyridin-1(2H)-yl)propanamide (Scheme 1), was shown to be the most potent and selective inhibitor with a Ki value of 21 nM against the sole β1i subunit. The mechanism of the non-covalent inhibition of this compound has been investigated using advanced molecular dynamics methods such as MD-binding (MDB) and Binding Pose MetaDynamics (BPMD). Moreover, an advanced docking method, Induced Fit Docking (IFD), was used for further investigation [24]. These approaches allowed for analyzing the binding mechanism, gaining insights into the ligand entrance pathway, and elucidating the dynamic behavior of the ligand within the binding cavity, providing a dynamic point of view for the definition of the pharmacophore features.
In this study, we employed these collected outcomes to obtain a dynamic pharmacophore model for β1i inhibitors. Additionally, to find new inhibitors of the β5i subunit, we developed a structure-based study using docking and pharmacophore models approaches, starting from the β5i subunit co-crystallized with PR-957 inhibitor (PDB ID: 3UNF).
Both pharmacophore models were used for a virtual screening campaign on commercially available compound collection of 2 million compounds (hereby called Scheme 1. Structure of the selective β1i inhibitor 1. The mechanism of the non-covalent inhibition of this compound has been investigated using advanced molecular dynamics methods such as MD-binding (MDB) and Binding Pose MetaDynamics (BPMD). Moreover, an advanced docking method, Induced Fit Docking (IFD), was used for further investigation [24]. These approaches allowed for analyzing the binding mechanism, gaining insights into the ligand entrance pathway, and elucidating the dynamic behavior of the ligand within the binding cavity, providing a dynamic point of view for the definition of the pharmacophore features.
In this study, we employed these collected outcomes to obtain a dynamic pharmacophore model for β1i inhibitors. Additionally, to find new inhibitors of the β5i subunit, we developed a structure-based study using docking and pharmacophore models approaches, starting from the β5i subunit co-crystallized with PR-957 inhibitor (PDB ID: 3UNF).
Both pharmacophore models were used for a virtual screening campaign on commercially available compound collection of 2 million compounds (hereby called compound library) that allowed the selection of thirty-four potential inhibitors from three commercial databases for inhibition assays against β1i/β5i subunits.

β1i Pharmacophore Modeling
The previous MD study [24] allowed identifying three representative poses (pose 1, pose 2, and pose 3) that elucidated the dynamic behavior and stability of compound 1 within the binding cavity. The interactions pattern of inhibition included two H-bonds between the amide group and Ser21 and Gly47, a π-stacking between the benzyl group and Phe31, and a cation-pi stacking interaction between the 2-pyridone moiety and the epsilon amino group of Lys33. These pose models were used to carry out a pharmacophore approach for the β1i subunit.
For each pose, pharmacophore models were generated using LigandScout [25] PHASE [26]. In LigandScout, pose1 showed three features: a hydrogen bond acceptor (HBA), a hydrogen bond donor (HBD), and a hydrophobic interaction (Figure 1a). Pose2 showed a hydrogen bond acceptor (HBA) and a hydrophobic interaction (Figure 1b). Pose3 was characterized by one hydrogen bond acceptor (HBA), one hydrogen bond donor (HBD), and two hydrophobic interactions (Figure 1c). A merged pharmacophore model [27,28] including all features observed for each pose was constructed and improved by removing the redundant features. The final model is characterized by five features: two hydrogen bond acceptors (HBA), one hydrogen bond donor (HBD), and two hydrophobic interactions ( Figure 1d). The model was retrospectively validated using the active and decoy datasets and performing the receiver operating characteristic (ROC) curve analysis. The ROC curves plot the sensitivity (true positive rate-TPR) against the specificity (false positive rate-FPR). The area under the curve (AUC) is used to assess the accuracy of the model. The grey dashed line represents the random line (AUC = 0.5); the blue curve represents the AUC of the model (AUC must be greater than 0.5). The enrichment factor is the number of active compounds found compared to their presence in the entire database. The model built with LigandScout showed values of AUC 1% = 1 and EF 1% = 15.3 ( Figure 2a). compound library) that allowed the selection of thirty-four potential inhibitors from three commercial databases for inhibition assays against β1i/β5i subunits.

β1i Pharmacophore Modeling
The previous MD study [24] allowed identifying three representative poses (pose 1, pose 2, and pose 3) that elucidated the dynamic behavior and stability of compound 1 within the binding cavity. The interactions pa ern of inhibition included two H-bonds between the amide group and Ser21 and Gly47, a π-stacking between the benzyl group and Phe31, and a cation-pi stacking interaction between the 2-pyridone moiety and the epsilon amino group of Lys33. These pose models were used to carry out a pharmacophore approach for the β1i subunit.
For each pose, pharmacophore models were generated using LigandScout [25] PHASE [26]. In LigandScout, pose1 showed three features: a hydrogen bond acceptor (HBA), a hydrogen bond donor (HBD), and a hydrophobic interaction (Figure 1a). Pose2 showed a hydrogen bond acceptor (HBA) and a hydrophobic interaction (Figure 1b). Pose3 was characterized by one hydrogen bond acceptor (HBA), one hydrogen bond donor (HBD), and two hydrophobic interactions (Figure 1c). A merged pharmacophore model [27,28]    The same procedure was carried out with PHASE. Pose1 showed three features: two hydrogen bond acceptors, one hydrogen bond donor, and two aromatic rings ( Figure 3a). In Pose2, one aromatic ring was recovered ( Figure 3b). Pose3 was characterized by one hydrogen bond acceptor, one hydrogen bond donor, and one aromatic ring ( Figure 3c). As previously, the PHASE-merged pharmacophore model was then obtained with all features observed in the three poses. Two hydrogen bond acceptors, one hydrogen bond donor, and two aromatic rings (Figure 3d). The same procedure was carried out with PHASE. Pose1 showed three features: two hydrogen bond acceptors, one hydrogen bond donor, and two aromatic rings ( Figure 3a). In Pose2, one aromatic ring was recovered ( Figure 3b). Pose3 was characterized by one hydrogen bond acceptor, one hydrogen bond donor, and one aromatic ring ( Figure 3c).    As previously, the PHASE-merged pharmacophore model was then obtained with all features observed in the three poses. Two hydrogen bond acceptors, one hydrogen bond donor, and two aromatic rings (Figure 3d).

β5i Docking and Pharmacophore Modeling
For the β5i subunit, a structure-based study was carried out on the experimental structure of the immunoproteasome (PDB ID: 3UNF). In this experimental structure of immunoproteasome, the β5i subunit is complexed with the PR-957 inhibitor. This ligand was covalently bound to Thr1 of the β5i subunit. The PR-957 inhibitor was structurally modified by breaking the covalent bond and filling in the open valence [29]. The portion of the epoxy ring responsible for the covalent attack on Thr1 was opened (Figure 4a,b). This model was also retrospectively validated and showed values of AUC1% = 0.60 and EF1% = 4.9 ( Figure 2b).

β5i Docking and Pharmacophore Modeling
For the β5i subunit, a structure-based study was carried out on the experimental structure of the immunoproteasome (PDB ID: 3UNF). In this experimental structure of immunoproteasome, the β5i subunit is complexed with the PR-957 inhibitor. This ligand was covalently bound to Thr1 of the β5i subunit. The PR-957 inhibitor was structurally modified by breaking the covalent bond and filling in the open valence [29]. The portion of the epoxy ring responsible for the covalent a ack on Thr1 was opened (Figure 4a  A re-docking protocol was performed using Glide with XP precision to validate the ability to reproduce the experimental position of the co-crystallized ligand. The experimental and docking positions differ only in the epoxy moiety responsible for the covalent interactions. The root means square deviation (RMSD) between the re-docked and crystal conformation of PR-957 was 1.77 Å, less than 2 Å, indicating that Glide XP has a reliable scoring function for this study. The docking protocol was further validated using a β5i validation database of active and decoy compounds. The model showed values of AUC1% = 0.90 and EF1% = 27 (Figure 4c,d).
Moreover, the docking method employed reports more than simple contact scoring by combining structural and energetic information between PR-957 and the active site. The docking post-processing outcomes were used to develop an e-pharmacophore hypothesis. The pharmacophore features with an energetic value <0.5 kcal/mol were A re-docking protocol was performed using Glide with XP precision to validate the ability to reproduce the experimental position of the co-crystallized ligand. The experimental and docking positions differ only in the epoxy moiety responsible for the covalent interactions. The root means square deviation (RMSD) between the re-docked and crystal conformation of PR-957 was 1.77 Å, less than 2 Å, indicating that Glide XP has a reliable scoring function for this study. The docking protocol was further validated using a β5i validation database of active and decoy compounds. The model showed values of AUC 1% = 0.90 and EF 1% = 27 (Figure 4c,d).
Moreover, the docking method employed reports more than simple contact scoring by combining structural and energetic information between PR-957 and the active site. The docking post-processing outcomes were used to develop an e-pharmacophore hypothesis. The pharmacophore features with an energetic value < 0.5 kcal/mol were retained and used to compose the pharmacophore hypothesis. The e-pharmacophore hypothesis includes one hydrogen bond acceptor (A7), two hydrogen bond donors (D8 and D11), and two aromatic rings (R17 and R18). The energy values for the favorable features in the hypothesis are A3 −0.63 kcal/mol, D8 −1.60 kcal/mol, D11 −0.70 kcal/mol, R17 −0.52 kcal/mol, and R18 −1.60 kcal/mol (Figure 5a). The higher absolute value of the feature energy indicates that the ligand atom mapping has a higher interaction energy with amino acids. This five-point pharmacophore model was retrospectively validated. It showed AUC 1% = 0.40 and EF 1% = 2, so we decided to recalculate these values by omitting at least one feature. Figure 5b,c shows the ROC curves for the e-pharmacophore hypothesis with features zero and one omitted. The AUC 1% improved up to 0.84. retained and used to compose the pharmacophore hypothesis. The e-pharmacophore hypothesis includes one hydrogen bond acceptor (A7), two hydrogen bond donors (D8 and D11), and two aromatic rings (R17 and R18). The energy values for the favorable features in the hypothesis are A3 −0.63 kcal/mol, D8 −1.60 kcal/mol, D11 −0.70 kcal/mol, R17 −0.52 kcal/mol, and R18 −1.60 kcal/mol (Figure 5a). The higher absolute value of the feature energy indicates that the ligand atom mapping has a higher interaction energy with amino acids. This five-point pharmacophore model was retrospectively validated. It showed AUC1% = 0.40 and EF1% = 2, so we decided to recalculate these values by omi ing at least one feature. Figure 5b,c shows the ROC curves for the e-pharmacophore hypothesis with features zero and one omi ed. The AUC1% improved up to 0.84.

β1i and β5i Virtual Screening
The pharmacophore models obtained from the β1i dynamic studies were used to perform a virtual screening on the compound library (VSβ1i).
The β5i pharmacophore models were used to perform the virtual screening campaign on the same compound library (VSβ5i) together with a docking analysis. The VSβ5i started with HTVS, SP, and XP docking screening. 1% of the prioritized ligands were retained and then submi ed to pharmacophore matching.

β1i and β5i Virtual Screening
The pharmacophore models obtained from the β1i dynamic studies were used to perform a virtual screening on the compound library (VSβ1i).
The β5i pharmacophore models were used to perform the virtual screening campaign on the same compound library (VSβ5i) together with a docking analysis. The VSβ5i started with HTVS, SP, and XP docking screening. 1% of the prioritized ligands were retained and then submitted to pharmacophore matching.
The hits obtained from VSβ1i and VSβ5i were ranked according to the consensus score described in the Section 3. The consensus score improved the selection and retrieval of molecules shared between the pharmacophore models for β1 and β5. Finally, the hit molecules retrieved after consensus scoring were further evaluated with Induced Fit Docking (IFD) on β1i and β5i subunits. IFD confers more flexibility to the protein side chains, allowing the ligands to adjust and optimize binding interactions within the active sites. Visual inspection of the IFD results helped to enrich the final selection. Thirty-four molecules were selected and purchased for in vitro screening against β1i and β5i subunits using an enzymatic assay.
These hits were evaluated for their ability to inhibit both subunits of i20S by measuring the rate of hydrolysis of the corresponding fluorogenic substrate Ac-Pro-Ala-Leu-AMC for β1i and Suc-Leu-Leu-Val-Tyr-AMC for β5i. MG-132, a reversible inhibitor of the immunoproteasome, was used as a positive control for β1i and β5i.
Initially, compounds were screened against each proteolytic subunit at 100 µM. Only compounds capable of inhibiting the enzymatic activity by more than 50% were characterized in detail using Tian continuous assays. We determined the dissociation constant, i.e., K i , of the non-covalent complex [E I], where K i is also the measure of the binding affinity, since the obtained linear progress curve exhibited a reversible, time-independent inhibition for all the active compounds.
To determine the K i values (Table 1), seven different concentrations were selected for each compound that passed the initial screening, ranging from those that minimally inhibited to those that completely inhibited the activity of the immunoproteasome subunit.

N-(2-(3-(1-phenylethyl)ureido)ethyl)-4-(trifluoromethyl)benzamide
24% 23%    The thirty-four identified hits are characterized by a central amide moiety except for compounds 17 and 31 which are, respectively, a purine derivative and a carboxylate derivative. Compounds 15, 25, and 35 have a central urea moiety. All the hits showed aromatic substituents on both sides of the amide moiety except for compounds 5, 7, 10, 27, and 29, which have just aromatic substituents at one side of the central core. Compounds bearing pyrrole, pyrazole, and/or imidazole moieties showed no inhibition or low inhibition on β1i and β5i subunits, respectively.

N-(3-(benzylamino)-3-oxopropyl)-2-fluoro-6-
As shown in Table 1, some of the tested compounds were active on the immunoproteasome, with K i values in the micromolar range (2, 3, 13, 14, 17, 30, and 34). Interestingly, as observed in the in silico studies, compound 2 inhibits both subunits (K i = 12.53 ± 0.18 µM for β1i and 31.95 ± 0.81 µM for β5i). In Figure 7, the binding modes of compound 2, obtained by IFD, are reported. The central amide moiety establishes an H-bond with the key residue Thr13 in both subunits. The furane ring forms a π-stacking interaction with Phe31 and Tyr169 in β1i and β5i subunits, respectively. An additional H-bond is established with Thr21 and the quinolinone moiety in β5i.

MG-132
0.13 ± 0.01 0.072 ± 0.014 The thirty-four identified hits are characterized by a central amide moiety except for compounds 17 and 31 which are, respectively, a purine derivative and a carboxylate derivative. Compounds 15, 25, and 35 have a central urea moiety. All the hits showed aromatic substituents on both sides of the amide moiety except for compounds 5, 7, 10, 27, and 29, which have just aromatic substituents at one side of the central core. Compounds bearing pyrrole, pyrazole, and/or imidazole moieties showed no inhibition or low inhibition on β1i and β5i subunits, respectively.
As shown in Table 1, some of the tested compounds were active on the immunoproteasome, with Ki values in the micromolar range (2, 3, 13, 14, 17, 30, and 34). Interestingly, as observed in the in silico studies, compound 2 inhibits both subunits (Ki = 12.53 ± 0.18 µM for β1i and 31.95 ± 0.81 µM for β5i). In Figure 7, the binding modes of compound 2, obtained by IFD, are reported. The central amide moiety establishes an Hbond with the key residue Thr13 in both subunits. The furane ring forms a π-stacking interaction with Phe31 and Tyr169 in β1i and β5i subunits, respectively. An additional Hbond is established with Thr21 and the quinolinone moiety in β5i.  Compounds 3 and 30 inhibit β1i subunits. Structurally, compound 3 shows an isosteric substitution of an oxygen atom in place of a sulfur atom and a naphthyl ring in place of a quinolinone moiety. This moiety is replaced by a triazolopyrimidine in compound 30, which is bound to the amide moiety without the thioether or ether chain as in the previous inhibitor. Compounds 13, 14, 17, and 34 inhibit β5i. Interestingly, compound 17 is the only inhibitor that does not bear an amide moiety in the central core.
These newly identified hits could be interesting for further optimization. The amide moiety remains a key pharmacophore that is surrounded by the right combination of aromatic moieties and could lead to more active non-peptide and non-covalent inhibitors (Scheme 2). compound 30, which is bound to the amide moiety without the thioether or ether chain as in the previous inhibitor. Compounds 13, 14, 17, and 34 inhibit β5i. Interestingly, compound 17 is the only inhibitor that does not bear an amide moiety in the central core.
These newly identified hits could be interesting for further optimization. The amide moiety remains a key pharmacophore that is surrounded by the right combination of aromatic moieties and could lead to more active non-peptide and non-covalent inhibitors (Scheme 2).

System and Ligand Preparations
For this study, we selected the amino acid chain K from the 20S mouse immunoproteasome crystallographic structure at 2.90 Å resolution (PDB ID: 3UNF) [30], corresponding to the β5i subunit with the co-crystallized ligand PR-957.
As reported in the previous paper, the reactive residue at the catalytic site, in this case, Thr1, was reconstituted after the removal of the covalent inhibitor by breaking the covalent bond and filling the open valence. Then, the structure and PR-957 were refined and optimized using the Protein Preparation Wizard (Schrödinger New York, 2021-1) and LigPrep (Schrödinger New York, 2021-1). OPLS2005 [31] was used as the force field and Epik (Schrödinger New York, 2021-1) was selected as the ionization tool at pH 7.2 ± 0.2.
A validation dataset containing active and decoy compounds was prepared for virtual screening. Active compounds on β1 and β5 subunits were obtained from the literature [23,[32][33][34][35] for both subunits, while decoys were obtained from the DUD-E database [36] and filtered to remove duplicates. For each active compound in the dataset, 50 decoys were generated. An in-house and two commercial libraries (Asinex and Bioscient) containing >300,000 compounds were prepared in SDF format. The libraries were pre-filtered to remove groups classified as Tox-Alerts, Pan Assays INterference compoundS (PAINS), and Rapid Elimination Of Swills (REOS), which can give false positives due to non-specific reactions during testing [37,38]. As for PR-957, validation, inhouse, and commercial libraries were prepared and optimized using LigPrep.

System and Ligand Preparations
For this study, we selected the amino acid chain K from the 20S mouse immunoproteasome crystallographic structure at 2.90 Å resolution (PDB ID: 3UNF) [30], corresponding to the β5i subunit with the co-crystallized ligand PR-957.
As reported in the previous paper, the reactive residue at the catalytic site, in this case, Thr1, was reconstituted after the removal of the covalent inhibitor by breaking the covalent bond and filling the open valence. Then, the structure and PR-957 were refined and optimized using the Protein Preparation Wizard (Schrödinger New York, 2021-1) and LigPrep (Schrödinger New York, 2021-1). OPLS2005 [31] was used as the force field and Epik (Schrödinger New York, 2021-1) was selected as the ionization tool at pH 7.2 ± 0.2.
A validation dataset containing active and decoy compounds was prepared for virtual screening. Active compounds on β1 and β5 subunits were obtained from the literature [23,[32][33][34][35] for both subunits, while decoys were obtained from the DUD-E database [36] and filtered to remove duplicates. For each active compound in the dataset, 50 decoys were generated. An in-house and two commercial libraries (Asinex and Bioscient) containing >300,000 compounds were prepared in SDF format. The libraries were pre-filtered to remove groups classified as Tox-Alerts, Pan Assays INterference compoundS (PAINS), and Rapid Elimination Of Swills (REOS), which can give false positives due to non-specific reactions during testing [37,38]. As for PR-957, validation, in-house, and commercial libraries were prepared and optimized using LigPrep.

Pharmacophore Modeling Generation
LigandScout and PHASE (Schrödinger, LLC, New York, NY, USA, 2021) were used to generate dynamic pharmacophore models for β1i and structure-based pharmacophore model based on the experimental structure of the β5i subunit.
The models were generated using the following merged pharmacophore creation settings: Feature tolerance scale factor = 1.0; Maximum number of resulting pharmacophores = 10; Number of omitted features for merged pharmacophore = 4.
The structure-based pharmacophore model for the β5i subunit was generated, starting from the PDB coordinate set of β5i subunit-PR-957. The resulting models were validated for their performance in discriminating between the active and decoy molecules using specific databases for β1i and β5i.
A consensus score was calculated. It combined the normalized scores obtained from the Pharmacophore-Fit Score from LigandScout and the Phase-Fit Score from Phase, as shown in Equation (1).

Molecular Docking
The docking study was performed using the Glide docking tool v. 9.0. [39].
The grid box was set on the PR-957 ligand coordinates. The Van der Waals radius was 1 Å, and the partial charge cut-off was 0.15 Å with flexible ligand sampling. Bias sampling torsion penalization for amides with non-planar conformation and Epik state penalties were added to the docking score. To validate the docking protocol and evaluate the ability to reproduce the experimental pose of the co-crystallized ligand, the re-docking test was performed with XP-level precision. The validation of the docking is carried out through the calculation of the RMSD of the cognate ligand, which gives an indication of the deviation between the experimental pose of the ligand and that reproduced by the docking algorithm. The reference value for RMSD is <2.0 for the scoring function.
The Virtual Screening was performed in three different docking steps: High Throughput Virtual Screening (HTVS), Standard Precision (SP), and Extra Precision (XP). Each docking step was run unconstrained.

In Vitro 20S Immunoproteasome Inhibition Assays
Human 20S immunoproteasome, obtained from the human spleen, was purchased from Enzo Life Science. The hydrolysis of the appropriate peptidyl 7-amino-4-methyl coumarin (AMC) substrate was monitored to measure the proteolytic activities of the immunoproteasome. The substrates Suc-Leu-Leu-Val-Tyr-AMC (Bachem, Bubendorf, Switzerland) for β5i and Ac-Pro-Ala-Leu-AMC (Biomol GmbH, Hamburg, Germany) for β1i subunits were used at 50 µM. Fluorescence of the hydrolysis product AMC was measured at 30 • C with a 380 nm excitation filter and a 460 nm emission filter using an Infinite 200 PRO microplate reader (Tecan, Männedorf, Switzerland). A preliminary screening at 100 µM inhibitor concentrations was carried out on the two proteolytic activities of the immunoproteasome.
An equivalent amount of DMSO as a negative control and compound 1 and MG-132 (a reversible immunoproteasome inhibitor) as a positive control were used for β1i and β5i, respectively. Compounds showing at least 50% inhibition at the screening concentration were then progressed into detailed assays. Continuous assays were performed at seven different concentrations ranging from 100 µM to 5 µM and by calculating the dissociation constants K i of the enzyme-inhibitor complex using the Cheng-Prusoff equation, K i = IC 50 /(1 + [S]K m −1 ). Inhibitor solutions were prepared from stocks in DMSO. Each independent assay was duplicated in 96-well plates with a total volume of 200 mL. For the assay on β1i and β5i subunits, a human 20S immunoproteasome was incubated at 30 • C at a final concentration of 4 µg/mL with the inhibitor at seven different concentrations. The reaction buffer comprised Tris·HCl (50 mM, pH 7.4), KCl (25 mM), NaCl (10 mM), MgCl 2 (1 mM), and 0.03% SDS. AMC released from substrate hydrolysis was monitored in a kinetic cycle for 10 min.

Conclusions
In this work, we report a virtual screening strategy built on a dynamic pharmacophore modeling approach onto the β1i subunit and a pharmacophore/docking approach onto the β5i subunit of the immunoproteasome. The outcomes of the computational approaches led to the identification of thirty-four hit compounds selected from three libraries. These hit compounds were tested in in vitro enzymatic assays and seven (2, 3, 13, 15, 17, 30, and 34) were active on the immunoproteasome, with K i values in the micromolar range. Compound 3 and compound 17 are the most active against the β1i subunit and β5i subunit, respectively. Interestingly, Compound 2 showed inhibitory activity on both subunits, and induced fit docking analysis revealed a binding pattern onto β1i with Thr1 and Phe31 that represent new key residues as reported in our previous work. Onto β5i subunit, it interacts with the key residues Thr1, Thr21, and Tyr169. The hit compounds identified represent an interesting starting point for further optimization. Compound 2 showed a dual inhibition and some key pharmacophoric features that could be exploited to identify more potent inhibitors. Moreover, we are planning to investigate the anti-inflammatory profile of β5i-inhibitors and the anticancer activity of β1i-inhibitors, since an up-regulation of the corresponding βi subunits has been detected in inflammatory diseases or hematological malignancies, respectively, in future studies.

Conflicts of Interest:
The authors declare no conflict of interest.