Synthesis, biochemical, pharmacological characterization and in silico pro ﬁ le modelling of highly potent opioid orvinol and thevinol derivatives

Morphine and its derivatives play inevitably important role in the m -opioid receptor (MOR) targeted antinociception. A structure-activity relationship study is presented for novel and known orvinol and thevinol derivatives with varying 3- O , 6- O , 17- N and 20-alkyl substitutions starting from agonists, an- tagonists and partial agonists. In vitro competition binding experiments with [ 3 H]DAMGO showed low subnanomolar af ﬁ nity to MOR. Generally, 6- O -demethylation increased the af ﬁ nity toward MOR and decreased the ef ﬁ cacy changing the pharmacological pro ﬁ le in some cases. In vivo tests in osteoarthritis in ﬂ ammation model showed signi ﬁ cant antiallodynic effects of thevinol derivatives while orvinol de- rivatives did not. The pharmacological character was modelled by computational docking to both active and inactive state models of MOR. Docking energy difference for the two states separates agonists and antagonists well while partial agonists overlapped with them. An interaction pattern of the ligands, involving the interacting receptor atoms, showed more ef ﬁ cient separation of the pharmacological pro ﬁ les. In rats, thevinol derivatives showed antiallodynic effect in vivo . The orvinol derivatives, except for 6-O-desmethyl-dihydroetor ﬁ n ( 2c ), did not show antiallodynic effect. © 2020 Elsevier Masson SAS. All rights reserved.


Introduction
Pain modulation is mainly regulated through the activation of the three classical types of opioid receptors, i.e. the m-, dand kopioid receptors (MOR, DOR and KOR, respectively) expressed in the neurons of the central and peripheral nervous system. The opioid receptors are members of the G-protein coupled receptors (GPCRs), the largest receptor family in the human genome, sharing the distinctive seven helical hydrophobic transmembrane helix domain [1e5]. Their activation leads to the inhibition of adenylyl cyclase which results in hyperpolarisation and inhibits neurotransmitter release [6,7]. The main target of the antinociceptive drugs in the treatment of pain is MOR.
Endogenous opioid peptides such as Met-and Leu-enkephalin [8], b-endorphin [9] and dynorphin-A [10] are produced in the brain. Two endomorphin tetrapeptides, endomorphin-1 and endomorphin-2 [11] were found to be highly selective endogenous agonists for MOR. Morphine is a prototype opioid agonist binding to MOR and is still the most frequently used drug in pain medication. Beside pain relief and analgesia, it has serious side effects including decreased respiratory effort, low blood pressure and it also has a high potential for addiction and abuse [12e14].
A number of structure-activity relationship studies dealing with thevinol and orvinol derivates are available [26,27], but the biochemical and pharmacological properties of our target compounds have not been reported except for 8 [25], The aim of present study was to compare the receptor binding properties and the MOR, DOR and KOR selectivity of some Bentley compounds in rat and guinea pig brain membrane preparations. The ligands were also investigated in [ 35 S]GTPgS functional binding assays to examine G-protein activation via opioid receptors. The effect of the investigated derivates was observed in vivo nociceptive tests.
The presence or absence of specific functional groups in the orvinol and thevinol derivatives can not be straightforwardly related to their pharmacological profiles. As an example, the 17-Nsubstituent serves as an acknowledged pharmacological switch between agonists and antagonists being methyl or cyclopropylmethyl, respectively. However, it is highly ambiguous within this class of opioids, regarding that 17-N-cyclopropylmethyl derivative can be full agonist as well [28,29] which may be a consequence of the bigger size of these opiates resulting in a more complex interaction pattern with the receptor. According to this, it seems plausible to investigate the interacting residues or atoms of the receptor leading to the specific response, i.e. pharmacological feature.
Former computational studies attempted to identify the interacting residues in MOR responsible for different pharmacological actions. However, no distinguishable interaction pattern was found for a structurally highly diverse set of agonist, partial agonist and antagonist to predict pharmacological activities using the inactive receptor state [30].
Compound 1b with a mixed agonist-antagonist (partial m-agonist/ k-antagonist) profile is clinicaly used as analgesic in the treatment of postoperative and/or cancer-related pain and in the substitution therapy of opioid dependent humans. 1c and 1d are nonselective opioid receptor agonists.
In the present study we performed a selective 6-O-demethylation of several compounds in order to study the structure-activity relationships. Binding affinities of 6-O-desmethyl-orvinols to opioid receptors were not yet investigated and their pharmacological/biochemical characterisation is currently not available in the scientific literature. In contrast to 3-O-demethylation the 6-Odemethylation of Bentley compounds is less explored. 6,14-Ethenomorphinans with a free tertiary hydroxyl group in position-6 were inaccessible before 1986. The first selective 6-Odemethylation of 7a-aminomethyland 7a-aldoxime-type 6,14ethenomorphinan derivatives was reported by Kopcho & Schaeffer [37]. Lithium aluminium hydride in tetrahydrofuran containing a halogenated co-solvent (CCl 4 ) was utilized as demethylating system. A six membered ring aluminum complex was hypothesized to play an important role in this unusual O-dealkylation. Subsequently, Lever et al. [16,17] extended this special 6-O-demethylation method for other 6,14-ethenomorphinans with  tert-alcohol functions in the position-20. Another method for the 6-O-demethylation of 6,14-ethenomorphinans has been developed by Breeden et al. [38], in 1999. respectively, in rat and guinea pig brain homogenates. All derivatives exhibited higher binding affinity in m-opioid receptor system than the selective peptide ligand DAMGO (Fig. 1A), some of them had extremely low K i values ( Table 1). For DOR the ligands showed comparable binding affinities than the selective DOR agonist Ile 5,6 -deltorphin II peptide ligand (Fig. 1B) except 8 (K i > 3000 nM). In the KOR binding assays, performed in guinea pig brain membranes, the analogues still displayed nanomolar affinities (Fig. 1C).

Functional GTPgS binding stimulation assay
The effect of the ligands on receptor-mediated G-protein activation was investigated in [ 35 S]GTPgS binding assays in rat brain membranes (Fig. 2). The highest stimulations were observed with 2d, 4, 5, 7 and 9, therefore they can be considered as full agonists. 1a, 1f, 2a and 2b did not produced dose-dependent increases, so they behave as neutral or pure antagonists ( Table 2). The remaining compounds (1b, 1e, 2c and 8) exhibiting intermediate levels of Gprotein activation are partial agonist ligands in this in vitro system.
The pure opioid antagonist 1a successfully reversed the efficacy of almost all compounds to basal activity with the exception of 2d which showed some remaining activation in the presence of equimolar 1a. Maximal stimulation produced by the ligands was mostly not elevated further when the full agonist 9 was coadministered ( Fig. 3A). However, in the case of 1f and 2b the copresence of 9 was able to effectively stimulate G-protein activation (Table 3).
Increasing concentrations of the partial agonists were also investigated in the presence of 10 mM 9 producing maximal stimulation ( Fig. 3B). All four compounds were able to inhibit the activation mediated by 9, although with relatively low efficacy and potency. This weak antagonizing effect in the presence of a full agonist validates that 1b, 1e, 2c and 8 are indeed partial agonist ligands for opioid receptors (Table 3).

In vivo studies
The basal withdrawal threshold of the non-inflamed side was 45 ± 0.5 g, and MIA caused significant decrease in paw withdrawal threshold on the injected side (24 ± 0.6 g). Only the largest dose of 9 treatments caused significant enhancement in the pain threshold on the non-inflamed side, therefore, results were analysed only on the inflamed paws. The different ligands showed different potencies, therefore, they were compared to distilled water (as negative control) in the ANOVA analysis, but the curve for 9 was also demonstrated as a positive control group with the lowest ED 30 value ( Table 4).
All of the thevinol derivatives showed significant antiallodynic effects (Fig. 4); however, the regression analysis revealed a lower in vivo potency compared to 9 as indicated by the ED 30 values, even it could not be calculated for 7 (Table 4). Regarding 4, the treatment was close to significant (Table 4). Time and their interaction showed significant effects, and the post-hoc analysis showed decreased allodynia in several time-points compared to the control group (Fig. 4A). ANOVA for repeated measurements showed Table 2 The maximal G-protein efficacy (E max ) of the morphine analogues in [ 35 S]GTPgS binding assays in rat brain membrane homogenates. The values were calculated according to dose-response binding curves in Fig. 2   Regarding 2a and 2b treatments, there were no significant effects and ED 30 values could not be calculated ( Fig. 5A and B, Table 4). However, 2c administration resulted in significant effects of treatment, time, and their interaction, with a relatively low ED 30 value ( Table 4). The post hoc comparison revealed that the largest dose of 2c caused similar degree of antiallodynic effect as did 9; even more prolonged effect was observed (Fig. 5C). 1e, 1f and 2d did not produce significant antinociceptive effects (Fig. 5D, E, F), and ED 30 values could not be calculated either.

Docking
Almost all of the investigated compounds showed higher binding affinity to MOR and thus they are expected to exert Gprotein activation through MOR. According to this, their pharmacological behaviour was modelled on the MOR crystal structures.
The target compounds are members of three pharmacological types showing agonism, antagonism or partial agonism at MOR. Thus the target compounds and several known agonists, antagonists and partial agonists, composing a set of 48 compounds (Tables S3e1,  Tables S3e2) [25,26,28,29,39], were docked to both the active and inactive receptor models to reveal whether their characteristics in docking experiments can be related to that of ligands with known pharmacological character. However, a simple visual inspection of the docked positions did not reveal specific features to explain the pharmacological diversity neither at the active nor at the inactive receptor state, therefore, the ligands were further characterized with their docking energies and the contacting receptor atoms. It is noteworthy that among the lowest energy docking poses only 16 out of 96 (48 ligands docked to both receptor states) originated from energy minimized conformers (see Experimental section 4.4.3.) suggesting the role of the flexibility of aliphatic rings.

Analysis of the docking energies and ligand efficiencies obtained for the active and inactive receptor states
Three kinds of docking energy measures were investigated: docking energies calculated by AutoDock Vina (E), ligand efficiency (docking energy divided by the number of non-hydrogen atoms of the ligand, LE) [40] and a similar value obtained from the docking energy divided by the count of the interacting atoms of the ligand Table 3 The maximal G-protein efficacy (E max ) the morphine analogues in the absence or presence of the opioid antagonist and agonist, 1a and 9, respectively in rat brain membrane homogenates. The values were calculated according to dose-response binding curves in Fig. 3 152.15 ± 0.85 100.50 ± 3.00 e e Table 4 The applied drugs and the cumulative dose procedure, the number of the animals and ANOVA results in each group and their in vivo potency as ED 30 with confidence interval (CI).
Drug (LEIAC). The effect of the receptor state on the docking energy measures was analysed by two-sided paired t-tests (Table 5.). According to this, agonists and partial agonists clearly differentiate between the receptor states by docking energies and LE while antagonists do not. In case of LEIAC agonists could not distinguish the receptor states. The energetic preference of the different ligand types for the receptor states, i.e. the difference between the docking energies obtained for the active and inactive states, was also investigated (Table S4). The energy difference would be negative for agonists showing their physically feasible preference for the active receptor and it should be the opposite (positive) for antagonists. Although this expectation was only partially fulfilled within the series of compounds investigated here, agonists were well separated from antagonists and partial agonists by two-sided unpaired t-tests using docking energy and LE values (Table 6.). Antagonists and partial agonists did not differ significantly. This is an interesting result however, because the geometric differencies between the receptor states do not seem significant (Table S5, Fig. S23), nevertheless the ligands distinguished them by binding energy.

Multivariate statistical analysis of docking energy values
Because both the receptor preference of the ligand types (Table 5.) and the difference between the types (Table 6.) were partially fulfilled, principal component analysis (PCA) was performed with five principal components to reveal the relationships between docking energy values and pharmacological features. Input data were the docking energies (E), LE and LEIAC values and their differences for the two receptor (Table S4). Additionally, the energy contributions decomposed for the specific interacting atom pair types, calculated by BINANA, were also involved. The different energy measures, however, did not contribute equally well in PCA to separate the experimentally determined pharmacological types. The efficiency of clustering was assessed by considering hierarchical clustering, J 2 statistics [41] and pairwise cluster overlap statistics using the program package "pca-utils" [42] for the five PCA components. Hierarchical clustering statistics and pairwise cluster overlap probabilities conformed each other showing that the separation of agonists from antagonists and partial agonists performed well while distinguishing antagonists from partial agonists highly depended on the energy measures used in PCA (Table 7. J 2 statistics). J 2 measures the fuzziness of the clusters, less J 2 means more compact cluster. It is interesting that using the decomposed energy contributions (ES) resulted in more efficient cluster separation but less compact clusters (Table 7. J 2 statistics, cluster overlapping statistics). Considering all three quality matrix, the best classification was obtained by the use of the docking energies and their active-inactive receptor differences (Ea, Ei, Ea-Ei, LEa-LEi, LEIACa-LEIACi). Results for the first two PCA dimensions and the distance matrix for five PCA variables are shown in (Fig. 6). The results show that the ligands can be classified to their known pharmacological groups using docking energies and related measures, albeit with significant overlap (Fig. 6A). However, a closer look on 3D representation of the first three PCA dimensions (Fig. 7) revealed better separation of the ligand types. There is a little overlap between agonists and antagonists while partial agonists overlap with both. According to PCA results 2c, a newly synthesized partial agonist, is among antagonists and close to established partial agonists (BU08028 (11), BU61 (14)). However, 1f, despite being antagonist experimentally, showed up among agonists in all kinds of PCA calculations. Furthermore, compound 30, considered to be a full agonist despite its N-cyclopropylmethyl substituent [29], was among antagonists in this model. It is worth to note that the use of energy differences alone gave almost the same classification efficiency.

Characterization of the ligands by the interacting receptor atoms
Because the ligands exert their effects through interactions with the receptor, the details of these interactions (interaction pattern)   Table 6 Comparison of docking energy differences between ligand types.

Ligand types F-test t-test
agonist -antagonist 0.020 2.00 Â 10 À5 * agonist -partial agonist 0.023 2.98 Â 10 À4 * partial agonist -antagonist 0.845 0.371 may reflect the pharmacological features. The interaction pattern of a ligand is the list of the interacting receptor atoms depicted in both receptor states. Additionally, each atom is marked whether it was involved in stabilizing (attractive) and/or destabilizing (repulsive) interaction. The interacting receptor atoms were extracted from the output of BINANA [43] analysis of the receptor-ligand complex listing up all particular interacting ligand-receptor atom pairs using the default parameters of BINANA. The stabilizing/destabilizing nature of an interaction between the ligand and receptor depends on the type of the interacting atoms. There are compatible and incompatible atom types (or atom classes) [44] resulting in stabilizing and destabilizing interactions, respectively. The atom class scheme was adopted to the AutoDock atomtypes from Sobolev et al. [44]. The classification of the ligands using their interaction patterns was performed by multiple correspondence analysis (MCA). A variety of interaction patterns, holding different information of the interactions, were created to find the most useful one for pharmacological classification of the ligands: i) both receptors, individual atoms, stabilization flags, ii) active receptor, individual atoms, stabilization flags, iii) inactive receptor, individual atoms, stabilization flags, iv) both receptors, individual atoms, without stabilization flags, v) both receptors, residues only, stabilization flags. As it was expected, the best classification was obtained with the most information rich input, i.e. with both receptors, individual atoms, stabilization flags (i) and the worst separation was obtained using only the interacting residues and stabilization flags (Table 7 J2). The MCA results are shown in Fig. 8 and the 3D representation in Fig. 9. Comparing the energy based PCA and the interaction pattern based MCA results suggests that the interaction patterns resulted in better pharmacological classification of the ligands, although the antagonist 1f was an outliner in both cases being close to the agonist group. It is worth to note that the cluster separation statistics was better for all kinds of MCA results compared with PCA, although the fuzziness of the clusters increased (Table 7).
Biological properties of our thirteen synthetic morphine analogues were investigated by in vitro biochemical and in vivo pharmacological experiments. In radioligand binding assays performed using MOR-, DOR-or KOR-selective [ 3 H]labelled primary ligands, all analogues exhibited excellent affinities for the multiple types of opioid receptors. At the MOR binding sites 2c, 4 and 5 were the most potent ligands, but all of the remaining compounds represented quite high affinities. They were also potent competitors in the DOR-selective receptor binding assays, displaying still    showed the highest affinities at KOR ligand binding sites, although all other analogues produced low nanomolar binding affinities. Taken together, the thevinol and orvinol derivatives are very high affinity opioid ligands, with a general preference for MOR > KOR > DOR binding sites.
Receptor mediated G-protein activation experiments were conducted in vitro using rat brain membrane preparations and [ 35 S] GTPgS binding stimulation assays. Transmembrane signalling properties of the compounds were variable, and the ligands used can be divided into three biochemical pharmacological groups based on their stimulation features. Full agonists, such as 2d, 4, 5, 7 and 9 produced the highest efficacy (E max values) and they can be characterised by culminating sigmoid stimulation curves with a plateau. Partial agonist ligands, 1b, 1e, 2c and 8, exhibited submaximal efficacies and rather decreased potencies in activating G i /G o regulatory proteins. Another previous term for partial agonists has been mixed agonist-antagonist ligands. When both a full agonist and partial agonist are present, the partial agonist actually acts as a competitive antagonist, competing with the full agonist for receptor occupancy and producing a net decrease in the receptor activation observed with the full agonist alone (Fig. 6B). The third cluster of ligands are pure or neutral opiate antagonists, such as 1a, 2a, 2b and 1f. They are characterised by horizontally linear dose-response curves indicating no changes in the basal G-protein activity.
6-O-demethylation usually increased the binding affinity to MOR (1a vs. 2a, 1b vs. 2b, 2c vs. 9, 4 vs. 5). In accordance with this finding, increasing the size of the 6-O-substituent (4 vs. 7) decreased the binding affinity to MOR ( Table 1). The opposite trend was observed for the efficacy in some cases, resulting in a shift of the pharmacological profiles, i.e. from partial agonist to antagonist (1b vs. 2b), from agonist to partial agonists (2c vs. 9). In the phenethyl thevinol series (2d, 4, 5, 7) 6-O-demethylation did not decrease the efficacy. On the contrary, the bigger 6-O-substituent in 7 slightly decreased the efficacy, however, it still can be considered as a full agonist ( Table 3). The changes in the pharmacological features caused by the different substituents are graphically summarized in Fig. S37.
Behavioural nociceptive properties of the ligands were studied in vivo, using Wistar rats in a chronic osteoarthritic inflammatory pain model. 9 at cumulative doses of 0.1-0.3e1.0 nmol/kg remained the most effective compound in decreasing inflammatory pain. 2c, 4 and 5, which exhibited the highest ligand binding affinities at the MOR, produced less antinociception than 9 did. The other full (2d) or partial antagonists (1e, 7), as well as the two neutral antagonists 2a and 2b examined in this test were practically not effective even in ten times higher (1-3-10 nmol) cumulative doses.
Docking the ligands to both the active and inactive receptor states makes possible reasonable pharmacological classification of the ligands by docking energies. The positive effect of the differences of the docking energy measures on the pharmacological classification of the ligands also emphasize that, despite the moderate geometric difference between the active and inactive receptor states within the binding pocket, the docking energies can predict the pharmacological features of the ligands. Pharmacological classification was also attempted by the interaction pattern of the docked ligands, i.e. by the interacting receptor atoms and the type of the interactions (stabilizing or distabilizing). More information (i.e. more kinds of docking energy measures, two receptor states instead of single one) resulted in better classification in both cases.
Agonists and antagonist are separated quite well while partial agonists overlap with the others which is in accordance with molecular dynamics results comparing interactions of agonists, antagonists and partial agonists; nevertheless, pharmacological classification by the docking algorithm in the present paper is still closer to the high throughput methodology.

Conclusions
All investigated compounds showed subnanomolar binding affinity to MOR and a preference to MOR over DOR and KOR. The pharmacological effects of the compounds involved agonism, partial agonism and antagonism. Neither binding affinities nor pharmacological features could be directly related to particular organic functional groups.
The in vitro pharmacological effects and the in vivo antiallodynic effects were in accordance except the full agonist 2d. As the only exception, the newly synthesized compound, 2c despite its partial agonist feature, showed antiallodynic effect equal tot hat of the full agonist 9.
Due to the harmony between the in vitro and in vivo results, the in silico calculations were expected to explain the pharmacological profiles of the compounds. The unsupervised multivariate classification methods (using either docking energies or interaction fingerprints) applied to the docking results, obtained from active and inactive receptor states, were able to separate the agonists from antagonists with a good accuracy. Additionally, the third group, partial agonists, were partially differentiated from the other two groups. Due to the effectivity of the multivariates classification, their further improvement seems to be promising. Differentiating between ORs needs accurate docking calculations, however, if it is accurate enough it can differentiate between the receptor states as well. If it is so, the first step in the modelling scenario should be the pharmacophoric featuring (receptor state selection) for binding affinity prediction. Potassium hydroxide (3.2 g, 57 mmol) was dissolved in diethylene glycol (20 mL) at 110 C. The solution was allowed to cool to 70 C and the corresponding 3-O-methyl-thevinol derivative (2 mmol) was added. The reaction mixture was stirred under argon atmosphere at 210 C (internal temperature) for 90 min. The brownish product mixture was allowed to cool to room temperature and poured into a saturated ammonium chloride solution (40 mL). The suspension was extracted with diethyl ether (4 Â 45 mL). The combined organic layer was extracted successively with 5% sodium hydrogen sulfite solution (2 Â 30 mL) and water (2 Â 30 mL). The organic phase was dried (Na 2 SO 4 ) and the solvent was evaporated under reduced pressure. The residue was purified by means of column chromatography on silica gel using the appropriate solvent system.

Chemicals
MgCl 2 x 6H 2 O, EGTA, Tris-HCl, NaCl, GDP, the GTP analogue GTPgS, were purchased from Sigma-Aldrich (Budapest, Hungary). The highly selective MOR agonist enkephalin analogue DAMGO was obtained from Bachem Holding AG (Bubendorf, Switzerland). The highly selective KOR agonist diphenethylamine derivative, HS665 [45] were kindly offered by Dr. Helmut Schmidhammer (University of Innsbruck, Austria). The morphine analogues were provided by ABX GmbH (Radeberg, Germany). The highly selective DOR agonist Ile 5,6 -deltorphin II was synthesized in the Laboratory of Chemical Biology group of the Biological Research Centre (Szeged, Hungary). DAMGO, [Ile 5,6 ]-deltorphin II and HS665 were dissolved in water, morphine analogues were dissolved in ethanol and were stored in 1 mM stock solution at À20 C. The radiolabeled GTP analogue, [ 35

Animals
For membrane homogenate preparations male and female Wistar rats and guinea pigs were used. Animals were housed in the local animal house of BRC (Szeged, Hungary). All the animals were kept in a temperature controlled room (21e24 C) under a 12:12 h light/dark cycle and allowed free access to food and water. All housing and experiments were conducted according to the European Communities Council Directives (86/606/ECC) and the Hungarian Act for the Protection of Animals in Research (XXVIII.tv. 32.x). The total number of animals as well as their suffering was minimized.

Preparation
A crude membrane fraction of Wistar rat and guinea pig brains were prepared for ligand binding experiments according to Ref. [48] with changes. After decapitation, brains were rapidly removed and homogenised in 30 vol of ice-cold 50 mM TriseHCl (pH 7.4) buffer with a teflon-glass homogeniser. After centrifugation at 40000g for 20 min at 4 C, the resulting pellet were suspended in 30 vol of the same buffer and incubated for 30 min at 37 C to remove endogenous opioids. Centrifugation was then repeated as described above. The final pellet were suspended in 5 vol of 50 mM TriseHCl (pH 7.4) buffer and stored at À80 C. Membranes were thawed, diluted with buffer and the resulting pellet were taken up in appropriate fresh buffer and immediately used. For the [ 35 S]GTPgS binding experiments the final pellet of rat and guinea pig brain membrane homogenates were suspended in 5 vol of ice-cold TEM (Tris-HCl, EGTA, MgCl 2 ) and stored at À80 C for further use. as previously described [49]. Total binding was measured in the absence of test compounds, non-specific binding was determined in the presence of 10 mM unlabeled GTPgS and subtracted from total binding. The difference between total and non-specific binding values represents basal activity. The reaction was terminated by rapid filtration under vacuum (Brandel M24R Cell Harvester), and Whatman GF/B glass fiber filters were washed three times with 5 mL ice-cold 50 mM Tris-HCl (pH 7.4) buffer. The radioactivity of the dried filters was detected in Ultima GoldTM MV aqueous scintillation cocktail with Packard TriCarb 2300 TR liquid scintillation counter. [ 35 S]GTPgS binding experiments were performed in triplicates and repeated three times.

Data analysis
Experimental data were presented as means ± S.E.M. Points were fitted with the professional curve fitting program, GraphPad Prism 5.0 (GraphPad Prism Software Inc., San Diego, CA), using nonlinear regression analysis. In the [ 35 S]GTPgS binding assays the 'Sigmoid dose-response' fitting was used to establish the maximal efficacy (E max ) of the receptors' G-protein and the ligand potency (EC 50 ). Stimulation was given as percent of the specific [ 35 S]GTPgS binding observed over the basal activity, which was settled as 100%.
In the competition binding assays the 'One site competition' fitting was used to establish the equilibrium binding affinity (K i value). Inhibition was given as percent of the specific binding observed.

In vivo experiments
After institutional ethical approval had been obtained (Institutional Animal Care Committee of the Faculty of Medicine at the University of Szeged), male Wistar rats (Charles River strain, Bioplan, Budapest, Hungary; 378 ± 5.1 g; n ¼ 6e8/group) were used in the experiments. The animals were group-housed (4 animals/cage) with free access to food and water, and with a 12:12 h light/dark cycle. Animal suffering and the number of animals per group were kept at a minimum; therefore, the drugs were administered in cumulative doses in 30-min intervals (at 0th, 30th and 60th min), and injections were repeated 7 days apart for the same animal, as in our previous study [50].
Osteoarthritis was induced by injecting MIA (Sigma-Aldrich Ltd. Budapest, Hungary; 1 mg/30 mL) into the tibiotarsal joint of the right hind leg on two consecutive days. MIA treatments were given to gently restrained conscious animals, using a 27-gauge needle, without anesthesia so as to exclude any drug interaction. These injections did not elicit signs of major distress. Within 14 days MIA had consistently been shown to cause severe end-stage cartilage destruction resulting in osteoarthritis-like joint pain accompanied with moderate edema [50,52,53]. The observer was blind to the drug treatment administered.

Behavioral nociceptive testing
The threshold for withdrawal from mechanical stimulation to the plantar aspect of the hindpaws was assessed using a dynamic plantar aesthesiometer (Ugo Basile, Comerio, Italy), which consists of an elevated wire mesh platform to allow access to the ventral surface of the hindpaws. Prior to baseline testing, each rat was habituated to a testing box for at least 20 min. Measurements were done with a straight metal needle (diameter 0.5 mm) that exerts an increasing upward force at a constant rate (6.25 g/s) with a maximum cut-off force of 50 g over an 8 s period. The measurement was stopped when the paw was withdrawn, and the results were expressed as paw withdrawal thresholds in grams.

Experimental paradigm
After MIA injections, baseline pain thresholds (2 times with 15 min interval) were determined 14 days later and their means provided the baseline value.
Cumulative-dosing procedure was applied ( Table 4). The control group received distilled water. In the positive control group, animals were treated with 9 (0.1, 0.3, and 1.0 nmol/kg). The higher doses of 9 were also applied in a preliminary experiment (1.0 and 3.0 nmol/kg) but 3.0 nmol/kg dose induced catatonia and respiratory depression in 100% of the animals, thus we determined the maximum dose as 1 nmol/kg. Response latencies were measured at 15 min intervals and the increasing doses of the drugs (three doses) were administered following two recordings. After the highest dose injection the pain threshold was assessed in every 15 min for 1 h, then hourly at the second and third hour to determine the time course of drug effects. Although, the motor behavior and general status of the animals were not investigated and quantified systematically, altered behavior (excitation, flaccidity or motor weakness) or any signs of opioid overdose (catatonia or respiratory depression) were not observed.

Statistical analysis
Data are presented as means ± SEM. Data sets were examined by repeated measures of ANOVA. Post hoc comparisons were carried out with the Fisher LSD test. A p value lower than 0.05 was considered significant. The mean paw withdrawal thresholds obtained 15 and 30 min after the injections (calculated after the individual drug injections) were used for linear regression analysis to determine the effective dose 30 (ED 30 ) values with 95% confidence intervals. Mean of 50 g would mean the complete relief of hyperalgesia, while mean of control value (29 g) means the zero effects, thus the difference (21 g) is the possible maximal effect. ED 30 is equivalent to the dose that yielded 30% difference (7 g) in the paw withdrawal threshold compared to the baseline (pretreatment) values. Data analyses were performed with the STATISTICA (Statistica Inc., Tulsa, Oklahoma, USA) and GraphPad Prism 4.0 (GraphPad Softwre Inc. La Jolla, CA, US) softwares.

Molecular docking
The ligands were flexibly docked to the experimentally solved active and inactive states of MOR (http://www.rcsb.org, pdb codes 5c1m and 4dkl for the active and inactive receptor states, respectively). The receptor structures were kept rigid. Due to the known limitation of the docking algorithm i.e. the rings are kept rigid during the calculation, 10 different conformations were generated by molecular dynamics to emulate the conformational flexibility for aliphatic rings in accordance with experimental findings [54]. The 10 conformations and the initial energy minimized conformer were provided for docking for each ligand. Docking calculations were performed by the program PSOVINA [55], a variant of the formerly released VINA [56]. Coordinate files were converted to pdbqt format by AutoDock Tools [57] or PSOVINA. Docking box center was set to the geometric center of the ligands in the crystals. Due to the impact of the size of the docking box on the docking pose, namely on the set of interacting receptor atoms, the box size was specifically determined for each ligand [58] using the radius of giration (R g ) of the minimum energy conformer calculated by Open Babel. PSOVINA was used with default parameters, albeit docking was repeated 5 times instead of increased exhaustiveness parameter [59] and the lowest energy docking poses was kept for each ligand.

Analysis of docking results
The selected poses were analysed by the program BINANA [43] to extract the interacting ligand-receptor atom pairs, the docking energies and the energy contributions belonging to the specific atom pairs. Stabilizing or destabilizing assignments for the ligandreceptor contacts based on the method Soboljev et al. [44]. Analysis of the docking energies and the interacting receptor atoms were performed by principal component analysis (PCA) and multiple correspondence analysis (MCA), respectively, using program packages of the R programming environment v. 3.5.1 [60]. Calculations were performed by the PCA and MCA moduls of the R package "FactoMineR" [61]. Results were visualized using package "factoextra" [62]. For 3D visualization packages "rgl" and "pca3d" were used [63,64]. Docking calculations were performed on a Linux cluster running Rocks 6.2 cluster management program and CentOS 6.6 operation system and on the PRACE/NIIF cluster, and the analyses in Linux Mint 18.1 using in house bash and R scripts.

Generation of the ligand conformations
Open Babel v. 2.4.1 program package [65] was used in a two-step manner: 1. conformers were generated with the "-score energy" option to avoid high energy structures, 2. conformers were further transformed with "-rings" option to allow ring flexibility using MMFF94s force field.

Declaration of competing interest
The authors declare no conflict/competing financial of interest.