# The Effect of CYP2B6, CYP2D6, and CYP3A4 Alleles on Methadone Binding: A Molecular Docking Study

Authors: Nik Nur Syazana Bt Nik Mohamed Kamal; Theam Soon Lim; Gee Jun Tye; Rusli Ismail; Yee Siew Choong
Journal: Journal of Chemistry
 (2013)

---

## Abstract

Current methadone maintenance therapy (MMT) is yet to ensure 100% successful treatment as the optimum dosage has yet to be determined. Overdose leads to death while lower dose causes the opioid withdrawal effect. Single-nucleotide polymorphisms (SNP) in cytochrome P450s (CYPs), the methadone metabolizers, have been showen to be the main factor for the interindividual variability of methadone clinical effects. In this study, we investigated the effect of SNPs in three major methadone metabolizers (CYP2B6, CYP2D6, and CYP3A4) on methadone binding affinity. Results showed thatCYP2B6*11, CYP2B6*12, CYP2B6*18, and CYP3A4*12 have significantly higher binding affinity to R-methadone compared to wild type. S-methadone has higher binding affinity in CYP3A4*3, CYP3A4*11, and CYP3A4*12 compared to wild type. R-methadone was shown to be the active form of methadone; thus individuals with CYP alleles that binds better to R-methadone will have higher methadone metabolism rate. Therefore, a higher dosage of methadone is necessary to obtain the opiate effect compared to a normal individual and vice versa. These results provide an initial prediction on methadone metabolism rate for individuals with mutant type CYP which enables prescription of optimum methadone dosage for individuals with CYP alleles.


---

## Body

## 1. Introduction

Methadone is a synthetic opioid that was first synthesized in the 1940s for analgesia. Besides that, it has also a longer half-life and is cheaper and able to “pharmacologically block” heroin like euphoria effect [1]. Therefore, methadone is widely used in the methadone maintenance therapy (MMT) in the treatment of patients with opioid dependency [2, 3]. These in turn improve the health and social productivity of the patients. However, overdose of methadone can lead to severe side effects, for example, coma, convulsions, and death. Insufficient dosage on the other hand will result in opioid withdrawal symptoms [4].Although methadone is generally used to combat with both illicit heroin addiction and HIV infection, optimization of the dose is still an ongoing process. Marketed methadone is usually a racemic form of two enantiomers, theR- and S-methadones at the ratio of 50 : 50 [5]. R-methadone accounts for the opioid effects [6] as it is reported to have clinically significant μ-receptor agonist activity [7]. Studies showed that improper dosage of methadone may lead to many undesirable effects like severe respiratory depression, QTc interval prolongations and “torsades de pointes” which can also cause sudden death. Thus, any differences, for example, interindividual variability in metabolism rate by CYPs, that lead to this variation need to be well studied [8].Interindividual variability in methadone metabolism is highly influenced by the genetic polymorphisms of cytochrome P450s (CYPs) [12] and  metalloproteins with a heme group as the catalytic center. CYPs are known to be important protein as they are involved in 70–80% of metabolism of xenobiotics (including methadone), converting some of them into less-toxic products or inactive form [13, 14]. The metabolizing activity by CYPs is substrate specific. Single-nucleotide polymorphisms (SNPs) in CYPs may contribute to changes in drug efficacy by leading to different effects in maximal plasma concentration, half life, and clearance of the drug from the body [12, 15]. CYP2B6, CYP2D6, and CYP3A4 are three major isoforms that are involved in methadone metabolism [6].Here we reported the methadone binding affinity of CYP2B6, CYP2D6, and CYP3A4 alleles (http://www.cypalleles.ki.se/) (accessed January 2012) as predicted by molecular docking simulation. Molecular docking simulation is an in silico method that calculates binding affinity, defining how favourable the binding between a given acceptor/ligand (e.g., methadone) in a receptor (e.g., CYP) is. Higher binding affinity (in terms of free energy of binding) shows that the ligand (e.g., methadone) is a better substrate for the receptor, thus indicating higher substrate metabolism activity [16–19]. Hence, free energy binding derived from docking calculation could help to predict the effects of SNPs in methadone metabolism rate [20]. Results showed that molecular docking simulation is able to distinguish the CYP alleles in methadone binding affinity which leads to the prediction of methadone metabolism rate. Thus, these data can help to shed some light on future methadone dosing for individual with CYP alleles, towards a better MMT management.


## 2. Materials and Methods

### 2.1. Methadone

Cartesian coordinates for neutral form ofR- and S-methadone were generated using Hyperchem 7.0 (Hyperchem 2001). AutoDockTool (ADT) [10, 11] was used to set the torsion number (a total of 7 torsions) and to calculate the Gasteiger PEOE partial charges [21–23].


### 2.2. Wild-Type (WT) CYP

Starting Cartesian coordinates for WT CYP2B6, CYP2D6, and CYP3A4 were obtained from RCSB Protein Data Bank with PDB id 3IBD [24], 2F9Q [25], and 3NXU [26], respectively. Introduced mutations of Y226H and K262R in 3IBD and L230D and L231R in 2F9Q were reversed back to the WT using MODELLER9v8 [27] mutate-script. All residues in CYP3A4 crystal remained unchanged as no mutations were introduced during the crystallization process [26].


### 2.3. CYP Alleles

A total of 10, 14, and 12 SNPs (mutant type; MT) were identified from the lists in Human Cytochrome P450 (CYP) Allele Nomenclature Committee website (http://www.cypalleles.ki.se/) to represent alleles for CYP2B6, CYP2D6, and CYP3A4, respectively. MODELLER9v8 mutate script was used to generate Cartesian coordinates for the respective alleles.


### 2.4. Docking Simulation

Docking simulation was performed by Autodock 4.0 [22, 28]. All WT and MT CYPs were added polar hydrogen atoms using the program protonate, and charges were loaded using the kollua. amber option. As the selection of the grid center is depending on the active site residues [29, 30], therefore the grid box and grid center for each CYP were individually optimized (Table 1). Residues within the methadone binding site were set to be flexible. These residues were Leu363, Val367, Thr302, Thr303, Arg98, Ile114, His369, Arg434, and Thr306 in CYP2B6; Arg101, Phe120, Leu302, Thr309, Thr310, Ile369, Val370, Met374, and His376 in CYP2D6; Arg105, Val118, Ser119, Phe302, Thr309, Thr310, Ile369, Leu373, and Arg375 in CYP3A4. A total of 21, 20, and 21 torsions were assigned to the flexible residues for CYP2B6, CYP2D6, and CYP3A4, respectively. Grid spacing was set to 0.375 Å for all CYPs and generated by AutoGrid. A total of 100 docking runs were performed by employing Lamarckian genetic algorithm (LGA) with pseudo-Solis and Wets local search. Other docking parameters were population size of 300; energy evaluations of 25,000,000; maximum generations of 27,000; translational step of 0.2 Å; quaternion step of 50°; torsional step of 50.0°; clustering tolerance of 2.0; crossover rate of 0.80; mutation rate of 0.02; elitism of 1; local search rate of 0.06; 300 iterations per local search with termination value of 0.01 and consecutive successes or failures before doubling or reducing local search step size of 4.Table 1
Grid box and grid center for CYP2B6, CYP2D6, and CYP3A4.
Isoform
Grid box dimension
Grid center
X
Y
Z
CYP2B6
40
×
70
×
50
22.00
15.00
25.00
CYP2D6
70
×
50
×
70
−34.19
28.88
−47.00
CYP3A4
50
×
70
×
50
34.00
−18.00
25.00


### 2.5. Analysis

Interactions (e.g., formation of hydrogen bond,π-π interaction, and π-δ interaction) of methadone in CYP alleles were performed using Accelrys Discovery Studio 3.0 (DS; Accelrys Inc.) and visualized by Visual Molecular Dynamic 1.8.7 (VMD) [9].


## 2.1. Methadone

Cartesian coordinates for neutral form ofR- and S-methadone were generated using Hyperchem 7.0 (Hyperchem 2001). AutoDockTool (ADT) [10, 11] was used to set the torsion number (a total of 7 torsions) and to calculate the Gasteiger PEOE partial charges [21–23].


## 2.2. Wild-Type (WT) CYP

Starting Cartesian coordinates for WT CYP2B6, CYP2D6, and CYP3A4 were obtained from RCSB Protein Data Bank with PDB id 3IBD [24], 2F9Q [25], and 3NXU [26], respectively. Introduced mutations of Y226H and K262R in 3IBD and L230D and L231R in 2F9Q were reversed back to the WT using MODELLER9v8 [27] mutate-script. All residues in CYP3A4 crystal remained unchanged as no mutations were introduced during the crystallization process [26].


## 2.3. CYP Alleles

A total of 10, 14, and 12 SNPs (mutant type; MT) were identified from the lists in Human Cytochrome P450 (CYP) Allele Nomenclature Committee website (http://www.cypalleles.ki.se/) to represent alleles for CYP2B6, CYP2D6, and CYP3A4, respectively. MODELLER9v8 mutate script was used to generate Cartesian coordinates for the respective alleles.


## 2.4. Docking Simulation

Docking simulation was performed by Autodock 4.0 [22, 28]. All WT and MT CYPs were added polar hydrogen atoms using the program protonate, and charges were loaded using the kollua. amber option. As the selection of the grid center is depending on the active site residues [29, 30], therefore the grid box and grid center for each CYP were individually optimized (Table 1). Residues within the methadone binding site were set to be flexible. These residues were Leu363, Val367, Thr302, Thr303, Arg98, Ile114, His369, Arg434, and Thr306 in CYP2B6; Arg101, Phe120, Leu302, Thr309, Thr310, Ile369, Val370, Met374, and His376 in CYP2D6; Arg105, Val118, Ser119, Phe302, Thr309, Thr310, Ile369, Leu373, and Arg375 in CYP3A4. A total of 21, 20, and 21 torsions were assigned to the flexible residues for CYP2B6, CYP2D6, and CYP3A4, respectively. Grid spacing was set to 0.375 Å for all CYPs and generated by AutoGrid. A total of 100 docking runs were performed by employing Lamarckian genetic algorithm (LGA) with pseudo-Solis and Wets local search. Other docking parameters were population size of 300; energy evaluations of 25,000,000; maximum generations of 27,000; translational step of 0.2 Å; quaternion step of 50°; torsional step of 50.0°; clustering tolerance of 2.0; crossover rate of 0.80; mutation rate of 0.02; elitism of 1; local search rate of 0.06; 300 iterations per local search with termination value of 0.01 and consecutive successes or failures before doubling or reducing local search step size of 4.Table 1
Grid box and grid center for CYP2B6, CYP2D6, and CYP3A4.
Isoform
Grid box dimension
Grid center
X
Y
Z
CYP2B6
40
×
70
×
50
22.00
15.00
25.00
CYP2D6
70
×
50
×
70
−34.19
28.88
−47.00
CYP3A4
50
×
70
×
50
34.00
−18.00
25.00


## 2.5. Analysis

Interactions (e.g., formation of hydrogen bond,π-π interaction, and π-δ interaction) of methadone in CYP alleles were performed using Accelrys Discovery Studio 3.0 (DS; Accelrys Inc.) and visualized by Visual Molecular Dynamic 1.8.7 (VMD) [9].


## 3. Results and Discussion

Metabolism of methadone to its metabolites, for example, L-a-acetylmethadol (LAAM) is mediated by certain steps; for example, N-demethylation process occurred in CYP [31]. Substrate binding (e.g., methadone) requires reduction of CYPs; a process that involves heme iron [32]. Thus, studies recorded that active sites for CYPs are area around the heme [33, 34]. Preliminary docking using AutoDock3.05 [35, 36] showed R- and S-methadones were docked in the active site of CYP3A4 (see Figure S1 in the Supplementary Material available online at http://dx.doi.org/10.1155/2013/249642). However, active site superimposition of CYP3A4 with CYP2B6 and CYP2D6 showed that Val367 and Met374 in both CYP2B6 and CYP2D6 were overlapped with R- and S-methadones (see Figure S2). Thus, resulted in R- and S-methadone was docked outside the active site of CYP2B6 and CYP2D6 (see Figure S2). Therefore, flexible docking were performed by AutoDock4.0 [22, 28] to study the conformation of methadone in the CYPs active site.Visualization of CYP2B6, CYP2D6, and CYP3A4 revealed the differences in their overall structure, size, and shape of active site (Figure1). Thus, grid box and grid center for each CYP were individually optimized (Table 1). Residues (including Val367 and Met374) within 10 Å from the heme were set to be flexible (Figure 2). The upper part of heme was set as grid center as CYP crystal structures deposited in PDB showed that it is the binding site for many substrates. The selection of grid center was also supported by Protein Binding Site Detection (ProBis) [37, 38] as this upper of heme is a highly conserved region and predicted as the binding site (data not shown). Recent study by Moors and coworkers also reported that the upper part of heme is the active site for CYP2D6 [39].Active site of CYPs. Surface presents the active site of (a) CYP2B6, (b) CYP2D6, and (c) CYP3A4. Heme is in blue stick representation. Figure was generated using VMD1.8.7 [9].
(a)
(b)
(c)Grid box and flexible residues of (a) CYP2B6, (b) CYP2D6, and (c) CYP3A4, in ribbon representation. CPK represents SNP; blue CPKs represent mutation that lead to the lowest/highest binding affinity forR-/S-methadone in CYP. Heme is in stick presentation. Red square resembles the grid box. Figure was generated using ADT4.0 [10, 11].
(a)
(b)
(c)Docking simulation showed that bothR- and S-methadone are docked in the same active site (Figure 3) with similar binding affinity (Table 2) for WT. The binding of methadone with CYP2B6 and CYP2D6 alleles was more “centered” compared to “wider” distribution in CYP3A4 (Figure 3). A total of 1, 5, and 3 alleles of CYP2B6, CYP2D6, and CYP3A4, respectively, are located within 5 Å from the heme group.Table 2
Estimated free energy of binding (FEB) and the number of interactions in CYPs alleles from the docking simulation ofR- and S-methadone.
CYP
Allele
Mutation
Distance from heme (Å)
Free energy of binding
Number of interaction
(kcal/mol)
R
-methadone
S
-methadone
R
-methadone
S
-methadone
Hb
π-π
π-δ
Hb
π-π
π-δ
2B6
*1A
WT
—
−22.20
−20.69
1
0
0
0
0
0
*9
Q172H
11.0
−22.52
−19.64
0
0
0
0
0
0
*11
M46V
20.0
−19.64
−19.11
0
0
0
0
0
1
*12
G99E
8.0
−18.58
−22.58
0
0
2
0
0
0
*14
R140Q
15.0
−20.49
−20.82
0
0
1
0
0
0
*15
I391N
10.0
−20.86
−20.78
0
2
0
0
0
0
*18
I328T
19.0
−18.66
−21.49
0
0
0
0
1
0
*21
P428T
4.0
−21.72
−21.45
0
0
0
0
0
0
*24
G476D
12.0
−20.92
−21.62
0
0
1
0
0
1
*25
Q485L
21.0
−21.99
−19.90
0
0
0
0
0
1
*27
M198T
15.0
−20.36
−20.45
0
0
1
0
1
2
2D6
*1A
WT
—
−20.93
−20.81
0
1
0
0
0
0
*7
H324P
19.0
−23.21
−19.31
0
1
0
0
0
1
*23
A85V
14.0
−20.37
−21.72
0
0
1
0
0
0
*24
I297L
11.0
−20.71
−20.73
0
2
1
0
0
0
*26
I369T
5.0
−20.79
−21.69
0
1
0
0
1
1
*33
A237S
23.0
−21.17
−20.08
0
1
0
0
0
0
*34
R296C
12.0
−22.79
−20.89
0
0
1
0
0
0
*39
S486T
12.0
−20.85
−20.91
0
0
0
0
1
0
*48
A90V
10.0
−20.88
−21.03
0
0
0
0
1
0
*50
E156A
13.0
−20.88
−21.33
0
0
0
0
1
0
*62
R441C
2.0
−23.16
−21.79
0
0
0
0
2
0
*74
L91M
5.0
−20.79
−21.19
0
1
0
0
1
0
*75
R441H
2.0
−21.09
−20.25
0
1
0
0
0
0
—
R123L
8.0
−23.15
−20.95
0
1
0
0
0
1
—
G445E
4.0
−21.77
−21.22
0
0
0
0
1
1
3A4
*1A
WT
—
−20.97
−20.47
0
0
0
0
0
0
*2
S222P
23.0
−19.45
−21.64
1
1
0
1
1
0
*3
M445T
6.0
−21.32
−23.92
1
1
0
1
1
0
*4
I118V
4.0
−20.01
−22.72
0
0
0
0
1
2
*8
R130Q
3.0
−20.60
−20.34
0
0
0
0
0
0
*10
D174H
16.0
−20.39
−19.30
1
1
1
0
0
0
*11
T363M
10.0
−19.41
−24.05
0
0
0
1
0
0
*12
L373F
4.0
−24.46
−24.64
0
1
0
0
0
0
*16
T185S
9.0
−21.18
−20.61
1
1
0
0
0
0
*17
F189S
10.0
−22.01
−20.91
0
1
0
1
1
0
*18
L293P
14.0
−20.58
−22.65
0
1
0
0
0
0
—
I431T
11.0
−21.45
−22.71
0
0
0
0
0
0
—
F176V
12.0
−21.66
−21.26
1
1
1
0
1
0
Hb: Hydrogen bond.Docked conformations ofR- and S-methadone (stick presentation) with the highest binding affinity in the alleles of (a) CYP2B6, (b) CYP2D6, and (c) CYP3A4, in the active site (red surface presentation). Pink sticks are R-methadone; orange sticks are S-methadone; blue stick is R-methadone in wild type (WT) and green stick is S-methadone in WT. Heme is in blue stick presentation. Figure was generated using VMD 1.8.7 [9].
(a)
(b)
(c)Table2 summarized the interactions of methadone with CYPs. Most alleles consist of one π-π interaction and/or with additional π-δ interaction. Only a few CYP3A4 alleles formed hydrogen bond with R- or S-methadones while one (WT) was found in CYP2B6 and none in CYP2D6 alleles. Although there is no hydrogen bonding, π-π and/or π-δ interactions were observed for some R- or S-methadone conformations in the alleles, their free energy of binding (FEB) (Table 2) might be contributed by van der Waals and/or electrostatic interactions (indirectly measured by the number of residues appeared around the conformation) (see Figure S3).Estimated free energy of binding (FEB), dissociation constant(Kd), and binding constant (Kb) of methadone in CYP in can be further described by the function below:
(1)ΔG=RTlnKd,(2)Kd=1Kb,
where ΔG is FEB, R is gas constant of 1.985 cal/K/mol and T is absolute temperature which is assumed to be the room temperature of 298.15 K. More negative value of FEB will derive lower Kd value and thus presume higher Kb with higher metabolism rate of a substrate. Calculated FEB for most R- and S-methadones is in the range of ~−18.0 to −24.0 kcal/mol (Table 2). By correcting for the standard error by AutoDock4 (2.5 kcal/mol), FEB obtained by R-methadone with CYP2B6*11, CYP2B6*12, CYP2B6*18, and CYP3A4*12 was significantly different with that of WT (Table 2). As for S-methadone, the FEB obtained with CYP3A4*3, CYP3A4*11, and CYP3A4*12 can be distinguished with that of WT. For all three CYPs, CYP2B6, CYP2D6, and CYP3A4, binding affinity obtained by both R- and S-methadones with the WT was not significant (<2.5 kcal/mol). FEB for CYP2D6 alleles also did not show significant difference with that of WT. The observed less-significant differences in FEB for most R- and S-methadone in CYP2B6, CYP2D6, and CYP3A4 (Table 2) might be due to the SNPs located further from the active site (>5 Å from heme and beyond π-π and/or π-δ interactions). When a mutation occurred beyond the active/binding site, it possesses little or no effect on the methadone binding. There are reports suggesting that distanced mutation(s) from active site affect the metabolism rates which could be due to the open and close conformation of the channel and changes in hydrophobicity in the protein [40, 41]. However, the abovementioned effects were not studied here because of the current docking procedure limitations in handling protein folding.Binding affinity is inversely proportional toKd value.  Kd value can be further related to another function of Kd=koff/kon where the kinetics and thermodynamics of a simple one-step binding and one-step dissociation mechanism can be predicted.: When a conformation with lower Kd value, corresponding to stronger binding, the koff value thus is lower, and the occupancy time is longer [42]. In general, when the binding of a substrate towards a receptor is favorable, the substrate will be metabolised faster compared with that of unfavourable binding. In CYP-methadone perspective, when methadone is bound to CYP, it will be metabolised and thus the lower concentration of methadone will reach opioid receptors to generate the actual opioid effects [7]. Hence, in terms of MMT, the lower binding affinity of methadone in CYP is preferred, as lower dosage of methadone is needed for therapeutic effect.Data generated (Table2) especially for CYP2B6*11, CYP2B6*12, CYP2B6*18, CYP3A4*3, CYP3A4*11, and CYP3A4*12 cannot be compared directly to other studies or experimental data as these studies (http://www.cypalleles.ki.se/) were either done using other drug or/and focus on protein expression but not the study on CYP-methadone interaction. It is expected that the data generated by molecular modelling studies, for example, molecular docking simulation, are preliminary and experimental enzymatic assay which need to be conducted to further prove the concordant of the simulation data. Besides, the complexity in the involvement of more than one CYP isoform in methadone metabolism [6, 43], as well as the possibility that more than one substrate (other than methadone) can bind simultaneously within the same binding site [44, 45] will also need to be considered in overall in vivo methadone metabolism rate.


## 4. Conclusions

R- and S-methadones were docked into a similar active site for all CYPs and their alleles. The calculated free energy of binding was able to differentiate the effect of SNPs within 5 Å from heme towards methadone binding. R- and/or S-methadones with the higher binding affinity with CYP alleles were predicted to have higher methadone metabolism rate compared with that of WT, and vice versa. These preliminary predictions may possibly give some insights on optimum methadone dosing for individuals with SNPs. Ultimately, determination of optimum methadone dosage is important to ensure MMT can be continuously used as a potent corrective treatment for heroin addiction. This disease associated with heroin injection in hopes can be reduced.



---

*Source: 249642-2013-09-09.xml*