# Insights into molecular mechanisms of drug metabolism dysfunction of human CYP2C9*30

Maxime Louet, Céline M. Labbé, Charline Fagnen, Cassiano M Aono, Paula

Homem-De-Mello, Bruno O Villoutreix, Maria A Miteva

► To cite this version:

Maxime Louet, Céline M. Labbé, Charline Fagnen, Cassiano M Aono, Paula Homem-De-Mello, et al. Insights into molecular mechanisms of drug metabolism dysfunction of human CYP2C9*30. PLoS ONE, 2018, 13 (5), pp.e0197249. 10.1371/journal.pone.0197249. hal-01801373

HAL Id: hal-01801373

https://hal.sorbonne-universite.fr/hal-01801373

Submitted on 28 May 2018

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
RESEARCH ARTICLE

# Insights into molecular mechanisms of drug metabolism dysfunction of human CYP2C9*30

Maxime Louet¹,², Céline M. Labbé¹,², Charline Fagnen¹,²,³, Cassiano M. Aono⁴, Paula Homem-de-Mello⁴, Bruno O. Villoutreix¹,², Maria A. Miteva¹,²*

1 Université Paris Diderot, Sorbonne Paris Cité, Inserm UMR-S 973, Molécules Thérapeutiques *In silico*, Paris, France, 2 INSERM, U973, Paris, France, 3 Université Pierre et Marie Curie, Sorbonne Universités, UMR 7590, Institut de Minéralogie, de Physique des Matériaux et de Cosmochimie, Paris, France, 4 Universidade Federal do ABC, Centro de Ciências Naturais e Humanas, Santo André, Brazil

* maria.mitev@inserm.fr

<figure><img src="image_2.png" /></figure>

## Abstract

Cytochrome P450 2C9 (CYP2C9) metabolizes about 15% of clinically administrated drugs. The allelic variant *CYP2C9*\*\*30 (A477T) is associated to diminished response to the antihypertensive effects of the prodrug losartan and affected metabolism of other drugs. Here, we investigated molecular mechanisms involved in the functional consequences of this amino-acid substitution. Molecular dynamics (MD) simulations performed for the active species of the enzyme (heme in the Compound I state), in the apo or substrate-bound state, and binding energy analyses gave insights into altered protein structure and dynamics involved in the defective drug metabolism of human CYP2C9.30. Our data revealed an increased rigidity of the key Substrate Recognition Sites SRS1 and SRS5 and shifting of the β turn 4 of SRS6 toward the helix F in CYP2C9.30. Channel and binding substrate dynamics analyses showed altered substrate channel access and active site accommodation. These conformational and dynamic changes are believed to be involved in the governing mechanism of the reduced catalytic activity. An ensemble of representative conformations of the WT and A477T mutant properly accommodating drug substrates were identified, those structures can be used for prediction of new CYP2C9 and CYP2C9.30 substrates and drug-drug interactions.

<figure><img src="image_3.png" /></figure>

### OPEN ACCESS

**Citation:** Louet M, Labbé CM, Fagnen C, Aono CM, Homem-de-Mello P, Villoutreix BO, et al. (2018) Insights into molecular mechanisms of drug metabolism dysfunction of human CYP2C9*30. *PLoS ONE* 13(5): e0197249. https://doi.org/10.1371/journal.pone.0197249

**Editor:** Xuhui Huang, Hong Kong University of Science and Technology, HONG KONG

**Received:** January 26, 2018

**Accepted:** April 30, 2018

**Published:** May 10, 2018

**Copyright:** © 2018 Louet et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

**Data Availability Statement:** All relevant data are within the paper and its Supporting Information files.

**Funding:** This work was supported by Inserm, the University Paris Diderot, ANR ToxME (MAM), FAPESP (CMA), CNPq (PHM). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

**Competing interests:** The authors have declared that no competing interests exist.

## Introduction

Cytochrome P450 2C9 (CYP2C9) is the most expressed member of the human CYP2C family and metabolizes more than 15% of clinically administrated drugs including hypoglycemic agents, anticonvulsants, anticoagulants, nonsteroidal anti-inflammatory drugs (NSAIDs), anti-hypertensives, diuretic drugs [1–3] and several endogenous compounds [4]. The genetic variations resulting mainly from single nucleotide polymorphisms (SNPs) explain a large part of the inter-individual variability in the activity of CYP enzymes [5]. More than 50 variants of CYP2C9 have shown decreased enzymatic activity [3, 6–8]. Polymorphisms of CYP2C9 can affect the clinical response of drugs metabolized by CYP2C9, especially those with a narrow therapeutic index and can induce adverse drug reactions (ADR). The rare allelic variant
CYP2C9*30 found in Japanese subjects has been associated to diminished response to the anti-hypertensive effects of the prodrug losartan [9]. Moreover, in vitro experiments have shown that the intrinsic clearance of losartan, diclofenac and glimepiride was markedly decreased for CYP2C9.30 (mutation A477T) (by > 70%) compared to that of the wild-type (WT) enzyme [10, 11]. It has been previously proposed that conformational changes occurring nearby F476 could play a role in the substrate recognition [12]. Recently the crystal structures of human CYP2C9 and its variant CYP2C9.30 co-crystalized with losartan have been reported [13]. However, the structural overlay of these two complexes did not reveal significant differences between the losartan-bound WT and A477T mutant crystal structural (Root Mean Square Deviation (RMSD) of ~0.2 Å). Three losartan molecules have been found to bind in the same positions in the co-crystallized WT and A477 mutant structures, neither of them was oriented with the losartan hydroxylation site toward the heme iron required for the hydroxylation reaction. Thus, the molecular mechanisms responsible for the mutation-induced dysfunction of the metabolism of CYP2C9.30 remain unclear.

Here we provide insights into the dysfunctional metabolism mechanism in human CYP2C9.30 based on molecular dynamics (MD), quantum mechanics and docking simulations. Molecular dynamics (MD) simulations combined with energetic analyses have been demonstrated to be a pertinent approach to understand molecular mechanisms of amino-acid substitutions affecting drug metabolism or regioselectivity of CYP [14–18]. In this work, we explored a large region of the conformational space of the WT of CYP2C9 and its mutant A477T, in the apo or substrate-bound states, for two chemically different substrates, diclofenac and losartan, with the site of hydroxylation oriented toward the heme iron consistent with the metabolic reaction. MD simulations of 1.5 µs were run in total. Structural and dynamics changes caused by the A447T mutation were revealed and analyzed to explain the decrease of the enzymatic activity of CYP2C9.30. The A477T mutant demonstrated an increased rigidity of key segments involved in substrate recognition that could affect the substrate entrance and accommodation. Binding energy analyses via molecular docking of five drugs substrates of CYP2C9 (shown in Fig 1) suggested molecular mechanisms implicated in the defective drug metabolism in carriers of CYP2C9*30. Our results highlighted conformational and dynamic changes believed to be involved in the governing mechanism of the catalytic activity reduction and confirmed that protein dynamics affected by missense mutations could play a crucial role in affected regulation of substrate entry, recognition and metabolism of CYP2C9.

## Results and discussion

### Conformational and dynamics behavior of CYP2C9 WT and A477T

The A477T mutation present in CYP2C9.30 variant is located in the Substrate Recognition Site 6 (SRS6) (Supplementary S1 Fig). The comparison of the two crystal structures of losartan-bound CYP2C9 WT and A477T mutant did not show important differences between the two complexes except a rotation of the side chain of Q214 forming thus a hydrogen bond (HB) with the site chain of T477 in the mutant protein (PDB IDs 5XXI and 5X23 [13]). Three losartan molecules were found to bind equally in the crystal structures of CYP2C9 WT and A477T mutant, two substrate molecules were present in the binding pocket and one of them was located in a peripheral site of the protein. One of the co-crystalized losartan molecule was bound far from the heme but nearby Q214 and the mutation residue 477 and involved in tight aromatic contacts with F476. The presence of the losartan molecule at this place does not facilitate understanding the role of the mutated residue 477 and its interactions with Q214 and F476. The residue 477 shows a low sequence conservation based on a sequence alignment of all CYP2 family members performed with the ConSurf web server [19]. A recent analysis of
<figure><img src="image_4.png" /><figcaption>**Fig 1.** Structures of five drug substrates of CYP2C9. The arrows show the sites of metabolism mediated by CYP2C9.</figcaption></figure>

https://doi.org/10.1371/journal.pone.0197249.g001

this mutation [3] carried out with well-established online computational tools predicting
mutation effects [20–23] did not suggested a hypothesis about its impact on the enzyme struc-
ture or function. To elucidate the molecular mechanisms involved in the defected metabolism
of CYP2C9.30 we performed molecular dynamics simulations of the active species of the WT
of CYP2C9 and the mutant A477T, both containing the heme in the Compound I state (Cpd I)
[24]. To better understand the structural and dynamics bases for the metabolism alteration we
have studied six CYP2C9 systems (1) WT apo, (2) A477T apo, (3) WT diclofenac-bound, (4)
A477T diclofenac-bound, (5) WT losartan-bound and (6) A477T losartan-bound. The initial
positions of the bound substrates were obtained by docking and corresponded to the Sites of
Metabolism (SOM) oriented toward the oxygen atom of Cpd I. It is to note that the positions
of the co-crystalized losartan in the WT and A477T mutant structures are not consistent with
the catalytic reaction because the site of hydroxylation is positioned away from the heme iron
[13]. Thus, these crystal binding modes are not appropriate for our MD simulations. For each
of the six systems, 5 MD trajectories of 50 ns were run. The MD simulations of 250 ns per stud-
ied system have shown stable protein structures of the CYP2C9 WT and its A477T mutant.
The all-heavy-atom RMSD over the merged five MD trajectories for each system showed maxi-
mal mean value of 1.59 Å for all systems. The backbone RMSD data (see S2 Fig) did not show
large conformational changes of the A477 mutant compared to the WT suggesting that subtle molecular mechanisms could be responsible for the altered metabolism.

The Root Mean Square Fluctuations (RMSF) results confirmed the stability of the WT and mutant A477T structures during the MD simulations. The mean RMSF values of the merged MD trajectories for each of the six systems varied between 1.08 Å and 1.18 Å. However, the comparison of the atomic fluctuations per residue (see Fig 2 and S3 Fig) revealed several discrepancies between the WT apo and the mutant apo proteins in key regions involved in substrate binding. The A447T mutant in apo form was found to be more rigid than the WT apo CYP2C9, in particular the residues from 367 to 369 of SRS5, and the residues from 271 to 276 of helix H (see Fig 2A and 2B). It is important to note that in the presence of a bound substrate (see Fig 2C, 2D, 2E and 2F), the WT and the mutant proteins showed similar flexibility probably due to the stabilization of the protein structures triggered by the interactions with the ligand. Different regions were influenced when CYP2C9 was bound to diclofenac or losartan that was likely due to the different locations of the two substrates in the binding pocket during the MD simulations (see below on the substrates dynamics).

In order to find noticeable conformational differences of the substrate binding pocket of CYP2C9 between the WT and the A477T mutant, we performed structural clustering of the five merged MD trajectories for each of the studied systems based on the binding site residues (see S4 Fig). The structural analysis of the centroid structures for the 15 most populated clusters of the six studied systems covering at least 55% of the binding site conformational space explored during all MD simulations showed structural differences of the β turn 4 including the mutation site. For the WT apo state, the most representative centroid structure corresponding to 11.3% of the total 250 ns MD was very similar to the initial apo crystal structure of CYP2C9 PDB ID 1OG2 [25]), (RMSD of the Cα atoms of 1.01 Å) with a conserved β turn 4 conformation and unchanged F100 and F114 side-chain orientations. The other centroids showed different side-chain orientations for F100 and F114 belonging to the BC loop, a protein segment known to be a part of previously identified substrate access channels (2a, 2b, 2c, 2ac, and 2e) [26, 27] (see below for channel dynamics). The 15 centroid structures of the A477T apo mutant showed different binding-site conformations for F476 belonging to the β turn 4. F476 has already been pointed out to be a key residue involved in the channel access to the catalytic pocket [25, 27]. To study the conformational changes of the β turn 4 we further analyzed the HB interactions of the β turn 4 during all MD simulations. The most important changes were detected for the residue S209 located on helix F with the residues 470–490 of the β turn 4. Two residues were identified to be involved in HBs with S209 (Fig 3A). We did not observe a HB bond between Q214 and T477 during the mutant apo MD simulations. In the apo mutant protein we detected possible HB breaks between the side chain of Q214 and the main chains of G475 (a HB conserved in 10.3% in the apo WT vs 5.6% in the apo mutant MD simulations) and F476 (a HB conserved in 26.3% in the apo WT vs 2.1% in the apo mutant MD simulations). Changes of the HB interactions of Q214 due to the A477 mutation were also observed in the case of bound losartan to CYP2C9 [13]. Therefore, significant interaction modifications occur in the mutant protein compared to the WT one corresponding to a movement of the β turn 4 in regard to the helix F of about 5 Å (see Fig 3B). This shift may alter the well-known channel access including F476 (or so-called warfarin-binding subpocket found in the co-crys talized CYP2C9–warfarin structure, PDB ID 1OG5 [25]) and may destabilize the CYP2C9-substrate interactions in this zone.

Hydrophobic interactions between F476 and other hydrophobic residues of the CYP2C9 binding pocket were monitored within a distance of 5 Å during the MD simulations (see Fig 4 and S5 Fig). We observed that in the presence of the A477T mutation, F476 increased its contacts with L361 and F476 remained in stable contacts with L366 in the cases of ligand-bound
<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>
**Fig 2.** Root Mean Square Fluctuations (RMSF) of the Cα atoms during the MD simulations of CYP2C9. The RMSF computed for five merged MD trajectories of (A) WT apo, (B) A477T apo, (C) WT diclofenac-bound, (D) A477T diclofenac-bound, (E) WT losartan-bound and (F) A477T losartan-bound are mapped on the CYP2C9 crystallographic structure (PDB ID: 1OG5) [25]. The color code of RMSF is ranging from 0.7 Å (blue) to 1.6 Å (red). The black arrows show noticeable differences between the WT and the mutant RMSF.

https://doi.org/10.1371/journal.pone.0197249.g002

states. Overall, F476 showed more stable contacts with L361, L362 and L366 belonging to SRS5 (shown in orange in Fig 4) in the mutant protein than in the WT one. These contacts explain the increased rigidity of SRS5 in the mutant compared to the WT (S3 Fig). When focusing on the apo state, we observed that the shift motion of the β turn 4 in the mutant protein modified the hydrophobic-interaction patterns of F476 as well. It is important to note that new contacts of F476 and L208 and L362 appeared during the MD simulations, which were not seen in the X-ray structure. The contacts of F476 with L208 strongly increased in the mutant structures confirming the movement of F476 toward L208 (helix F shown in red in Fig 4). This is likely related to the partial closing of the solvent channel S1 in the presence of the mutation A477T (see details on the channels dynamics below). Interestingly, F476 remained in stable contacts with the key F100 gate residue in the apo WT and mutant protein. This interaction is also conserved in the two X-ray structures (warfarin-bound PDB ID 1OG5 and flurbiprofen-bound PDB ID 1R9O).

## Substrate dynamics during the MD simulations

Substrate positions are crucial for the oxidation reaction. Fig 5 shows the diclofenac and losartan dynamics during the MD simulations for the WT and the A477T mutant. The diclofenac moved quickly out of its initial position and explored the entire binding pocket during three out of the five simulations for both WT and its A477T mutant (see Fig 5C). Despite of its mobility, the diclofenac always remained in contact with one of the SRS1, SRS2, SRS3, and SRS4 regions during both WT and A477T mutant simulations. However, the diclofenac lost completely its contacts with SRS5 during the mutant MD simulations (see Fig 5A).

On the contrary, the losartan reached an equilibrium position (called here the preferential conformation), remaining in the same subpocket in contact with the heme similarly to its

<figure><img src="image_11.png" /><figcaption>**Fig 3.** Interaction patterns between the helix F and the β-turn 4. (A) Percentage of the HBs between S209 and N474 or F476 during the MD simulations. (B) Representation of the shift motion of the β-turn 4 in regard to helix F. The interactions of S209 with N474 (blue) and F476 (red) are noted.</figcaption></figure>

https://doi.org/10.1371/journal.pone.0197249.g003
<figure><img src="image_12.png" /><figcaption>**Fig 4.** Hydrophobic contacts of F476 and other residues of the CYP2C9 binding pocket monitored over the MD simulations. Residues found in contact with F476 are depicted in sticks. Three distinct hydrophobic clusters are found and colored in blue, red and orange. The purple sphere represents the location of the mutation.</figcaption></figure>

https://doi.org/10.1371/journal.pone.0197249.g004

initial docking position in both WT and A477T mutant proteins (see Fig 5C and S6 Fig).
Along the 250 ns MD simulations for the WT, 1.5%, 86% and 11.5% of all generated conformations corresponded to the initial structure, the preferential structure and other diverse
<figure><img src="image_13.png" /><figcaption>**Fig 5.** Substrate dynamics during the MD simulations for the WT and the A477T mutant. Volumetric map of the pocket space explored by diclofenac (A) and losartan (B) during the entire MD simulations are shown. The initial positions of the ligands and Cpd I are shown as sticks. The protein is represented as a white cartoon. SRS regions are colored in green. The purple sphere represents the location of the mutation. The blue (WT) and red (A477T) surfaces represent regions</figcaption></figure>
explored by the ligands over 95% of all simulation snapshots. (C) Distance between the Site Of Metabolism (SOM) of diclofenac (top) and losartan (bottom) and the reactive oxygen atom of Cpd I along the MD simulations.

https://doi.org/10.1371/journal.pone.0197249.g005

structures, respectively. For the A477T mutant, 8.5%, 52.7% and 38.8% of all MD conformations corresponded to the initial structure, the preferential structure and other structures, respectively. These observations suggested that the preferential losartan position in the catalytic site consistent with the catalytic reaction is more stable in the WT than in the A477T mutant.

The RMSD of Cpd I along the MD simulations (S7 Fig) showed that Cpd I remained slightly less stable in the A477T mutant than in the WT. The Cpd I instability could also alter the catalytic activity of the A477T mutant. Not surprisingly, the interactions between Cpd I and the bound substrates stabilized Cpd I in the holo states in comparison with the apo states.

## Channels dynamics of CYP2C9 WT and A477T

We explored the structure and dynamics of the substrate access channels for the WT and the A477T mutant. Significant changes of the substrate access channels were found to occur in the presence of the mutation A477T. The open channels detected along all MD simulations for the WT and the mutant proteins are summarized in Table 1. CYP2C9 channels were previously studied [26, 27] and our predicted channels were labeled accordingly (see Fig 6A and 6B). We identified similar channels as the previously described ones, except for the so-called solvent channel S, which was found here to be divided into two independent channels named S1 and S2. In all cases, the solvent channel S1 remained more open than S2.

We observed that several channels were open in more than 20% of the MD simulations of the diclofenac-bound state (2a, 2b, 2e 2f, 5, S1 and S2). This observation is in accordance with the large movements of diclofenac observed during the MD simulations suggesting a dynamic behavior of the channels. Similar results have been reported in a previous study [27] for CYP2C9 bound to flurbiprofen, a NSAID similar to diclofenac. On the contrary, the large losartan molecule showing a very stable position in the binding pocket during the MD simulations triggered a very different behavior regarding the opening and closing of the channels. In fact, only three channels (2c, S1 and S2) were open during the MD simulations of losartan-bound CYP2C9 in accordance with the two structures of the WT and A477T mutant co-crystallized with losartan, PDB ID 5XXI and 5X23 [13], respectively. Interestingly, channel 2c was open in the losartan-bound structures, but closed in all other studied CYP2C9 states. Recently, it has been suggested that the peripheral site situated nearby channel 4 accommodating losartan and palmitic acid molecules in CYP2C9 and CYP2C8, respectively, could play some role in the substrate recognition [13]. However, our analyses indicated that channel 4 was closed in the MD simulations for all WT and mutant protein structures, either in apo or substrate-bound states. Thus, it is not likely this area to be a preferred access channel site.

<table><caption>Table 1. Percentage of CYP2C9 channels opened during the MD simulations of the WT and the A477T mutant in the apo, diclofenac- and losartan-bound states.</caption><thead><tr><th>Simulation system</th><th>2a</th><th>2b</th><th>2c</th><th>2ac</th><th>2e</th><th>2f</th><th>4</th><th>5</th><th>S1</th><th>S2</th></tr></thead><tbody><tr><td>WT apo</td><td>17.0</td><td>32.1</td><td>5.9</td><td>0.6</td><td>34.9</td><td>14.2</td><td>0.7</td><td>13.8</td><td>79.4</td><td>31.3</td></tr><tr><td>A477T apo</td><td>20.4</td><td>21.7</td><td>6.8</td><td>0.1</td><td>33.6</td><td>9.1</td><td>0.3</td><td>8.0</td><td>67.8</td><td>30.1</td></tr><tr><td>WT diclofenac</td><td>40.7</td><td>25.0</td><td>10.9</td><td>1.9</td><td>23.6</td><td>17.4</td><td>2.3</td><td>21.2</td><td>67.5</td><td>33.1</td></tr><tr><td>A477T diclofenac</td><td>24.5</td><td>23.8</td><td>3.7</td><td>14.6</td><td>20.6</td><td>37.7</td><td>0.3</td><td>15.5</td><td>65.8</td><td>28.5</td></tr><tr><td>WT losartan</td><td>6.2</td><td>5.9</td><td>45.6</td><td>0.0</td><td>12.6</td><td>6.1</td><td>1.6</td><td>8.0</td><td>71.7</td><td>12.1</td></tr><tr><td>A477T losartan</td><td>16.3</td><td>21.3</td><td>49.7</td><td>0.0</td><td>5.7</td><td>12.5</td><td>0.1</td><td>10.3</td><td>60.5</td><td>22.6</td></tr></tbody></table>

https://doi.org/10.1371/journal.pone.0197249.t001
<figure><img src="image_14.png" /><figcaption>**Fig 6.** CYP2C9 channels identified during the MD simulations. The channels are labelled according to Cojocaru et al. [26]. The protein is represented in white cartoon and the heme is represented in sticks.</figcaption></figure>

https://doi.org/10.1371/journal.pone.0197249.g006

The comparison between the apo WT and mutant proteins showed several important changes due to the mutation. Most of the substrate and solvent access channels were more closed in the A477T mutant as compared to the WT. More precisely, the channels 2b, 2f, 5 and S1 were more closed in the apo state of the mutant than of the WT (see Table 1). The channel 5 is partly formed by the β turn 4, where the mutation is located. As described above, this protein region changed its conformation in presence of the mutation that may lead to the closing of channel 5. The channel 2b is located far from the mutation residue 477 (> 10 Å). However, the movement of F476 toward L361, L362 and L366 belonging to SRS5 in the A477T mutant structure, as discussed above, likely alters the 2b formation depending on the SRS5 conformation. In addition, the shift of F476 toward L208, as explained above, favors the closing of the channel S1 in presence of the mutation A477T. All together the channels accessibility results suggest that the F476 shift in the apo A477T mutant plays a critical role in the channel dynamics.

## Substrate docking into CYP2C9 WT and A477T structures

We carried out molecular docking calculations to identify different conformations of the WT and A477T mutant of CYP2C9 properly accommodating the substrates for the catalytic reaction. In order to explore the binding spectrum of CYP2C9, we studied five drugs metabolized by CYP2C9 by molecular docking: diclofenac, flurbiprofen, warfarin, losartan and glimepiride. It has been shown that the metabolism of diclofenac, losartan and glimepiride is strongly affected in the A477T mutant [9-11]. To take into account the plasticity of the CYP2C9 binding pocket we studied multiple conformations of the WT and mutant protein generated by the MD simulations. To avoid redundancy in docking experiments we focused on the centroids structures generated by structural clustering of the binding site over the MD simulations. The most representative 40 centroids (corresponding to clusters covering at least 85% of the MD simulations, S4 Fig) generated for each of the six studied systems (WT apo, WT diclofenac-bound, WT losartan-bound, A477T apo, A477T diclofenac-bound, A477T losartan-bound)
were used for docking. In total, the five drug substrates were docked into 5 crystal and 240 CYP2C9 WT/mutant conformations (see S8 Fig).

Ten independent docking calculations were performed for each protein centroid structure. The three top-scored substrate poses were taken for the analysis. This resulted in 360 substrate poses in MD generated WT CYP2C9 conformations and 360 poses in MD generated A477T mutant conformations for each of the five drugs. The two small NSAID flurbiprofen and diclofenac were found to bind in different positions in the CYP2C9 binding pocket, while the other ligands showed preferable binding zones. The comparison of the WT and the mutant showed differences for the diclofenac and flurbiprofen docking in the so called "flurbiprofen-binding" subpocket corresponding to the catalytic site subpocket close to the heme (Fig 7A and 7B). These two ligands were accommodated closer to the helix I (SRS4) in the mutant structures compared to the WT ones. The large glimepiride and losartan molecules were oriented toward the Cpd I oxygen atom and no important differences were detected for their docking positions when comparing the WT and the mutant structures (Fig 7C and 7D). Interestingly, the warfarin was found to bind preferably in the "warfarin-binding" sub-pocket corresponding to the access channel and to an intermediate binding state before reaching a position for the catalysis (Fig 7E and 7F). The warfarin also changed its docking positions in the mutant structures moving toward helix I as it was observed for diclofenac and flurbiprofen.

Then, we focused on the SOM positions regarding the Cpd I oxygen atom, which was included in the docking calculations. We considered that the docking poses showing a distance between the SOM (Fig 1) and the Cpd I catalytic oxygen atom < 6 Å corresponded to a substrate binding mode appropriate for the catalysis (see S8 Fig and Table 2). As seen in S8 Fig, docking of substrates only into the crystal structures of CYP2C9 is not sufficient to explore the binding behavior of the substrates. We did not observe significant differences between the WT and the mutant regarding the number of acceptable poses. However, noticeable differences between the WT and the A477T mutant were found (see Table 2) when focusing on docking poses obtained in the highly representative centroids of the MD simulations (representativity > 6% of the MD simulations). A larger number of acceptable poses for the catalytic reaction were obtained for the WT as compared to the mutant for the most representative centroids. These results indicated that some conformations of the A477T mutant were still able to accommodate the substrates, however, mutant conformations that properly bound the substrates were visited for a shorter period of time of the MD simulations compared to the WT. The detailed analysis showed that the 3 top-scored docking poses with distance < 6 Å between the SOM and the Cpd I oxygen were obtained in: 3 centroids (cluster numbers 1, 2 and 3) generated by WT apo MD simulations; 2 centroids (cluster numbers 0 and 1) generated by WT diclofenac-bound MD simulations; 3 centroids (cluster numbers 0, 1 and 2) generated by WT losartan-bound MD simulations, 1 centroid (cluster number 0) generated by A477T apo MD simulations; 2 centroids (cluster numbers 0 and 1) generated by A477T diclofenac-bound MD simulations; 3 centroids (cluster numbers 0, 1 and 2) generated by A477T losartan-bound MD simulations. They represented 25.76%, 19.64%, 24.98%, 8.44%, 20.31% and 23.36% of the entire 250 ns MD simulations for each of the six studied systems, respectively (see S4A Fig and S1 Table). We should note that only 1 centroid structure suitable for substrate docking was identified from the A477T apo mutant MD simulations, thus confirming the altering effect of the mutation on the binding site conformation. As observed from the MD RMSF analysis, the centroid conformations generated by the MD simulations in presence of a bound substrate were less affected by the mutation regarding the substrate docking. In total, 8 WT and 6 mutant centroids (noted here as "best" centroids) were identified to properly accommodate the substrates with distance < 6 Å between the SOM and the Cpd I oxygen atom required for the catalytic reaction.
<figure><img src="image_15.png" /></figure>

<figure><img src="image_16.png" /></figure>

<figure><img src="image_17.png" /></figure>

<figure><img src="image_18.png" /></figure>

<figure><img src="image_19.png" /></figure>

<figure><img src="image_20.png" /></figure>
**Fig 7.** The most representative three docking positions obtained by structural clustering with RMSD of 4 Å over each **substrate**. The ligands are represented as blue sticks for WT and red sticks for A477T mutant. (A) diclofenac, (B) flurbiprofen, (C) glimepiride, (D) losartan and (E) warfarin. (F) The experimentally known positions of flurbiprofen and warfarin, both represented in sticks. CYP2C9 is represented in white cartoon, Cpd I in sticks and A/T477 residue is represented by a purple sphere.

https://doi.org/10.1371/journal.pone.0197249.g007

In order to analyze other possible positions of the losartan in the CYP2C9 active site we compared the docking poses of losartan with the recently reported co-crystalized structures, PDB ID 5XXI (WT–losartan bound) and 5X23 (A477T mutant–losartan bound). The most similar positions were found when the losartan was docked into the active site of the centroid number 4 of CYP2C9 WT apo representing 5.76% of the entire MD simulations for WT apo (RMSD = 5.41 Å) and of the centroid number 0 of A477T mutant–losartan bound representing 6.59% of the entire MD simulations for the A477T apo mutant–losartan bound (RMSD = 4.36 Å). These docking poses suggest that losartan can also adopt positions similar to the co-crystalized ones in CYP2C9 WT and CYP2C9 A477T mutant–losartan bound. Although RMSD are not very low, the losartan is similarly oriented in the docked and experimental structures (see S9 Fig). However, such positions are not consistent with the hydroxylation catalytic reaction.

In the three best WT apo MD derived centroids the residues F476 and L362 adopted different conformations (Fig 8). The residue L366 kept its hydrophobic contacts with the heme as in the apo X-ray structures. On the contrary, L366 was observed to shift far from the heme in the best apo mutant centroid. In the best centroids derived from the WT diclofenac-bound state MD simulations L361, L362, F100 and F114 show different orientations in the binding pocket due to the large motions of the diclofenac during the MD simulations. Compared to the WT diclofenac-bound MD based centroids, the mutant diclofenac-bound MD based centroids showed very important shift of F100 and the β turn 4 containing F476 (see Fig 8). In the centroids generated from the MD simulations of the WT and A477T mutant losartan-bound states, L362 adopted different conformations. Overall the above cited residues showed large mobility and adopted different conformations depending on the presence of the mutation A477T or a bound substrate during the MD simulations.

## Conclusions

We investigated molecular mechanisms implicated in the altered substrate metabolism in CYP2C9.30 exhibiting < 50% of the WT enzyme activity. MD simulations performed for the active species of the enzyme, Compound I, in apo or substrate-bound states, and binding energy analyses for several drugs revealed molecular mechanisms involved in the defective drug metabolism in carriers of CYP2C9*30. Our data demonstrated that the A477T mutation present in CYP2C9.30 led to increased rigidity of the SRS1 and SRS5 regions and to shifting of

<table><caption>**Table 2.** Number of the top 3 scored docking poses within a distance of 6 Å between the SOM and the Cpd I oxygen docked in 120 WT and 120 A477T structures corresponding to the 40 MD centroids for each MD case. The centroids showing > 6% representativity correspond to a cluster covering at least 6% of the conformers generated along the 250 ns MD simulations.</caption><thead><tr><th>Substrate</th><th>Nº of poses into 120 WT centroids</th><th>Nº of poses into 120 A477T centroids</th><th>Nº of poses into WT centroids showing > 6% representativity</th><th>Nº of poses into A477T centroids showing > 6% representativity</th></tr></thead><tbody><tr><td>Diclofenac</td><td>85</td><td>83</td><td>8</td><td>5</td></tr><tr><td>Flurbiprofen</td><td>92</td><td>102</td><td>12</td><td>7</td></tr><tr><td>Glimepiride</td><td>63</td><td>72</td><td>7</td><td>2</td></tr><tr><td>Losartan</td><td>83</td><td>89</td><td>7</td><td>3</td></tr><tr><td>Warfarin</td><td>86</td><td>77</td><td>10</td><td>4</td></tr></tbody></table>

https://doi.org/10.1371/journal.pone.0197249.t002
<figure><img src="image_21.png" /><figcaption>Fig 8. Key residues of the active site of the best identified CYP2C9 structures for substrate docking generated from MD simulations of diclofenac-bound CYP2C9. The side chains of the diclofenac-bound WT structures are colored in grey (centroid 0), cyan (centroid 1), light green (centroid 2). The side chains of the diclofenac-bound A477T structures are colored in red (centroid 0) and violet (centroid 1).</figcaption></figure>

https://doi.org/10.1371/journal.pone.0197249.g008

the β turn 4 of SRS6 toward helix F thus decreasing the substrate access to some protein channels. We observed different behaviors of the substrates diclofenac and losartan during the MD simulations for the WT and the A477T mutant, suggesting that A477T can alter differently the metabolism of various substrates. The overall stabilization of the WT and the mutant structures in the cases of a bound substrate observed for the performed here MD simulations as well as for the available crystal structures demonstrated the importance to study structural and dynamics properties of the apo protein in order to detect molecular mechanisms due to the mutation without screening effects of a bound ligand. Molecular docking of five drug substrates of CYP2C9 into WT and A477T mutant structures generated by MD simulations indicated that some conformations of the A477T mutant were able to correctly accommodate the substrates, however, their number was very low compared to the WT MD conformations. All together these conformational and dynamic changes are believed to be involved in the governing mechanism of the catalytic activity dysfunction. Finally, an ensemble of representative conformations of the WT and A477T mutant properly accommodating the substrates were
identified, these structures can be used for future prediction of unknown CYP2C9.1 and CYP2C9.30 substrates and drug-drug interactions.

## Materials and methods

### System preparation

Six crystallographic structures of CYP2C9 were taken from the Protein Data Bank (PDB) and analysed: the native sequence flurbiprofen-bound one (PDB ID 1R9O) [28], losartan-bound one (PDB ID 5XXI [13]), losartan-bound CYP2C9.30 (containing the mutation A477T) (PDB ID 5X23 [13]), and three structures of CYP2C9 containing the mutations K206E, I215V, C216Y, S220P, P221A, I222L, I223L: the apo CYP2C9 one (PDB ID 1OG2 [25]), the warfarin-bound one (PDB ID 1OG5 [25] and the inhibitor 2QJ-bound one (PDB ID 4NZ2 [29]. We decided to construct our initial model of CYP2C9 based on the only one available apo crystal structure of CYP2C9 (PDB ID 1OG2) in order to get away possible effects of bound ligands on the initial enzyme structure. We restored the 7 original CYP2C9 residues mutated in the crystallized enzyme. In addition, the side-chain orientation of R108 was rotated in order to permit its hydrogen bond formation with bound substrates. The pKa values of the titratable groups of the CYP2C9 WT and A477T mutant proteins were calculated with the FDPB approach using the web server PCE [30]. The hydrogen atoms were then added using the DockPrep utility of Chimera [31]. In order to check the effects of the 7 mutations present in the crystallized apo 1OG2 structure on our starting model conformation of the CYP2C9 WT, we performed initial molecular dynamics (MD) simulations of the crystallized apo 1OG2 structure containing these mutations and the apo 1OG2 structure with restored native CYP2C9 residues. The same MD protocol of 250 ns used for the other MD simulations and described in the next section was employed. The analysis has shown that the mutations present in 1OG2 and 1OG5 do not critically affect the overall protein dynamics and structure neither the dynamic behaviour of the BC loop and the β turn 4 shown in this study to be affected by the A477T mutation in CYP2C9.30.

The initial models of diclofenac-bound and losartan-bound WT and A477T mutant structures were generated using docking with Autodock 4.2 program [32]. The best substrate conformations showing the lowest binding energies and a correct position of SOM (see Fig 1) (i.e. SOM within distance of 6 Å to the heme catalytic oxygen) were chosen as starting conformations for MD simulations.

### Molecular dynamics simulations

The six initial CYP2C9 models (WT apo, A477T apo, WT diclofenac-bound, A477T diclofenac-bound, WT losartan-bound, A477T losartan-bound) were fully solvated with explicit TIP3P water molecules [33] and neutralized resulting in systems of about 60,000 atoms. The CHARMM36 force field was used [34]. The CGenFF force field [35] was used to parameterize diclofenac with acceptable penalty score [36]. In contrary, very high penalty score (> 205) were obtained for losartan using CGenFF. Thus, we performed quantum-mechanical calculations by means of density functional theory (DFT) to obtain the losartan parameters. Geometry optimization was performed by using B3LYP functional and 6-31g(d) basis set. The obtained losartan structure was located at a minimum of the potential energy surface through frequency calculations at the same level of energy (see S10 Fig). Atomic charges derived from the electrostatic potential were calculated using the ChelpG scheme. All calculations were performed with the Gaussian 09 software (version C.01, Gaussian, Inc., Wallingford CT, 2010). Atomic parameters of losartan obtained with quantum mechanical calculations are given in S2 and S3 Tables. The force field for the heme Cpd I was taken from [37].
The NAMD2.6 program [38] was used for minimization and MD simulations. All systems containing water and ion molecules were firstly minimized with 5000 steps of conjugate gradient to remove poor contacts and relax the system. The systems were then submitted to a 1 ns MD simulation in order to equilibrate water and ion molecules located in the vicinity of the protein while the protein atoms were kept fixed. In a second stage, all constraints were removed and 5 MD production runs, each of 50 ns, were executed for each of the six CYP2C9 systems. Velocities were randomly assigned with different seeds for each trajectory. The particle mesh Ewald [39] method was used for electrostatic interactions, and the van der Waals interactions were computed with a switching function applied in the range 9–10 Å. The NPT ensemble was used (1.01325 bar and 300 K) with Langevin dynamics and a Nosé-Hoover-Langevin piston pressure control [40]. The integration step was set to 1 fs. The conformations were saved every 5 ps, resulting in 10,000 conformations per MD trajectory. The MD simulations were performed for 250 ns for each of the six studied systems. To probe the impact of longer MD simulations, we run an additional single MD simulation of 200 ns for the A477T mutant apo. Interestingly, the single MD simulation of 200 ns explored much less conformational space than the combined 5x50ns MD simulations and the conformational space explored by the single 200 ns run was covered by the 5 x50 ns MD simulations for the CYP2C9 A477T mutant (see S4 Fig).

Representative structures of the substrate binding pocket

To find the most representative structures of the substrate binding pocket for each of the six systems we merged the five MD trajectories. For the structural clustering, we took 25,000 protein conformations (1 conformation every 10 ps) for each system and used the quality threshold algorithm [41] of the VMD software [42]. Following our previous protocols for multiple receptor ensemble generation [43, 44], a Root Mean Square Deviation (RMSD) distance of at least 1.0Å for all atoms of the binding site was chosen as the clustering criteria. The following residues were considered as part of the enzyme binding pocket: I99, F100, L102, A103, A106, N107, R108, V113, F114, L201, N204, I205, I207, L208, S209, K232, L233, L234, A291, V292, D293, F295, G296, A297, G298, T299, E300, T301, L366, F476 and A(T)477. We took the centroid structure of each cluster in order to define a representative set of protein conformations for subsequent analysis.

Substrate docking

Five drugs known to be CYP2C9 substrates were docked into the crystal CYP2C9 structures PDB IDs: 1OG5 and 4NZ2 with restored native residues, 1R9O, 5XXI, 5X23, as well as into the initial models of the WT or the mutant A477T and into 240 representative conformations selected after the MD simulations for the six studied systems. The heme was kept in the intermediate Cpd I state. Drug protonation states were calculated at pH 7.4 with the major microspecies option of the ChemAxon package (ChemAxon www.chemaxon.com). Interestingly, ChemAxon predicted two different pKa values for losartan, 4.29 or 8.0, starting from two different smile structures downloaded from BindingDB and DrugBank, respectively. We analyzed possible protonation of the losartan molecules co-crystallized with CYP2C9 in the PDB structures ID 5X23 and 5XXI and decided to consider the protonated state of the losartan, which corresponds better to the losartan environment in the binding sites. Gasteiger atom charges were added using the AutoDockTools package. Molecular dockings were carried out with the Autodock 4.2 program [32]. A grid spacing of 0.375 Å was used. The search space was centered on the protein binding sites with a box size of 24x25x30 Å. Ten independent docking runs were performed using the Lamarckian genetic algorithm [45] with an initial population
of 150 random individuals and a maximum number of 25,000,000 energy evaluations. During docking, ligand molecules were flexible while the protein was kept rigid.

## Channel dynamics determining

The channel dynamics for the WT and the A477T mutant were investigated with the CAVER software [46]. All CYP2C9 conformations generated from the MD simulations were analyzed using CAVER (50,000 for each system). A probe radius of 1.2 Å was used for channel prediction. The starting point was set at the center of the binding pocket. In case of a bound substrate during the MD simulations, the ligand was removed prior to channel searching.

## Supporting information

<caption>S1 Fig. 3D structure of CYP2C9.</caption>

(PDF)

<caption>S2 Fig. Root Mean Square Deviations of all heavy atoms for each MD simulation.</caption>

(PDF)

<caption>S3 Fig. RMSF differences of Ca atoms between the WT and A477T proteins for the apo, diclofenac bound and losartan bound states.</caption>

(PDF)

<caption>S4 Fig. Structural clustering of the MD generated conformers.</caption>

(PDF)

<caption>S5 Fig. Hydrophobic contacts between F476 and other residues of the CYP2C9 binding pocket monitored over the MD simulations.</caption>

(PDF)

<caption>S6 Fig. Initial and preferential positions of losartan in the CYP2C9 active site during the MD simulations of the WT and A477T variant losartan-bound systems.</caption>

(PDF)

<caption>S7 Fig. RMSD of Compound I during the MD simulations.</caption>

(PDF)

<caption>S8 Fig. Docking scores of docking poses as a function of distance between the catalytic oxygen of Cpd I and SOM.</caption>

(PDF)

<caption>S9 Fig. Docking poses of losartan similar to the co-crystallized structures.</caption>

(PDF)

<caption>S10 Fig. Molecular structure of losartan as obtained with quantum-mechanical geometry optimization.</caption>

(PDF)

<caption>S1 Table. Centroid structures with representativity ≥ 6% along the MD simulations allowing to dock the substrates with distances between SOM and Cpd I catalytic oxygen within 6 Å.</caption>

(PDF)

<caption>S2 Table. Atomic charges of losartan used as parameters for the molecular dynamics simulations.</caption>

(PDF)
S3 Table. Atomic parameters of losartan obtained with quantum mechanical calculations.
(PDF)

## Acknowledgments

Supports from the French National Research Institute Inserm, the University Paris Diderot and ANR (ToxME) are greatly appreciated. PHM and CMA thank Brazilian agencies FAPESP and CNPq. We thank Pr. O. Taboureau for helpful discussions.

## Author Contributions

**Conceptualization:** Bruno O. Villoutreix, Maria A. Miteva.

**Data curation:** Céline M. Labbé, Charline Fagnen.

**Formal analysis:** Céline M. Labbé, Charline Fagnen, Cassiano M. Aono, Paula Homem-de-Mello.

**Investigation:** Maxime Louet, Céline M. Labbé, Charline Fagnen, Cassiano M. Aono, Paula Homem-de-Mello, Maria A. Miteva.

**Methodology:** Maxime Louet, Charline Fagnen, Cassiano M. Aono, Paula Homem-de-Mello, Bruno O. Villoutreix, Maria A. Miteva.

**Project administration:** Maria A. Miteva.

**Supervision:** Maria A. Miteva.

**Validation:** Maxime Louet.

**Visualization:** Céline M. Labbé.

**Writing – original draft:** Maxime Louet, Maria A. Miteva.

**Writing – review & editing:** Maxime Louet, Bruno O. Villoutreix, Maria A. Miteva.

## References

1. Hirota T, Eguchi S, leir il. Impact of genetic polymorphisms in CYP2C9 and CYP2C19 on the pharmacokinetics of clinically used drugs. Drug Metab Pharmacokinet. 2013; 28:28–37. PMID: 23165865

2. Moroy G, Martiny VY, Vayer P, Villoutreix BO, Miteva MA. Toward in silico structure-based ADMET prediction in drug discovery. Drug Discov Today. 2012; 17(1–2):44–55. Epub 2011/11/08. S1359-6446(11) 00362-X [pii] https://doi.org/10.1016/j.drudis.2011.10.023 PMID: 22056716.

3. Isvoran A, Louet M, Vladoiu DL, Craciun D, Loriot MA, Villoutreix BO, et al. Pharmacogenomics of the cytochrome P450 2C family: impacts of amino acid variations on drug metabolism. Drug Discov Today. 2017; 22:366–76. https://doi.org/10.1016/j.drudis.2016.09.015 PMID: 27693711.

4. Dai D, Zeldin DC, Blaisdell JA, Chanas B, Coulter SJ, Ghanayem BI, et al. Polymorphisms in human CYP2C8 decrease metabolism of the anticancer drug paclitaxel and arachidonic acid. Pharmacogenetics. 2001; 11(7):597–607. Epub 2001/10/23. PMID: 11668219.

5. Fujikura K, Ingelman-Sundberg M, Lauschke VM. Genetic variation in the human cytochrome P450 supergene family. Pharmacogenet Genomics. 2015; 25(12):584–94. Epub 2015/09/05. https://doi.org/10.1097/FPC.0000000000000172 PMID: 26340336.

6. Wang B, Wang J, Huang SQ, Su HH, Zhou SF. Genetic polymorphism of the human cytochrome P450 2C9 gene and its clinical significance. Curr Drug Metab. 2009; 10(7):781–834. Epub 2009/11/21. PMID: 19925388.

7. Hiratsuka M. In vitro assessment of the allelic variants of cytochrome P450. Drug Metab Pharmacokinet. 2012; 27(1):68–84. Epub 2011/11/02. JST.JSTAGE/dmpk/DMPK-11-RV-090 [pii]. PMID: 22041138.
8. Martiny VY, Miteva MA. Advances in molecular modeling of human cytochrome P450 polymorphism. J Mol Biol. 2013; 425(21):3978-92. Epub 2013/07/17. S0022-2836(13)00442-7 [pii] https://doi.org/10.1016/j.jmb.2013.07.010 PMID: 23856621.

9. Yin T, Maekawa K, Kamide K, Saito Y, Hanada H, Miyashita K, et al. Genetic variations of CYP2C9 in 724 Japanese individuals and their impact on the antihypertensive effects of losartan. Hypertens Res. 2008; 31(8):1549-57. https://doi.org/10.1291/hypres.31.1549 PMID: 18971529.

10. Maekawa K, Fukushima-Uesaka H, Tohkin M, Hasegawa R, Kajio H, Kuzuya N, et al. Four novel defective alleles and comprehensive haplotype analysis of CYP2C9 in Japanese. Pharmacogenet Genomics. 2006; 16(7):497-514. Epub 2006/06/22. https://doi.org/10.1097/01.fpc.0000215069.14095.c6 PMID: 16788382.

11. Maekawa K, Harakawa N, Sugiyama E, Tohkin M, Kim SR, Kaniwa N, et al. Substrate-dependent functional alterations of seven CYP2C9 variants found in Japanese subjects. Drug Metab Dispos. 2009; 37 (9):1895-903. https://doi.org/10.1124/dmd.109.027003 PMID: 19541829.

12. Melet A, Assrir N, Jean P, Pilar Lopez-Garcia M, Marques-Soares C, Jaouen M, et al. Substrate selectivity of human cytochrome P450 2C9: importance of residues 476, 365, and 114 in recognition of diclofenac and sulfaphenazole and in mechanism-based inactivation by tienilic acid. Arch Biochem Biophys. 2003; 409(1):80-91. PMID: 12464247.

13. Maekawa K, Adachi M, Matsuzawa Y, Zhang Q, Kuroki R, Saito Y, et al. Structural Basis of Single-Nucleotide Polymorphisms in Cytochrome P450 2C9. Biochemistry. 2017; 56:5476-80 https://doi.org/10.1021/acs.biochem.7b00795 PMID: 28972767

14. Sano E, Li W, Yuki H, Liu X, Furihata T, Kobayashi K, et al. Mechanism of the decrease in catalytic activity of human cytochrome P450 2C9 polymorphic variants investigated by computational analysis. J Comput Chem. 2010; 31(15):2746-58. Epub 2010/09/15. https://doi.org/10.1002/jcc.21568 PMID: 20839301.

15. Roberts AG, Cheesman MJ, Primak A, Bowman MK, Atkins WM, Rettie AE. Intramolecular heme ligation of the cytochrome P450 2C9 R108H mutant demonstrates pronounced conformational flexibility of the B-C loop region: implications for substrate binding. Biochemistry. 2010; 49(40):8700-8. https://doi.org/10.1021/bi100911q PMID: 20815369; PubMed Central PMCID: PMC3005909.

16. Lee MY, Borgiani P, Johansson I, Oteri F, Mkrtchian S, Falconi M, et al. High warfarin sensitivity in carriers of CYP2C9*35 is determined by the impaired interaction with P450 oxidoreductase. Pharmacogenomics J. 2014; 14:343-9. https://doi.org/10.1038/tpj.2013.41 PMID: 24322786

17. Dodani SC, Kiss G, Cahn JK, Su Y, Pande VS, Arnold FH. Discovery of a regioselectivity switch in nitrating P450s guided by molecular dynamics simulations and Markov models. Nat Chem. 2016; 8 (5):419-25. https://doi.org/10.1038/nchem.2474 PMID: 27102675; PubMed Central PMCID: PMC4843824.

18. Wang B, Li C, Dubey KD, Shaik S. Quantum mechanical/molecular mechanical calculated reactivity networks reveal how cytochrome P450cam and Its T252A mutant select their oxidation pathways. J Am Chem Soc. 2015; 137:7379-90. https://doi.org/10.1021/jacs.5b02800 PMID: 26011529

19. Ashkenazy H, Abadi S, Martz E, Chay O, Mayrose I, Pupko T, et al. ConSurf 2016: an improved methodology to estimate and visualize evolutionary conservation in macromolecules. Nucleic Acids Res. 2016; 44(W1):W344-50. https://doi.org/10.1093/nar/gkw408 PMID: 27166375; PubMed Central PMCID: PMCPMC4987940.

20. Adzhubei I, Jordan DM, Sunyaev SR. Predicting functional effect of human missense mutations using PolyPhen-2. Curr Protoc Hum Genet. 2013; Chapter 7:Unit7.20. https://doi.org/10.1002/0471142905. hg0720s76 PMID: 23315928; PubMed Central PMCID: PMCPMC4480630.

21. Dehouck Y, Kwasigroch JM, Gilis D, Rooman M. PoPMuSiC 2.1: a web server for the estimation of protein stability changes upon mutation and sequence optimality. BMC Bioinformatics. 2011; 12:151. https://doi.org/10.1186/1471-2105-12-151 PMID: 21569468; PubMed Central PMCID: PMCPMC3113940.

22. Frappier V, Chartier M, Najmanovich RJ. ENCoM server: exploring protein conformational space and the effect of mutations on protein function and stability. Nucleic Acids Res. 2015; 43(W1):W395-400. https://doi.org/10.1093/nar/gkv343 PMID: 25883149; PubMed Central PMCID: PMCPMC4489264.

23. Schymkowitz J, Borg J, Stricher F, Nys R, Rousseau F, Serrano L. The FoldX web server: an online force field. Nucleic Acids Res. 2005; 33(Web Server issue):W382-8. https://doi.org/10.1093/nar/gki387 PMID: 15980494; PubMed Central PMCID: PMCPMC1160148.

24. Rittle J, Green MT. Cytochrome P450 compound I: capture, characterization, and C-H bond activation kinetics. Science. 2010; 330(6006):933-7. https://doi.org/10.1126/science.1193478 PMID: 21071661.

25. Williams PA, Cosme J, Ward A, Angove HC, Matak Vinkovic D, Jhoti H. Crystal structure of human cytochrome P450 2C9 with bound warfarin. Nature. 2003; 424(6947):464-8. https://doi.org/10.1038/nature01862 PMID: 12861225.
26. Cojocaru V, Winn PJ, Wade RC. The ins and outs of cytochrome P450s. Biochim Biophys Acta. 2007; 1770(3):390-401. https://doi.org/10.1016/j.bbagen.2006.07.005 PMID: 16920266.

27. Cojocaru V, Balali-Mood K, Sansom MS, Wade RC. Structure and dynamics of the membrane-bound cytochrome P450 2C9. PLoS Comput Biol. 2011; 7(8):e1002152. https://doi.org/10.1371/journal.pcbi.1002152 PMID: 21852944; PubMed Central PMCID: PMCPMC3154944.

28. Wester MR, Yano JK, Schoch GA, Yang C, Griffin KJ, Stout CD, et al. The structure of human cytochrome P450 2C9 complexed with flurbiprofen at 2.0-A resolution. J Biol Chem. 2004; 279(34):35630-7. https://doi.org/10.1074/jbc.M405427200 PMID: 15181000.

29. Branden G, Sjogren T, Schnecke V, Xue Y. Structure-based ligand design to overcome CYP inhibition in drug discovery projects. Drug Discov Today 2014; 19:905-11. https://doi.org/10.1016/j.drudis.2014.03.012 PMID: 24642031

30. Miteva MA, Tuffery P, Villoutreix BO. PCE: web tools to compute protein continuum electrostatics. Nucleic Acids Res. 2005; 33(Web Server issue):W372-5. https://doi.org/10.1093/nar/gki365 PMID: 15980492; PubMed Central PMCID: PMC1160126.

31. Pettersen EF, Goddard TD, Huang CC, Couch GS, Greenblatt DM, Meng EC, et al. UCSF Chimera—a visualization system for exploratory research and analysis. J Comput Chem. 2004; 25(13):1605-12. https://doi.org/10.1002/jcc.20084 PMID: 15264254.

32. Morris GM, Huey R, Lindstrom W, Sanner MF, Belew RK, Goodsell DS, et al. AutoDock4 and Auto-DockTools4: Automated docking with selective receptor flexibility. J Comput Chem. 2009; 30(16):2785-91. https://doi.org/10.1002/jcc.21256 PMID: 19399780; PubMed Central PMCID: PMC2760638.

33. Jorgensen WL, Chandrasekhar J, Madura JD, Impey RW, Klein ML. Comparison of simple potential functions for simulating liquid water. The Journal of Chemical Physics. 1983; 79:926-35. https://doi.org/10.1063/1.445869

34. Best RB, Zhu X, Shim J, Lopes PEM, Mittal J, Feig M, et al. Optimization of the Additive CHARMM All-Atom Protein Force Field Targeting Improved Sampling of the Backbone φ, ψ and Side-Chain χ1 and χ2 Dihedral Angles. Journal of Chemical Theory and Computation. 2012; 8:3257-73. https://doi.org/10.1021/ct300400x PMID: 23341755

35. Vanommeslaeghe K, Hatcher E, Acharya C, Kundu S, Zhong S, Shim J, et al. CHARMM general force field: A force field for drug-like molecules compatible with the CHARMM all-atom additive biological force fields. Journal of Computational Chemistry. 2010; 31:671-90. https://doi.org/10.1002/jcc.21367 PMID: 19575467

36. Vanommeslaeghe K, Raman EP, MacKerell AD Jr. Automation of the CHARMM General Force Field (CGenFF) II: assignment of bonded parameters and partial atomic charges. J Chem Inf Model. 2012; 52(12):3155-68. https://doi.org/10.1021/ci3003649 PMID: 23145473; PubMed Central PMCID: PMCPMC3528813.

37. Bathelt CM, Zurek J, Mulholland AJ, Harvey JN. Electronic structure of compound I in human isoforms of cytochrome P450 from QM/MM modeling. J Am Chem Soc. 2005; 127(37):12900-8. https://doi.org/10.1021/ja0520924 PMID: 16159284.

38. Phillips JC, Braun R, Wang W, Gumbart J, Tajkhorshid E, Villa E, et al. Scalable molecular dynamics with NAMD. Journal of Computational Chemistry. 2005; 26:1781-802. https://doi.org/10.1002/jcc.20289 PMID: 16222654

39. Darden T, York D, Pedersen L. Particle mesh Ewald: An N·log(N) method for Ewald sums in large systems. The Journal of Chemical Physics. 1993; 98:10089-92. https://doi.org/10.1063/1.464397

40. Feller SE, Zhang Y, Pastor RW, Brooks BR. Constant pressure molecular dynamics simulation: The Langevin piston method. The Journal of Chemical Physics. 1995; 103:4613-21. https://doi.org/10.1063/1.470648

41. Heyer LJ, Kruglyak S, Yooseph S. Exploring expression data: identification and analysis of coexpressed genes. Genome Res. 1999; 9(11):1106-15. PMID: 10568750; PubMed Central PMCID: PMC310826.

42. Humphrey W, Dalke A, Schulten K. VMD: visual molecular dynamics. Journal of Molecular Graphics. 1996; 14:33-8, 27-8. PMID: 8744570

43. Martiny VY, Carbonell P, Chevillard F, Moroy G, Nicot AB, Vayer P, et al. Integrated structure- and ligand-based in silico approach to predict inhibition of cytochrome P450 2D6. Bioinformatics. 2015; 31(24):3930-7. https://doi.org/10.1093/bioinformatics/btv486 PMID: 26315915.

44. Moroy G, Sperandio O, Rielland S, Khemka S, Druart K, Goyal D, et al. Sampling of conformational ensemble for virtual screening using molecular dynamics simulations and normal mode analysis. Future Med Chem. 2015; 7(17):2317-31. https://doi.org/10.4155/fmc.15.150 PMID: 26599419.
45. Morris GM, Goodsell DS, Halliday RS, Huey R, Hart WE, Belew RK, et al. Automated Docking Using a Lamarckian Genetic Algorithm and an Empirical Binding Free Energy Function. Journal of Computational Chemistry. 1998; 19(14):1639–62.

46. Chovancova E, Pavelka A, Benes P, Strnad O, Brezovsky J, Kozlikova B, et al. CAVER 3.0: a tool for the analysis of transport pathways in dynamic protein structures. PLoS Comput Biol. 2012; 8(10): e1002708. https://doi.org/10.1371/journal.pcbi.1002708 PMID: 23093919; PubMed Central PMCID: PMC3475669.