Both Reactivity and Accessibility Are Important in Cytochrome P450

Article Cite This: J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

pubs.acs.org/jcim

Both Reactivity and Accessibility Are Important in Cytochrome P450 Metabolism: A Combined DFT and MD Study of Fenamic Acids in BM3 Mutants Rasmus Leth,† Bogac Ercig,†,‡ Lars Olsen,†,§ and Flemming Steen Jørgensen*,† †

Department of Drug Design and Pharmacology, Faculty of Health and Medical Sciences, University of Copenhagen, Universitetsparken 2, DK-2100 Copenhagen, Denmark

J. Chem. Inf. Model. Downloaded from pubs.acs.org by WEBSTER UNIV on 02/14/19. For personal use only.

S Supporting Information *

ABSTRACT: Cytochrome P450 102A1 from Bacillus megaterium (BM3) is a fatty acid hydroxylase that has one of the highest turnover rates of any mono-oxygenase. Recent studies have shown how mutants of BM3 can produce metabolites of known drug compounds similar to those observed in humans. Single-point mutations in the binding pocket change the regioselective metabolism of fenamic acids from aromatic hydroxylation to aliphatic hydroxylation. This study is concerned with the individual contribution from accessibility and reactivity for drug metabolism with a future goal to develop fast methods for prediction. For a BM3 M11 mutant as well as the M11 V87F and M11 V87I mutants, we studied the metabolism of the nonsteroidal anti-inflammatory drugs (NSAIDs) mefenamic acid, meclofenamic acid, tolfenamic acid, and diclofenac. Density functional theory (DFT; B3LYP and B3LYP-D3) calculations for all possible reactions were performed using a porphyrin model reacting with the four substrates. Molecular dynamics (MD) simulations were used to determine the potential sites of metabolism that are accessible. Finally, we combine reactivity and accessibility for each potential site to interpret the experimentally determined metabolism. Generally, the 3 and 5 positions (on the ring containing the acidic substituent) and the 2′, 3′, and 4′ positions are most reactive, whereas 4, 5, 3′, and 4′ are most accessible. Combining reactivity and accessibility show that the 5, 3′, and 4′ positions are predicted to be most prone to be metabolized, in agreement with experimentally observed data. Reactivity seems to be the dominant factor in the CYP-mediated metabolism of these NSAIDs, which is consistent with previously published methods based solely on reactivity.

■

INTRODUCTION Cytochrome P450 (CYP) are a family of monooxygenase enzymes involved in the metabolism of endogenous and exogenous drugs, e.g., drug compounds.1 The metabolism of drug compounds can lead to the formation of reactive metabolites and inhibition of cytochromes.2 It is therefore important to get a detailed understanding of the mechanisms of drug metabolism to accurately predict the fate of drug compounds when they enter the human body. CYPs are also present in plants and many bacteria, where they synthesize numerous molecules, metabolites, which are important for plant development, defense, and for bacterial growth.3,4 Cytochrome P450 102A1 (BM3) from Bacillus megaterium is a fatty acid hydroxylase. In nature the BM3 P450 system is fused with the NADPH-P450 reductase domain into a single polypeptide and displays the highest catalytic rate of any mono-oxygenase.5 Mutants of the bacterial cytochrome BM3 are able to produce human metabolites of some drug compounds.6−8 Among these compounds are the fenamic acid analogs mefenamic acid, meclofenamic acid, tolfenamic acid, and diclofenac (Figure 1), commonly used nonsteroidal anti-inflammatory drugs (NSAIDs) because of their ability to inhibit the COX enzymes.9,10 Metabolism of these compounds © XXXX American Chemical Society

is highly relevant as it can lead to the generation of reactive quinone imines.11 Interestingly, point mutations in the binding pocket of the BM3 mutant called M11 can change the regioselectivity and the product distribution of fenamic acid metabolites.12 The three major metabolites for fenamic acids in BM3 are the monohydroxylated products 3′-OH-methyl-, 4′OH-, or 5-OH-fenamic acid. Generally, for all BM3 mutants aliphatic hydroxylation at the 3′ position and aromatic hydroxylation in the 4′ position are favored over aromatic hydroxylation in the 5 position. The BM3 M11 mutant slightly favors the 4′ aromatic reaction, and the BM3 M11 V87F mutant greatly favors the 4′ aromatic reaction. In contrast, the BM3 M11 V87I mutant favors 3′-OH aliphatic hydroxylation. There are many different in silico methods for predicting P450 metabolite generation.13,14 Docking and molecular dynamics (MD) studies are typically used to model the binding mode of the substrate and thus identify which sites are accessible for the oxo-ferryl of the heme group. Density functional theory (DFT), on the other hand, has been used to study how likely the subsequent reaction is to occur by Received: October 26, 2018

A

DOI: 10.1021/acs.jcim.8b00750 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling

Figure 1. Three fenamic acid analogues (mefenamic acid, meclofenamic acid, and tolfenamic acid) and diclofenac studied in this paper. Numbering scheme is shown for the mefenamic acid. Product distributions (in percentages) and site-of-metabolism data for the three fenamic acid analogues are from Venkataraman et al.12 for the M11/M11 V87F/M11 V87I mutants, respectively.

corrected for zero-point vibrations at the bs1 level. All transition state barriers are given relative to the separated substrate and Cpd I model. All DFT calculations were performed with Turbomole.28 Docking and MD Simulations. MD simulations were based on a crystal structure of BM3 M11 (PDB code 5E9Z). BM3 M11 V87F and BM3 M11 V87F mutants were generated with Maestro in the Schrödinger package version 9.3.29 All protein structures were prepared with Protein Preparation Wizard.30 The ligands mefenamic acid, meclofenamic acid, tolfenamic acid, and diclofenac were built and prepared with Ligprep.30 Pyramidalization of the amines was treated in accordance with published crystal structures.31 Both neutral and charged ligands were generated and subsequently docked into each protein structure using Glide version 5.8 with default settings (van der Waals scaling factor and partial charge cutoff 0.80 and 0.15, respectively; standard precision and flexible ligand sampling; no constraints).32 We chose to continue with charged ligands for further study based on published pKa values33 and consistency with previously published studies on mefenamic acid in BM3.12,34 The tleap program of the Amber12 package was used to create prmtop and prmcrd, topology and coordinate files of the proteins, respectively, with the ff99SB force field.35,36 The HEME was prepared using parameters from Shahrokh et al.37 The system was neutralized with counterions. The tool ACPYPE was used to convert prmtop and prmcrd to GROMACS format, and GROMACS 5.0.2 was used for periodic box generation, solvation (approximately ∼13 500 solvent molecules), equilibration, and MD simulations.38−41 All protein structures were solvated in a truncated octahedron box with TIP3P water. Berendsen thermo- and barostats were used to maintain the temperature and pressure of the system.42 A time step of 2 fs was used while van der Waals and shortrange electrostatic interactions were evaluated explicitly within a 9 Å cutoff. Long-range electrostatic interactions was included using the particle mesh Ewald method.43 Equilibration of the systems includes minimization of the solvent, gradual heating from 0 to 300 K in the NVT ensemble, and 1 ns simulation at 300 K and 1.01325 bar in the NPT ensemble. Finally, MD production runs of 10 ns were carried out. The MD simulations were analyzed with VMD 1.9.1.44

determining the activation barriers. For the human CYP enzymes, DFT calculations have been used for prediction of likely drug metabolites.14−18 We have previously studied the metabolism of mefenamic acid by the three previously mentioned BM3 M11 mutants by subsequent application of docking, MD simulations, and free energy calculations. It turned out that it was necessary to refine the initial docking poses by MD simulations in order to obtain binding modes, which explained the formation of experimentally observed metabolites. The free energy calculations enabled a correct quantification of the metabolites. In this study, we have extended and refined the approach by considering not only mefenamic acid but also three related compounds. DFT calculations and MD simulations were used to study the regioselective metabolism of the three fenamic acids and the structurally related diclofenac by three different BM3 mutants. DFT, at the B3LYP and B3LYP-D3 level, was used to characterize all possible oxidation reactions for the four substrates to identify the most reactive sites. Subsequently, docking was used to identify which binding modes are possible, and MD was used on the top three binding modes. The MD trajectories were used to determine the accessibility of each potential site of metabolism for each ligand in each of the BM3 mutants. Finally, we combined the reactivity calculated with DFT and accessibility from the MD simulations and compared with the experimental product distributions in order to determine the relative importance of reactivity and accessibility in CYP-mediated metabolism.

■

METHODS DFT Calculations. Compound I (Cpd I) was modeled as porphyrin without side chains and with CH3S− and O2− as axial ligands. All calculations were done at the DFT/B3LYP level of theory. Both doublet and quartet spin states were calculated for all complexes. Transition state structures were confirmed by frequency calculations. For calculations with dispersion the doublet spin state was calculated with B3LYPD3.19,20 The geometry optimizations and frequency calculations were performed with the Ahlrichs VDZ basis set of Schäfer et al.,21 enhanced with a p function on the Fe atom and the 6-31G(d) basis set on the other atoms (bs1).22−24 Final energies were determined with a larger 6-311++G(2d,2p) basis set, except Fe, for which we used the Ahlrichs VDZ basis set of Schäfer et al.,25,26 enhanced with s, p, d, and f functions (bs2).27 Energies were corrected for zero-point vibrational energy. Unless otherwise stated, we report the bs2 energies

■

RESULTS DFT-Determined Activation Barriers. The activation barrier for oxidation of all possible sites of metabolism for B


Article


Table 1. B3LYP Activation Barriers for Sites of Metabolism of Mefenamic Acid, Meclofenamic Acid, Tolfenamic Acid, and Diclofenaca reaction type

position

mefenamic acid (kJ/mol)

meclofenamic acid (kJ/mol)

tolfenamic acid (kJ/mol)

diclofenac (kJ/mol)

H-abstraction rebound aromatic aromatic aromatic aromatic aromatic aromatic aromatic aliphatic aliphatic/aromatic

amine amine 3 4 5 6 4′ 5′ 6′ Me-2′ Me-3′/3′

57 59 75 91 67 81 72 81 81 62 66

69 86 87 94 73 80 86 100

59 76 80 96 70 97 79 86 80 68

41 65 84 87 73 82 81 91b

70

73

a

Calculated as transition state energy relative to the energy of the reactant complex. bThe lone-pair of the N atom is on the same side as the attack of the oxo-ferryl group.

turn lowers the barrier for aromatic hydroxylation in the 4′ position. Interestingly, the barrier for the 3′ position in diclofenac is also low (63 kJ/mol), which is most likely caused by the presence of a chlorine in both the ortho and the para positions. The 5′ position has a higher activation barrier (91 kJ/mol), although it also has two chlorine atoms in ortho and para positions. However, the lone pair of the N atom is oriented differently and is positioned on the same side of the aromatic ring as the attack from the oxo-ferryl group, which increases the barrier. In tolfenamic acid the chloride substitution is in the 3′ position and thus ortho to the reaction in the 4′ position, where it has both a stabilizing electronic and a destabilizing steric effect, resulting in a small increase in activation barrier compared to mefenamic acid. We found a large activation barrier for tolfenamic acid in the 5′ position, because it has the chlorine and the amine in meta position. Generally, the reaction in the 6′ position is sterically hindered in the transition state. For the 6′ position to become available for the Cpd I oxygen, the two aromatic rings have to rotate, which reduces the conjugation across the amine and the electronic effect is lost. In the other ring, the lowest barriers are found for the 3 and 5 positions, again due to the stabilizing effect of the amine. These barriers are slightly higher for diclofenac, most likely because the presence of the carboxylic acid, now one carbon atom further from the ring, leads to a less pronounced stabilizing effect. We have calculated activation barriers for the most reactive positions (Table 2) with B3LYP-D3, a functional containing a correction for dispersive forces. The dispersive forces favor van der Waals (vdW) interactions and thus stabilize the transition state complex more than the separated species.47 Accordingly, we observe very low activation barriers with B3LYP-D3 as experienced previously.48 Generally, the D3 functional decreases the activation barrier for attack on the aromatic carbon atom compared to the aliphatic hydroxylation reaction, thus decreasing the difference in energy between these two reactions. Because the substrate is closer to Cpd I in the transition state for aromatic hydroxylation, the favorable dispersive forces increase further. This phenomenon is more pronounced for mefenamic acid compared to meclofenamic acid, which could be due to more attractive dispersive forces of the methyl group (Me-2′) compared to the chlorine in the 2′ position (see Table 2). To summarize the DFT and DFT-D3 data, the most reactive sites in the four substrates studied are

mefenamic acid, meclofenamic acid, tolfenamic acid, and diclofenac were determined using DFT calculations at the B3LYP level (Table 1 and Figure 2 and Figures S1−S3, Supporting Information). The activation energies for aliphatic hydroxylation in positions 2′ and 3′ are between 50 and 63 kJ/mol (Table S1, Supporting Information), which is similar to what has been determined previously for, e.g., toluene.15 Thus, the methyl and amine substituents do not affect the hydrogen abstraction energies much. For oxidation of the amines, both direct oxidation and hydrogen abstraction from the amine were investigated. Hydrogen abstraction from the amine is more favorable, and in most cases the rebound reaction is the ratedetermining step, which is consistent with previous findings.45 The activation barriers for the hydrogen abstraction and the subsequent rebound are higher than what has been reported for other amines. For example, for N-methylaniline the TS energy relative to the reactant state is 11 and 48 kJ/mol for amine hydrogen abstraction and the subsequent rebound, respectively.45 For mefenamic acid, the activation barriers are 57 and 59 kJ/mol (relative to the separated species, cf. Table S1 in Supporting Information). This is probably due to limited accessibility of the NH group in these substrates. Additionally, the internal hydrogen bond between the carboxylic acid group and the NH group is broken in the TS, thus further increasing the difference in energy between the TS and the reactant complex and the substrate (Figure 2). The aromatic sites in the 5 and 4′ positions, para position to the amine, have the lowest activation barriers of the aromatic hydroxylation reactions. The barriers of about 75 kJ/mol correspond to previously determined barriers and are significantly lower than an unsubstituted benzene ring, which has a barrier of 87 kJ/mol.16 This effect was previously studied and showed that an amine in the para position decreases the activation barrier.46 Electronically the 3 and 6′ positions, which are ortho to the amine, could also be expected to have low activation barriers, although interactions with the neighboring substituent can have an effect as well. Indeed, the aromatic 3 position also has a low activation barrier in mefenamic and tolfenamic acid. The addition of chloride in the 2′ and 6′ positions (i.e., meta to 4′) only has a marginal effect on the activation barriers for aromatic hydroxylation in the 4′ position, in agreement with published structure−activity relationships for aromatic hydroxylation.46 This effect is less evident for diclofenac because it has no 3′ substituent, which in C


Article


Figure 2. Transition state structures for the N-hydroxylation (amine), aromatic hydroxylation (3, 4, 5, 6, 4′, 5′, and 6′), and aliphatic hydroxylation (Me-2′ and Me-3′) reactions of mefenamic acid. Labels below the structures refer to the site of metabolism, cf. numbering in Figure 1.

Table 2. Gas-Phase B3LYP-D3 Activation Barriers for the Sites of Metabolism of Mefenamic Acid, Meclofenamic Acid, Tolfenamic Acid, and Diclofenaca reaction type

position

mefenamic acid (kJ/mol)

meclofenamic acid (kJ/mol)

tolfenamic acid (kJ/mol)

diclofenac (kJ/mol)

aliphatic aromatic aromatic

Me-3′ 4′ 5

65 62 51

70 63 65

64 75

71 66

a

Calculated as transition state energy relative to the energy of the reactant complex. The reported energies are doublet spin state energies calculated with bs1 and correction for zero-point vibration.

ligand in the proteins. From these studies we identified three possible binding modes corresponding to hydroxylation at the aliphatic sites in the 3′ position for mefenamic and meclofenamic acids (color coded blue in figures and tables) and aromatic hydroxylations in positions 4′ and 5′ (color

the methyl groups in the 2′ and 3′ positions and the aromatic carbon atoms in the 3, 5 and 4′ positions. Criteria Used in the Analysis of the MD Simulations. The substrates were initially docked into BM3 mutants M11, M11 V87F, and M11 V87I, to study the orientation of the D


Article


Figure 3. Docking poses of mefenamic acid in BM3 M11 representing aliphatic hydroxylation (left) and aromatic hydroxylation on both aromatic rings (center and right). Poses correspond to the initial structures for the subsequent MD simulations. Docking poses for meclofenamic acid, diclofenac, and tolfenamic acid are shown in Figure S4, Supporting Information. See text for the blue, red, and green color coding.

coded red) and 4, 5, 6 (color coded green) for all four substrates (see Figure 3). These binding modes were among the top three ranking poses obtained and match previously published binding modes for mefenamic acid in BM3.34 No binding modes representing amine oxidation were identified, because the binding pocket of BM3 is too narrow and makes it impossible for the amine to reach the heme group. Because metabolism of each NSAID in each of the BM3 mutants is experimentally documented,12 it was assumed that each ligand binds and thus ligand channeling was not studied. For each of the binding modes, short molecular dynamics simulations were performed to evaluate the accessibility of each site of metabolism. Equilibrated sampling was confirmed by measured RMSD values (see Table S4). Geometrical information of the DFT structures was used to assess whether a reaction is possible for the given binding mode in the MD simulations. The O−C distances in the reactant complexes from the DFT calculations are on average 2.6 ± 0.1 and 3.4 ± 0.1 Å for the aliphatic and aromatic reactions, respectively (Table S2). In our analysis of the MD trajectories we used a cutoff value for the O−C distance of 3.5 Å, as well as extended criteria of 4.5 and 5.5 Å, as it was often observed that the O−C distance was longer during the MD simulations. The observation of extended distances is not surprising, as atoms in a classical force field cannot overlap in a complex that will lead to a chemical reaction. The values have been used in the literature. Lonsdale et al.49 used a Fe−C distance criterion of 4 Å for modeling of aromatic hydroxylation of diclofenac, while de Graaf et al.50 used a distance criterion of 6 Å to the heme Fe atom to select for potential sites of metabolism. A range of criteria has therefore been applied in this study to investigate how these affect the outcome. The orientation of the C atoms corresponding to those observed in the TS structures was used to ensure that a reaction is possible. The Fe−O−C angles are 120 ± 1° for the aliphatic hydroxylation and 130 ± 7° for the aromatic hydroxylation (Table S3). The angular dependence of these two Fe−O−C angles was investigated for both types of reactions (Table 3). From the measured angles in the DFT calculations we developed angle criteria for analysis of the reactive conformations in the MD simulations. The aliphatic hydroxylation reaction is more dependent on the value of

Table 3. Gas-Phase Energies Associated with Increasing the Angle between Cpd I Iron, Cpd I Oxygen, and the Site of Metabolism Carbon (Fe−O−C) in Mefenamic Acid for the Aliphatic and Aromatic Hydroxylation Reactions in the Transition State Complexa aliphatic hydroxylation Fe−O−C angle (deg) 120 (freely optimized) 130 140 150

relative energy (kJ/mol) 0 3 15 18

aromatic hydroxylation Fe−O−C angle (deg) 129 (freely optimized) 140 150 160

relative energy (kJ/mol) 0 3 6 6

a

The energies are calculated relative to the transition state complex with freely optimized angles. The reported energies are bs1 energies with no corrections.

this angle, and therefore, a more narrow interval was applied (115−130°), while a larger interval was applied for the aromatic hydroxylation, since it seems less affected by a change in this angle (115−160° was applied). We did not use an angle criteria lower than 115°, and because of the narrow binding pocket we rarely found any frames in which this occurs. From three DFT calculations a Fe−O−N angle of 140° in the transition state geometries of amine oxidation was observed. This angle is larger than what has previously been described,51 probably due to limited accessibility of the amine of fenamic acids. The orientation of the methyl group in the 3′ position in the docking (Figure 3, left) indicates that the hydrogen abstraction occurs from the H atom perpendicular to the phenyl ring, whereas the H atom in the plane of the phenyl ring seems to be abstracted in the aromatic binding mode (Figure 3, center). To quantify this difference in reactivity between geometries of the DFT calculations and the docking or MD simulations, we performed a scan of the orientation of the CH2 group on a mefenamic acid radical. The frozen radical method has previously been shown to correlate well with TS energies for hydrogen abstraction reactions.15 We find that abstraction of the H atom in the plane of the phenyl ring is associated with an energy increase on 38 kJ/mol, which most likely would result in a nonproductive binding mode (Figure 4). Another geometrical effect is the orientation of the nitrogen lone pair with respect to the aromatic hydroxylation reaction. E


Article


Figure 4. Molecular representation of a frozen mefenamic acid radical after 3′-aliphatic hydroxylation. The graph shows the stability of the frozen radical with the CH2 group constrained at different torsional angles (dotted line). Energies are bs1 energies with no corrections.

Figure 5. Molecular representation of the model system used for evaluating the effect of varying the torsional angle around the amine (dotted line). The graph shows the energies relative to the freely optimized complex. Energies are bs1 energies with no corrections.

Figure 6. Histograms from the MD simulations of mefenamic acid, meclofenamic acid, tolfenamic acid, and diclofenac in three BM3 mutants, showing the percentages of frames in which the given atom is fulfilling both distance (3.5 Å) and angle criteria. For the histograms for the distance criteria on 4.5 and 5.5 Å, see Figure S8 in Supporting Information. For each ligand the top plot corresponds to BM3 M11, middle plot to BM3 M11 V87F, and bottom plot to BM3 M11 V87I. Different colors refer to data from simulations started from the binding poses corresponding to the structures shown in Figure 3.

F


Article


Figure 7. Snapshots from the MD simulations of mefenamic acid in the three different binding modes for the three BM3 mutants. Corresponding representations of meclofenamic acid, tolfenamic acid, and diclofenac are shown in Figures S5−S7, Supporting Information. See text for the blue, red, and green color coding.

the accessibilities if a O−C distance criterion of 4.5−5.5 Å is used instead. Meclofenamic acid displays accessibilities comparable to mefenamic acid (Figure 6) with M11 and M11 V87I being similar and different from M11 V87F. For both mefenamic and meclofenamic acid, the aromatic binding mode has the Me-3′ group pointing toward the FeO group in an orientation with the H atom that is not favorable for the hydrogen abstraction, i.e., in the plane of the phenyl ring of the substrate. This is in contrast to the aliphatic binding mode, in which the H atom being abstracted is perpendicular to the phenyl ring that typically would have a lower activation barrier (Figures 3 and 4). By relating the calculated activation barriers shown in Figure 4 and the binding modes observed in docking and MD simulations, we can quantify which reactions are more likely in the given conformation. Our findings suggest that the aliphatic binding mode will most likely be associated with hydroxylation at the 3′ position, while the aromatic binding mode will most likely be associated with hydroxylation at the 4′ and 5 positions. The orientation of the ligand for aliphatic hydroxylation is different for the three BM3 mutants (Table 4). For mefenamic acid, the orientation in the V87F mutant favors aliphatic hydroxylation. These results are contradicting the experimental results in which aliphatic hydroxylation is more favorable in the V87I mutant. This tendency is reversed for the simulations with meclofenamic acid, in which the V87I mutation imposes a favorable conformation of the ligand compared with the other two mutations. This is consistent with experimental results and

For example, the difference in the 3′ and 5′ hydroxlation reaction for diclofenac (Table 1) is a consequence of this effect. Figure 5 shows that rotation of the torsional angle around the amine is associated with an energy increase of 13.7 kJ/mol. Thus, although this effect is smaller than the effect from the orientation of the H atom to be abstracted (Figure 4), it is not negligible. Determination of Accessibility. For mefenamic acid the shortest O−C distance criteria of 3.5 Å show that the Me-3′ and 4′ positions are easily accessible and that the 3, 4, 5, and 6 positions, at least to a certain extent, fulfill the distance and angle criteria (Figure 6). There is a tendency that the 4 position is more accessible for the FeO than the other C atoms in that ring. Interestingly, the V87F mutant generally seems to have fewer structures that fulfill the distance and angle criteria than the two other mutants. For example, there are very few structures fulfilling the aromatic criteria for 4′ hydroxylation. As shown in Figure 7, π−π stacking is observed between the aromatic ring of the ligand and the phenylalanine side chain of residue 87. This could prevent the substrate from moving closer to the FeO group, but sterically the reaction is not hindered. This finding is in agreement with the V87F mutant being a back mutation, toward the wild type, which does not bind drug-like compounds.6 The number of productive binding modes of mefenamic acid, of course, increases if distance criteria of 4.5 or 5.5 Å are used (cf. Figure S8). The majority of frames still have the Me-3′ and 4′ fulfilling the criteria, whereas in the other ring, the 3, 4, 5, and 6 positions are all accessible. There is not a dramatic change in G


Article

Journal of Chemical Information and Modeling Table 4. Aliphatic Hydroxylation Ratiosa ligand

BM3 M11

BM3 M11 V87F

BM3 M11 V87I

mefenamic acid meclofenamic acid

12 (40) 81 (36)

50 (21) 63 (9)

22 (75) 99 (63)

seems that the aliphatic reaction is easily accessible. However, the torsional angle between the methyl group of the ligand and the reactive center reveals that there are only a few favorable orientations for the hydrogen abstraction. For the aromatic binding mode there are, as noted previously, also favorable orientations of the methyl group for the hydrogen abstraction reaction. However, the orientation of the methyl group follows the same trend in both binding modes for mefenamic acid (Table 4). It is also observed that the orientation is favorable most of the time in the V87F mutant. The binding energy of the aliphatic and aromatic binding modes of mefenamic acid in the BM3 M11, M11 V87F, and M11 V87I mutants has recently been studied by free energy calculations.34 Here it was shown by thermodynamic integration that the aromatic binding mode is favored over the aliphatic by roughly 10 and 13 kJ/mol in the M11 and V87F mutants, respectively. The aliphatic binding mode was shown to be favored over the aromatic by roughly 9 kJ/mol in the V87I mutant. From these data it is evident that the binding energy, at least in part, determines whether the aliphatic or aromatic reaction is preferred. The second aromatic binding mode corresponding to hydroxylation in the 5 position was not investigated in this study. The 5 position of all ligands has a low activation barrier and is easily accessible in the molecular dynamics simulations. However, the hydroxylated product in this position ranks very low in the experimental product formation ratios. One reason for this is that we have no measure for the frequency of occurrence for each binding mode. The fenamic acid substrates might bind with the carboxylic acid moiety pointing away from the reactive center of the enzyme more frequently, e.g., due to the interaction with Ser72 or because it is stabilized by the solvent. Lonsdale et al. previously discussed the binding energy for different poses.49 Comparison of the three mutants reveals that the most notable difference in the MD accessibilities is that the V87F mutant seems to bind the substrates at a longer distance from the reactive center. This is especially pronounced for the 4′ position corresponding to aromatic hydroxylation. Experimentally, the V87F mutant favors the aromatic hydroxylation over the aliphatic hydroxylation reaction.12 The binding modes (Figure 7, center and right) from the MD simulations of the V87F mutant show the π−π stacking between the aromatic ring of the ligand and the phenylalanine side chain of residue 87. The performance of DFT and MD has been evaluated quantitatively by two different approaches. The area under the curve (AUC) of receiver operating characteristics (ROC) plots identifying any correctly predicted sites-of-metabolism (SOMs) has been calculated considering all hydroxylation sites. The AUC is 0.88 for DFT and 0.83, 0.86, and 0.84 for the three mutants using MD. The ability to identify the experimentally observed SOMs within top-three (mefenamic acid and meclofenamic acid) or top-two (tolfenamic acid and diclofenac) guesses was evaluated showing that 21 of 30 sites and 19 of 30 sites were correctly predicted by DFT and MD, respectively. These numbers indicate that although both methods perform quite well, reactivity seems to be more relevant to consider in this case. Interestingly, we observe that some experimental SOMs are only marginally accessible in the MD simulations. However, all SOMs are characterized by activation energy lower than 74 kJ/mol (see Table S1). In future work, a combination of reactivity and accessibility is desirable. Potentially, the activation barriers in Table 3 and

a

Percentages for which the torsional angle between atom number 4′, 3′, the carbon in the 3′ methyl group of the ligand, and the oxygen of Cpd I during a MD simulation is between 45 and 135°. Numbers in parentheses are percentage of experimentally observed aliphatic hydroxylation from Venkataraman et al.12

probably due to the chlorine substituents in meclofenamic acid. For tolfenamic acid, the strictest O−C distance criterion of 3.5 Å seems to eliminate most of the structures from the trajectory for the M11 and M11 V87F mutants (Figure 6). However, the V87I mutant binds tolfenamic acid closer to the FeO group. Using 4.5 and 5.5 Å as distance criteria does not change the accessibility significantly for the M11 and M11 V87F mutants (cf. Figure S8). Comparing the binding mode of tolfenamic acid and mefenamic acid in M11 V87F, we observe that they bind in a similar manner, and one might speculate that it is the Cl substituent in the 3′ position that prevents the 4′ position to move closer. For diclofenac, using the distance criterion of 3.5 Å, both the 5 and the 4′ positions are accessible in all three mutants, although M11 and M11 V87I show a preference for the 5 position and M11 V87F for the 4′ position (Figure 6). At longer O−C distances, all unsubstituted aromatic C atoms are accessible.

■

DISCUSSION Although it is generally accepted that metabolism only takes place if the site is both reactive and accessible, only a limited number of methods take both effects into account and evaluate their mutual importance.13,14 In a previous computational study on aflatoxin B1 (AFB1) metabolized by CYP 1A2 and 3A4,52 it was shown that AFB1 binds to the CYP in a conformation where the atom in the ligand is accessible for the heme group. The reactivity of the ligand atom should also be sufficiently high to allow the reaction to take place. For the rigid compound AFB1 in CYP 1A2 and 3A4, the reactivity was the most important factor, but whether this conclusion can be extended to other, more flexible, ligands and others CYPs is not evident. In this study the B3LYP functional was chosen for the DFT calculations. Although additional functionals and QM methods have been developed since B3LYP, comparisons suggest that B3LYP is well suited and provides consistency with previous studies.53,54 The most reactive sites from the DFT calculations on the four compounds studied in this paper are the aromatic hydroxylation reactions in 5 and 4′ positions which, including the dispersion, are as reactive as the Me-2′ and Me-3′ hydroxylation. From the MD simulations, most of the aromatic sites are generally equally accessible. Thus, it can be concluded that the reactions occur for the most reactive sites. The Me-2′ site is not accessible and is therefore not observed experimentally. This is also the case for the N-hydroxylation reactions that are not observed experimentally either. Interestingly, the aliphatic hydrogen abstraction reactions for mefenamic acid and meclofenamic acid are less favored in the M11 and V87F mutants as opposed to the V87I mutant that favors the aliphatic reaction. From the MD simulations it H


Article

Journal of Chemical Information and Modeling the torsional angle variations shown in Figures 4 and 5 could be combined. In this way, one could calculate reactivity as a function of geometrical parameters from the MD simulations to provide a time-averaged reactivity without having to implement time-consuming methods such as QM/MM calculations.49,55

three different binding modes for the three BM3 mutants; snapshots from the MD simulations of tolfenamic acid in the two different binding modes for the three BM3 mutants; snapshots from the MD simulations of diclofenac in the two different binding modes for the three BM3 mutants; histograms of accessibilities from the molecular dynamics simulations using O−C distances of 4.5 and 5.5 Å as criterion; ligand RMSD values for MD simulation of meclofenamic acid in the three BM3 mutants; ligand RMSD values for MD simulation of tolfenamic acid in the three BM3 mutants; ligand RMSD values for MD simulation of diclofenac in the three BM3 mutants; B3LYP or B3LYP-D3 activation energies for the hydroxylation of mefenamic acid, meclofenamic acid, tolfenamic acid, and diclofenac calculated as the difference between the transition state energy and the energy of the separated species; distances and angles in the reactant state for mefenamic acid, meclofenamic acid, tolfenamic acid, and diclofenac; distances and angles in the transition state for mefenamic acid, meclofenamic acid, tolfenamic acid, and diclofenac; average RMSD values for the MD simulations of the three BM3 mutants with the four ligands (PDF)

■

CONCLUSION We investigated the reactivity of the possible SOMs in four structurally related compounds, three fenamic acid analogues and diclofenac, by DFT calculations. The reactivity studies were combined with an evaluation of the accessibility determined by docking and subsequent MD simulations. From the DFT B3LYP calculations we observe that the aliphatic sites in the 2′ and 3′ positions are the most reactive. This effect is dampened in the B3LYP-D3 calculations in which the aromatic reaction are stabilized by dispersive forces and become equally favorable to the aliphatic reactions. From the DFT calculation we conclude that the most reactive sites are the 2′, 3′, 4′, 3, and 5. On the basis of the TS and reactant structure from DFT calculations accessibility criteria were determined. Angle criteria of 115−130° and 130−160° were used for aliphatic and aromatic hydroxylation, respectively. O− C distance criteria of 3.5, 4.5, and 5.5 Å were investigated. While 4.5 Å is sensitive enough to define the 4′ position accessible in the V87F mutant, only 5.5 Å is sensitive enough for tolfenamic acid. Generally sites 3, 4, 5, 6, 3′, 4′, and 5′ are accessible for the substrates in all BM3 mutants. Combining these data with the DFT calculations leaves the 3, 5, 3′, and 4′ as the top sites of metabolism which is consistent with the 5-, 3′-, and 4′substituted products being the experimentally observed metabolites.12 In the drug discovery process identification of the most likely sites of metabolism is important, and accordingly, the most predictive methods have focused on predicting site of metabolism. Our study suggests that the combined approach presented here is sufficient to both capture the binding and explain the metabolism of conformationally flexible compounds. Although accessibility should not be ignored, these findings expand the field of predictive drug metabolism and provide confidence in previous work, indicating that reactivity may be more important than accessibility.

■

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Flemming Steen Jørgensen: 0000-0001-8040-2998 Present Addresses ‡

Department of Plasma Proteins, Sanquin-AMC Landsteiner Laboratory, Amsterdam, The Netherlands. § Protein Engineering, Novozymes A/S, Krogshøjvej 36, DK2880 Bagsvaerd, Denmark. Author Contributions

The manuscript was written through contributions of all authors. All authors have given approval to the final version of the manuscript. Notes

The authors declare no competing financial interest.

■

ASSOCIATED CONTENT

S Supporting Information *

REFERENCES

(1) Ortiz de Montellano, P. R. Cytochrome P450: Structure, Mechanism, and Biochemistry, 4th ed.; Springer International Publishing, 2015. (2) Orr, S. T. M.; Ripp, S. L.; Ballard, T. E.; Henderson, J. L.; Scott, D. O.; Obach, R. S.; Sun, H.; Kalgutkar, A. S. Mechanism-Based Inactivation (Mbi) of Cytochrome P450 Enzymes: Structure−Activity Relationships and Discovery Strategies to Mitigate Drug−Drug Interaction Risks. J. Med. Chem. 2012, 55, 4896−4933. (3) Xu, J.; Wang, X.-y.; Guo, W.-z. The Cytochrome P450 Superfamily: Key Players in Plant Development and Defense. J. Integr. Agric. 2015, 14, 1673−1686. (4) Kelly, S. L.; Kelly, D. E. Microbial Cytochromes P450: Biodiversity and Biotechnology. Where Do Cytochromes P450 Come from, What Do They Do and What Can They Do for Us? Philos. Trans. R. Soc., B 2013, 368, 20120476. (5) Munro, A. W.; Leys, D. G.; McLean, K. J.; Marshall, K. R.; Ost, T. W. B.; Daff, S.; Miles, C. S.; Chapman, S. K.; Lysek, D. A.; Moser, C. C.; Page, C. C.; Dutton, P. L. P450 Bm3: The Very Model of a Modern Flavocytochrome. Trends Biochem. Sci. 2002, 27, 250−257.

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jcim.8b00750. Transition state structures for the N-hydroxylation (amine), aromatic hydroxylation (3, 4, 5, 6, 4′, and 5'), and aliphatic hydroxylation (Me-3′) reactions of meclofenamic acid; transition state structures for the Nhydroxylation (amine), aromatic hydroxylation (3, 4, 5, 6, 4', 5′, and 6′), and aliphatic hydroxylation (Me-2′) reactions of tolfenamic acid; transition state structures for the N-hydroxylation (amine) and aromatic hydroxylation (3, 4, 5, 6, 3′, 4′, and 5′) reactions of diclofenac; docking poses of meclofenamic acid, diclofenac, and tolfenamic acid in BM3 M11 representing aliphatic hydroxylation (only for meclofenamic acid) and aromatic hydroxylation on both aromatic ring; snapshots from the MD simulations of meclofenamic acid in the I


Article

Journal of Chemical Information and Modeling (6) van Vugt-Lussenburg, B. M.; Damsten, M. C.; Maasdijk, D. M.; Vermeulen, N. P.; Commandeur, J. N. Heterotropic and Homotropic Cooperativity by a Drug-Metabolising Mutant of Cytochrome P450 Bm3. Biochem. Biophys. Res. Commun. 2006, 346, 810−818. (7) van Vugt-Lussenburg, B. M.; Stjernschantz, E.; Lastdrager, J.; Oostenbrink, C.; Vermeulen, N. P.; Commandeur, J. N. Identification of Critical Residues in Novel Drug Metabolizing Mutants of Cytochrome P450 Bm3 Using Random Mutagenesis. J. Med. Chem. 2007, 50, 455−461. (8) Damsten, M. C.; van Vugt-Lussenburg, B. M.; Zeldenthuis, T.; de Vlieger, J. S.; Commandeur, J. N.; Vermeulen, N. P. Application of Drug Metabolising Mutants of Cytochrome P450 Bm3 (Cyp102a1) as Biocatalysts for the Generation of Reactive Metabolites. Chem.-Biol. Interact. 2008, 171, 96−107. (9) Cryer, B.; Feldman, M. Cyclooxygenase-1 and Cyclooxygenase-2 Selectivity of Widely Used Nonsteroidal Anti-Inflammatory Drugs. Am. J. Med. 1998, 104, 413−421. (10) Mancy, A.; Antignac, M.; Minoletti, C.; Dijols, S.; Mouries, V.; Ha Duong, N.-T.; Battioni, P.; Dansette, P. M.; Mansuy, D. Diclofenac and Its Derivatives as Tools for Studying Human Cytochromes P450 Active Sites: Particular Efficiency and Regioselectivity of P450 2cs. Biochemistry 1999, 38, 14264−14270. (11) Madsen, K. G.; Skonberg, C.; Jurva, U.; Cornett, C.; Hansen, S. H.; Johansen, T. N.; Olsen, J. Bioactivation of Diclofenac in Vitro and in Vivo: Correlation to Electrochemical Studies. Chem. Res. Toxicol. 2008, 21, 1107−1119. (12) Venkataraman, H.; Verkade-Vreeker, M. C.; Capoferri, L.; Geerke, D. P.; Vermeulen, N. P.; Commandeur, J. N. Application of Engineered Cytochrome P450 Mutants as Biocatalysts for the Synthesis of Benzylic and Aromatic Metabolites of Fenamic Acid Nsaids. Bioorg. Med. Chem. 2014, 22, 5613−5620. (13) Olsen, L.; Oostenbrink, C.; Jørgensen, F. S. Prediction of Cytochrome P450 Mediated Metabolism. Adv. Drug Delivery Rev. 2015, 86, 61−71. (14) Kirchmair, J.; Williamson, M. J.; Tyzack, J. D.; Tan, L.; Bond, P. J.; Bender, A.; Glen, R. C. Computational Prediction of Metabolism: Sites, Products, Sar, P450 Enzyme Dynamics, and Mechanisms. J. Chem. Inf. Model. 2012, 52, 617−648. (15) Olsen, L.; Rydberg, P.; Rod, T. H.; Ryde, U. Prediction of Activation Energies for Hydrogen Abstraction by Cytochrome P450. J. Med. Chem. 2006, 49, 6489−6499. (16) Rydberg, P.; Ryde, U.; Olsen, L. Prediction of Activation Energies for Aromatic Oxidation by Cytochrome P450. J. Phys. Chem. A 2008, 112, 13058−13065. (17) Shaik, S.; Kumar, D.; de Visser, S. P.; Altun, A.; Thiel, W. Theoretical Perspective on the Structure and Mechanism of Cytochrome P450 Enzymes. Chem. Rev. 2005, 105, 2279−2328. (18) Meunier, B.; de Visser, S. P.; Shaik, S. Mechanism of Oxidation Reactions Catalyzed by Cytochrome P450 Enzymes. Chem. Rev. 2004, 104, 3947−3980. (19) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (Dft-D) for the 94 Elements H-Pu. J. Chem. Phys. 2010, 132, 154104. (20) Grimme, S.; Ehrlich, S.; Goerigk, L. Effect of the Damping Function in Dispersion Corrected Density Functional Theory. J. Comput. Chem. 2011, 32, 1456−1465. (21) Schafer, A.; Horn, H.; Ahlrichs, R. Fully Optimized Contracted Gaussian-Basis Sets for Atoms Li to Kr. J. Chem. Phys. 1992, 97, 2571−2577. (22) Hariharan, P. C.; Pople, J. A. Influence of Polarization Functions on Molecular-Orbital Hydrogenation Energies. Theor. Chim. Acta 1973, 28, 213−222. (23) Hehre, W. J.; Ditchfield, R.; Pople, J. A. Self-Consistent Molecular-Orbital Methods. 12. Further Extensions of Gaussian-Type Basis Sets for Use in Molecular-Orbital Studies of Organic Molecules. J. Chem. Phys. 1972, 56, 2257−2260. (24) Francl, M. M.; Pietro, W. J.; Hehre, W. J.; Binkley, J. S.; Gordon, M. S.; Defrees, D. J.; Pople, J. A. Self-Consistent Molecular-

Orbital Methods. 23. A Polarization-Type Basis Set for 2nd-Row Elements. J. Chem. Phys. 1982, 77, 3654−3665. (25) Krishnan, R.; Binkley, J. S.; Seeger, R.; Pople, J. A. SelfConsistent Molecular-Orbital Methods. 20. A Basis Set for Correlated Wave-Functions. J. Chem. Phys. 1980, 72, 650−654. (26) McLean, A. D.; Chandler, G. S. Contracted Gaussian-Basis Sets for Molecular Calculations. 1. 2nd Row Atoms, Z = 11−18. J. Chem. Phys. 1980, 72, 5639−5648. (27) Rulisek, L.; Jensen, K. P.; Lundgren, K.; Ryde, U. The Reaction Mechanism of Iron and Manganese Superoxide Dismutases Studied by Theoretical Calculations. J. Comput. Chem. 2006, 27, 1398−1414. (28) Turbomole3, Version 6; University of Karlsruhe and Forschungszentrum Karlsruhe Gmbh, 1989−2007. (29) Maestro, Version 9.8; Schrödinger, LLC., New York, 2014. (30) Madhavi Sastry, G. M.; Adzhigirey, M.; Day, T.; Annabhimoju, R.; Sherman, W. Protein and Ligand Preparation: Parameters, Protocols, and Influence on Virtual Screening Enrichments. J. Comput.-Aided Mol. Des. 2013, 27, 221−234. (31) Wester, M. R.; Johnson, E. F.; Marques-Soares, C.; Dijols, S.; Dansette, P. M.; Mansuy, D.; Stout, C. D. Structure of Mammalian Cytochrome P450 2c5 Complexed with Diclofenac at 2.1 Å Resolution: Evidence for an Induced Fit Model of Substrate Binding. Biochemistry 2003, 42, 9335−9345. (32) Halgren, T. A.; Murphy, R. B.; Friesner, R. A.; Beard, H. S.; Frye, L. L.; Pollard, W. T.; Banks, J. L. Glide: A New Approach for Rapid, Accurate Docking and Scoring. 2. Enrichment Factors in Database Screening. J. Med. Chem. 2004, 47, 1750−1759. (33) Domańska, U.; Pobudkowska, A.; Pelczarska, A.; WiniarskaTusznio, M.; Gierycz, P. Solubility and Pka of Select Pharmaceuticals in Water, Ethanol, and 1-Octanol. J. Chem. Thermodyn. 2010, 42, 1465−1472. (34) Capoferri, L.; Leth, R.; ter Haar, E.; Mohanty, A. K.; Grootenhuis, P. D.; Vottero, E.; Commandeur, J. N.; Vermeulen, N. P.; Jørgensen, F. S.; Olsen, L.; Geerke, D. P. Insights into Regioselective Metabolism of Mefenamic Acid by Cytochrome P450 Bm3Mutants through Crystallography, Docking, Molecular Dynamics, and Free Energy Calculations. Proteins: Struct., Funct., Genet. 2016, 84, 383−396. (35) Hornak, V.; Abel, R.; Okur, A.; Strockbine, B.; Roitberg, A.; Simmerling, C. Comparison of Multiple Amber Force Fields and Development of Improved Protein Backbone Parameters. Proteins: Struct., Funct., Genet. 2006, 65, 712−725. (36) Case, D. A.; Darden, T. A.; Cheatham, T. E.; Simmerling, C. L.; Wang, J.; Duke, R. E.; Luo, R.; Walker, R. C.; Zhang, W.; Merz, K. M.; Roberts, B.; Hayik, S.; Roitberg, A.; Seabra, G.; Swails, J.; Goetz, A. W.; Kolossváry, I.; Wong, K. F.; Paesani, F.; Vanicek, J.; Wolf, R. M.; Liu, J.; Wu, X.; Brozell, S. R.; Steinbrecher, T.; Gohlke, H.; Cai, Q.; Ye, X.; Wang, J.; Hsieh, M. J.; Cui, G.; Roe, D. R.; Mathews, D. H.; Seetin, M. G.; Salomon-Ferrer, R.; Sagui, C.; Babin, V.; Luchko, T.; Gusarov, S.; Kovalenko, A.; Kollman, P. A. Amber; University of California: San Francisco, CA, 2012. (37) Shahrokh, K.; Orendt, A.; Yost, G. S.; Cheatham, T. E. Quantum Mechanically Derived Amber-Compatible Heme Parameters for Various States of the Cytochrome P450 Catalytic Cycle. J. Comput. Chem. 2012, 33, 119−133. (38) Berendsen, H. J. C.; van der Spoel, D.; van Drunen, R. Gromacs: A Message-Passing Parallel Molecular Dynamics Implementation. Comput. Phys. Commun. 1995, 91, 43−56. (39) Wang, J.; Wang, W.; Kollman, P. A.; Case, D. A. Automatic Atom Type and Bond Type Perception in Molecular Mechanical Calculations. J. Mol. Graphics Modell. 2006, 25, 247−260. (40) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. Development and Testing of a General Amber Force Field. J. Comput. Chem. 2004, 25, 1157−1174. (41) Sousa da Silva, A. W.; Vranken, W. F. Acpype - Antechamber Python Parser Interface. BMC Res. Notes 2012, 5, 367. (42) 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. J


Article

Journal of Chemical Information and Modeling (43) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103, 8577−8593. (44) Humphrey, W.; Dalke, A.; Schulten, K. Vmd: Visual Molecular Dynamics. J. Mol. Graphics 1996, 14, 33−38. (45) Seger, S. T.; Rydberg, P.; Olsen, L. Mechanism of the NHydroxylation of Primary and Secondary Amines by Cytochrome P450. Chem. Res. Toxicol. 2015, 28, 597−603. (46) Bathelt, C. M.; Ridder, L.; Mulholland, A. J.; Harvey, J. N. Mechanism and Structure-Reactivity Relationships for Aromatic Hydroxylation by Cytochrome P450. Org. Biomol. Chem. 2004, 2, 2998−3005. (47) Lonsdale, R.; Harvey, J. N.; Mulholland, A. J. Effects of Dispersion in Density Functional Based Quantum Mechanical/ Molecular Mechanical Calculations on Cytochrome P450 Catalyzed Reactions. J. Chem. Theory Comput. 2012, 8, 4637−4645. (48) Rydberg, P.; Lonsdale, R.; Harvey, J. N.; Mulholland, A. J.; Olsen, L. Trends in Predicted Chemoselectivity of Cytochrome P450 Oxidation: B3lyp Barrier Heights for Epoxidation and Hydroxylation Reactions. J. Mol. Graphics Modell. 2014, 52, 30−35. (49) Lonsdale, R.; Houghton, K. T.; Zurek, J.; Bathelt, C. M.; Foloppe, N.; de Groot, M. J.; Harvey, J. N.; Mulholland, A. J. Quantum Mechanics/Molecular Mechanics Modeling of Regioselectivity of Drug Metabolism in Cytochrome P450 2c9. J. Am. Chem. Soc. 2013, 135, 8001−8015. (50) de Graaf, C.; Oostenbrink, C.; Keizers, P. H.; van der Wijst, T.; Jongejan, A.; Vermeulen, N. P. Catalytic Site Prediction and Virtual Screening of Cytochrome P450 2d6 Substrates by Consideration of Water and Rescoring in Automated Docking. J. Med. Chem. 2006, 49, 2417−2430. (51) Rydberg, P.; Olsen, L. Do Two Different Reaction Mechanisms Contribute to the Hydroxylation of Primary Amines by Cytochrome P450? J. Chem. Theory Comput. 2011, 7, 3399−3404. (52) Bonomo, S.; Jorgensen, F. S.; Olsen, L. Dissecting the Cytochrome P450 1a2- and 3a4-Mediated Metabolism of Aflatoxin B1 in Ligand and Protein Contributions. Chem. - Eur. J. 2017, 23, 2884−2893. (53) Leth, R.; Rydberg, P.; Jørgensen, F. S.; Olsen, L. Density Functional Theory Study on the Formation of Reactive Benzoquinone Imines by Hydrogen Abstraction. J. Chem. Inf. Model. 2015, 55, 660− 666. (54) Tentscher, P. R.; Arey, J. S. Geometries and Vibrational Frequencies of Small Radicals: Performance of Coupled Cluster and More Approximate Methods. J. Chem. Theory Comput. 2012, 8, 2165−2179. (55) Lonsdale, R.; Rouse, S. L.; Sansom, M. S.; Mulholland, A. J. A Multiscale Approach to Modelling Drug Metabolism by MembraneBound Cytochrome P450 Enzymes. PLoS Comput. Biol. 2014, 10, e1003714.

K


Both Reactivity and Accessibility Are Important in Cytochrome P450

Recommend Documents