Identification of PRMT5 inhibitors with novel scaffold structures through virtual screening and biological evaluations

Protein arginine methyltransferase 5 (PRMT5), an important member in PRMT family, has been validated as a promising anticancer target. In this study, through the combination of virtual screening and biological experiments, we have identified two PRMT5 inhibitors with novel scaffold structures. Among them, compound Y2431 showed moderate activity with IC50 value of 10.09 μM and displayed good selectivity against other methyltransferases. The molecular docking analysis and molecular dynamics (MD) simulations suggested that the compound occupied the substrate-arginine binding site. Furthermore, Y2431 exhibited anti-proliferative activity to leukemia cells by inducing cell cycle arrest. Overall, the hit compound could provide a novel scaffold for further optimization of small-molecule PRMT5 inhibitors.


Introduction
Protein arginine methylation is an important posttranslational modification catalyzed by protein arginine methyltransferases (PRMTs). By using S-adenosyl-L-methionine (SAM) as methyl donor, PRMTs catalyze the methylation of arginine residues in a variety of proteins [1]. To date, there are nine members of PRMTs that have been reported in mammalian cells. Based on the difference in the catalytic products, they are divided into three forms: Type I (PRMT1, 2, 3, 4, 6, and 8), Type II (PRMT5, 9) and Type III (PRMT7) [2]. Both of type I PRMTS and type II PRMTs could catalyze the formation of monomethylarginines (MMA), type I PRMTs could also catalyze the formation of asymmetric dimethylarginines (ADMA), while type II PRMTs catalyze the production of symmetric dimethylarginines (SDMA). As the only III type PRMTs, PRMT7 solely generates MMA [3].
Due to the important role of PRMT5 in cancer physiology and pathology, it is urgent to develop effective inhibitors. Up to now, great efforts have been made to develop PRMT5 inhibitors (Fig. 1). According to the structure characteristics, they can be classified into SAM analogues (e.g., Sinefungin [12], DS-437 [13] and LLY-283 [14]) and non-SAM analogues (e.g., DC_Y134 [15], GSK-3326595 [16], DC-C01 [17] and CMP5 [18] As one of computer aided drug design (CADD) methods, virtual screening has been widely used in academia and industry to discover lead compounds [19]. Previously, we identified compound DC_Y134 as a PRMT5 inhibitor by virtual screening targeting the binding site composed of the SAM-binding pocket and the substrate-binding pocket [15]. Here, to obtain compounds with better selectivity, we fine-tuned our strategies to focus on the arginine-binding pocket, and discovered a new PRMT5 inhibitor which showed good selectivity against other methyltransferases and anti-proliferative activity in several PRMT5-associated cancer cell lines.

Protein structure preparation
In this study, the crystal structure of the human PRMT5:MEP50 in complex with substrate H4 peptide was chosen as the receptor (PDB ID: 4GQB) [12]. After deleting the water molecules and partner MEP50 protein, the Protein Preparation Wizard Workflow integrated in the Maestro 9.0 graphical user interface was used to prepare the structure of PRMT5. The pH value was defined as 7.0 ± 2.0 and residues within 6 Å around substrate arginine in H4 peptide were defined as binding sites. Other parameters were set as the default.

Ligand database preparation
Two hundred twelve thousand two hundred fifty-five small molecule compounds from Specs database (https:// specs. net/) was chosen for virtual screening. To select drug candidates that could be better used by organisms, Lipinski's rule of five [20] was applied to filter the database, and "pan-assay interference compounds (PAINS)" [21] were removed with Pipeline Pilot, version 7.5 (Accelrys Software Inc., San Diego, CA, USA). The generating 18,0000 compounds were optimized to generate quality 3D structures with the precise chirality adopting the default settings by LigPrep provided in the Schrödinger suite [22] and the protonation state was generated by Epik with a pH of 7.0 ± 2.0 [23].

Virtual screening
The crystal structure of PRMT5 with H4 peptide (PDB ID: 4GQB) was adopted as the target for the following virtual screening procedures, and the flowchart is shown in Fig. 2a. First of all, the generating 18,000 compounds were screened through Glide [24] SP and XP modes, leading to the topranked 300 compounds. Then, to ensure the structural diversity of candidates, the remaining compounds were classified to 30 groups by "Clustering Molecules" protocol in Pipeline Pilot, version 7.5. Considering the interaction of molecules with key residues such as Q309, S310, F327, E444 or F580, one to three compounds were selected from each group. Finally, 80 compounds were selected for subsequent biological evaluations.

In vitro PRMT5 enzyme inhibition and selectivity evaluation
In vitro enzyme inhibitory activities and selectivity assay of compounds were tested by Shanghai Chempartner Co., Ltd. using 3 H-labeled radioactive methylation and AlphaLISA assay, which were adopted by the same protocol as described previously [17,25,26]. The PRMT5 protein was purchased from BPS (Cat. No. 51045) for enzyme activity inhibition test and adding SAH only is control assays. The following formula (Eq. (1)) [26] was used to calculate the inhibition rate: (1) % inhibition = (max signal − compound signal)∕(max signal − min signal) × 100 The reaction value of enzyme and substrate is the maximum signal; the reaction value of substrate only is the min signal.
The IC 50 values were analyzed in GraphPad Prism 8.0 software with Eq. (2) [26] and the compound with lowest IC 50 was chosen for selective evaluating: In selectivity assay, five other different types of methyltransferases were chosen for further investigation, including representatives of type I arginine methyltransferase (PRMT1, PRMT4), DNA methyltransferase 1 (DNMT1), enhancer of zeste homolog 2 (EZH2), histone methyltransferase G9a and lysine-specific demethylase 1 (LSD1). 80 compounds were diluted to 50 μM and 100 μM for selectivity evaluation. Except for DNMT1 and PRMT5, which were detected by radioactive methylation, all the enzymes were analyzed by AlphaLISA and Eq. (1) was carried out to calculate the inhibition value.

Molecular dynamics simulation
Molecular dynamics simulations (MD) was performed using Amber14 software (University of California, San Francisco, USA) [27,28] package. The small molecules were parametrized using the generalized AMBER force field (GAFF) [29] and AM1-BCC charge model [30] by antechamber [31]. The tleap module was employed to generate the topology of each system. Before MD simulations, all the systems were solvated into the TIP3P [32] water box with a 10 Å buffer distance between the box wall and the nearest solute atoms. Na + ions were added to neutralize the simulation systems. The three systems were subsequently minimized with the Amber ff14SB force field. After minimization, each system was heated from 0 to 300 K within 100 ps, then 1 ns equilibration was carried out. Finally, 100 ns MD simulations on each system were then carried out using the pmemd. cuda module with constant temperature (NPT), and periodic boundary condition. The CPPTRAJ [33] program was applied for all MD trajectory analyses and the VMD software was used to perform visual observation [34].

MM/GBSA calculations
The MM/GBSA method [35] in Amber14 was implemented to calculate binding free energies of PRMT5 complex. From the molecular dynamics trajectory, 100 snapshots of the last equil 10 ns were carried out on the MD trajectory. The equation (Eq. TΔS is entropic contribution to ΔG at absolute temperature T.

Free energy decomposition of inhibitor-residue pairs
Determining the energy contribution of inhibitor-residue pairs to estimate the components in the interactions between the residue and ligand is important [36]. In the current study, residues within the range of PTMT5 substrate active site 12 Å are selected for calculation, and pairwise decomp with 1-4 terms (1-4 EEL and 1-4 VDW) were added to internal potential terms.

Anti-proliferative assay
According to the results of enzymatic activity tests, we discovered that Y2431 possesses predominant inhibitory activity against PRMT5 in vitro. To further explore whether it would affect the proliferation of human Leukemia cell line, MV4-11 and K562 were adopted to assess the inhibitory effect of Y2431 on cancer cells. The cells survival rate was ΔG solv = ΔG pol + ΔG nonpol measured using the WST-1 Cell Proliferation Assay Kit. Two cells were seeded in 96-well plates and adding RPMI 1640 medium with 10% fetal bovine serum and then treated with Y2431 at different concentrations or DMSO control for 8 days in a cell culture incubator at 37 °C with a supply of 5% CO 2 . WST-1 solution was added according to the user guide. Three repetition and replicates were done in each assay. The data of IC 50 were obtained by GraphPad Prism 8.0.

Cell-cycle assays
K562 cell line was plated in 6-well plates and then cultured with Y2431 at a diverse mix of concentrations or DMSO. After 48 h, cells collected were re-suspended in 70% ethanol and fixed overnight at 4 °C. After washing with PBS, the cells were incubated with Propidium Iodide/RNaseA at room temperature for 30 min and detected using flow cytometer (ACEA NovoCyte).

Virtual screening against substrate-binding pocket
To date, a number of crystal structures of human PRMT5 have been determined, which indicate that PRMT5 possess two ligand binding sites of cofactor SAM and substrate arginine. Previously, we identified several PRMT5 inhibitors by virtual screening targeting the pocket that contains the two ligand binding sites [15]. In contrast to the binding site of SAM, the protein sequence of substrate-binding site varies in different PRMTs. Consequently, to obtain more PRMT5 small-molecule inhibitors with good selectivity, an optimized approach that focused on substrate-binding site was carried out. Here, the crystal structure of PRMT5 with H4 peptide (PDB ID: 4GQB) was adopted as the target for virtual screening procedures, and the flowchart is shown in Fig. 2a. After a series of virtual calculation and screening, 80 compounds were selected for the following biochemical evaluation.

Enzyme inhibition and selectivity assay
The 80 compounds selected by virtual screening were tested for PRMT5 inhibitory activity by 3 H-labeled radioactive methylation assay. Among these candidates, two compounds (Y2431, Y0102) inhibited PRMT5 activity with IC 50 values of 10.09 µM and 31.53 µM, respectively (Fig. 2b). Subsequently, we tested the selectivity of the compound Y2431 against other methyltransferases including PRMT1, PRMT4, DNMT1 EZH2, G9a and LSD1. As shown in Table 1, Y2431 exhibited much lower inhibition against other methyltransferase, demonstrating that this compound is a selective inhibitor for PRMT5.

Molecular dynamics simulation
To investigate the dynamics and energy properties of the two PRMT5/SAM/inhibitor complexes and PRMT5/SAM system, 100 ns molecular dynamics simulations were carried out. The Fig. 3(a-b) showed dynamic stability of the three systems during simulations. We monitored root mean square deviation (RMSD) values of protein backbone and inhibitors relative to the initial structure along the entire MD trajectories, which suggested that all of the three systems reached equilibrium states after 80 ns. And the RMSD of inhibitor Y2431 retain lower value than that of Y0102, which indicated that Y2431 is more stable in the pocket. Although its benzene ring was reversed during the 50-80 ns  RMSD (a, b), RMSF (c) and Rg (d) analysis over the entire 100 ns MD trajectory of PRMT5 (red) and PRMT5 complexes with compound Y2431 (orange lines) and Y0102 (blue lines), respectively  period, resulting in a large change in RMSD, it was still in a stable binding state and returned to its initial conformation after 80 ns. To investigate the dynamic features, root mean square fluctuation (RMSF) values were also calculated, which revealed that the fluctuations of residues are similar in the three systems, especially for the residues near the substrate binding site, indicating that binding of the two inhibitors has little effect on the flexibility of residues in PRMT5 (Fig. 3b). Since the initial structures for MD simulation didn't contain the partner MEP50 protein, the residue fluctuations of the two regions binding to MEP50 (residues 54-74 and 158-180) significantly increased [37]. Radius of gyration (Rg) is a parameter representing the compactness of protein.
To estimate the effects of Y2431 and Y0102 bindings on protein, Rg values for the three systems were calculated (Fig. 3c), which showed that the difference between the three protein systems was only 0.5 Å after 90 ns, indicating that the looseness of the protein does not change much and the conformation is quite stable, which is consistent with the RMSD analysis results.

Binding free energy analysis
Based on the results of MD simulations, MM/GBSA method was utilized to calculate the binding free energy of the two compounds with PRMT5. As shown in Table 2, compound Y2431 binds to PRMT5 with a ΔG value of − 44.03 kcal/ mol, while Y0102 only got a ΔG value of − 27.57 kcal/mol, indicating that Y2431 has a higher binding affinity with PRMT5, which is consistent with the experimental data of bioassays. Besides, the energy component of van der Waals (ΔE vdw ) is favorable for the binding of both compounds, while the polar solvation is not beneficial for the binding.
To obtain a detailed view of the protein-inhibitor binding, the binding energy of the residues in the substrate binding pocket was decomposed to investigate the partial energy contributions of each residue upon substrate binding. The residues that have a contribution of < − 1 kcal/mol were summarized in Table 3, suggesting that some residues that play important roles in substrate H4 peptide binding [12] are also critical for inhibitors binding, especially for Y2431.
Subsequently, we compared the energy contribution components of the key residues that are important for binding of H4 peptide [12]. As shown in Fig. 4a, compared to Y0102, key residues contain Y304, Y307, F327 and F580 formed stronger van der Waals interactions with Y2431. Meanwhile, the electrostatic interactions between these key residues and Y2431 are also stronger, especially for the residue E435 and E444 that form the catalytic active site of PRMT5 (Fig. 4b), providing a molecular basis to explain the inhibitory activity of Y2431 towards PRMT5.

Binding mode analysis of PRMT5/Y2431 complex
In order to disclose the binding modes of Y2431 with PRMT5, we extracted the structure of last frame (100 ns) from the MD trajectory to analyze molecular interaction pattern. As shown in Fig. 5a, Y2431 occupies the binging site of residue Arg3 in the substrate H4 peptide, while the methyl group of Y0102 almost deviates from the arginine pocket, providing further proof for the results of biochemical experiments. Furthermore, Y2431 fits into the hydrophobic pocket which is made up of residues W579, F580, T323, Q309, V307, V326, Q322, S578 (Fig. 5b-c). In particular, hydrogen bond is formed between the guanidine group and residue E444 and E435, whose occupancy in simulation is 87.47% and 68.63%. In addition, E435 and E444 also form

Anticancer activity and cycle arrest analysis
Since PRMT5 is a potential target of leukemia, two human leukemia cell lines, MV4-11 and K562 were adopted to test the inhibitory effect of Y2431 on cancer cells. As shown in Fig. 6a-b, Y2431 could inhibit the proliferation of MV4-11 and K562 cells with IC 50 values of 10.27 μM and 10.12 μM at 8 days, respectively, suggesting modest inhibitory effect in cells of this compound. To further investigate the antiproliferation mechanism of Y2431 towards K562 cells, flow cytometry was carried out to assess the effect on cell cycle. Figure 6c indicated that Y2431 arrested the cell cycle in the G2/M phase in a concentration-dependent manner upon treatment for 48 h, which revealed the cytotoxic activity of Y2431 against PRMT5-releated leukemia cells.

Conclusion
In conclusion, as one of the most promising biological anticancer targets in the PRMT family, PRMT5 has attracted more and more attention. In this study, through virtual screening and biological evaluations, compound Y2431 has been discovered as a novel PRMT5 inhibitor with IC 50 values of 10.09 μM. Notably, Y2431 displayed good selectivity against other methyltransferases. Furthermore, MD simulations and binding free energy analysis revealed that Y2431 has a stronger binding affinity with PRMT5 than compound Y0102, which was consistent with inhibitory data. The binding free energy decompositions and binding mode analysis suggested that Y2431 occupied the substrate-binding pocket. Notably, the electrostatic interactions between positively charged guanidine group and the charged negatively binding pocket play a significant role in the binding of Y2431. E444, which is one of methyltransferase catalytic residues, interacts with Y2431 through electrostatic attraction. These results provided a molecular basis for explaining the inhibitory activity of Y2431 toward PRMT5. In addition, Y2431 also has good cell membrane permeability and can inhibit the proliferation of several PRMT5-related cancer cells. In summary, this study provided a reliable structure-based virtual screening method for identifying novel PRMT5 inhibitors, as well as a novel chemical scaffold for further optimization.