Metabolic regio- and stereoselectivity of cytochrome P450 2D6 towards 3,4-methylenedioxy-N-alkylamphetamines: In silico predictions and experimental validation

Article abstract of DOI:10.1021/jm050338+

A series of 3,4-methylenedioxy-N-alkylamphetamines (MDAAs) were automatically docked and subjected to molecular dynamics (MD) simulations in a cytochrome P450 2D6 (CYP2D6) protein model. The predicted substrate binding orientations, sites of oxidation, and relative reactivities were compared to the experimental data of wild-type and Phe¹²⁰Ala mutant CYP2D6. Automated docking results were not sufficient to accurately rationalize experimental binding orientations of 3,4-methylenedioxy-N-methylamphetamine (MDMA) in the two enzymes as measured with spin lattice relaxation NMR. Nevertheless, the docking results could be used as starting structures for MD simulations. Predicted binding orientations of MDMA and sites of oxidation of the MDAAs derived from MD simulations matched well with the experimental data. It appeared the experimental results were best described in MD simulations considering the nitrogen atoms of the MDAAs in neutral form. Differences in regioselectivity and stereoselectivity in the oxidative metabolism of the MDAAs by the Phe¹²⁰Ala mutant CYP2D6 were correctly predicted, and the effects of the Phe¹²⁰Ala mutation could be rationalized as well.

Full text of DOI:10.1021/jm050338+

J. Med. Chem. 2005, 48, 6117-6127
6117
Metabolic Regio- and Stereoselectivity of Cytochrome P450 2D6 towards
3,4-Methylenedioxy-N-alkylamphetamines: in Silico Predictions and
Experimental Validation
Peter H. J. Keizers,†,‡ Chris de Graaf,†,‡ Frans J. J. de Kanter,§ Chris Oostenbrink,‡ K. Anton Feenstra,‡
Jan N. M. Commandeur,‡ and Nico P. E. Vermeulen*,‡
Leiden Amsterdam Center for Drug Research (LACDR)/Division of Molecular Toxicology and Division of Organic and
Inorganic Chemistry, Department of Chemistry and Pharmacochemistry, Vrije Universiteit, De Boelelaan 1083,
1081 HV Amsterdam, The Netherlands
Received April 12, 2005
A series of 3,4-methylenedioxy-N-alkylamphetamines (MDAAs) were automatically docked and
subjected to molecular dynamics (MD) simulations in a cytochrome P450 2D6 (CYP2D6) protein
model. The predicted substrate binding orientations, sites of oxidation, and relative reactivities
were compared to the experimental data of wild-type and Phe120Ala mutant CYP2D6. Automated
docking results were not sufficient to accurately rationalize experimental binding orientations
of 3,4-methylenedioxy-N-methylamphetamine (MDMA) in the two enzymes as measured with
spin lattice relaxation NMR. Nevertheless, the docking results could be used as starting
structures for MD simulations. Predicted binding orientations of MDMA and sites of oxidation
of the MDAAs derived from MD simulations matched well with the experimental data. It
appeared the experimental results were best described in MD simulations considering the
nitrogen atoms of the MDAAs in neutral form. Differences in regioselectivity and stereoselec-
tivity in the oxidative metabolism of the MDAAs by the Phe120Ala mutant CYP2D6 were
correctly predicted, and the effects of the Phe120Ala mutation could be rationalized as well.
Introduction
refined, and validated experimentally. The homology
models of CYP2D6 suggest that substrates interact with
two or three aromatic/hydrophobic residues, e.g., Phe120
and Phe483, and a carboxylic acid residue, Glu216 or
Cytochromes P450 (CYPs) are heme containing en-
zymes which can be found in virtually all organisms.
This large family of enzymes is capable of oxidizing and
reducing a broad range of endogenous and exogenous
substrates, such as steroids, carcinogens, and drugs.1,2
In humans, one of the most relevant drug metabolizing
CYPs is CYP2D6. Although the expression levels of
CYP2D6 are only 2% of all hepatic CYPs, following
CYP3A4, it is the second most important drug metabo-
lizing enzyme, involved in the metabolism of about 30%
of the currently marketed drugs, including â-blockers,
neuroleptics, antidepressants, and antiarythmics.3-5
Large interindividual differences exist in CYP2D6
activity, due to gene multiplicity and polymorphisms,
thus further increasing its clinical importance.6,7 The
rationalization and prediction of potential CYP2D6
substrates is, therefore, advantageous in the discovery
and development of new drugs. Nowadays, the first
crystal structures of mammalian CYPs are becoming
available.8-10 However, so far, a crystal structure of
CYP2D6 remains to be resolved, and any structural
information on this enzyme still depends on homology
modeling, mutagenesis, and spectroscopic studies.
Asp301 12,13
The relevance of these residues was sup-
.
ported by site-directed mutagenesis studies,14-18 and
these interactions are consistent as well with those
derived from pharmacophore models of inhibitors and
substrates.19,20 Molecular modeling can provide infor-
mation on active site characteristics and the importance
of specific amino acid residues in enzyme-substrate
interactions, but it can also be used to rationalize and
predict regio- and stereoselectivity in metabolism by
CYPs. It has recently also been shown that automated
docking can successfully be applied to predict sites of
oxidation in substrates using CYP crystal structures.21
Molecular dynamics (MD) simulations in addition can
account for distributions of multiple binding conforma-
tions and, thus, give a more comprehensive description
of multiple sites of oxidation in substrates catalyzed by
CYPs.22-25
In this study, a molecular modeling approach has
been used to study the binding orientation and the sites
of oxidation of a series of 3,4-methylenedioxy-N-alkyl-
amphetamines (MDAAs), or XTC analogues, by CYP2D6.
The primary aim was to evaluate an integrated molec-
ular modeling approach to rationalize and predict
substrate binding and metabolism by CYP2D6. When
regioselectivity and stereoselectivity in drug metabolism
and the effect of active site mutations can be rational-
ized for a series of closely related compounds, such an
approach can be considered to be accurate and reliable.
CYP2D6 is one of the CYP isoforms studied most
extensively using molecular modeling.11 Several homol-
ogy model structures of CYP2D6 have been built,
* Author for correspondence: e-mail npe.vermeulen@few.vu.nl;
phone 0031 20 5987590; fax 0031 20 5987610.
† Both authors contributed equally.
‡ Division of Molecular Toxicology.
§ Division of Organic and Inorganic Chemistry.
10.1021/jm050338+ CCC: $30.25 © 2005 American Chemical Society
Published on Web 08/19/2005

6118 Journal of Medicinal Chemistry, 2005, Vol. 48, No. 19
Keizers et al.
Table 1. Estimated Dissociation Constants, Ks, (µM) of the
Substrates for Wild-Type and Phe120Ala Mutant CYP2D6a
The molecular modeling predictions were based on a
homology model of CYP2D6. Reactive substrate poses
were generated using automated docking and were
obtained from MD simulations of wild-type and mutant
CYP2D6 protein containing all three MDAAs. Two
different protonation states of the basic nitrogen atom
were studied, as it is known that this can drastically
influence the oxidation by CYP2D6.26 The computa-
tional predictions were validated experimentally. The
effects of mutating Phe120 into alanine on the binding
modes and metabolic products were studied for five
MDAAs, i.e., 3,4-methylenedioxy-N-methylamphetamine
(MDMA), its stereoisomers, 3,4-methylenedioxy-N-ethyl-
amphetamine (MDEA), and 3,4-methylenedioxy-N-prop-
ylamphetamine (MDPA). For the recently described
Phe120Ala mutant, large effects on the regioselective
oxidation of MDMA were demonstrated.18 Dissociation
constants for the three MDAAs were determined by
spectral titration. The metabolic products were ana-
lyzed, and furthermore, the binding orientations of
MDMA in the wild-type and mutant active sites were
explored experimentally by NMR spin lattice relaxation
rate measurements, a technique previously used to
determine hydrogen atom to heme distances in
CYPs.27-30
substrate
Ks, wild-type
Ks, Phe120Ala
R/S-MDMA
R-MDMA
S-MDMA
MDEA
28 ( 3
43 ( 6
40 ( 5
19 ( 2
18 ( 1
28 ( 2
53 ( 2
38 ( 4
15 ( 4
12 ( 2
MDPA
a All values are the means ( standard deviation (SD) of at least
two independent experiments, as described in the Materials and
Methods section. R/S-MDMA refers to the racemic mixture of 3,4-
methylenedioxy-N-methylamphetamine (MDMA), and R-MDMA
and S-MDMA refer to the pure R- and S-isomers of MDMA,
respectively. For the corresponding N-ethyl- and N-propylamphet-
amines (MDEA and MDPA), racemic mixtures were used.
Results
Figure 1. Scheme of the oxidation of 3,4-methylenedioxy-N-
alkylamphetamines (MDAAs) by Phe120Ala mutant CYP2D6.
The MDAAs are O-demethylenated to 3,4,-dihydroxy-N-alkyl-
amphetamine (1), N-dealkylated to 3,4-methylenedioxy-
amphetamine (2), N-hydroxylated to 3,4-methylenedioxy-N-
hydroxy-N-alkylamphetamine (3), or ω/ω-1-hydroxylated to
3,4-methylenedioxy-ω/ω-1-hydroxy-N-alkylamphetamine (4).
The C3-methyl is indicated, and R is methyl, ethyl, and propyl
for MDMA, MDEA, and MDPA, respectively.
Expression and Purification of Wild-Type and
Phe120Ala Mutant CYP2D6. The expression levels of
both wild-type and Phe120Ala mutant CYP2D6 were
similar as described previously.18 The enzymes which
contained a C-terminal 6 × histidine tag were purified
using Ni-NTA-agarose. Eluting the enzymes from the
Ni-NTA-agarose using imidazole led to a large amount
of enzyme in the high spin state, and even after
overnight dialysis, high spin remained the predominant
state. Therefore, 0.2 mM L-histidine was preferred as
eluting agent over imidazole, because it did not cause a
change in spin state after dialysis yielding pure enzymes
which were stable for at least 4 h at 24 °C in the
presence of 5% glycerol and 100 µM of substrate, as
assessed by CO-difference spectroscopy. Hardly any
P420 (i.e. inactive CYP detected as an absorbance
maximum at 420 nm in the CO-difference spectrum)
was present in the batches of purified enzyme.
Spectral Titration of Substrates to the Enzymes.
The addition of all MDAAs to purified wild-type or
Phe120Ala mutant CYP2D6 led to a type I change in the
visible absorbance spectrum. By titration of the com-
pounds to the solutions of enzyme, the spectral dissocia-
tion constants, Ks, could be estimated (using eq 1, data
shown in Table 1). While there were minor differences
in affinity between the wild-type and the mutant, for
both enzymes, a slightly lower Ks was found for MDEA
and MDPA than for MDMA. The enantiomers of MDMA
had about equal affinity for the wild-type, but a slight
preference toward S-MDMA was observed for the mu-
tant enzyme.
both enzymes. MDEA was O-demethylenated by wild-
type and mutant CYP2D6, and additionally, the mutant
N-dealkylated and N-hydroxylated MDEA. A fourth
oxidation product of MDEA by the mutant was detected
at a tR of 18.5 min, with a m/z of 222, fragmenting to
m/z ) 163. This product could be the aldehyde or keto
form of ω- or ω-1-hydroxylation, respectively.
MDPA was O-demethylenated by wild-type and the
mutant CYP2D6; the latter enzyme also N-dealkylated
and N-hydroxylated this substrate. A peak with a m/z
of 254 was found at 19.2 min with the mutant, indicat-
ing a double hydroxylation. No fragmentation took
place, so the positions of the hydroxyl groups remain
unknown.
Product formation from all substrates was linear for
at least 10 min for both enzymes under the conditions
used, so the kinetic parameters were determined for the
major products (Table 2). The overall rate of oxidation
by the Phe120Ala mutant was higher than that of wild-
type CYP2D6 for all MDAAs. The rates for O-de-
methylenation were two to five times higher for the
mutant, and the fact that more products are formed
indicates an even higher increase in oxidation rates. To
quantify the turnover rates, a synthetic reference
compound was only available for MDA; for the other
products, turnover rates were only compared in relative
terms. Because the fluorescence peak areas of equal
concentrations of MDA, MDMA, MDEA, and MDPA are
equal, it is likely that the N-alkyl chain has no influence
Metabolism of MDAAs. As reported previously,18
MDMA is oxidized by the Phe120Ala mutant to MDA and
N-OH-MDMA, in addition to 3,4-OH-MA, the only
product formed by the wild-type enzyme (Figure 1).
Elongation of the N-methyl chain to an ethyl or propyl
yielded the corresponding products after incubation with

CYP2D6 Catalytic Selectivity
Journal of Medicinal Chemistry, 2005, Vol. 48, No. 19 6119
Table 2. Product Formation of the Oxidation of MDAAs by Wild-Type and Phe120Ala Mutant CYP2D6d
wild-type
Phe120Ala
a
substrate
product
structure
Km
Vmax
Km
Vmax
R/S-MDMA
3,4-OH-MA
MDA
N-OH-MDMA
1
2
3
1.9 ( 0.5
c
c
1.2 ( 0.4
c
c
11.8 ( 3.8
14.5 ( 2.4
11.2 ( 1.4
2.2 ( 0.0a
3.5 ( 0.1b
7.1 ( 0.5a
R-MDMA
S-MDMA
MDEA
3,4-OH-MA
MDA
1
2
3
2.7 ( 0.4
1.2 (0.2
13.2 ( 0.6
26.1 ( 3.3
20.6 ( 2.5
2.4 ( 0.1a
c
c
c
c
5.9 ( 0.7b
N-OH-MDMA
10.7 ( 1.3a
2.4 ( 0.0a
3.8 ( 0.1b
8.4 ( 0.4a
2.9 ( 0.2a
4.2 ( 0.1b
1.4 ( 0.1a
3,4-OH-MA
MDA
N-OH-MDMA
1
2
3
3.3 ( 0.5
c
c
1.2 ( 0.0
c
c
4.9 ( 0.1
10.3 ( 0.6
6.9 ( 0.6
3,4-OH-EA
MDA
N-OH-MDEA
1
2
3
1.1 ( 0.1
c
c
1.2 ( 0.2
c
c
3.3 ( 0.0
4.5 ( 0.8
5.7 ( 0.3
MDPA
3,4-OH-PA
MDA
N-OH-MDPA
1
2
3
1.0 ( 0.1
c
c
1.3 ( 0.1
c
c
5.2 ( 0.6
9.1 ( 0.5
28.0 ( 6.5
7.5 ( 0.6a
7.3 ( 1.8b
3.0 ( 0.2a
b
a Vmax expressed in 1 × 105 fluorescence units min-1 nmol-1 CYP. Vmax expressed in min-1
.
c Not observed. d All values are the means
of at least two independent experiments ( SD, as described in the Materials and Methods section. Structure numbers refer to products
drawn in Figure 1; Km is expressed in µM.
Figure 3. Ratios of S-MDMA catalytic efficiency (Vmax/Km)
over R-MDMA catalytic efficiency (Vmax/Km) by wild-type (gray)
and by Phe120Ala mutant CYP2D6 (black) for the products 3,4,-
dihydroxy-N-alkylamphetamine (3,4-OH-MA), 3,4-methylene-
dioxyamphetamine (MDA), and 3,4-methylenedioxy-N-hydroxy-
N-alkylamphetamine (N-OH-MDMA).
enantiomers did not differ (Table 2). Eventually, these
changes lead to a different product ratio for the enan-
Figure 2. Vmax/Km product ratios (min-1 µM-1 for N-dealky-
lation, AU for the other products) for the O-demethylenation
tiomers: O-demethylenation is more prominent for
by wild-type (A), and O-demethylenation (B), N-dealkylation
S-MDMA than for R-MDMA, followed by N-hydroxyla-
(C), and N-hydroxylation (D) by Phe120Ala mutant CYP2D6 of
tion and N-dealkylation (Figure 3).
the different MDAAs.
NMR Spin Lattice Relaxation Rate Measure-
on the fluorescence yield of the compounds and, there-
fore, it is reasonable that this holds true for the
O-demethylenated and N-hydroxylated products of the
MDAAs as well.
ments. Enzyme bound MDMA hydrogen atom to heme
iron distances were determined by spin lattice relax-
ation NMR to validate the computationally predicted
distances. The 1H NMR spectrum of MDMA in the
presence of glycerol and each of the two enzymes was
well resolved; all signals of the different hydrogen atoms
were clearly visible (Figure 4). Because the doublet of
the C3-methyl (δ 1.18 ppm) and the singlets of the
N-methyl (δ 2.61 ppm) and methylene (δ 5.89 ppm)
hydrogen atoms could be best quantified and are the
most distal groups in MDMA, these were used to
determine the substrate orientation in the active site
of the two enzymes (inversion recovery spectra shown
in the Supporting Information). An increase in Tl of the
MDMA hydrogen atoms was found after reducing the
enzymes with dithionite, indicating that the paramag-
netic effect of the enzymes was diminished. From Tl,
the distances from the methylene, C3-methyl, and
N-methyl hydrogen atoms to the heme iron atom were
calculated (using eq 2) and tabulated (in Table 3). The
fast exchange condition (see the Materials and Methods
The turnover rates of MDMA, MDEA, and MDPA
were equal for wild-type CYP2D6 (Table 2). The cata-
lytic efficiency (Vmax/Km) indicated that MDPA is oxi-
dized most efficiently and MDMA is oxidized least
efficiently (Figure 2A). The mutant showed a higher
catalytic efficiency in O-demethylenating MDPA over
MDMA and MDEA, while the efficiency in N-hydroxy-
lation was the highest for MDMA (Figure 2B-D).
MDEA and MDPA were N-dealkylated more efficiently
than MDMA. So, elongation of the N-alkyl chain leads
to an altered catalytic regioselectivity by the mutant
enzyme.
Stereoselectivity in MDMA oxidation by wild-type
CYP2D6 was hardly found (Figure 3). However, the
mutant enzyme showed at least 2-fold lower Km values
for all three products of S-MDMA when compared to
R-MDMA, while the Vmax values for the products of both

6120 Journal of Medicinal Chemistry, 2005, Vol. 48, No. 19
Keizers et al.
Table 3. Hydrogen Atom to Heme Iron Distances (Å) of the MDMA Enantiomers to Wild-Type and Phe120Ala Mutant CYP2D6 as
Determined by NMR Spin Lattice Relaxation Rate Measurements and by MD Simulations of Neutral and Charged MDMAa
wild-type
Phe120Ala
position
method
NMR
MD (charged)
MD (neutral)
NMR
MD (charged)
MD (neutral)
NMR
MD (charged)
MD (neutral)
R-MDMA
S-MDMA
R-MDMA
S-MDMA
methylene (δ 5.89 ppm)
6.4
3.5
4.3
7.3
11.0
5.9
6.9
11.8
5.3
6.3
3.6
4.5
7.2
10.1
6.5
6.7
11.0
6.6
6.4
4.6
5.0
7.5
10.1
8.4
6.6
8.7
6.9
6.5
3.1
5.5
7.2
9.4
6.0
6.7
11.5
6.5
N-methyl (δ 2.61 ppm)
C3 (δ 1.18 ppm)
a NMR derived values are the means of at least two independent measurements with SD less than 5%, as described in the Materials
and Methods section. MD simulations derived values are averaged over time and for three individual runs starting from different docking
poses, according to eq 4.
starting conformations for MD calculations, and an
example of charged R-MDMA is shown (Figure 5). These
three distinct docked poses could be generally described
as one binding orientation corresponding to O-demeth-
ylenation, with the substrate nitrogen atom in close
contact with Glu216 (pose 1, generated with AutoDock),
another corresponding to O-demethylenation, with the
substrate nitrogen atom placed between Glu216 and
Asp301 (pose 2, generated with GOLD), and one corre-
sponding to N-dealkylation or N-hydroxylation (pose 3,
generated with AutoDock). In cases where pose 3 was
not observed for the substrate in the charged form (see
Figure 4. 1H NMR spectrum of 20 mM R-MDMA in deuter-
Table 4A), the corresponding pose of the neutral state
ated KPi-glycerol in the presence of 5 µM Phe120Ala mutant
was protonated and taken as a starting conformation
CYP2D6 showing the resonance assignments.
for MD simulations.
MD Simulations. Simulations of the enzymes with
substrates bound in the active site, starting from the
three different automatically docked poses, remained
stable during 10 ns of unrestrained MD at 300 K, with
final atom-positional root-mean-square deviation values
of 1.3-3.8 Å. Throughout the simulations, structure
determining hydrogen bonds were observed and the
majority of the orientations of the substrates in the
enzymes were within a reactive distance to the iron
atom (74% and higher, Table 4B). The substrates
appeared to have considerable freedom within the
binding site. Complete reorientation from the initial
coordinates was frequently observed, indicating that the
unrestrained MD simulations sampled the substrate
conformational space sufficiently. The distances be-
tween substrate hydrogen atoms and the heme iron
atom of wild-type and Phe120Ala mutant CYP2D6 were
derived for two different protonation states of the
MDMA enantiomers by averaging the distances over the
three simulations (using eq 4), and these were subse-
quently compared to experimentally determined dis-
tances (Table 3). Relatively large differences in averaged
hydrogen atom to iron distances were found within
individual MD simulations starting from different sub-
strate docking poses. However, averaged distances of
individual MD simulations of the same protonation state
and starting from the same docking pose were ap-
proximately equal. An exception was the difference in
averaged distance of the methylenedioxy hydrogen
atoms of both neutral MDMA enantiomers in simula-
tions starting from pose 3 (corresponding to N-dem-
ethylation/hydroxylation) in wild-type CYP2D6 (∼7.5 Å)
compared to the Phe120Ala mutant (∼11.0 Å). Also, the
averaged distances derived from all three simulations
section) was validated by measuring the temperature
dependence of RlP from 281 K to 307 K (see the
Supporting Information). For the three groups of hy-
drogen atoms studied, there was a linear increase of RlP
with an increase in the reciprocal temperature (1/T),
showing that the exchange between substrate molecules
in the enzyme active site and in solution is fast on the
NMR time scale. Only minor differences were observed
for the hydrogen atoms to iron distances between the
enantiomers of MDMA and each of the two enzymes.
The measured hydrogen atoms to iron distances indicate
that the methylene moiety is closest to the heme iron
at about 6.4 Å, the N-methyl is the most distal to the
heme at about 7.3 Å, and the C3 hydrogen atoms are at
a distance of about 6.7 Å.
Automated Docking Studies. All automatically
docked poses of the MDAAs, in two protonation states
and in both enzymes, were found to be predominantly
within a reactive distance of 6 Å from the heme iron
atom (69% and higher, Table 4A). Of these reactive
poses, the large majority corresponded to O-demethyl-
enation (62-100%). The relative occurrence of docking
poses corresponding to N-dealkylation (0-38%) and
N-hydroxylation (0-13%) of the MDAAs was found to
be enzyme, stereoisomer, and protonation state depend-
ent. In both enzymes, the neutral forms of the com-
pounds showed the widest range of reactive substrate
orientations, whereas the wild-type allowed for more
diversity. This observation is in contrast with the
experimental finding that wild-type CYP2D6 only cata-
lyzes O-demethylenation, while the Phe120Ala mutant
allows for alternative product formation.
Three distinct, energetically most favorable automati-
cally docked poses were selected of all MDAAs as

CYP2D6 Catalytic Selectivity
Journal of Medicinal Chemistry, 2005, Vol. 48, No. 19 6121
Table 4. Distributions of Substrate Orientations Corresponding to Specific Sites of Oxidation (%), According to Automated Docking
Studies (A) and MD Simulations (B)a
wild-type
NDM
Phe120Ala
NDM
substrate
MDMA
reactive
ODM
NOH
ωOH
reactive
ODM
NOH
ωOH
(A) Automated Docking Studies
R
S
R
S
R
S
charged
neutral
charged
neutral
charged
neutral
charged
neutral
charged
neutral
charged
neutral
100
88
92
99
91
100
69
99
70
92
74
82
100
86
0
10
0
14
1
9
1
0
4
0
-
-
-
-
0
0
0
0
0
0
0
0
100
93
100
80
100
100
100
100
99
100
90
100
99
99
98
99
93
99
97
62
73
0
3
0
7
0
1
0
0
0
0
0
0
0
0
-
-
-
-
0
0
0
0
0
0
0
0
100
73
0
13
0
0
MDEA
MDPA
99
1
91
0
2
99
0
1
78
22
21
38
0
0
7
79
0
1
62
0
100
95
100
3
100
99
0
38
27
1
0
(B) MD Simulations
MDMA
MDEA
MDPA
R
S
R
S
R
S
charged
neutral
charged
neutral
charged
neutral
charged
neutral
charged
neutral
charged
neutral
97
77
100
83
100
88
99
85
100
92
100
100
100
100
0
3
0
1
0
2
0
2
0
0
0
0
0
12
0
-
-
89
96
100
69
100
64
99
71
100
61
100
58
0
3
0
28
0
-
-
99
-
100
92
0
-
80
11
0
-
23
0
18
0
32
0
22
0
1
13
0
-
74
1
11
0
0
0
0
0
0
64
1
2
0
2
0
20
0
0
95
2
93
9
100
82
0
76
0
6
79
5
97
97
100
100
0
99
0
0
100
98
1
0
100
99
0
0
98
0
a The columns labeled “reactive” indicate the percentage of all binding orientations that displayed a site of oxidation within 6 Å of the
heme iron. The subsequent columns show the distribution of the reactive conformations for the different sites of oxidation, corresponding
to O-demethylenation (ODM), N-demethylation (NDM), N-hydroxylation (NOH), and “omega-hydroxylation” (ωOH). The automated docking
studies are based on 50 independent AutoDock and GOLD docking runs; the MD simulations are Boltzmann weighted averages over
three independent simulations starting from three different docking poses. Substrates were simulated in neutral and charged form.
of individual protonation states were approximately the
same for the different MDMA enantiomers in wild-type
and Phe120Ala mutant CYP2D6. The averaged hydrogen
atom to heme iron distances of the neutral MDMA
enantiomers agreed well with the distances determined
by NMR spin lattice relaxation rate measurements.
Apparently, considering the substrates to be neutral
describes the experimental situation best. With the
exception of the N-methyl and C3 hydrogen atoms of
R-MDMA in the Phe120Ala mutant, all hydrogen atom
to heme iron distances calculated from MD simulations
of neutral MDMA enantiomers were shorter than the
experimentally derived distances. The methylenedioxy
hydrogen atom to heme iron distances were slightly
lower than those of the N-methyl and C3 hydrogen
atoms. In contrast, the averaged N-methyl and C3
hydrogen atom to heme iron distances of the charged
MDMA enantiomers violated the experimental dis-
tances consistently.
For all MDAAs, differences between the averaged
enzyme-substrate interaction energies calculated using
SCORE were within 3 kJ mol-1 (∼kT) for all three MD
simulations of the same protonation states and within
5 kJ mol-1 for all six MD simulations. The relative
weights of the MD simulations of the three starting
orientations used in the Boltzmann averaging (eq 5, see
the Materials and Methods section) over the simulations
ranged between 0.11 and 0.54.
binding pocket, forming a hydrogen bond between its
protonated nitrogen and the carboxylate group of Glu216
(Table 4B and Figure 5). Relatively high probabilities
of O-demethylenation poses were found for neutral
MDAAs (>58%), although orientations corresponding to
N-demethylation/hydroxylation were also observed (ex-
cept for MDPA in wild-type CYP2D6). N-demethylation
or N-hydroxylation of neutral MDAAs had significantly
higher probabilities in the Phe120Ala mutant than in
wild-type CYP2D6. Furthermore, substrate dependent
catalytic regioselectivity was observed for the mutant:
elongation of the N-alkyl chain led to a decreased
N-hydroxylation efficiency with an optimal N-dealky-
lation efficiency for MDEA (Table 4B). These trends
were equal to those found experimentally and became
even more pronounced when the MDAAs were averaged
over enantiomers and protonation states and Boltzmann
weighted (Figure 6).
Significant stereoselective differences in the prob-
abilities of O-demethylenation vs N-demethylation/
hydroxylation were not observed for wild-type CYP2D6.
The Phe120Ala mutant, however, did in all cases dis-
criminate the enantiomers of the neutral MDAAs.
Differences between the probabilities of N-demethyla-
tion and N-hydroxylation of MDMA enantiomers were
found, and R-MDPA was found in orientations leading
to N-demethylation and ω-hydroxylation more often
than S-MDPA (Table 4B).
Already during the equilibration time or soon after
the start of the MD simulations with charged MDAAs,
reorientation of pose 3 into pose 1 (corresponding to
O-demethylenation) occurred. In the case of MDEA, this
led to nonreactive binding orientations outside of the
Enzyme-substrate interactions observed in the MD
simulations of neutral MDAAs showed a significant
decrease in van der Waals/aromatic interactions with
residue 120 in the Phe120Ala mutant compared to wild-
type CYP2D6 (Table 5). Furthermore, a small decrease

6122 Journal of Medicinal Chemistry, 2005, Vol. 48, No. 19
Keizers et al.
Figure 6. Probabilities (%) of reactive binding orientations
corresponding to O-demethylenation (gray), N-dealkylation
(white), N-hydroxylation (diagonal stripes), and ω-hydroxyla-
tion (black) of the MDAAs in wild-type (A) and in Phe120Ala
mutant (B) CYP2D6, as derived from 12 MD simulations.
Shown are the Boltzmann weighted averages of two different
protonation states and three different docking poses (examples
in Figure 5), assuming a 1:1 ratio between R- and S-enanti-
omers.
Table 5. Enzyme-Substrate Interactions (%) during MD
Simulationsa
wild-type
Phe120Ala
VdW/arom H-bond VdW/arom H-bond
substrate
MDMA
Phe120
Glu216
Ala120
Glu216
R
S
R
S
R
S
81
80
55
66
71
81
65
69
51
70
74
52
12
23
2
54
47
46
57
52
61
MDEA
MDPA
1
5
7
Figure 5. MD simulations starting orientations of charged
R-MDMA, generated by automated docking studies in wild-
type (A) and Phe120Ala mutant (B) CYP2D6 homology models.
Two orientations are corresponding to O-demethylenation,
with the substrate nitrogen atom either in close contact with
Glu216 (pose 1, in yellow) or placed between Glu216 and Asp301
(pose 2, in cyan). A third orientation is corresponding to
N-demethylation or N-hydroxylation (pose 3, in purple). Water
oxygen atoms, as predicted by GRID (see Materials and
Methods section), are depicted in blue.
a Enzyme-substrate interactions during the MD simulations
were monitored in terms of the occurrence of atom-atom distances
of 3.5 Å or less between the Phe/Ala120 side chain and the MDAA
ring system, representing aromatic/van der Waals interactions,
and hydrogen bonds between the Glu216 carboxylate oxygen atoms
and the MDAA nitrogen atom. Values of neutral and charged
substrates were averaged.
were found to form H-bonds with the substrate. At most,
one water molecule was observed in other regions of the
active site.
in hydrogen-bond interactions with the carboxylic acid
of Glu216 was observed for the Phe120Ala mutant as a
direct result of the relatively high probability of binding
orientations corresponding to N-dealkylation or N-
hydroxylation. In both wild-type and Phe120Ala mutant
CYP2D6, however, hydrogen-bond interactions were
found to be significantly decreased for neutral MDAAs
(up to 23% on average for each enantiomer) when
compared to charged MDAAs (up to 92% on average for
each enantiomer). During the MD simulations, water
molecules were primarily observed to fill the part of the
Discussion
The primary aim of this study was to develop an
integrated molecular modeling approach to analyze
substrate oxidation by CYP2D6. When regioselectivity
and stereoselectivity in metabolic oxidation and the
effect of enzyme active site mutations can be visualized
and quantified using a series of related compounds, such
an approach could be a useful tool in the rationalization
and prediction of metabolism. Using this approach, the
experimentally determined binding modes and the
product formation of a series of MDAAs by wild-type
and Phe120Ala mutant CYP2D6 could be rationalized.
Substrate Binding Orientations. In the present
study, molecular modeling was used to explain the
active site which is connected with a so-called “water
31
channel”,
which was already partially occupied by
predicted water molecules in the automated docking
starting structures (Figure 5). These water molecules

CYP2D6 Catalytic Selectivity
Journal of Medicinal Chemistry, 2005, Vol. 48, No. 19 6123
apparent discrepancy between experimentally observed
differences in regioselective oxidation and the stereo-
selective preference by the Phe120Ala mutant toward
MDMA (Figure 1) and the average substrate orientation
in the mutant and wild-type active sites as determined
by NMR spin lattice relaxation rate measurements
(Table 3). The average hydrogen atom to iron distances
calculated from the Tl of the MDMA enantiomers in both
enzymes were equal, with the methylene moiety in both
cases closest to the heme iron atom. Although the
distances between the hydrogen atoms at potential sites
of oxidation in MDMA and the heme iron atom were
somewhat larger than comparable distances found in
CYP X-ray structures (4-6 Å), they were within the
range of distances determined in other NMR studies.27-30
This difference can be explained by the fact that
experimentally determined distances are averages of all
possible substrate orientations during access to, resi-
dence in, and exit from the active site. As a result, they
contain more information than just that of a reactive
substrate orientation. NMR spin lattice relaxation rate
measurements have been used before to determine the
orientation of codeine in CYP2D6, where the measured
hydrogen atoms to iron atom matched the orientation
expected from codeine’s metabolic profile.28 Codeine is
a more rigid and bulky substrate than MDMA, which
may explain the good match. Caffeine, a small and low
affinity substrate of CYP1A2, was also a subject of NMR
spin lattice relaxation rate measurements, but the
orientation of caffeine found did not match the metabolic
profile.30 From other experimental studies, it is known
that substrates can be quite mobile in CYP active sites,32
and the existence of multiple substrate binding modes
and substrate mobility in CYP2D6 is clearly supported
by our molecular modeling studies. Averaged distances
derived from the present studies, considering three
distinct automatically docked substrate binding modes
of the MDAAs in combination with long MD simula-
tions, were in good agreement with the experimentally
determined metabolic profiles when considering the
MDAAs to be in their neutral form (Table 3). Still, this
does not explain the experimentally observed change in
regioselective oxidation and stereoselective preference
of the Phe120Ala mutant toward MDMA compared to the
wild-type. So, it is likely that a small substrate like
MDMA is very mobile in the CYP2D6 active site and
that the reactivity of possible sites of oxidation also
determines which products are being formed.
efficiency was that of MDEA. Stereoselectivity in oxida-
tion of MDMA by wild-type CYP2D6 was hardly found
experimentally or predicted by molecular modeling.
However, the experimentally observed preference for
S-MDMA by the Phe120Ala mutant could not be ratio-
nalized by molecular modeling. This may be due to the
fact that the stereoselectivity is predominantly caused
by a difference in Km and not in Vmax, while affinities of
the substrates to the enzymes were not predicted in a
rigorous manner in this study. On the other hand, the
small violations of the calculated to the experimentally
determined hydrogen atom to iron distances for neutral
R-MDMA bound to the Phe120Ala mutant do indicate
that this enantiomer experiences more difficulties in
finding its optimum orientation in the binding pocket.
Influence of Phe120Ala Mutation on Oxidation
Reactions. By removing the phenyl ring of Phe120, more
space is created in the active site of CYP2D6 and a
specific aromatic interaction point is eliminated.18 From
the MD simulations, it became clear that this space was
not occupied by the MDAAs or by water molecules, but
that critical substrate interactions were lost after
mutating Phe120 into an alanine (Table 5). This observa-
tion and the fact that the stereo- and regioselectivities
of the mutant enzyme differed from those of the wild-
type show that Phe120 is important for specific interac-
tions with the MDAAs in CYP2D6. There is not suffi-
cient space in the active sites of wild-type or the
Phe120Ala mutant CYP2D6 to accommodate two mol-
ecules of the MDAAs at the same time. However, it
could be possible that two substrate molecules bind
simultaneously to different parts of the protein, e.g., in
the active site and in substrate access channels31 or in
a peripheral binding site.10
Combination of Automated Docking and MD
Simulations. In earlier studies, it was shown that
automated docking can successfully be applied to predict
sites of oxidation by CYPs.21 In contrast, the current
study shows that automated docking alone is not
suitable for accurate determination of relative prob-
abilities of different substrate binding modes and for
discrimination between substrates and enzymes of high
similarity, such as regio/stereoisomers and mutants of
enzymes. To reproduce the experimental data, which are
an average over many different substrate and enzyme
orientations in time, a dynamic treatment of both
substrate and enzyme is required. Only when averaging
over different positions and orientations was taken into
account could experimentally determined substrate
hydrogen atoms to heme iron atom distances be repro-
duced accurately. The experimentally observed differ-
ences in substrate oxidation by wild-type and Phe120Ala
mutant CYP2D6 could only be explained from the
different binding modes observed in the MD simula-
tions. Multiple and extensive MD simulations are
needed to catch these subtle differences between sub-
strates, enzyme structures, and their dynamics.
Regioselectivity and Stereoselectivity in Oxida-
tion. Experimentally observed trends in regioselectivity
in the oxidation of the MDAAs by wild-type and Phe120
-
Ala mutant CYP2D6 (Table 2) were in good agreement
with the relative probabilities of different binding modes
of the MDAAs observed in the MD simulations (Table
4B). In the wild-type enzyme, only O-demethylenation
was observed experimentally, while in the Phe120Ala
mutant N-demethylation and N-hydroxylation also oc-
curred. The molecular modeling studies showed that
N-demethylation/hydroxylation probabilities for neutral
MDAAs were significantly higher in the Phe120Ala
mutant than in the wild-type enzyme (Table 4B). Both
the experimental and modeling studies showed a de-
creased N-hydroxylation efficiency with elongation of
the N-alkyl chain and that the highest MDA formation
Substrate Protonation States. The molecular mod-
eling was performed with neutral and charged sub-
strates because this difference markedly influenced the
binding and mobility of the MDAAs in the active site of
the two enzymes. Even though the MDAAs in physi-
ological solution will occur predominantly in their
charged form, the preferred protonation state within the

6124 Journal of Medicinal Chemistry, 2005, Vol. 48, No. 19
Keizers et al.
temperature on a Pharmacia Ultrospec 2000 spectrometer. In
short, 1 mL of 0.5 µM purified enzyme in a KPi-glycerol buffer
was divided over two cuvettes; to the first, 5 µL of a 1 mM
solution of compound in the same buffer was added, and to
the second, the same volume of just buffer was added.
Difference spectra from 350 to 450 nm were taken before and
after each addition of compound. The difference in absorbance
at 390 and 425 was plotted against the substrate concentration
to estimate the dissociation constant (Ks) using eq 1 in Graph
Pad Prism 4.0 software, where B is the absorbance difference
(390-425) and Bmax is the absorbance difference when [S] )
∞.
binding pocket can differ significantly; thus, reductions
of apparent pKa values for substrates upon binding to
CYP2D6 of up to 2 pH units have been reported.26 From
a direct comparison between the experimentally derived
hydrogen atom to heme iron distances and the compu-
tationally determined values, it became clear that the
MD simulations reproduced the experimental data best
when the MDAAs were represented in their neutral
form.
Conclusions
The approach of combining automated docking and
MD simulations in a protein model was a successful way
to rationalize the oxidative metabolism of a series of 3,4-
methylenedioxy-N-alkylamphetamines (MDAAs) by
CYP2D6. Differences in oxidative metabolism of these
closely related substrates could be predicted in details
such as regioselectivity and stereoselectivity. Further-
more, the effects of the active site Phe120Ala mutation
on the selectivity of the substrate could be rationalized.
The presented integrated modeling method is a promis-
ing tool in the prediction of the metabolic properties of
new (druglike) compounds with respect to interactions
with CYP2D6 and might be applicable to other drug
metabolizing enzymes as well.
Bmax[S]
B )
(1)
Ks + [S]
Metabolism of MDAAs. Metabolic reactions were carried
out in 200 µL of 50 mM KPi (pH 7.4) and 5 mM MgCl2
supplemented with 10 concentrations ranging from 0 to 200
µM of one of the MDAAs and E. coli membranes corresponding
to 25 nM wild-type or Phe120Ala mutant CYP2D6. After 5 min
of preincubation at 37 °C, the reactions were initiated with
an NADPH regenerating system, resulting in final concentra-
tions of 0.2 mM NADPH, 0.3 mM glucose-6-phosphate, and
0.4 units mL-1 glucose-6-phosphate dehydrogenase. The reac-
tion was allowed to proceed for 10 min before it was stopped
by the addition of 20 µL of 23% HClO4. After centrifugation
(10 min, 6800g), 30 µL aliquots of the supernatant were
analyzed by high-performance liquid chromatography (HPLC).
To measure O-demethylenation, analytes were separated
isocratically using a C18 column (Phenomenex Luna 5u 150
mm × 4.6 mm), with a mobile phase consisting of 22%
acetonitrile and 0.1% triethylamine, set to pH 3 with HClO4
at a flow rate of 0.6 mL min-1. Other oxidation products were
found using the same column with a gradient (A, 5% ACN with
20 mM ammonium acetate; B, 90% ACN with 10 mM am-
monium acetate). All oxidation products were detected by
fluorescence (λex ) 280 nm, λem ) 320 nm) and identified by
liquid chromatography-mass spectroscopy (LC-MS). Peak
areas of all products were quantified with Shimadzu Class VP
4.3 software and analyzed using nonlinear regression with one
site binding in Graph Pad Prism 3.0.
LC-MS. To identify the metabolic products of the MDAAs,
incubations were carried out for 15 min as described above
with 100 µM MDMA, MDEA, or MDPA and E. coli membranes
corresponding to 50 nM wild-type or Phe120Ala mutant CYP2D6.
Volumes of 100 µL of supernatant were injected and separated
using a C18 column (Phenomenex Luna 5u 150 mm × 4.6 mm)
with a flow rate of 0.6 mL min-1 and analyzed by MS. The
analytes were eluted using a gradient starting with a 5% ACN
eluens, supplemented with 20 mM ammonium acetate for 7
min, then, increasing linearly to 90% ACN with 10 mM
ammonium acetate in 14 min, and remaining there for 5 min.
APCI positive ionization was used on a LCQ Deca mass
spectrometer (Thermo Finnigan, Breda, The Netherlands) with
a vaporizer temperature of 450 °C, N2 as a sheath (40 psi) and
an auxiliary gas (10 psi), a needle voltage of 6000 V, and a
heated capillary at 150 °C.
Materials and Methods
Materials. The pSP19T7LT plasmid containing human 2D6
with a C-terminal His6-tag bicistronically coexpressed with
human NADPH-cytochrome P450 reductase was kindly pro-
vided by Prof. Dr. Ingelman-Sundberg. 3,4-Methylenedioxy-
N-alkylamphetamines (MDAAs) were synthesized as described
before.33,34 S- and R-MDMA‚HCl were obtained from the
Division of Neuroscience and Behavioral Research of the
National Institute on Drug Abuse (Bethesda, MD). Emulgen
911 was purchased at KAO Chemicals (Tokyo, Japan). Ni-
NTA-agarose was from Qiagen (Westburg, Leusden, The
Netherlands). All other chemicals were of analytical grade and
obtained from standard suppliers.
Expression and Purification of CYP2D6. Both the
Phe120Ala mutant and the wild-type pSP19T7LT plasmids
were transformed into Escherichia coli, strain JM109. Expres-
sion and membrane isolation were carried out as described.18
Membranes were resuspended in 0.4% of the original culture
volume of a KPi-glycerol buffer (50 mM potassium phosphate
buffer, pH 7.4, 10% glycerol), and CYP contents were deter-
mined by CO-difference spectra.35
Enzymes were solubilized from membranes by stirring in
KPi-glycerol supplemented with 0.5% Emulgen 911 for 2 h
at 4 °C. Insoluble parts were removed by centrifugation (60
min, 120000g at 4 °C). The supernatant was incubated, while
gently rocking, with Ni-NTA-agarose for 30 min at 4 °C. The
Ni-NTA-agarose was retained in a polypropylene tube with a
porous disk (Pierce, Perbio Science, Ettenleur, The Nether-
lands) and washed with a KPi-glycerol buffer containing 2
mM histidine. CYP2D6 was eluted with 0.2 M histidine. After
overnight dialysis in the KPi-glycerol buffer, the sample was
concentrated on a Vivaspin 20 filtration tube (10.000 MWCO
PES, Sartorius, Nieuwegein, The Netherlands). For NMR spin
lattice relaxation rate measurements, the buffer was ex-
changed for a deuterated KPi-glycerol buffer by repeated
rounds of adding a larger volume of deuterated buffer and
concentrating. The deuterated KPi-glycerol buffer was made
by repeated rounds of evaporating water or D2O from a KPi
buffer and dissolving the residue again in D2O. Then 5%
glycerol was added, and traces of iron were removed with
Chelex 100 (Biorad, Richmond, CA).
Products of MDMA metabolism: 3,4-OH-MA (tR 9.5 min, m/z
182, MS/MS m/z 151), MDA (tR 17.0 min, m/z 180, MS/MS m/z
163), and N-OH-MDMA (tR 20.4 min, m/z 210, MS/MS m/z
163). Products of MDEA metabolism: 3,4-OH-EA (tR 13.0 min,
m/z 196, MS/MS m/z 151), MDA (tR 17.0 min, m/z 180, MS/
MS m/z 163), N-OH-MDEA (tR 21.6 min, m/z 224, MS/MS m/z
163), and a fourth product for the mutant (tR 18.5 min, m/z
222, MS/MS m/z 163). Products of MDPA metabolism: 3,4-
OH-PA (tR 16.0 min, m/z 210, MS/MS m/z 151), MDA (tR 17.0
min, m/z 180, MS/MS m/z 163), N-OH-MDPA (tR 23.4 min, m/z
238, MS/MS m/z 163), and a fourth product for the mutant
(tR 19.2 min, m/z 254, no fragmentation observed).
NMR Spin Lattice Relaxation Rate Measurements.
Distances from substrate hydrogen atoms to the enzyme heme
iron atom were calculated from their longitudinal relaxation
rates determined using 1H NMR.27-30 1H NMR measurements
Optical Titrations. Dissociation constants of substrates
to the enzymes were determined by spectral titration, accord-
ing to the method of Jefcoate.36,37 Spectra were taken at room

CYP2D6 Catalytic Selectivity
Journal of Medicinal Chemistry, 2005, Vol. 48, No. 19 6125
were performed on a Bruker Avance 400 MHz spectrometer;
signals were referenced to HDO at 4.75 ppm. Longitudinal
relaxation times (Tl) of substrate hydrogen atoms were mea-
sured with a standard inversion recovery sequence of 180°-
τ-90°. The spectra were recorded with 20 values of τ, ranging
from 0.1 to 12.5 s, and a relaxation delay of 5-10 times Tl.
The temperature was kept at 300 K, except for in the variable
temperature experiments. Tl was determined by plotting peak
areas against delay times, τ, and by exponential fitting (It )
CYP2D6 homology model using AutoDock and GOLD. Active
site water molecules, predicted using a protocol based on the
program GRID48 and a ligand-based cut off, were included. The
preparation of enzyme and substrate structures in the auto-
mated docking studies was performed using default settings
and predicted water molecules, as described before.21 For each
substrate, 50 independent AutoDock and GOLD docking runs
were performed using default parameter values. Docking poses
generated by both docking algorithms were clustered (maxi-
mum root-mean-square differences between members of the
cluster of 1 Å) and rescored with the SCORE scoring function.49
The three energetically most favorable, distinct docking poses
of each MDAA species were selected as starting structures for
MD calculations, yielding in total 12 (3 docking poses times 2
enatiomers times 2 protonation states) simulations for each
substrate in both wild-type and Phe120Ala mutant CYP2D6.
MD Simulations. Substrate conformations selected from
the automated docking studies were used as starting struc-
tures for MD simulations, using the GROMACS molecular
simulation package50 and the GROMOS-96 force field, param-
eter set 43A1.51,52 Heme parameters were taken from this
parameter set, and additionally, a covalent bond was defined
between the heme iron atom and the sulfur atom of Cys443 with
an ideal bond distance of 0.240 nm. The enzyme, including
substrate and heme, was energy minimized in a vacuum using
the steepest descent method, first, with harmonic position
restraints using a force constant of 1000 kJ mol-1 nm-2 on all
non-hydrogen atoms and, then, without position restraints.
Subsequently, preequilibrated water was added in a dodeca-
hedral box with a minimum distance of 1 nm between the
enzyme and box edges, followed by minimization. For 1 ps,
the water was allowed to relax, while the whole enzyme was
position restrained. Then, first, only the side chains of the
residues not involved in substrate binding were released for 1
ps, and next, also the backbone of these residues was released
for 10 ps. Finally, only the backbone of the binding residues
was restrained for 100 ps. The equilibration scheme was
followed by an unrestrained 10 ns MD production run. MD
simulations were carried out with a time step of 2 fs, a twin-
range cutoff of 0.8/1.4 nm, and a relative dielectric constant
of 1.0. The LINCS algorithm was used to constrain the length
of all covalent bonds.53 Pressure was maintained at 1 bar and
temperature at 300 K by weak coupling to an external bath,54
with relaxation times of 0.1 and 1.0 ps, respectively. Coordi-
nates were saved at a 1.0 ps interval for subsequent analysis.
Enzyme-substrate interactions were estimated using the
SCORE scoring function and averaged over the trajectories.
Analysis of MD Simulations. Hydrogen atoms to iron
atom distances, rMD, for different protonation states of the
MDMA enantiomers were derived from MD simulations by
averaging distances (r) observed in the three independent
simulations with different starting positions as
A + Be-τ/T ) with Bruker XWINNMR software. The measured
l
Tl was transformed into the longitudinal paramagnetic relax-
3+
ation rate (RlP) according to eq 2. Tl,Fe is the relaxation time
-1
RlP ) (Tl,Fe3+
)
-1 - (Tl,Fe2+
)
(2)
+CO
of substrate hydrogen atoms in the presence of ferric enzymes,
2+
and Tl,Fe +CO is the relaxation time in the presence of ferrous
CO-bound enzymes. The hydrogen atoms to iron atom distance,
rNMR, can be calculated from RlP using the simplified Solomon
and Bloembergen equation27,38 (eq 3)
r
NMR ) C(RRlPF(τc))1/6
(3)
where C is a constant taking into account the spin state as
well as other nuclear and electronic factors, with a value of
813 for high spin CYP, 30 and R is the fraction of enzyme-bound
substrate, determined as the ratio of the concentration of
enzyme and the sum of the substrate concentration and its
spectral dissociation constant (Ks). F(τc) is a function of the
correlation time, τc, for which the value of 3 × 10-10 s is taken,
as determined for CYP2D6 by others.28 In the calculations of
hydrogen atoms to heme distances according to eqs 2 and 3,
the assumption has been made that on the NMR time scale
there is a fast exchange between substrate molecules in
solution and in the enzyme active site. This was verified by
measuring the temperature dependence of 1/Tl, which should
show a linear increase with the reciprocal value of the
temperature.29,30
In a typical experiment, 5 µM enzyme and 20 mM substrate
were used in a volume of 500 µL. Recording the spectra for
the determination of relaxation times with oxidized enzyme
took about an hour; then, sodium dithionite was added, and
the sample was treated for 30 s with CO to do the measure-
ments with the reduced CO-bound enzyme.
Homology Modeling. A protein homology model of CYP2D6
was constructed based on the crystal structures of dimethyl-
sulfophenazole derivative and diclofenac-bound rabbit CYP2C5
(PDB codes 1N6B and 1NR6, respectively).39,40 Homology
modeling, model refinement, and model validation were per-
formed according to the approach described previously.12 In
short, one hundred CYP2D6 models were generated with the
restraint-based comparative modeling program Modeler,41
using the same amino acid alignment as before.12 Subse-
quently, three models were selected with the best loop
conformations, as determined by visual inspection, stereo-
chemical parameters using PROCHECK,42 and side-chain
environment using Errat43 and Verify3D.44 Homology models
had approximately the same protein quality check scores as
the crystal structure templates. Finally, one single model was
selected which could accommodate codeine best in the experi-
mental binding orientation.28 This final model was validated
on its ability to reproduce substrate binding orientations
corresponding to metabolic profiles with AutoDock,45 FlexX,46
and GOLD47 automated docking studies (De Graaf, unpub-
lished results). The final model of wild-type CYP2D6 was used
as a template for modeling of the Phe120Ala mutant. The Phe
residue at position 120 was mutated to Ala using the homology
module of Insight II (Biosym, San Diego, CA), after which an
energy minimization and a 1 ps MD simulation with position
restraints on the protein backbone were carried out as
described in MD Simulations of the Results section.
-1/6
r
MD ) ( r-6
)
(4)
These distances could be directly compared to the experimen-
tally obtained distances, rNMR (eq 3).
Sites of oxidation were predicted based on enzyme-
substrate conformations generated by automated docking and
MD simulations. Distances of 6.0 Å or less between the heme
iron and an atom in the substrate (e.g., corresponding to
O-demethylenation, N-dealkylation, N-hydroxylation, or ω-hy-
droxylation) were considered as reactive in the sense that
activation of the C-H bond by the heme-Fe-O moiety is
possible.55 Site of oxidation predictions were performed to
determine both substrate reactivity (expressed in the number
of reactive MD configurations out of the total) and regioselec-
tivity. The catalytic regioselectivity was calculated by deter-
mining which site of oxidation was closest to the heme iron
atom for each reactive docking pose or MD configuration.
Enzyme-substrate interactions during the MD simulations
were monitored in terms of the occurrence of atom-atom
distances of 3.5 Å or less between the Phe/Ala120 side chain
and the MDAA ring system (representing aromatic/van der
Automated Docking Studies. The R- and S-enantiomers
of MDMA, MDEA, and MDPA in their neutral (basic) and
charged (acidic) forms were docked into the active site of the

6126 Journal of Medicinal Chemistry, 2005, Vol. 48, No. 19
Keizers et al.
Waals interactions) and hydrogen bonds between the Glu216
carboxylate oxygen atoms and the MDAA N-atom. Substrate
reactivity, regioselectivity, and enzyme-substrate interactions
were calculated for any given protonation state and stereo-
and regioisomer, as a Boltzmann weighted average over the
three independent MD simulations with different substrate
starting orientations. The average of an individual simulation
i was weighted by a relative probability, pi, which was
calculated from the average enzyme-substrate interaction
estimates as obtained from the program SCORE (EiSCORE).
(14) Ellis, S. W.; Hayhurst, G. P.; Smith, G.; Lightfoot, T.; Wong, M.
M.; Simula, A. P.; Ackland, M. J.; Sternberg, M. J.; Lennard,
M. S.; Tucker, G. T.; et al. Evidence that aspartic acid 301 is a
critical substrate-contact residue in the active site of cytochrome
P450 2D6. J. Biol. Chem. 1995, 270, 29055-29058.
(15) Guengerich, F. P.; Hanna, I. H.; Martin, M. V.; Gillam, E. M. J.
Role of glutamic acid 216 in cytochrome P450 2D6 substrate
binding and catalysis. Biochemistry 2003, 42, 1245-1253.
(16) Paine, M. J.; McLaughlin, L. A.; Flanagan, J. U.; Kemp, C. A.;
Sutcliffe, M. J.; Roberts, G. C.; Wolf, C. R. Residues glutamate
216 and aspartate 301 are key determinants of substrate
specificity and product regioselectivity in cytochrome P450 2D6.
J. Biol. Chem. 2003, 278, 4021-4027.
(17) Flanagan, J. U.; Marechal, J. D.; Ward, R.; Kemp, C. A.;
McLaughlin, L. A.; Sutcliffe, M. J.; Roberts, G. C.; Paine, M. J.;
Wolf, C. R. Phe120 contributes to the regiospecificity of cyto-
chrome P450 2D6: mutation leads to the formation of a novel
dextromethorphan metabolite. Biochem. J. 2004, 380, 353-360.
(18) Keizers, P. H. J.; Lussenburg, B. M. A.; de Graaf, C.; Mentink,
L. M.; Vermeulen, N. P. E.; Commandeur, J. N. M. Influence of
phenylalanine 120 on cytochrome P450 2D6 catalytic selectivity
and regiospecificity: crucial role in 7-methoxy-4-(aminomethyl)-
coumarin metabolism. Biochem. Pharmacol. 2004, 68, 2263-
2271.
(19) Ekins, S.; de Groot, M. J.; Jones, J. P. Pharmacophore and three-
dimensional quantitative structure activity relationship methods
for modeling cytochrome P450 active sites. Drug Metab. Dispos.
2001, 29, 936-944.
(20) Vermeulen, N. P. E. Prediction of drug metabolism: the case of
cytochrome P450 2D6. Curr. Top. Med. Chem. 2003, 3, 1227-
1239.
(21) de Graaf, C.; Pospisil, P.; Pos, W.; Folkers, G.; Vermeulen, N. P.
E. Binding mode prediction of cytochrome P450 and thymidine
kinase protein-ligand complexes by consideration of water and
rescoring in automated docking. J. Med. Chem. 2005, 48 (7),
2308-2318.
i
i
SCORE/kB
SCORE/kBT
pi ) (e-∆E
T)/( e-∆E
)
(5)
∑
i
In this equation, kB is the Boltzmann constant (8.3441 kJ mol-1
K-1), and T is the temperature (300 K).
Acknowledgment. We thank Prof. Dr. Magnus
Ingelman-Sundberg and Dr. Mats Hidestrand for pro-
viding the pSP19T7LT plasmid containing CYP2D6 and
the human NADPH-CYP reductase. We thank Dr. Hari
Sing from the NIDA-NIH for providing the MDMA
enantiomers. We thank Ed Groot and Jon de Vlieger
for technical assistance.
Supporting Information Available: Figures showing the
1H NMR spectra of racemic MDMA with varying values of τ
and the temperature dependence of T1 for specified hydrogen
atoms of MDMA in wild-type CYP2D6. This material is
available free of charge via the Internet at http://pubs.acs.org.
(22) Harris, D.; Loew, G. Prediction of regiospecific hydroxylation of
camphor analogs by cytochrome-P450(Cam). J. Am. Chem. Soc.
1995, 117, 2738-2746.
References
(23) Harris, D. L.; Loew, G. H. Investigation of the proton-assisted
pathway to formation of the catalytically active, ferryl species
of P450s by molecular dynamics studies of P450eryF. J. Am.
Chem. Soc. 1996, 118, 6377-6387.
(24) Keseru, G. M.; Kolossvary, I.; Bertok, B. Cytochrome P-450
catalyzed insecticide metabolism. Prediction of regio- and ste-
reoselectivity in the primer metabolism of carbofuran: A theo-
retical study. J. Am. Chem. Soc. 1997, 119, 5126-5131.
(25) Park, J. Y.; Harris, D. Construction and assessment of models
of CYP2E1: Predictions of metabolism from docking, molecular
dynamics, and density functional theoretical calculations. J.
Med. Chem. 2003, 46, 1645-1660.
(26) Miller, G. P.; Hanna, I. H.; Nishimura, Y.; Guengerich, F. P.
Oxidation of phenethylamine derivatives by cytochrome P450
2D6: the issue of substrate protonation in binding and catalysis.
Biochemistry 2001, 40, 14215-14223.
(27) van de Straat, R.; de Vries, J.; de Boer, H. J.; Vromans, R. M.;
Vermeulen, N. P. E. Relationship between paracetamol binding
to and its oxidation by two cytochromes P-450 isozymessa
proton nuclear magnetic resonance and spectrophotometric
study. Xenobiotica 1987, 17, 1-9.
(1) Goeptar, A. R.; Scheerens, H.; Vermeulen, N. P. E. Oxygen and
xenobiotic reductase activities of cytochrome P450. Crit. Rev.
Toxicol. 1995, 25, 25-65.
(2) Guengerich, F. P. Cytochrome P450: what have we learned and
what are the future issues? Drug Metab. Rev. 2004, 36, 159-
197.
(3) Shimada, T.; Yamazaki, H.; Mimura, M.; Inui, Y.; Guengerich,
F. P. Interindividual variations in human liver cytochrome P450
enzymes involved in the oxidation of drugs, carcinogens and toxic
chemicals: studies with liver microsomes of 30 Japanese and
30 Caucasians. J. Pharmacol. Exp. Ther. 1994, 270, 414-423.
(4) Bertilsson, L.; Dahl, M. L.; Dalen, P.; Al-Shurbaji, A. Molecular
genetics of CYP2D6: clinical relevance with focus on psycho-
tropic drugs. Br. J. Clin. Pharmacol. 2002, 53, 111-122.
(5) Zanger, U. M.; Raimundo, S.; Eichelbaum, M. Cytochrome P450
2D6: overview and update on pharmacology, genetics, biochem-
istry. Naunyn-Schmiedeberg’s Arch. Pharmacol. 2004, 369, 23-
37.
(6) Oscarson, M. Pharmacogenetics of drug metabolising enzymes:
importance for personalised medicine. Clin. Chem. Lab. Med.
2003, 41, 573-580.
(7) Ingelman-Sundberg, M. Pharmacogenetics of cytochrome P450
and its applications in drug therapy: the past, present and
future. Trends Pharmacol. Sci. 2004, 25, 193-200.
(8) Williams, P. A.; Cosme, J.; Sridhar, V.; Johnson, E. F.; McRee,
D. E. Mammalian microsomal cytochrome P450 monooxyge-
nase: structural adaptations for membrane binding and func-
tional diversity. Mol. Cells 2000, 5, 121-131.
(9) Williams, P. A.; Cosme, J.; Ward, A.; Angove, H. C.; Vinkovic,
D. M.; Jhoti, H. Crystal structure of human cytochrome P450
2C9 with bound warfarin. Nature 2003, 424, 464-468.
(10) Williams, P. A.; Cosme, J.; Vinkovic, D. M.; Ward, A.; Angove,
H. C.; Day, P. J.; Vonrhein, C.; Tickle, I. J.; Jhoti, H. Crystal
structures of human cytochrome P450 3A4 bound to metyrapone
and progesterone. Science 2004, 305, 683-686.
(11) de Graaf, C.; Vermeulen, N. P. E.; Feenstra, K. A. Cytochrome
P450 in silico: an integrative modeling approach. J. Med. Chem.
2005, 48, 2725-2755.
(12) Venhorst, J.; ter Laak, A. M.; Commandeur, J. N.; Funae, Y.;
Hiroi, T.; Vermeulen, N. P. Homology modeling of rat and human
cytochrome P450 2D (CYP2D) isoforms and computational
rationalization of experimental ligand-binding specificities. J.
Med. Chem. 2003, 46, 74-86.
(13) Kemp, C. A.; Flanagan, J. U.; van Eldik, A. J.; Marechal, J. D.;
Wolf, C. R.; Roberts, G. C.; Paine, M. J.; Sutcliffe, M. J.
Validation of model of cytochrome P450 2D6: an in silico tool
for predicting metabolism and inhibition. J. Med. Chem. 2004,
47, 5340-5346.
(28) Modi, S.; Paine, M. J.; Sutcliffe, M. J.; Lian, L. Y.; Primrose, W.
U.; Wolf, C. R.; Roberts, G. C. A model for human cytochrome
P450 2D6 based on homology modeling and NMR studies of
substrate binding. Biochemistry 1996, 35, 4540-4550.
(29) Poli-Scaife, S.; Attias, R.; Dansette, P. M.; Mansuy, D. The
substrate binding site of human liver cytochrome P450 2C9: an
NMR study. Biochemistry 1997, 36, 12672-12682.
(30) Regal, K. A.; Nelson, S. D. Orientation of caffeine within the
active site of human cytochrome P450 1A2 based on NMR
longitudinal (T1) relaxation measurements. Arch. Biochem.
Biophys. 2000, 384, 47-58.
(31) Wade, R. C.; Winn, P. J.; Schlichting, I.; Sudarko A survey of
active site access channels in cytochromes P450. J. Inorg.
Biochem. 2004, 98, 1175-1182.
(32) Prasad, S.; Mitra, S. Role of protein and substrate dynamics in
catalysis by Pseudomonas putida cytochrome P450cam. Bio-
chemistry 2002, 41, 14499-14508.
(33) Braun, U.; Shulgin, A. T.; Braun, G. Centrally active N-
substituted analogues of 3,4-methylenedioxyphenylisopropy-
lamine (3,4-methylenedioxyamphetamine). J. Pharm Sci. 1980,
69, 192-195.
(34) Onderwater, R. C.; Venhorst, J.; Commandeur, J. N. M.; Ver-
meulen, N. P. E. Design, synthesis, and characterization of
7-methoxy-4-(aminomethyl)coumarin as a novel and selective
cytochrome P450 2D6 substrate suitable for high-throughput
screening. Chem. Res. Toxicol. 1999, 12, 555-559.

CYP2D6 Catalytic Selectivity
Journal of Medicinal Chemistry, 2005, Vol. 48, No. 19 6127
(35) Omura, T.; Sato, R. The carbon monoxide-binding pigment of
liver microsomes. II. Solubilization, purification, and properties.
J. Biol. Chem. 1964, 239, 2379-2385.
(46) Rarey, M.; Kramer, B.; Lengauer, T.; Klebe, G. A fast flexible
docking method using an incremental construction algorithm.
J. Mol. Biol. 1996, 261, 470-489.
(47) Jones, G.; Willett, P.; Glen, R. C.; Leach, A. R.; Taylor, R.
Development and validation of a genetic algorithm for flexible
docking. J. Mol. Biol. 1997, 267, 727-748.
(48) Goodford, P. J. A computational procedure for determining
energetically favorable binding sites on biologically important
macromolecules. J. Med. Chem. 1985, 28, 849-857.
(36) Jefcoate, C. R. Measurement of substrate and inhibitor binding
to microsomal cytochrome P450 by optical-difference spectros-
copy. Methods Enzymol. 1978, 52, 258-279.
(37) Hanna, I. H.; Krauser, J. A.; Cai, H.; Kim, M. S.; Guengerich,
F. P. Diversity in mechanisms of substrate oxidation by cyto-
chrome P450 2D6. Lack of an allosteric role of NADPH-
cytochrome P450 reductase in catalytic regioselectivity. J. Biol.
Chem. 2001, 276, 39553-39561.
(38) Solomon, I.; Bloembergen, N. Nuclear magnetic interactions in
the HF molecule. J. Chem. Phys. 1956, 25, 261-266.
(39) Wester, M. R.; Johnson, E. F.; Marques-Soares, C.; Dansette,
P. M.; Mansuy, D.; Stout, C. D. Structure of a substrate complex
(49) Wang, R. X.; Liu, L.; Lai, L. H.; Tang, Y. Q. SCORE: A new
empirical method for estimating the binding affinity of
protein-ligand complex. J. Mol. Model. 1998, 4, 379-394.
(50) Lindahl, E.; Hess, B.; van der Spoel, D. GROMACS 3.0:
a
a
package for molecular simulation and trajectory analysis. J. Mol.
Model. 2001, 7, 306-317.
(51) Daura, X.; Mark, A. E.; van Gunsteren, W. F. Parametrization
of aliphatic CHn united atoms of GROMOS96 force field. J.
Comput. Chem. 1998, 19, 535-547.
of mammalian cytochrome P450 2C5 at 2.3
A resolution:
evidence for multiple substrate binding modes. Biochemistry
2003, 42, 6370-6379.
(40) Wester, M. R.; Johnson, E. F.; Marques-Soares, C.; Dijols, S.;
Dansette, P. M.; Mansuy, D.; Stout, C. D. Structure of mam-
malian cytochrome P450 2C5 complexed with diclofenac at 2.1
A resolution: evidence for an induced fit model of substrate
binding. Biochemistry 2003, 42, 9335-9345.
(52) Scott, W. R. P.; Hunenberger, P. H.; Tironi, I. G.; Mark, A. E.;
Billeter, S. R.; Fennen, J.; Torda, A. E.; Huber, T.; Kruger, P.;
van Gunsteren, W. F. The GROMOS biomolecular simulation
program package. J. Phys. Chem. A 1999, 103, 3596-3607.
(53) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M.
LINCS: A linear constraint solver for molecular simulations.
J. Comput. Chem. 1997, 18, 1463-1472.
(54) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.;
DiNola, A.; Haak, J. R. Molecular dynamics with coupling to an
external bath. J. Chem. Phys. 1984, 81, 3684-3690.
(55) Meunier, B.; de Visser, S. P.; Shaik, S. Mechanism of oxidation
reactions catalyzed by cytochrome p450 enzymes. Chem. Rev.
2004, 104, 3947-3980.
(41) Sali, A.; Blundell, T. L. Comparative protein modelling by
satisfaction of spatial restraints. J. Mol. Biol. 1993, 234, 779-
815.
(42) Laskowski, R. A.; Macarthur, M. W.; Moss, D. S.; Thornton, J.
M. Procheck - a program to check the stereochemical quality
of protein structures. J. Appl. Crystallogr. 1993, 26, 283-291.
(43) Colovos, C.; Yeates, T. O. Verification of protein structures -
patterns of nonbonded atomic interactions. Protein Sci. 1993,
2, 1511-1519.
(44) Luthy, R.; Bowie, J. U.; Eisenberg, D. Assessment of protein
models with 3-dimensional profiles. Nature 1992, 356, 83-85.
(45) Morris, G. M.; Goodsell, D. S.; Halliday, R. S.; Huey, R.; Hart,
W. E.; Belew, R. K.; Olson, A. J. Automated docking using a
Lamarckian genetic algorithm and an empirical binding free
energy function. J. Comput. Chem. 1998, 19, 1639-1662.
JM050338+