The human soluble epoxide hydrolase (hsEH, EC 3.3.2.10) is a bifunctional homodimeric enzyme weighing 62 kDa. It is a member of the epoxide hydrolase superfamily and is composed of two independently folded domains: N-terminal domain (NTD) with phosphatase activity and C-terminal domain (CTD) with hydrolase activity. hsEH is encoded by the EPHX2 gene and can be found in both cytosol and peroxisomes[1–3]. The enzyme shows the highest activity in the liver; however, it is distributed in the majority of tissues, including the kidney, intestine, adipocytes, and vessels[4–6]. The clinical relevance of hsEH is linked with a C-terminus motif which performs the hydrolysis of arachidonic acid epoxides and other natural epoxy-fatty acids, especially the transformation of four regioisomers of epoxyeicosatrienoic acids (EETs): 5,6-, 8,9-, 11,12-, and 14,15-EETs to the corresponding dihydroxyeicosatrienoic acids (DHETs)[7, 8]. The conversion of EETs into the corresponding diols performed by hsEH generates non-bioactive molecules, thus the enzyme inhibition would be expected to enhance the EETs bioavailability and their beneficial properties[9, 10]. EETs have been found to modulate the activities of many molecular targets and signalling pathways. They have been shown to possess vasodilator and anti-inflammation properties, as well as act as mediators of autocrine and paracrine[9, 11, 12]. These characteristics are why many studies implicate their significant meaning in several cardiovascular disorders, e.g. blood pressure, hypertension, and atherosclerosis[13, 14]. EETs have also been found in the brain, where they regulate cerebral blood flow and protect against ischemic brain injury[15, 16]. Several studies were carried out to investigate the effect of the sEH gene deletion[17–19]. Samokhvalov et al. reported that such a deletion provides a protective effect against lipopolysaccharide (LPS)-induced cardiotoxicity by maintaining mitochondrial function[20]. In their research associated with Parkinson's disease, Ren et al. reported that sEH gene deletion or enzyme inhibition may protect against MPTP (1-methyl-4-phenyl-1,2,3,6-tetrahydropyridine) induced neurotoxicity[21]. EETs from the cytochrome P450/soluble epoxide hydrolase pathway, together with the prostaglandins derived from the cyclooxygenase (COX) pathway, have also been found to be important eicosanoids that may regulate angiogenesis and tumorigenesis[22, 23]. Therefore, a great effort has been made to develop inhibition strategies towards sEH.
Some of the earliest sEH inhibitors were chalcone oxide, which can form a covalent intermediate with the enzyme. However, they were slowly hydrolysed and then released from the active site due to time-dependent inhibition and they were rather unstable, making in vivo usage limited[9, 24–26]. Urea, carbamate, and amide compounds with the hydrophobic group substitution were the first successful generation of hsEH inhibitors. They were identified as significantly more potent inhibitors and still represent the most popular class of hsEH inhibitors[26, 27]. These compounds can mimic epoxide binding and are not released from the active site. However, they have rather poor physical properties which include, low solubilities and high melting points, resulting in inappropriate bioavailability[26, 28, 29]. To improve their properties, while maintaining a high potency, a polar functional group at a proper distance from the urea or amide moiety was added. Soon after, the second generation of inhibitors, based on the 12-(3-adamantan-1-ylureido)dodecanoic acid (AUDA) and 1-adamantan-1-yl-3-{5-[2-(ethoxyethoxy)ethoxy]pentyl} urea (AEPU), were developed[9, 26, 28, 30–32]. But, these compounds remain flawed, due to their rapid in vivo metabolism and limited pharmacokinetics parameters.
Bruce Hammock’s group proposed the so-called “Pharmacophore Model”, which has been widely used to develop new sEH inhibitors[33]. The “Pharmacophore Model” is as follows: the “primary pharmacophore” (P1) may consist of carbamate, amide, or urea (preferably) and bear a bulky, hydrophobic substituent, whereas the “secondary pharmacophore” (P2) should consist of a polar group, located at a distance of five or six carbons (7.5 Å) from the urea carbonyl group. The “secondary pharmacophore” could be any group that is capable of serving as a proton acceptor (e.g., amide, ketone, alcohol, etc.); however, moieties that are too polar offered less potent inhibition. The model may also show an optional “tertiary pharmacophore” (P3), located at least 12 Å from the urea carbonyl group. The hydrophobic linkers, L1 and L2, should connect the polar moieties from the primary, secondary and tertiary pharmacophores[33–35].
Several pharmaceutical companies and academic institutions and pharmaceutical companies, such as the University of California-Davis, Columbia University, Chinese Academy of Sciences, Arête, Boehringer Ingelheim, Taisho, Takeda, Merck, Roche, etc., have entered the arena of developing sEH inhibitors. Thousands of various hsEH inhibitors were synthesised over the last three decades, but only individual compounds have been qualified into clinical trials. Such an example is the 1-(1-acetyl-piperidin-4-yl)-3-adamantan-1-yl-urea (AR9281) inhibitor from Arête Therapeutics[36]. The compound has completed phase I of the clinical trial successfully[37] but has failed the IIa clinical trial due to lack of efficiency[38]. Esters and salts of AUDA, commonly used as an experimental sEH reference inhibitor, are good inhibitors of sEH but their clinical use has been restricted due to metabolic instability and limited solubility in water and many organic solvents[33, 34, 39]. Also, GSK2256294 ((1R,3S)-N-(4-cyano-2-(trifluoromethyl)benzyl)-3-((4-methyl-6-(methylamino)-1,3,5-triazin-2-yl)amino)cyclohexanecarboxamide) inhibitor has shown promising effects in two phase I studies, however, there is no information on the results of phase II clinical trials[40].
Various research groups began using the fragment-based drug discovery with biophysical or biochemical assays and a high-concentration of fragments to speed up the research of new inhibitors. For instance, Amano et al. conducted fragment screening against sEH using X-ray crystallography together with biochemical assays[41, 42]. They were able to identify numerous fragments that inhibit sEH and reveal their binding modes. Using the optimisation techniques, they determined structures of fragment-enzyme complexes and discovered a novel scaffold for sEH inhibitors. In 2015, Öster et al. also presented a fragment-based crystallisation approach, where a fragment-based lead generation was executed with a limited high-throughput screening against the sEH[43]. The obtained complexes extended the structural knowledge about possible dual conformations or multiple binding sites. In 2016, Xue et al. used a similar approach to completely map the active site pocket to gain a deeper knowledge of ligand interactions and to guide the rational structure-based drug design[44].
In the last decade, research groups also applied various approaches to discover new hsEH inhibitors based primarily on in silico studies. In 2011, Tanaka et al. identified almost 70 new sEH inhibitors using a docking-based virtual screening[45]. They analysed different physicochemical properties of the hit compounds, e.g., size and lipophilicity, to find a compound with high ligand efficiency. Xing et al. used a similar approach to construct combinatorial libraries of potential inhibitors based on a benzoxazole template[46]. The libraries were further evaluated by docking and by analysing the capacity to make hydrogen bonds. Research was then expanded into exploring the structure-activity relationship by the amide bond formation. In 2012, Moser et al. developed a pharmacophore model based on co-crystallised inhibitors of hsEH[47]. They used a combination of Protein Ligand Interaction Fingerprints with the QueryGenerator tool to model the shape of the binding pocket indicated as excluded volumes. Then, they proposed a new pharmacophore model (with a special focus on the catalytic triad, as described by Morisseau and Hammock[26]) by examination of the clusters of receptor pharmacophore annotation points within the binding pocket. This pharmacophore model was later used by Moser et al. to find dual 5-lipoxygenase (5-LO)/soluble epoxide hydrolase inhibitors[48]. In 2016, Waltenberger et al. created multiple structure- and ligand-based pharmacophore models to cover different putative binding modes of the ligands and to cover the active space more thoroughly[49]. In 2018, Tripathi et al. applied a 3D-QSAR approach using sEH inhibitors from literature (with a wide range of activity and structural diversity) to design a pharmacophore model[50]. A new pharmacophore hypothesis was presented by Bhagwati et al. in recent research, based on the molecular modelling studies, which include ligand-based pharmacophore modelling, virtual screening, molecular docking, and all-atom molecular dynamics simulations[51].
The use of in silico methods has undoubtedly helped propose new hsEH inhibitors, however, the research was mainly focused on exploring the main binding pocket and developing new pharmacophore models based on the structures of known inhibitors binding in this area. Exploration of other parts of the hsEH enzyme, e.g., accessible after the conformation changes that could host additional functional groups, was omitted. For the exploration of new, potential binding sites in the protein, molecular dynamics simulation with (co)solvent molecules could be applied. Such an approach, called mixed-solvent molecular dynamics (MixMD) or cosolvent-based molecular dynamics (CMD) simulations, mixes small organic molecules with a water solvent to stimulate dynamic protein motions and induce partial conformational changes of binding pocket residues appropriate for the binding of diverse ligands[52–56]. Based on such abundant information, it is possible to identify the ligand binding modes, preferable binding regions (often called hot-spots), and present structural characteristics of a given place to build a pharmacophore model[57–59]. But, the MixMD has several issues. There is no universal procedure or methodology for performing this type of simulation. Various research groups select different types of cosolvent molecules and the concentration of cosolvent in the systems can vary significantly[52, 53, 67–69, 57, 60–66]. One problem is the selection of cosolvents to properly map all of the possible protein-ligand interactions and regions, such as hydrogen bonding (donors and acceptors), hydrophobic regions, and positive and negative ionic interactions[60, 70, 71]. Moreover, the water:cosolvent ratio could also affect the results, since a higher concentration of particular molecules could disturb the protein’s structure and the cosolvent molecules could aggregate[54, 72–74]. MixMD simulations were usually used for relatively small proteins, such as thermolysin[52, 63, 75–77], p53[76, 78, 79], HIV-1 protease[80–82], hen egg-white lysozyme[83], Hsp90[81], and several kinases[72, 84] to identify potential binding regions on the protein’s surface exclusively. Therefore, an analysis for proteins with buried active sites (as in the case of hsEH) remains a challenge. For a long time, no bioinformatics tools were available that would allow accurate tracking of the behaviour of solvent molecules inside the protein’s core and its further analysis from an intramolecular voids the point of view during the molecular dynamics simulations. The situation changed when new and easy-to-use AQUA-DUCT (AQ) software appeared[85]. This tool was primarily used to track the molecular trajectories of ligands (solvent molecules) throughout molecular dynamics simulations. The new version, released in 2019, also allows the identification of structurally important residues and regions of macromolecules responsible e.g., for ligand trapping, which can vary depending on the tracked ligand[86]. Thus, such functionality can be used for in silico drug design and for identifying new potential binding sites for hsEH inhibitors.
In this paper, we have addressed two main questions: i) is it possible to expand the targeted space of the hsEH interior to propose unique binding sites for new classes of inhibitors? and ii) is it possible to identify amino acids that, so far, were not involved with inhibitors’ interactions and could accept functional groups to increase the inhibitors solubility? To address these questions, we used the combination of small molecules tracking and classical and mixed-solvent molecular dynamics simulations to analyse the intramolecular voids of the hsEH. Such an approach allows us to analyse the local distribution of particular solvent molecules employed as molecular probes. The densest regions of water and cosolvent pathways (hot-spots), are used to predict particular regions of protein-ligand interactions and could be used for further pharmacophore design. Since hsEH has a relatively high number of crystallographic protein-inhibitor complexes deposited at the Protein Data Bank (PDB) database[87], we could use this abundant information to gain structural insight into the binding of known inhibitors and compare them with the obtained results. We believe that such an approach can be applied for any molecular target and can work as a universal pipeline in the drug discovery pathway.