First-in-Class Inhibitors of the Ribosomal Oxygenase MINA53

MINA53 is a JmjC domain 2-oxoglutarate-dependent oxygenase that catalyzes ribosomal hydroxylation and is a target of the oncogenic transcription factor c-MYC. Despite its anticancer target potential, no small-molecule MINA53 inhibitors are reported. Using ribosomal substrate fragments, we developed mass spectrometry assays for MINA53 and the related oxygenase NO66. These assays enabled the identification of 2-(aryl)alkylthio-3,4-dihydro-4-oxoypyrimidine-5-carboxylic acids as potent MINA53 inhibitors, with selectivity over NO66 and other JmjC oxygenases. Crystallographic studies with the JmjC demethylase KDM5B revealed active site binding but without direct metal chelation; however, molecular modeling investigations indicated that the inhibitors bind to MINA53 by directly interacting with the iron cofactor. The MINA53 inhibitors manifest evidence for target engagement and selectivity for MINA53 over KDM4–6. The MINA53 inhibitors show antiproliferative activity with solid cancer lines and sensitize cancer cells to conventional chemotherapy, suggesting that further work investigating their potential in combination therapies is warranted.

p. S14 Figure S7. Dose response inhibition of hydroxylation activity with 2OG inhibitors in MINA53 and NO66 MS hydroxylation assays.
p. S18 Table S6. Data collection and refinement statistics for the 8-KDM5B complex structure.
p. S21 Table S7. Docking assessment results obtained from random conformation cross-docking.
p. S22 Table S8. Total docking score obtained from the docking of the three tautomers with MINA53 and NO66.
p. S26 Figure S15. MD C1 analysis of ligand RMSD with respect to pose C (x axis) and pose A (y axis). Points are colored according to calculated KDE density.
p. S27 Table S10. MM/GBSA calculated binding free energies for selected poses of compound 8 with MINA53 and NO66.
p. S28 Table S11. MM/GBSA calculated binding free energies for selected poses of compound 10 with MINA53 and NO66.
p. S34 Figure S18. MINA53 substrate trapping with inhibitors 9 and 10 in cells. p. S35 Figure S19. Western blot analysis and relative densitometric analysis of the phosphorylated form of H2AX (ser139) on U-87MG glioma cells exposed for p. S36 S4 the indicated times to 1 (negative control) or 9 and 10 (MINA53 inhibitors) at 10 μM.
Purity control by HPLC for compounds 7-10 p. S37    Enzyme activity of MINA53 (M1-V464) or NO66 (Gln116-Asn641) is not inhibited by addition of 1% (v/v) aqueous DMSO, n = 1. NO66 (Gln116-Asn641) appears to be more sensitive to DMSO than MINA53, but tolerates 1% (v/v) DMSO well. Note that the data presented here came from an initial purification of both enzymes; subsequent optimized purifications led to more active material.  Table S4 for assay details. Values are means ± 95% confidence, n = 4.  Room Temperature a See methods for detailed assay procedures. S12 Figure S4. Optimization of MS hydroxylation assay conditions for NO66. Percentage substrate turnover for NO66 (Ala167-Asn641) after 60 min incubation with 50 mM MES, HEPES, Tris, Bis-Tris or phosphate buffers with a range of NaCl or KCl salt concentrations. NO66 was used at a concentration of 300 nM, with 5 µM 2OG, 5 µM peptide substrate, 50 µM ferrous ammonium sulphate (FAS) and 100 µM L-ascorbic acid (L-AA), N=1.      (Met1-Val464) and NO66 (Ser183-Asn641) for 2,4-PDCA (pyridine-2,4-dicarboxylate) and NOG (N-oxalylglycine), two established broad-spectrum inhibitors 2 of 2OG oxygenases, and IOX1. 3 The plot for GSK-J1 is not shown for clarity. Values are means ± 95% confidence, n = 3. Note, although these compounds are relatively broad spectrum 2OG oxygenase inhibitors, their potency against different 2OG oxygenases varies, as manifest in our data for MINA53 and NO66, where 2,4-PDCA and NOG are more potent than IOX-1, which is active against multiple JmjC KDMs. 3 Compound IC50 value NOG    Figure S8. Counter-screen of 2-(aryl)alkylthio-3,4-dihydro-4-oxopyrimidine-5-carboxylic acids for AlphaScreen signal interference. Inhibitors were pre-incubated with the biotinylated-product and a product-specific antibody and added to a streptavidin-coated donor and Protein A-coated acceptor bead mixture. Data are normalised to a DMSO control (100%), with no biotin-substrate as the baseline; n= 3 +/-STD.   , where F obs is the observed structure factor and F calc is the calculated structure factor. § R free is the same as R cryst except calculated with a subset (5%) of data that were excluded from the refinement calculations. † † Engh Huber (1991).

Molecular Modeling
Compounds 8 and 10 can adopt different tautomeric forms in which the exchangeable hydrogen can be either on the pyrimidine's nitrogens (Tauto-1 and Tauto-2, Figure S11), or on the oxygen atom at the pyrimidine C-4 position (Tauto-3, Figure S11). As no experimental evidence is available on the preferred tautomeric forms, quantum mechanical calculations were done to evaluate their relative stabilities. Structure-based (SB) studies using docking simulations were applied to investigate likely ligand-protein interactions. To this, a docking assessment 4 was initially undertaken to select the best scoring function/docking algorithm in Smina 5 and Plants, 6 which are freely available. The best performing program was used to investigate the likely tautomeric forms of 8 and 10 in complex with MINA53 and NO66 using reported crystal structures. Molecular dynamics (MD) simulations were then run on the proposed docked poses in the complexes to investigate the basis of the experimentally observed activities and selectivies.

Docking Assessment
X-ray structures of MINA53 and NO66 were used to evaluate the ability of a docking program to predict the correct geometry of protein-ligand complex (re-docking). The Plants and Smina software were investigated considering three scoring functions for each, Chemplp, Plp, Plp95 and Vina, Vinardo, Ad4 scoring, respectively. Experimental (EC) and randomized (RC) ligand conformation re-docking (RD) and cross-docking (CD) methods were used to assess the software in reproducing experimentally observed binding modes of co-crystallized inhibitors. The root mean square deviation (RMSD) 7 was calculated between the docked ligand conformation and the crystallographically observed one to investigate the docking accuracy (DA) 8,4 . Analysis of calculated docking accuracy percentage (DA%) revealed Plants with the scoring function Plp95 as the program showing the highest DA% value. (Table S7) S22

Molecular Docking of Compounds 8 and 10
To investigate the protein-tautomer interaction stabilities of 8 and 10, a sort docking consensus was carried out as follows: each tautomer was docked, using Plants with the Plp95 scoring function, into MINA53 and NO66 structures available in PDB; the first five docked conformations were evaluated with a docking score (Table S8).
Although the docking score values for the three 8 and 10 tautomers docked into MINA 53 and NO66 are in a narrow range, it is possible to propose a preference of binding of 8Tauto-1 for MINA53, as this combination gives the best docking score value of -117.95. For compound 10, the best selectivity is shown by 10Tauto-1 and MINA53, with a docking score of -120.75.
The best docked conformations for 8 and 10 with MINA53 were very similar ( Figure S13 A and C), with both being predicted to interact via π-π stacking with Tyr167 and chelation of the Mn 2+ ion (substituting for catalytic Fe 2+ ). By contrast with NO66, the best docked conformations were clearly different ( Figure S13 B and D). However, with MINA53 the phenyl ring of 8 and with NO66 the pyrimidine ring of 10 are involved in a π-π stacking with one of the two histidine residues involved in Mn 2+ ligation ( Figure S13).

MD Investigations on the 8/KDM5B System
As from crystallographic data it is not possible to establish which tautomer of derivative 8 preferentially binds to KDM5B, MD simulations were performed. Starting from the reported KDM5B crystal structure (PDB ID: 5FZI) the three tautomers (8Tauto-1, 8Tauto-2, and 8Tauto-3) were modeled using the experimentally observed binding mode as a template; these were labeled as conformations A and C (8Pose-A and 8Pose-C). A total of six different simulations (25 ns each) were performed with all combinations of the poses and tautomers: 8Pose A/Tauto-1, 8Pose A/Tauto-2, 8Pose A/Tauto-3, 8Pose C/Tauto-1, 8Pose C/Tauto-2, 8Pose C/Tauto-3. Crystallographic water molecules in the active site were retained in the MD simulations, as they are involved in establishing briding interactions between the ligand and the Mn +2 ion. Consistent with this, preliminary MD simulation runs without the crystallographic waters manifested large conformational changes in the inhibitor, or loss of binding during the equilibration phase. 8Pose A/Tauto-1 was the only combination showing some stability without the crystallographic waters; further in this case solvent waters were observed to move to occupy the space originally occupied by the deleted active site crystallographic waters (data not shown). Subsequent MD simulations thus were performed with the crystallographically observed water molecules.
The MD analyses of the 8 tautomers bound into KDM5B were accomplished by calculating the ligand root mean square deviation (RMSD) along the trajectories ( Figure S14A). All trajectories returned an RMSD below 2 Å, confirming the stability of the complexes. A kernel density estimation (KDE) analysis on RMSD curves enabled analysis of the distribution of RMSD values ( Figure S14B). The results indicate that 8Pose A/Tauto-1 is the most stable combination; this shows a narrow RMSD distribution centered at ~0.35 Å, the lowest value observed in all MD simulations. This analysis also indicates that for all 8 tautomers resolved conformation C is, in general, less stable than the A. Moreover, it should be noted that, for pose A, all the three tautomers could in principle bind KDM5B, but that 8Tauto-1 shows the best stability.
Trajectories were used to calculate ligand binding free energies using the molecular mechanics energies combined with the generalized Born and surface area continuum solvation (MM/GBSA) method. 9 The calculated free energies (Table S9) were in good agreement with the RMSD analysis: 8Tauto-1 is the 8 preferred tautomer, though with a small ΔΔG value. However, the binding free energy value calculated for 8Pose C/Tauto-1 was slightly lower than that for 8Pose A/Tauto-1. A visual inspection of the analyzed MD trajectory 8Pose C/Tauto-1 showed the ligand switching from alternative pose C to that of A during the first ns ( Figure S15).
The joint interpretation of QM and MD simulation results clearly shows that 8Tauto-1 is the preferred tautomer and that experimental alternative pose A is the preferred for binding with KDM5B.

ΔG (GB) (kcal/mol) Std ΔG (GB) (kcal/mol) Std ΔG (GB) (kcal/mol) Std
Pose  To investigate the inhibitory activity of 8 against MINA53 and NO66, molecular docking simulations coupled with MD simulations were performed as described above for KDM5B. The tautomers of 8 were docked with the MINA53 and NO66 crystal structures (see molecular docking section). Three binding poses were selected for MINA53 and NO66 for MD-based investigations: the best docked (BD, the pose characterized by the lowest score), and the poses presenting the lowest RMSD (LR) with respect to the crystalized alternative conformations (PDB ID: 5FZI) A (LRA) and C (LRC). A total of six 25 ns MD simulations were performed for 8.
Binding poses obtained from docking were complexed with their respective locks (MINA53 key/lock systems: BD/4BXF, LRA/2XDV, LRC/2XDV; NO66 key/lock systems: BD/4CCK, LRA/4CCK, LRC/4CCK) and subjected to MD simulation. The resulting MD trajectories were analyzed and MM/GBSA binding free energies were calculated. The calculated ΔGs (Table S10) are in agreement with the observed activity profiles.
The lowest energy binding modes for 8 on MINA53 (BD, 8Tauto-1) and NO66 (LRC, 8Tauto-3) were visually inspected using the most representative frame of the simulation. A KDE analysis was performed on the first two components of the PCA calculated from ligand cartesian coordinates collected along the MD simulation. The frame with highest density, taken as the most representative, was then visually inspected.

Molecular Dynamics Investigations of 10 into MINA53 and NO66
Calculations to investigate the activity of 10 versus MINA53 and NO66 used the same procedure as done for 8. 10 was docked into the MINA53 and NO66 structures. Three poses for each protein were selected for analysis, based on best docked (BD) and lowest RMSD (LR) poses with respect to the crystallographically observed (PDB ID: 5FZI) alternative conformations A (LRA) and C  (Table S11), i.e. predict a higher affinity for 10 for MINA53 compared to NO66. The best binding mode identified by calculations for 10 with NO66 was LRA (10Tauto-1), while the best with MINA53 was BD (10Tauto-1). The BD and LRA calculated binding free energies with MINA53 had large negative ΔG values principally due to Mn 2+ chelation.
As for 8, for 10 the most representative frame of the simulation was extracted through KDE analysis performed on the first two components of the PCA calculated from ligand cartesian coordinates collected along the MD simulation. The frame with the highest density, taken as the most representative one, was then visually inspected. A RMSD analysis was conducted on both protein backbone and ligand of complexes 8Tauto-1/MINA53 (BD/4BXF), 8Tauto-3/NO66 (LRC/4CCK), 10Tauto-1/MINA53 (BD/4BXF), and 10Tauto-1/NO66 (LRA/4CCK) trajectories ( Figure S16). While the backbone RMSDs displayed a comparable stability for all the investigated systems ( Figure S16A and B), the ligand RMSD timeseries ( Figure S16C) and KDE distributions ( Figure S16D) clearly evidenced the higher stability of compounds 8Tauto-1 and 10Tauto-1 binding modes to MINA53 over those proposed for NO66, in agreement with the MM/GBSA calculations.

Computational Methods
Crystal Data. X-ray crystal structures of ribosomal histidinyl hydroxylases (MINA53 and NO66) were downloaded from Protein Data Bank (PDB), i.e. two structures for MINA53, 2XDV (in complex with NOG) and 4BXF (in complex with 2OG), and one for NO66, 4CCK (in complex with NOG). 4BXF and 4CCK have 2 and 4 protein chains, respectively, each complexed with the proper ligand. These multi-chains PDB were loaded through UCSF Chimera v1.14, visually inspected, then split into a single chain and treated as single protein. Missing residues were built through homology modeling using Modeller v.10.1. Crystallographically observed solvent molecules were deleted, and missing hydrogens were added; residues side chains were ionized (protonated/deprotonated) at physiological pH (7.4). Geometric optimization was done by means the OpenMM tool using an in-house Python script. The minimized complexes were SB aligned on the alpha carbon atoms by means of UCSF Chimera MatchMaker module using the PDB: 2xdv structure as a reference.
Tautomeric structures. The three tautomeric forms of the 8 and 10 were built using Marvin Sketch software, then processed depending to the computational approach. For docking studies, the hydrogens at pH 7.4 were added and a random three-dimensional conformation was produced with Marvin sketch. A genetic algorithm-based conformational search using Obonformer was then used to optimize this 3-D conformation.

Molecular docking settings. Plants settings.
The docking of proteins with ligands was performed using Plants v1.2 version with three different scoring functions at the default speed (SPEED1). The docking tools generated 10 conformations for each docked ligand. The docking binding site was centered at the molecule's mean center and enlarged to a radius of 12 Å. Docking was performed using three different scoring functions: Chemplp, Plp and Plp95. Autodock Vina settings. Intermediary steps, such as pdbqt files for protein and ligand preparation, were completed using obabel software. The box was created by mean the Vina graphical interface module in UCSF Chimera. For each calculation, ten poses were obtained and ranked according to the scoringfunctions.

QM Calculations.
Starting from the three tautomeric forms of 8 and 10 represented as SMILES, a first three-dimensional geometry was obtained; the six structures were then subjected to a genetic algorithm-based conformational search returning the 30 lowest energy conformers from the last generation, optimized for RMSD diversity. This first manipulation was accomplished using OpenBabel (v. 2.4.0). 10 The resulting conformer structures were optimized at the semiempirical quantum mechanical level of theory PM6. 11 Each tautomer's lowest energy conformer was then subjected to further geometric optimization at a higher level of theory using the Minnesota DFT functional M06 12 equipped with the Pople split-valence triple-zeta basis set with polarization and diffuse functions on heavy atoms and hydrogens, 6-311++G** 13 Harmonic frequencies calculation were accomplished at the end of the geometry optimizations to evaluate zero-point vibrational energies, thermochemical corrections and to verify the nature of the stationary point. PCM solvation model was used for all QM calculations: PM6 and DFT ones. QM calculations were performed using the General Atomic and Molecular Electronic Structure System (GAMESS) software. 14 MD Simulations. Protein structures were prepared as reported in crystal structure section: solvent molecules and crystallization residues were removed; gaps were filled by means of Modeller 15 and residue protonation states were defined using PROPKA3. 16 Ligand general amber force field 2 (GAFF2) 17 parameters were calculated by means of antechamber 18 at the semiempirical AM1-BCC level of theory. 19 The ff14SB force field 20 was used for proteins, while Li/Merz parameters were used for ions. 2121, 22 The complexes were solvated in an orthorhombic box using the fourpoint OPC water model 23 setting to 1.5 nm the distance between the box boundaries and the solute. Na + and Clions were added to neutralize the system. Complete parameters and topology files were obtained using tLeap. All parametrization procedures were accomplished using the AmberTools 20 suite. 24 The MD integration timestep was set to 2 fs, the hydrogen bond length was constrained using the LINCS algorithm, 25 long range electrostatics were treated using the particle mesh Ewald (PME) 26 method, temperature was imposed using a velocity-rescale thermostat. 27 MD simulations started with a 50000 step minimization applying harmonic position restraints to backbone atoms (k=1000 kJ/mol/nm 2 ) followed by a second unconstrained minimization. Systems were then heated gradually to the target temperature (300.15 K) during 1 ns constrained (k=100 kJ/mol/nm 2 ) simulation in an NVT ensemble. A second constrained equilibration (k=10 kJ/mol/ nm 2 ) occurred in an NPT ensemble where pressure control (1 atm) is achieved by means of a Berendsen barostat. 28 Finally, 25 ns production simulations in an NPT ensemble using a Parrinello-Rahman barostat 29 to control pressure (1 atm). MD simulations were performed with GROMACS (v. 2019.6). 30 Each simulation analysis was performed after the system had reached equilibration (i.e. stabilization of the backbone RMSD), so the first 5 ns of the trajectories were discarded. Binding free energy calculations were performed with the endpoint generalized Born and surface area continuum solvation (MM/GBSA) method 9 by means of the MMPBSA.py 31 Python script delivered with the AmberTools 18 suite. The entropy contribution was neglected, since the benefits coming from its calculation remain controversial and normal model calculations are computationally expensive. 32 The GB method used for MM/GBSA calculations was GB-Neck2. 33 This GB variant has the same form as the GBn method, 34 but uses a different parameter set. MD analysis were performed using the MDTraj 35 Python library. Statistical analysis such as KDE and PCA were performed with scipy 36 and the sckit-learn 37 Python library, respectively.   . Anti-HA IPs were immunoblotted for RPL27A (17 kDa) and HA. β-Actin was used as loading control. Evidence for MINA53 substrate trapping is observed in all cell lines following incubation with 9 and 10, but not 1. Note that the H179A mutant served as a control likely lacking or having reduced substrate binding capability. Figure S19. Western blot analysis and relative densitometric analysis of the phosphorylated form of H2AX (ser139) on U-87MG glioma cells exposed for the indicated times to 1 (negative control), 9 and 10 (MINA53 inhibitors) at 10 μM. Doxorubicin was used as a positive control (0.5 μM, 8 h). GAPDH expression was used as a control loading. Results are expressed as fold changes over untreated cells.