Exhaustive computational search of ionic-charge clusters that mediate interactions between mammalian cytochrome P450 (CYP) and P450-oxidoreductase (POR) proteins

Exhaustive computational search of ionic-charge clusters that mediate interactions between mammalian cytochrome P450 (CYP) and P450-oxidoreductase (POR) proteins

Computational Biology and Chemistry 34 (2010) 42–52 Contents lists available at ScienceDirect Computational Biology and Chemistry journal homepage: ...

522KB Sizes 0 Downloads 55 Views

Computational Biology and Chemistry 34 (2010) 42–52

Contents lists available at ScienceDirect

Computational Biology and Chemistry journal homepage: www.elsevier.com/locate/compbiolchem

Research article

Exhaustive computational search of ionic-charge clusters that mediate interactions between mammalian cytochrome P450 (CYP) and P450-oxidoreductase (POR) proteins Alexander Zawaira a,∗ , Marco Gallotta b , Natasha Beeton-Kempen a , Lauren Coulson a , Patrick Marais b , Michelle Kuttel b , Jonathan Blackburn a a b

Division of Medical Biochemistry, Institute for Infectious Disease & Molecular Medicine, Faculty of Health Sciences, University of Cape Town, Observatory 7925, South Africa Department of Computer Science, University of Cape Town, Private Bag X3, Rondebosch 7701, South Africa

a r t i c l e

i n f o

Article history: Received 26 September 2009 Received in revised form 20 October 2009 Accepted 21 October 2009 Keywords: Cytochrome P450 Cytochrome P450 reductase P450-oxidoreductase Protein–protein interfaces Computational modelling Ionic interactions

a b s t r a c t In this work, a model for the interaction between CYP2B4 and the FMN domain of rat P450-oxidoreductase is built using as template the structure of a bacterial redox complex. Amino acid residues identified in the literature as cytochrome P450 (CYP)–redox partner interfacial residues map to the interface in our model. Our model supports the view that the bacterial template represents a specific electron transfer complex and moreover provides a structural framework for explaining previous experimental data. We have used our model in an exhaustive search for complementary pairs of mammalian CYP and P450-oxidoreductase (POR) charge clusters. We quantitatively show that among the previously defined basic clusters, the 433K–434R cluster is the most dominant (32.3% of interactions) and among the acidic clusters, the 207D–208D–209D cluster is the most dominant (29%). Our analysis also reveals the previously not described basic cluster 343R–345K (16.1% of interactions) and 373K (3.2%) and the acidic clusters 113D–115E–116E (25.8%), 92E–93E (12.9%), 101D (3.2%) and 179E (3.2%). Cluster pairings among the previously defined charge clusters include the pairing of cluster 421K–422R to cluster 207D–208D–209D. Moreover, 433K–434R and 207D–208D–209D, respectively the dominant positively and negatively charged clusters, are uncorrelated. Instead our analysis suggests that the newly identified cluster 113D–115E–116E is the main partner of the 433K–434R cluster while the newly described cluster 343R–345K is correlated to the cluster 207D–208D–209D. © 2009 Elsevier Ltd. All rights reserved.

1. Introduction The cytochromes P450 (CYP, EC 1.14.14.1) are a superfamily of heme-thiolate enzymes that catalyze the monooxygenation of hydrophobic endogenous and xenobiotic substrates (Bridges et al., 1998). CYP homologs have been sequenced from all lineages of life—including eukaryotes and bacteria (Sevrioukova et al., 1999b; Zawaira et al., 2008). The cytochrome P450 enzymatic cycle includes substrate binding, first electron transfer, oxygen binding, second electron transfer, substrate oxidation and finally, product dissociation. Hence the CYP reaction cycle involves two distinct

Abbreviations: CYP, cytochrome P450; POR, P450-oxidoreductase; FMN, flavin mononucleotide; BM3, Bacillus megaterium 3. ∗ Corresponding author at: Room N3.03, Wernher & Beit Building North, Department of Medical Biochemistry, Institute for Infectious Disease & Molecular Medicine, Faculty of Health Sciences, University of Cape Town, Anzio Road, Observatory 7925, South Africa. Tel.: +27 21 406 6453. E-mail addresses: [email protected], [email protected] (A. Zawaira). 1476-9271/$ – see front matter © 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.compbiolchem.2009.10.006

electron transfer steps (Bridges et al., 1998; Kuznetsov et al., 2006). Biologically relevant CYP-catalyzed monooxygenation reactions include steroid metabolism and biogenesis, xenobiotic metabolism, fatty acid metabolism and antibiotic metabolism (Sevrioukova et al., 1999b). Bacterial CYPs are soluble proteins that have high substrate specificity while their mammalian cousins are endoplasmic reticulum-bound CYPs mostly involved in xenobiotic metabolism and have broad substrate specificity (Nelson and Strobel, 1987). Despite these differences, surveys of known CYP structures from across the different lineages of life show that the CYP structure/fold is highly conserved (Hasemann et al., 1995; Sansen et al., 2007; Schoch et al., 2004; Wester et al., 2004; Yano et al., 2004). A notable feature of this highly conserved CYP fold architecture is that a proximal surface, where the heme cofactor comes closest to the protein surface, can readily be discerned from the distal surface where the heme group is farthest from the protein surface (Hasemann et al., 1995; Poulos et al., 1987; Sansen et al., 2007; Schoch et al., 2004; Wester et al., 2004; Yano et al., 2004). The cytochrome P450 enzyme sources electrons from redox partner systems. Bacterial and mitochondrial CYPs obtain catalytic

A. Zawaira et al. / Computational Biology and Chemistry 34 (2010) 42–52

cycle electrons from small soluble iron–sulfur electron transport proteins such as adrenodoxin (also known as ferredoxin-1, FDX1) and putidaredoxin (Sevrioukova et al., 1999b). On the other hand, mammalian microsomal CYPs source their electrons from NADPH dependent P450-oxidoreductase (POR). POR contains the flavin cofactors FAD and FMN (Sevrioukova et al., 1999b). POR is organised into four domains which are, from the protein N-terminus to the C-terminus, the FMN binding domain, the connecting domain, the FAD binding domain and the NADPH-binding domain (Wang et al., 1997). The FAD cofactor accepts electrons from NADPH and transfers them to the FMN cofactor, which in turn, transfers electrons to the heme iron of substrate-bound CYP (Sevrioukova et al., 1999b; Wei et al., 2007). An important question in CYP research is the elucidation of the protein–protein interactions involved in the process of electron transfer from redox partners to substrate-bound CYP. Following the determination of the first structure of a CYP (CYP450Cam from Pseudomonas putida) in 1985 (Poulos et al., 1985), the most parsimonious answer to that question would have been: CYP is likely to interact with redox partners at the proximal side where the heme centre is most easily accessible. Early models for the CYP–POR complex supported the parsimonious view (Stayton et al., 1989). The current mainstream views on CYP–POR interaction surfaces ultimately trace their origins to the studies on the reduction of cytochrome c by cytochrome b5 (Stonehuerner et al., 1979). Stonehuerner et al. (1979) showed that the interaction between cytochrome c and cytochrome b5 is mediated by complementary electrostatic interactions between specific lysine residues on cytochrome c and specific acidic residues (aspartate or glutamate) on cytochrome b5. Similar experiments have been performed demonstrating the role of electrostatic interactions in mammalian CYP–POR (Bernhardt et al., 1988; Davydov et al., 2000; Kelley et al., 2005). Extensions from the cytochrome c–cytochrome b5 model system to CYP target systems by Poulos and co-workers (1989) and by Davydov et al. (1992) provided the first information about CYP–redox partner interfacial residues. Poulos and co-workers’ extensions to the bacterial putidaredoxin–CYP450Cam system suggested that the P450Cam residues Arg112, Lys344 and Arg364 formed salt bridges with putidaredoxin residues Glu48, Glu44 and Asp60 respectively (Stayton et al., 1989). The P450Cam residues Arg112, Lys344 and Arg364 are homologous to CYP2B4 residues Arg125, Lys421 and Arg443 (Bridges et al., 1998). Davydov et al. (1992) extensions to the CYP2B4–P450-oxidoreductase system suggested that CYP2B4 positions 121–145 are involved in POR recognition. Several other workers have contributed to the identification of CYP–redox partner interfacial residues. These include Juvonen et al. (1992) who showed the involvement of Arg129 in the binding of CYP2A5 to cytochrome b5. Bridges et al. (1998) delimited the following CYP2B4 residues as CYP2B4–POR interfacial residues: R122, R126, R133, F135, K433, R422 and R443. Shen and Kasper (1995) searched for surface exposed acidic clusters within the FMN-binding domain of POR and identified two interfacial residue clusters comprising 207D, 208D and 209D and 213E, 214E and 215D. Zhao et al. (1999) identified the cluster comprising 142D, 144D and 147D. Of the bacterial redox systems that have been used as model systems to address the problem of redox partner recognition in mammalian systems, the cytochrome P450-BM3 bears the closest resemblance to mammalian systems. Flavocytochrome P450-BM3 is a 119 kDa self-sufficient fatty-acid monooxygenase from Bacillus megaterium. The protein consists of a heme domain (BMP) and a FMN/FAD-containing P450-oxidoreductase domain linked together in a single polypeptide chain (Sevrioukova et al., 1999b). Sevrioukova et al. (1999a,b) have determined the structure

43

of a P450-BM3 construct comprising the heme/FMN-containing domain and it has been debated whether this structure is a plausible prototypical representation of a specific electron transfer complex. Here we investigate the relationship between CYP2B4–POR interfacial residues one can delimit using the P450-BM3 structure as template and those which have been delimited to date by the many various experimental approaches listed above. We have also used our model of the CYP2B4–POR complex to implement an exhaustive search for oppositely charged clusters mediating the ionic interactions between mammalian CYP and POR. Comparison of our pairings with Karyakin’s experimentally delimited pairings (Karyakin et al., 2007) has been used to validate our assignments. 2. Materials and methods 2.1. Sequences and structures Sequences were retrieved from the Universal Protein Resource Knowledgebase (UniProtKB/Swiss-Prot) database via the European Bioinformatics Institute website (http://www.ebi.ac.uk). The accession numbers were according to the UniProtKB database unless otherwise stated. Structures were retrieved from the RCSB Protein Data Bank (www.pdb.org). The sequence accession numbers used herein are as follows. CYP1A1: P04798 (Homo sapiens), CYP1A2: P05177 (H. sapiens), CYP2B1: P00176 (Rattus norvegicus) and CYP2B4: P00178 (Oryctolagus cuniculus). The structure PDB codes are as follows. CYP450BM-3: 1BVY (B. megaterium), CYP450Cam: 3CPP (P. putida), CYP1A2: 2HI4 (H. sapiens), CYP2B4: 2BDM (O. cuniculus), CYP2C9: 1OG5 (H. sapiens) and CYP450reductase: 1AMO (R. norvegicus). 2.2. Sequence and structure alignments Protein sequences were aligned using the program ClustalW (Jeanmougin et al., 1998; Thompson et al., 1997, 1994). Structure visualisations and alignments were done using the program Pymol (www.pymol.org). Structure alignments of CYP catalytic domains were evaluated by applying a panel of visual criteria such as the simultaneous achievement of stacking of heme cofactors in parallel planes with less than 3 Å distance separation between corresponding atoms and complete superimposition of I-helices. General good alignment of F and G helices was also considered in accepting or rejecting alignments of structures. Structure alignments of FMN binding domains were evaluated by considering the achievement of stacking of FMN cofactors in parallel planes with corresponding atoms less than 1 Å apart. The root mean square deviation (RMSD) of the alignments was also considered. 2.3. Generation of a structural model for the CYP2B4–P450-oxidoreductase FMN domain complex A structural model for the CYP2B4–P450-oxidoreductase FMN binding domain complex was built by swapping the CYP450BM3 heme domain (1BVY.pdb, chain A) with CYP2B4 (2BDM.pdb) and the CYP450BM-3 FMN binding domain (1BVY.pdb, chain F) with the FMN binding domain of P450-oxidoreductase (1AMO.pdb, chain A, residues 1–235) using the program FATCAT (Ye and Godzik, 2004). In order to facilitate FATCAT chain swapping, Chain names in PDB files were changed using the search and replace facilities of common text editors such as Microsoft Word. 2.4. Treatment of side-chain conformations and preliminary searches for paired charge clusters A structural model for the CYP2B4–POR FMN domain complex is derived as described above making use of separately deter-

44

A. Zawaira et al. / Computational Biology and Chemistry 34 (2010) 42–52

mined structures of CYP2B4 (2BDM.pdb) and POR–FMN domain (1AMO.pdb). The side chains in 2BDM.pdb and 1AMO.pdb are not expected to be in the optimal orientations that are related to the role of ionic interactions in the recognition and proper relative orientation of CYP and POR. Optimisation of side chains to capture re-organisations which occur upon complex formation remains a major and largely inadequately addressed problem in computational biology (Schueler-Furman et al., 2005; Wang et al., 2005). A preliminary search for ionic-charge clusters mediating the CYP–POR interaction began with the use of the program SCWRL (Wang et al., 2008) to optimise side-chain interactions between CYP2B4 and POR–FMN domain in our structural model. The program SCWRL was run with default settings on an IBM PC. We also applied to the initial FATCAT-derived model, other tools for side-chain optimisation such as normal mode analysis-based algorithms for the refinement of side-chain implemented by the Normal Mode Analysis, Deformation and Refinement Server (NOMAD-Ref, http://lorentz.immstr.pasteur.fr/nomad-ref.php) by Lindahl et al. (2006). We subsequently delimited ionic interactions from these models using the Protein Interactions Calculator Server (http://crick.mbu.iisc.ernet.in/∼PIC/) by Tina et al. (2007) and the ESBRI Server (http://bioinformatica.isa.cnr.it/ESBRI/) by Costantini et al. (2008). 2.5. Calculation of protein surfaces and interfacial residues Protein solvent accessible surfaces (Totrov and Abagyan, 1996) were calculated from PDB coordinates by applying algorithms that achieve the computational equivalent of rolling a spherical probe (representing a water molecule) of radius 1.4 Å on the Van der Waals surface of the protein molecules. The marching tetrahedra algorithm (Chan and Purisima, 1998) was used to calculate solvent accessible surfaces. This algorithm is related to the marching cubes algorithm (Lorensen and Cline, 1987) and like marching cubes, requires a three dimensional grid of distances from the surface. Each entry in this grid is the distance from the location of the point it represents to the surface. Entries are positive for points exterior to the surface and negative for points enveloped inside the surface. The marching tetrahedra algorithm uses this grid of distances as a contour map by approximating the surface where the distance is zero. This generates the solvent accessible surface as a set of triangular panels, i.e. the method returns a triangulated numerical surface. In this work, we use Kd-trees (Goodman et al., 2004; Lee and Wong, 1977) to calculate the grid of distances that form the basis of the marching tetrahedra calculations. A Kd-tree containing the coordinates of all N atoms in the protein is constructed in O (N log2 N) time. The distances from each point in the grid are calculated by performing a nearest neighbour query on the Kd-tree in O (log N) time. The distance to the solvent accessible surface is then calculated as the distance to this nearest atom’s centre plus the radius of the atom plus the radius of the water molecule used as a probe. This gives an overall efficiency of O (N log2 N + M log N) to compute an entire grid of M points. Paul Bourke’s polygonisation libraries (http://local.wasp.uwa.edu.au/∼pbourke/geometry/polygonise/) were used together with the Kd-tree functions we defined to achieve an implementation of the marching tetrahedra algorithm. Protein–protein interfaces were calculated from protein binary complexes by considering the region of the solvent accessible surface of one partner that lies within a given Euclidean distance cut-off (i.e. is buried) of the other binding partner. A distance cutoff of 4.5 Å was used as the default as this value is commonly used as the default in such calculations (Jackson et al., 1998; Lu et al., 2003; Ofran and Rost, 2003; Pulim et al., 2008; Tovchigrechko and Vakser, 2005, 2006). We also considered interface residue delim-

itations at greater distance cut-offs such as 6, 6.5 and 7 Å. In our calculation of the interface from a given protein–protein complex, all triangles of the solvent accessible surface of one binding partner (as calculated using the marching tetrahedra algorithm) whose centroids are within the selected distance cut-off (4.5, 6, 6.5 or 7 Å) from the centre of an atom in the other binding partner are considered as part of the binding interface. The rest of the triangles are deleted. In order to improve the production of continuous interfaces, morphological operations of dilation and erosion (Giardina and Dougherty, 1988) are applied on the current interface. In our implementation, the dilation step adds all triangles within a specified distance from the current interface, closing all gaps under a certain size. The erosion step then removes all triangles in the current interface that are within this same distance from the complement of the interface—thereby removing all triangles that were added in the dilation step but did not close any holes. A distance cut-off of 1.4 was chosen as the default for the dilation step after optimisation studies. Given a calculated interface, amino acid residues that are considered to be part of this interface (i.e. making a contribution to this interface) were calculated as follows. A normal is calculated from each triangle in the interface. The triangle is assigned to the atom with smallest vertical distance from its centre to this normal. Next, the triangle is assigned to the amino acid residue that contains this atom. Finally, amino acids with at least one triangle assigned to them are assigned to the interface.

2.6. Projection of literature CYP–redox partner interfacial residues onto the CYP2B4–P450-oxidoreductase FMN domain complex The literature was searched for experimentally delimited CYP–redox partner interfacial residues. The goal was to project these CYP–redox partner interfacial residues from the literature onto the CYP450BM-3 (1BVY.pdb)-derived structural model for the CYP2B4–POR FMN binding domain complex. This would permit direct visual assessment of the agreement between the CYP450BM-3-derived model for the CYP2B4–POR FMN binding domain complex and the experimentally delimited interfacial residues. CYP2B4 residues homologous to experimentally delimited CYP–POR interfacial residues from other members of the CYP2 family were identified by pair-wise alignment using the program ClustalW (Thompson et al., 1997; Thompson et al., 1994). This was done for Arg129 of CYP2A5 from Mus musculus, Arg144 of CYP2C9 from H. sapiens and residues 116–134 of CYP2B1 from R. norvegicus. Wherever it was possible, our assignments were checked against assignments in the literature. CYP2B4 residues homologous to CYP–redox partner interfacial residues delimited in CYPs (bacterial and eukaryotic) other than the CYP2 family were identified by aligning the structure of CYP2B4 (2BDM.pdb) with the structure of the CYP for which the interfacial residues have been experimentally delimited. In the absence of a structure for that CYP with experimentally delimited interfacial residues, a structure of a member from the same family was used and residue positions were then transferred via pair-wise alignment using ClustalW (Jeanmougin et al., 1998; Thompson et al., 1994). This strategy was used in transferring Shimizu’s interfacial residue delimitations in R. norvegicus CYP1A2 (Mayuzumi et al., 1993; Shimizu et al., 1991) to CYP2B4 via H. sapiens CYP1A2 and also to transfer Cvrk’s CYP1A1 interfacial residues (Cvrk and Strobel, 2001) to CYP2B4 via H. sapiens CYP1A2. In this work, no assignments to CYP2B4 of experimentally delimited interfacial residue positions from sequences outside the CYP2 family were performed without the mediation of structure alignments. The CYP2B4 residues identified as analogous to experimentally delimited CYP–redox partner interfacial residues

A. Zawaira et al. / Computational Biology and Chemistry 34 (2010) 42–52

were viewed in the CYP2B4–POR model using the program Pymol (www.pymol.org). 2.7. Exhaustive pairing of CYP2B4 basic clusters to P450-oxidoreductase acidic clusters A more exhaustive yet simple strategy for the pairing of positively charged surface residues on CYP2B4 with negatively charged surface residues on CYP2B4–POR FMN binding domain that seeks to circumvent the limitations of contemporary tools for side-chain optimisation begins with the analysis of the CYP2B4–POR FMN domain structural model with the software ESBRI. ESBRI is a web server that evaluates salt bridges in protein structures (Costantini et al., 2008). Salt bridges are delimited as interactions between positively charged (i.e. protonated) Histidines, Lysines or Arginines and negatively charged Aspartates or Glutamates. In this work, interactions involving Histidine residues are ignored. Furthermore, salt bridges were searched at 4, 10, 15, 20, 25 and 30 Å distance cut-offs. This was done because the side chains in the 1BVY.pdb-derived structural model for the CYP2B4–P450-oxidoreductase complex are not necessarily expected to be orientated for salt bridge interaction and application of the simple scoring scheme described below to the ESBRI-delimited salt bridges allows us to take into account these sub-optimal orientations of side chains in salt bridge delimitation. The results from the ESBRI runs at the distance cut-offs described above were analyzed to establish the smallest distance cut-off (among them) that allows the inclusion of all previously described charge clusters on the CYP2B4 surface and on the POR surface. The structures of CYP2B4 (2BDM.pdb) and POR FMN binding domain (1AMO.pdb, residues 1–235) were searched for stretchedout Lysine, Arginine, Aspartate and Glutamate residues. The distances (in Angstroms) from the alpha-carbon atom (CA) to the protonated terminal nitrogens (NZ in Lysine and NH1 and NH2 in Arginine) and from the CA to the carboxylic acid side chain oxygen atoms (OD1 and OD2 in Aspartate and OE1 and OE2 in Glutamate) were measured in Pymol and recorded. These values were compared to those obtained from building stretched-out Arg-Arg-Arg, Lys-Lys-Lys, Asp-Asp-Asp and Glu-Glu-Glu tripeptides and were found to be similar. The extended peptides were built using PEPBUILD (http://www.imtech.res.in/bvs/pepbuild/team.html) and distances were measured in Pymol. The set of distances found from searching the CYP2B4 structure (2BDM.pdb) and the POR FMN binding domain structure (1AMO.pdb, residues 1–235) was used in all subsequent calculations. These distances are named as follows: distance from Lysine CA to Lysine NZ is ˛, from Arginine CA to Arginine NH1 is ˇ, from Arginine CA to Arginine NH2 is , from Aspartate CA to Aspartate OD1 is ı, from Aspartate CA to Aspartate OD2 is ε, from Glutamate CA to Glutamate OE1 is  and from Glutamate CA to Glutamate OE2 is . The larger of ˇ and  (i.e. max(ˇ, )) was used to represent the distance from the Arginine CA to the Arginine protonated salt-bridge nitrogen in subsequent calculations. Similarly, the larger of ı and ε (max(ı, ε)) was used in subsequent calculations to represent the distance from Aspartate CA to the salt bridge-forming carboxylate side chain and the larger of  and  (max(, )) was similarly used in Glutamate. The ESBRI-identified pairs were then evaluated by applying a simple test to check if there exists the theoretical possibility of the side chains of the identified residues coming close enough to form a salt bridge. Here, the distance between the CA atoms of the ESBRI-identified amino acid residue pairs are calculated in Pymol. This CA–CA distance is referred-to as ω for a given pair of residues identified by ESBRI to be salt bridge interaction partners. Next, the closest distance of approach between the side chains of the ion-pair residues is calculated by subtracting (from the CA–CA distance, ω): the sum of ˛ and max(ı, ε) for Lysine–Aspartate pairs, the sum of

45

˛ and max(, ) for Lysine–Glutamate pairs, the sum of max(ˇ, ) and max(ı, ε) for Arginine–Aspartate pairs and the sum of max(ˇ, ) and max(, ) for Arginine–Glutamate pairs. Differences strictly below the cut-off of 4 Å are accepted as indicative of salt-bridges involving the associated pair of amino acid residues. This distance cut-off for salt-bridge interactions is the same as the default distance cut-off in ESBRI (Costantini et al., 2008). 2.8. Mathematical scoring of the degree of co-relatedness of pairs of oppositely charged clusters Given a cluster (let us call it X) containing x positively charged residues (Lysines and/or Arginines) and a cluster (let us call it Y) containing y negatively charged residues (Aspartates and/or Glutamates), we can use combinatorics to calculate the maximum number of pair-wise interactions between residues in X and residues in Y. Let us call this number max(X, Y). When we compare max(X, Y) with the number of interactions delimited on the 1BVY.pdb-based structural model using the approaches described above, we can achieve a mathematical description of the degree of co-relatedness of the two clusters. This comparison is achieved by expressing the number of interactions delimited from the 1BVY.pdb structure using the schemes described above as a fraction of max(X, Y). Observe that this score is simply a mathematical score and makes no considerations of the physics of interactions (mechanics, energetics, etc.). From the combinatorics result that the number of distinct pairwise interactions in a group of N objects (e.g. number of distinct handshakes in a group of people) is (Zawaira and Hitchcock, 2009) (N)(N − 1) 2

(Formula 1)

it follows that, max(X, Y ) =

(x + y)(x + y − 1) (x)(x − 1) (y)(y − 1) − − 2 2 2

(1)

That is, apply (Formula 1) on the total collection of objects in X and Y (i.e. x + y) then subtract all the non-permissible interactions between objects in the same set (these are repulsive interactions, the numbers are obtained by applying (Formula 1) on set X (contains x objects) and on set Y (contains y objects)). 2.9. Determination of charged residues in CYP2B4 and P450-oxidoreductase analogous to experimentally delimited salt bridge-forming residue pairs in the CYP450BM-3 system Karyakin et al. (2007) have recently inferred salt bridge interactions that occur upon BMP–FMN complex formation from Fourier transform infrared (FTIR) spectroscopy experiments. The saltbridge interactions delimited by Karyakin and co-workers are as follows: FMN binding domain Arg 498–BMP Asp 242, FMN binding domain Asp 607–BMP Lys 306 and FMN binding domain Asp 542–BMP Lys 59. We are interested in the determination of CYP2B4–POR residue pairs analogous to Karyakin’s pairs and the comparison of these to the ones inferred from the approach of applying ESBRI (Costantini et al., 2008) and our scoring scheme onto the 1BVY.pdb-derived structural model of the CYP2B4–POR complex. The identification of CYP2B4 charged residues (Arginines and Lysines) analogous to basic residues in Karyakin’s pairs was done via structural alignment of the 2BDM.pdb structure of CYP2B4 and the BMP domain in 1BVY.pdb using the program Pymol. The FMN-binding domain of R. norvegicus POR (residues 1–235 in 1AMO.pdb) was aligned to the FMN-binding domain of CYP450BM3 (in 1BVY.pdb) using Pymol to identify POR residues analogous

46

A. Zawaira et al. / Computational Biology and Chemistry 34 (2010) 42–52

Table 1 CYP2B4–P450-oxidoreductase interfacial residues calculated from the 1BVY.pdb-derived structural model at 6 Å distance cut-off. P450-oxidoreductase

CYP2B4 residues

64V, 65K, 84Y, 86S, 87Q, 88T, 89G, 90T, 92E, 93E, 96N, 97R, 100K, 111S, 113D, 115E, 116E, 140Y, 141G, 145P, 146T, 147D, 148N, 150Q, 151D, 175N, 177T, 178Y, 208D, 209D, 210G, 211N, 213E and 214E

85R, 88L, 89V, 90D, 91Q, 92A, 93E, 94A, 96S, 129L, 130A, 132M, 133R, 134D, 136G, M137, E143, L340, D341, R343, A344, K345, P347, D350, H354, T372, K373, K384, F413, D415, A419, L420, K421, R422, N423, E424, F426, M427, F429, S430, L431, G432, K433, I435, C436, L437, E439, G440 and R443

Residues that have been described in the literature as interfacial residues are underlined. The references are given in Table S2 in the supplementary material.

to the acidic residues in Karayakin’s pairs. The alignments were assessed using the criteria described above. 3. Results 3.1. Structural model for the CYP2B4–P450-oxidoreductase FMN-binding domain complex The structural model for the CYP2B4–POR FMN-binding domain complex is shown in Fig. 1. CYP–redox partner residues described in the literature are projected onto this model. This generates a (qualitative) visual representation of the relationship between the interface delimited by adopting the structure by Sevrioukova et al. (1999b) (1BVY.pdb) as the prototypical CYP–redox partner interface and the set of interfacial residues delimited by various experimental approaches described in the literature. As can be seen in Fig. 1, there is in general good concordance between the interfacial residues implied in 1BVY.pdb and literature interfacial residues.

A more quantitative evaluation of the relationship between the 1BVY.pdb-derived CYP–redox partner interface and the interfacial residue in the literature has been achieved by calculating the interface as described in Section 2. The interface defined on the surface of R. norvegicus POR FMN-binding domain (1AMO.pdb, residues 1–235) calculated at the 6 Å distance cut-off (residues listed in the left column of Table 1) shows that 100% of literature interfacial residues (listed in Table S2) map to the interface calculated from our model. In this calculation a literature residue in Table S2 for which the residue before or after (in sequence) is listed in Table 1 is also included as an interfacial residue predicted by our model. Conversely, the results in Table 1 and Table S2 suggest that 29.4% of interfacial residues derived from our model (listed in Table 1) have been described previously in the literature. As before, a calculated interfacial residue in Table 1 that occupies the residue position before or after a residue described in the literature (i.e. listed in Table S2) is also included in the calculation. When we consider the interface defined on the surface of CYP2B4 calculated at the 6 Å distance cut-off (residues listed in the right column of Table 1), and using the same rules as above for inclusion/exclusion of residues in calculations, we see that 55.3% of literature interfacial residues (listed in Table S2) map to the interface calculated from our model (listed in Table 1). Conversely, it may also be seen that 36.7% of interfacial residues derived from our model (listed in Table 1) have been described previously in the literature (listed in Table S2). Table 1 suggests that 14.5% of the 235 residues in the R. norvegicus POR FMN-binding domain are involved in CYP-binding as interfacial residues. It also suggests that 10% of the 491 residues in CYP2B4 are involved in POR-binding as interfacial residues. The list of residues delimited at the 4.5 Å distance cut-off is shown in Table S1a that delimited at distance cut-off 6.5 Å in Table S1b and finally, that delimited at distance cut-off 7 Å in Table S1c.

3.2. Treatment of side-chain conformations and preliminary searches for paired charge clusters

Fig. 1. Visualisation of the relationship between the 1BVY.pdb-derived CYP2B4–P450-oxidoreductase complex and experimentally delimited interfacial residues. CYP2B4 is shown in green and the FMN-binding domain of POR is shown in yellow. CYP–redox partner interfacial residues were sourced from the literature and mapped onto the 1BVY.pdb-derived structure of the CYP2B4–POR complex. The interfacial residues sourced from the literature are listed in Table S2 in the supplementary material. The experimentally delimited CYP2B4 interfacial residues are shown in red and the experimentally delimited POR interfacial residues are shown in blue. The residues on the top far right of the figure marked a (232R), b (251K) and c (253R) are generally far removed from the interface and do not conform to the general location of the rest of the residues in Table S2. The residue 262R is also far removed from the interface and, relative to the shown view, is located at the far back of CYP2B4. The model for the CYP–POR complex was rendered in Pymol (www.pymol.org).

Analysis of the initial FATCAT-derived model for the CYP2B4–POR FMN domain complex (no optimisation of sidechain conformations) using the Protein Interaction Calculator Server (6 Å default distance cut-off for the delimitation of ionic interactions) and the ESBRI Server (4 Å default distance cut-off for the delimitation of ionic interactions) reveals that Aspartate 113 and Glutamate 115 of POR FMN domain form a cluster that makes ionic interactions with Lysine 433 of CYP2B4. Analysis of the model whose side-chain interactions have been optimised using SCWRL identifies Aspartate 113, Glutamate 115 and Glutamate 116 of POR FMN domain as forming a negatively charged cluster of residues that make ionic interactions with Lysine 433 of CYP2B4. Furthermore, the SCWRL optimised model reveals that Glutamate 92 of POR FMN domain forms ionic interactions with Arginine 133 of CYP2B4. While the preliminary analysis of both the models (with and without SCWRL side-chain optimisation) affords the identification of the previously not described POR FMN domain acidic charge clusters Aspartate 113, Glutamate 115 and Glutamate 116 and Glutamate 92, the failure to identify all other

A. Zawaira et al. / Computational Biology and Chemistry 34 (2010) 42–52

47

Table 2 Correlated charge clusters in the 1BVY.pdb-derived structural model for the CYP2B4–P450-oxidoreductase FMN-binding domain complex. Correlated pairs (CYP2B4:POR)

Number of distinct pair-wise interactions made by residues CYP2B4 residues (residue:number)

POR residues (residue:number)

(85R:142E), (85R:144D), (85R:147D), (125R:113D), (133R:92E), (133R:93E), (133R:101D), (343R:207D), (343R:208D), (343R:209D), (343R:213E), (345K:209D), (373K:116E), (421K:179E), (421K:209D), (422R:144D), (422R:207D), (422R:208D), (422R:209D), (433K:92E), (433K:113D), (433K:115E), (433K:116E), (433K:147D), (434R:92E), (434R:113D), (434R:115E), (434R:116E), (434R:147D), (443R:208D) and (443R:213E). Total number of pair-wise interactions = 31

(85R:3), (125R:1), (133R:3), (343R:4), (345K:1), (373K:1), (421K:2), (422R:4), (433K:5), (434R:5) and (443R:2). Total number of pair-wise interactions = 31

(92E:3), (93E:1), (101D:1), (113D:3), (115E:2), (116E:3), (142E:1), (144D:2), (147D:3), (179E:1), (207D:2), (208D:3), (209D:4) and (213E:2). Total number of pair-wise interactions = 31

In the left column are shown ordered (CYP2B4, POR) ion-pairs. Pairs in which at least one member has not been described in the literature (and is therefore a charge-charge interaction residue first identified in this work) are shown in bold. The previously not described residue in the pair is shown underlined. In the right-hand column is given a list of CYP2B4 and POR residues and the number of distinct pair-wise interactions in which they are involved, i.e. given in the format (residue, number of distinct pair-wise interactions involving residue). Pairs in which the associated residue has not been previously described in the literature are given in bold type (these are the bold and underlined residues in the left column). Observe that POR 214E and 215D are not correlated to any positively charged residues.

Table 3 Pairing of previously defined CYP2B4 and P450-oxidoreductase charged clusters. Analysis of individual clusters

Analysis of paired CYP–POR clusters

CYP clusters

POR clusters

Positive cluster–negative cluster combination

Number of contacts as percentage of total

Cluster a: 433K–434R (10/31 = 32.3%)

e

Cluster b–cluster g (421K–422R)–(207D–208D–209D) *(4/6 = 66.7% correlated) Cluster a–cluster h (433K–434R)–(142E–144D–147D) *(2/6 = 33.3% correlated) Cluster a–cluster g (433K–434R)–(207D–208D–209D) *(0/6 = 0% correlated) Cluster d–cluster h (85R)–(142E–144D–147D) *(3/3 = 100% correlated)

(4/31 = 12.9%)

Cluster b: 421K–422R (6/31 = 19.4%)

Cluster c: 133R (3/31 = 9.7%)

Cluster g: 207D–208D–209D (9/31 = 29%)

f

Cluster h: 142E–144D–147D (6/31 = 19.4%)

(2/31 = 6.5%)

(0/31 = 0%)

(3/31 = 9.7%)

g

Cluster i: 213E–214E–215D (2/31 = 6.5%)

Cluster d: 85R (3/31 = 9.7%) Cluster e: 443R (2/31 = 6.5%) Cluster f: 125R (1/31 = 3.2%) Analysis of clusters whose residues have been previously described in the literature is presented in the left column. The number of salt-bridge interactions a given individual charge cluster is involved-in is enumerated from Table 2 by summing the number of interactions of the individual residues in the cluster and is expressed a fraction and percentage of the total (a total of 31 interactions are delimited in Table 2). The right-hand column gives correlations of some of the clusters. The number of interactions accounted-for by a given pair of clusters is calculated by summing all distinct ion-pair interactions (listed in Table 2) made by residues from the two clusters under consideration. The proportion of interactions accounted-for by the given cluster-pair is calculated by expressing the number of interactions it accounts-for as a fraction of the total number of pair-wise interactions. The percentage (and fractional) contribution is reported in the “Number of contacts as a percentage of total” column of the table. The degree of co-relatedness of cluster combinations considered here is shown for each cluster pair in the brackets beginning with asterisk (*). The numerator for this term is enumerated from Table 2 and the denominator is determined as described in Section 2. The purpose of this score is also described in Section 2. The cluster combination in bold type (cluster a–cluster g) emphasizes the non-correlation between two of the biggest negative and positive clusters described in the literature. This is in harmony with our discovery of two large, previously not described oppositely charged clusters that correlate (separately) with cluster a and cluster g. Notice that the cluster h correlates better with cluster d than it does with cluster a. e Designated as cluster 1 by Shen and Kasper (1995) and also investigated by Zhao et al. (1999). f Designated as cluster 3 by Zhao et al. (1999). g Designated as cluster 2 by Shen and Kasper and also investigated by Zhao et al. (1999).

known charge clusters highlights the need for a more powerful and exhaustive search strategy that circumvents the limitations of contemporary side-chain optimisation tools. This point is also emphasized by the relatively marginal differences in prediction outcomes between the SCWRL optimised and the non-optimised models. Side-chain optimization with the NOMAD-Ref Server (Lindahl et al., 2006) gives a more extensive delimitation of ionic interaction pairs (ordered as (FMN domain residue:CYP2B4 residue)): (Glutamate 92:Arginine 133), (Aspartate 113:Lysine 433), (Glutamate 115:Lysine 433), (Aspartate 208:Arginine 343), (Aspartate 208:Arginine 422), (Aspartate 208:Arginine 443), (Aspartate 209:Arginine 343), (Aspartate 209:Arginine 422). While ionic interaction pair-delimitations from the NOMAD-Ref simulations include more basic clusters (including the previously not

described Arginine 343 residue), the achieved coverage of known clusters still highlights the need for a more robust strategy for searching for ionic interaction pairs. 3.3. Exhaustive pairing of CYP2B4 basic clusters to P450-oxidoreductase acidic clusters The values of ˛, ˇ, , ı, ε,  and  (see Section 2 for definition and significance of these values) are reported in Table S3 of the supplementary material. The ESBRI (Costantini et al., 2008) run at the 15 Å distance cut-off includes all the previously described charge clusters on the surface of CYP2B4 and POR. Furthermore, our scoring system shows that no additional ion-pair interactions are included at the higher

48

A. Zawaira et al. / Computational Biology and Chemistry 34 (2010) 42–52

Table 4 Analysis of CYP2B4 and P450-oxidoreductase charge clusters newly identified in this work. Analysis of individual clusters

Analysis of paired CYP–POR clusters

CYP clusters

POR clusters

Positive cluster–negative cluster combination

Number of contacts as percentage of total

Cluster j: 343R–345K (5/31 = 16.1%)

Cluster q: 113D–115E–116E (8/31 = 25.8%)

(6/31 = 19.4%)

Cluster p: 373K (1/31 = 3.2%)

Cluster r: 92E–93E (4/31 = 12.9%)

Cluster a–cluster q (433K–434R)–(113D–115E–116E) *(6/6 = 100% correlated) Cluster j–cluster g (343R–345K)–(207D–208D–209D) *(4/6 = 66.7% correlated) Cluster j–cluster h (343R–345K)–(142E–144D–147D) *(0/2 = 0%) Cluster c–cluster r (133R)–(92E–93E) *(2/2 = 100% correlated)

Cluster s: 101D (1/31 = 3.2%) Cluster t: 179E (1/31 = 3.2%)

(4/31 = 12.9%)

(0/31 = 0%) (2/31 = 6.5%)

Analysis of clusters whose residues have not been previously described in the literature (i.e. identified in this work) is presented in the left column. The number of salt-bridge interactions a given individual charge cluster is involved-in is enumerated from Table 2 by summing the number of interactions of the individual residues in the cluster and is expressed a fraction and percentage of the total (a total of 31 interactions are delimited in Table 2). The right-hand column gives correlations of some of the clusters. The number of interactions accounted-for by a given pair of clusters is calculated by summing all distinct ion-pair interactions (listed in Table 2) made by residues from the two clusters under consideration. The proportion of interactions accounted-for by the given cluster-pair is calculated by expressing the number of interactions it accounts-for as a fraction of the total number of pair-wise interactions. The percentage (and fractional) contribution is reported in the “Number of contacts as a percentage of total” column of the table. The degree of co-relatedness of cluster combinations considered here is shown for each cluster pair in the brackets beginning with asterisk (*). The numerator for this term is enumerated from Table 2 and the denominator is determined as described in Section 2. The purpose of this score is also described in Section 2. We analysed the 1AMO.pdb structure using the program Ligplot (Wallace et al., 1995) and found that none of the residues in cluster q (113D–115E–116E), cluster r (92E–93E), cluster s (101D) or cluster t (179E) make either hydrophobic or hydrogen-bond contacts with the FMN cofactor. However, it is noteworthy that residue 178Y (preceding cluster t) makes hydrophobic contact with FMN and residues 90T and 91A (preceding cluster r) make H-bond contacts with FMN. These findings are shown in Fig. S1 in the supplementary material. The cluster combination in bold type (cluster j–cluster h) emphasizes the non correlation between the largest positively charge cluster we have identified in this work and one of the largest negatively charge cluster on the surface of POR. The analysis of the cluster j–cluster h pair in this table together with the analysis of the cluster a–cluster g pair in Table 3 highlights the presence of large but well-separated charge clusters on CYP and POR surfaces.

distance cut-offs of 20, 25 and 30 Å (data not shown). Ion-pair interactions in this work are therefore delimited from the ESBRI run performed at the 15 Å distance cut-off. As described above, our scoring system is a simple scheme to allow for sub-optimal side-chain orientations in the CYP2B4–POR structural model. Salt bridge pairings that score favourably (score < 4 Å) are reported in Table 2 and in Table S4 of the supplementary material. Calculation of the scores of the ESBRI-delimited interactions is shown in Table S4 of the supplementary material. In the literature, some of the residues in Table 2 have been grouped into clusters. We have used our model to delimit interactions between clusters and we have also used the simple mathematical scoring system described in Section 2 to achieve an indication of the extent of co-relatedness of these charge clusters.

Tables 2 and 3 give results from the first ever attempt to directly match the previously described positive and negative charge clusters. Table 4 gives pairings that involve clusters we have identified in this work (i.e. previously not described in the literature). 3.4. Determination of charged residues in CYP2B4 and P450-oxidoreductase analogous to experimentally delimited salt bridge-forming residue pairs in the CYP450BM-3 system The CYP2B4 and POR residues analogous to CYP450BM-3 heme domain and CYP450BM-3 FMN-binding domain residues in Karyakin’s pairs (Karyakin et al., 2007) were determined by pairwise alignments of CYP2B4 and CYP450BM-3 heme domain and POR FMN domain and CYP450BM-3 FMN domain structures as

Table 5 Comparison of Karyakin‘s salt bridge pairs (Karyakin et al., 2007) with 1BVY.pdb derived charge cluster pairings. Karyakin CYP450BM-3 salt bridge pair (CYP450BM-3 heme domain residue–CYP450BM-3 FMN-binding domain residue)

Analogous salt bridge pair in mammalian CYP–POR system (CYP2B4 residue–Rattus norvegicus POR residue)

(242D–498R)

a

(133R − 92E + 93E)

(306K–607D)

b

(343R + 345K − 207D + 208D + 209D)

(59K–542D)

c

(433K + 434R − 142E + 144D + 147D)

Comments

Notice that the location of charge types is reversed in this instance. The negative charge is located on the CYP450BM-3 heme domain and the positive charge is located on the CYP450BM-3 FMN-binding domain. However, structure-based assignments of charge centres are still possible. This assignment is in complete agreement with the cluster c–cluster r assignment in Table 4. This assignment is in complete agreement with the cluster j–cluster g assignment in Table 4. This assignment is in complete agreement with the cluster a–cluster h assignment in Table 3.

Notice that two of Karyakin’s assignments are in agreement with predictions from Table 4. Table 4 predictions include charge clusters that have not been described explicitly in the literature. Since the transfer of Karyakin’s residue pairs to mammalian CYP–POR is not trivial and has been first described in this work, our work gives the first direct description of the following charge clusters: 92E + 93E and 343R + 345K. a The residue 498R is proximal to both 92E and 93E and it is not obvious which of the two is “more analogous” to it. 92E and 93E are therefore lumped-up as a charge centre and collectively assigned to 498R. b Similarly, the residue 306K is proximal to both 343R and 345K and the residue 607D is proximal to 207D, 208D and 209D. c The residue 59K is proximal to 433K and 434R and is assigned to both using the same logic as above. The residue 542D is proximal to 142E, 144D and 147D.

A. Zawaira et al. / Computational Biology and Chemistry 34 (2010) 42–52

described in Section 2 and the results are shown in Table 5. A summary of the alignment runs in Pymol is given in the supplementary materials (TableS5).

4. Discussion Sevrioukova and co-workers have debated the biological significance of the CYP450BM-3 heme domain–CYP450BM-3 FMNdomain binding complex they determined (1BVY.pdb). While electrophoretic analysis of dissolved crystals revealed that the linker between the heme and the FMN-binding domain was proteolyzed (resulting in an asymmetric unit consisting of two heme domain molecules and only one FMN-binding domain), there are aspects of the structure of the complex that made it plausible it could represent a specific electron-transfer complex. Such aspects include the fact that in the structure, the methyl groups of the FMN cofactor of the FMN-binding domain are proximal to and pointed towards the heme-binding loop of the heme domain of P450BM-3 (Sevrioukova et al., 1999a,b). Starting from the observation that CYP secondary and tertiary structure is in general well conserved across different lineages of life (Hasemann et al., 1995; Sansen et al., 2007; Schoch et al., 2004; Wester et al., 2004; Yano et al., 2004), we have swapped the CYP450BM-3 heme domain with CYP2B4 and the CYP450BM3 FMN-binding domain with POR in 1BVY.pdb to build a model for the structure of the CYP2B4–POR complex. We have collected CYP–redox partner interfacial residues described in the literature. These residues have been determined separately on the CYP surface and on the CYP–redox partner surface. In cases where interfacial residues were determined on CYP–redox partner systems other than the CYP2B4–POR system, we have determined the corresponding residues in the latter system via structure alignments or combination of sequence alignment and structure alignment. Projection of these interfacial residues described in the literature onto our 1BVY.pdb-derived model for the CYP2B4–POR complex reveals that a significant proportion of the literature interfacial residues fall in the protein–protein interface region of our model. This agreement strongly argues for the biological relevance of the 1BVY.pdb structure and models extracted from it. Sevrioukova and co-workers followed-up structure determination work with mutagenesis studies of interfacial residues deduced from their P450BM-3 heme-FMN-binding domain structure. Mutation (to Cys) and subsequent attachment of dansylcystine to the proximal residue positions 104 and 387 in CYP450BM-3 heme domain resulted in a 10–20-fold reduction inter-domain electron transfer for the position 104 mutant in both the reconstituted system (i.e. the heme and FMN domains expressed as separate polypeptides) and the intact CYP450BM-3 protein. The position 387 mutant resulted in a 2-fold decrease in inter-domain electron transfer in the reconstituted system but had no effect in the intact CYP450BM-3 protein. This work supports the notion that the interface deduced from their structure of the P450BM3 heme-FMN-binding domain complex (i.e. the heme domain proximal side) represents a specific electron transfer interaction (Sevrioukova et al., 1999a,b). In our alignment of the CYP450BM3 heme domain and CYP2B4, CYP2B4 residues that are in the 4 Å neighbourhood of CYP450BM-3 heme domain residue position 104 are 129L, 130A, 131T, 132M, 134D, 135F, 136G and 137M. All these residues are listed as interfacial residues in Table 1. Furthermore, CYP2B4 residues that are in the 4 Å neighbourhood of CYP450BM3 heme domain residue position 387 are 81T, 82D, 84I, 85R, 86E, 87A, 88L, 89V, 90D, 424E, 426F, 427M and 431L—all but the first three residues are listed as interfacial residues in Table 1. This highlights the internal consistency of our model and its high accord with experimental data.

49

We went-on to analyze our CYP2B4–POR model by using the ESBRI web-server (Costantini et al., 2008) to search for possible saltbridge interactions and subsequently scored the ESBRI-calculated pairs of interacting negatively charged (Glutamate or Aspartate) and positively charged residues (Lysine or Arginine) by using a distance-based system that strictly rejects ESBRI pairs whose distance of closest approach of side chains is greater than 4 Å and strictly accepts pairs whose distance of closest approach is less than or equal to 4 Å. The parameters we used were chosen so as to maximise the likelihood of the distance of closest approach falling below 4 Å. Hence distance parameters from the alpha carbon to the salt-bridge forming nitrogen (NZ in Lysine and NH1 and NH2 in Arginine) and from the alpha carbon to the salt bridge forming oxygen (OD1 and OD2 in Aspartate and OE1 and OE2 in Glutamate) were extracted from stretched-out side chains. Furthermore, the larger of the distances from CA to OD1 and from CA to OD2, was used as the distance from CA to the salt-bridge forming oxygen in Aspartate. Parameters for Glutamate and Arginine were similarly extracted. It is noteworthy that our method for scoring salt-bridges implicitly makes the assumption that the model represents the closest approach between the interacting CYP2B4 and POR molecules. This assumption may result in misclassifications particularly in the case of interactions that score just above the 4 Å distance cut-off. This includes such interactions as 421K–142E (Entry 40 in Table S2, score = 4.2 Å), 373K–115E (Entry 31 in Table S2, score = 4.25 Å), 343R–215D (Entry 29 in Table S2, score = 4.24) and 345K–207D (Entry 5 in Table S2, score = 5.38). It could be argued that the CYP450BM-3 FMN domain and the CYP450BM-3 heme domain may dock still closer in solution than is shown in the 1BVY.pdb structure—for example, the CA atoms of the two closest residues between the two domains, 105P of the heme domain and 495G of the FMN domain, are 4.73 Å apart. Only slight movement (which could be realised in solution) would be required for interactions such as 421K–142E to score as salt-bridge interactions. In short, there is scope for interactions that just fail to score as salt-bridge interactions in our strict regime to exist in solution. Despite the limitations imposed by the simplistic criteria for acceptance or rejection pair-wise interactions—particularly the underlying assumption that the model represents the terminal closest approach between CYP2B4 and POR prior to electron transfer, our approach has permitted us to create a network of CYP2B4–POR pair-wise interactions that has substantial accord/agreement with experimental observations. Shen and Kapser’s (1995) work shows that mutations in the cluster 213E–214E–215D (cluster i, Table 3) do not affect PORsupported CYP-dependent benzphetamine N-demethylase activity while mutation of 208D in the 207D–208D–209D cluster (cluster g, Table 3) reduces activity by 63%. Our model predicts that, all other parameters fixed, complete disruption of cluster g residues’ pair-wise interactions (9 pairs of interactions) could affect the CYP–POR interaction more profoundly than the complete disruption of cluster i residues’ pair-wise interactions (2 pairs of interactions). Shen and Kapser’s work provides a very close experimental realisation of comparison of the impact of complete disruption of interactions in cluster g and cluster i (D207N–D208N–209D vs E213Q–E214Q–D215N). In this instance predictions may be made from our model with confidence (i.e. is unlikely to be confounded by details of the energetics of the pair-wise interactions such as the individual contributions of ionpairs to the total binding energy) because of the overwhelming majority of interactions in cluster g versus cluster i. When differences in counts of pair-wise interactions between or within clusters are marginal (≤2), such predictions are expected to be harder to build as the contributions of ion-pair types to overall protein stability are also dependent on the environmental context

50

A. Zawaira et al. / Computational Biology and Chemistry 34 (2010) 42–52

within which they occur (Horovitz et al., 1990; Makhatadze et al., 2003; Serrano et al., 1990; Waldburger et al., 1995). For example, our model predicts that mutation of residue 209D (involved in 4 pair-wise interactions) would have a more profound impact than mutation of residue 208D (involved in 3 pair-wise interactions) or 207D (involved in 2 pair-wise interactions) whereas Shen and Kasper (1995) observe that mutation of 208D has a far greater impact on catalytic activity than mutation of 207D or 209D. Zhao and co-workers’ observations of the impact of mutation of cluster g and i residues on POR-supported CYP2D6 dependent codeine O-demethylation activity are in close agreement with Shen and Kapser’s observations described above (Shen and Kasper, 1995; Zhao et al., 1999) and hence relate to our model in the same way as described above for Shen and Kasper’s data. Zhao and coworkers identified cluster h (Table 3) of acidic residues on the FMN-binding domain of POR, i.e. the residues 142E, 144D and 147D. They compared the impact on catalytic activity of mutation of residues in this cluster to mutation of residues in cluster g and cluster i. Our model predicts that, all other parameters fixed, complete disruption of cluster g residues’ pair-wise interactions (9 pairs of interactions) could affect the CYP–POR interaction more profoundly than the complete disruption of cluster h residues’ pairwise interactions (6 pairs of interactions) which, in turn, could affect the CYP–POR interaction more profoundly than complete disruption of cluster i residues’ pair-wise interactions (2 pairs of interactions). While Zhao and co-workers do not build and compare complete disruptions of interactions across clusters, their findings are in general agreement with predictions that can be made from our model. For example, the greatest reductions in activity that can be attained by single mutation in each cluster are ordered in the same way as the above prediction from our model—cluster g gives the largest reduction in activity (70% reduction), cluster h comes second (32%) and finally, cluster i (4%) (Zhao et al., 1999). In this work, we make the first direct identification of the acidic cluster 113D–115E–116E (we call it cluster q (Table 4)). While we predict confidently that complete disruption of interactions in this cluster would have a significant, measurable impact on POR supported CYP-dependent activity, the marginal differences in counts of pair-wise interactions involving cluster q (8 pairs of interactions) and cluster g (9 pairs of interactions) makes it unclear which cluster would give a greater change after complete disruption. Furthermore, it is also difficult to predict if mutation of cluster q would affect POR-supported CYP-dependent activity more than cluster h (6 interactions). Shimizu et al. (1991) have used site-directed mutagenesis in combination with measurement of kinetic and thermodynamic parameters of mutants to probe the role of Lysines and Arginines in the catalytic function of CYP1A2. We have identified (see Table S2) CYP2B4 amino acid residues that correspond to the CYP1A2 residues studied by Shimizu and co-workers. We find that Shimizu’s data of the dissociation constants (Kd ) of POR from the CYP mutants makes a better overall semi-quantitative fit than does Shimizu’s turnover data (turnover towards 7-ethoxycoumarin) with the counts of pair-wise interactions deduced from our analysis. When we perform a simple binary classification of Shimizu’s residues into two sets—those involved in high numbers of pairwise interactions (>3) according to our structure model and those involved in low numbers of pair-wise interactions (≤3) (thresholds chosen to achieve the lowest misclassification), we build the following set (set H) of high interaction number CYP1A2 residues: H (K99, K440, K453, R455). We also build the following set of low interaction number residues (set L): L (K94, K105, K463). We define the set of mutants with low Kd values (set LKd) from Shimizu’s data as mutants whose Kd values are lower than or

equal to 108 nM. On the other hand we define the set of mutants whose Kd values are higher than 108 nM as set HKd. This choice of thresholds is inspired by the classification of mutant types by Shimizu where it is noted that the CYP1A2 sequences (wild type and variants) can be divided into two groups where the Kd values of members of one group are roughly twice in magnitude those of members of the other group of sequences. The set LKd is as follows: LKd (K94, K105, R455, K463) and the set HKd is as follows: HKd (K99, K440, K453). Given that all other parameters are fixed, it would be expected that disruption of residues accounting for larger numbers of pair-wise interactions affects the formation and/or stabilization of CYP–POR complexes more than the disruption of residues accounting for smaller numbers of such pair-wise interactions. Hence under these assumptions, it is expected that set H (K99, K440, K453, R455) coincides with set HKd (K99, K440, K453) and that set L (K94, K105, K463) coincides with set LKd (K94, K105, R455, K463). We see that these sets match closely and the only misclassification is residue R455. We speculate that this misclassification might be because R455 might be the least important residue (making the fewest interactions) in the CYP1A2 contiguous charge cluster 455R–456R–457R. Its disruption would therefore not result in severe alteration of the formation and/or stabilization of CYP1A2–POR complex. We are also able to rationalise several aspects of Shimizu’s data using our counts of pair-wise interactions of the residues. For example, although the mutants at the CYP1A2 residue positions 94K and 105K have the same turnover number towards 7-ethoxycoumarin (0.18 nmol/min/nmol P450), they have hugely different Kd values (103 and <28 nM respectively) (Shimizu et al., 1991). We rationalise this as follows: CYP1A2 residue 94K corresponds to CYP2B4 residue 85R which is involved in 3 pair-wise interactions. On the other hand, residue 105K corresponds to residue 98R which is involved in 0 (zero) pair-wise interactions. Everything else constant, abrogation of residue 94K is therefore expected to have a more profound impact on Kd value (higher Kd value) than abrogation of 105K. We observe that Bridges et al. (1998) data of the dissociation constants (Kd ) of POR from the CYP2B4 mutants makes a good semiquantitative fit with the pair-wise interactions delimited from our model—much like the fit described above for Shimizu’s data over the set of charged residues common to Bridges’s list and our list in Table 2. We achieve a good fit with Bridges and co-workers’ data as follows. We let set H be the collection of positively charged residues (from among the set of residues that appear both in Table 2 of this work and in Bridges et al., 1998) which are involved in high numbers of pair-wise interactions. In this instance a residue involved in ≥3 pair-wise interactions is assigned to H . The set H is therefore as follows H (422R, 433R, 133R). We let set L be the collection of positively charged residue involved in a low (<3) number of pair-wise interactions. The set L is therefore L (421K, 443R). Next we define the set H Kd of residues that have a “high” Kd value. The cut-off we choose here for a “high” Kd value is any value greater than 10 times the Kd value for the dissociation of POR from the wild-type CYP2B4 protein (i.e. >0.2 ␮M—standard errors are discarded). We also define the set L Kd as the set of residues with “low” Kd values (≤0.2 ␮M). The set H Kd is therefore comprises the following H Kd (422R, 443R, 433K, 133R) and the set L Kd comprises L Kd (K421). Following arguments similar to those discussed above when considering Shimizu’s data, it is expected that the set H will overlap with set H Kd and the set L with set L Kd. Comparison of the two pairs of sets reveals that only the residue 443R is misclassified by our assignments—highlighting the great harmony that can be built between our model and experimental data. Karyakin et al. (2007) used Fourier transform infrared (FTIR) spectroscopy to determine residues in the CYP450BM-3 system forming salt-bridges upon complex formation. We have built pairs

A. Zawaira et al. / Computational Biology and Chemistry 34 (2010) 42–52

of residues in the CYP2B4–POR system that are analogous to Karyakin’s pairs via structure alignments and the residue pairs we obtained in this way are shown in Table 5. It can be observed that the pairing of the CYP2B4 cluster 133R to the POR FMN domain cluster 92E + 93E in Table 5 is identical to the cluster c–cluster r pairing in Table 4. Furthermore, the pairing of the CYP2B4 cluster 343R + 345K to the POR FMN domain cluster 207D + 208D + 209D in Table 5 is identical to the cluster j–cluster h pairing in Table 4. Finally, it can be seen that the pairing of the CYP2B4 cluster 433K + 434R to the POR FMN domain 142E + 144D + 147D in Table 5 is identical to the cluster a–cluster h pairing in Table 3. Hence our cluster pairings based on our 1BVY.pdb-derived model for the CYP2B4–POR FMN domain complex agree with all three of Karyakin’s pairs. We make the first direct description of the CYP2B4 clusters 343R + 345K and 373K and also of the POR FMN domain clusters 113D + 115E + 116E, 92E + 93E, 101D and 179E and show how they are paired in interactions of the CYP2B4–POR complex. We also develop a simple method that shows that the existence of cluster 113D + 115E + 116E, cluster 92E + 93E and cluster 343R + 345K can be inferred from Karyakin’s work. These relationships are highlighted for the first time in this work. In light of the good agreement between experimental data and the 1BVY.pdb-derived model of the complex between CYP2B4 and POR FMN domain, the fact that 100% of interfacial residues on POR surface described in the literature appear in our model yet only 29.5% of interfacial residues on the POR surface described in our model are also described in the literature could indicate that our model gives a more complete description of the mammalian CYP–POR interface. This notion is further reinforced by the fact that 55.3% of the interfacial residues on the CYP2B4 surface described in the literature appear in our model while 36.7% of interfacial residues on the CYP2B4 surface described in our model are also described in the literature. Intriguing ion-pair interaction patterns between CYP2B4 and POR are discerned from our model. For example, when one considers the interactions made by clusters that comprise more than one amino acid residue and are involved in at least 5 pair-wise interactions, the interaction pattern shown in Fig. S2 (supplementary material) emerges. This highlights the presence of three large positive clusters – 433K + 434R (cluster a), 421K + 422R (cluster b) and 343R + 345K (cluster j) and three large negative clusters – 113D + 115E + 116E (cluster q), 142E + 144D + 147D (cluster h) and 207D + 208D + 209D (cluster g). Residues in these clusters are involved in 87.1% of the pair-wise interactions delimited in this work (31 ion-pair interactions) and these clusters are therefore likely to dominate electrostatic interactions involved in the initial recognition (mediated by long-range electrostatic interactions) and proper relative orientation (mediated by short-range electrostatic interactions) of CYP2B4 and POR. We observe that among this collection of large clusters, the negative cluster q and the positive cluster j seem to be dedicated to one interaction partner (cluster a and cluster g respectively) while the rest of the clusters make interactions that “stitch” across more than one cluster. In summary, the assumption that the structure of the bacterial CYP450BM-3 heme domain–CYP450BM-3 FMN domain complex by Sevrioukova et al. (1999b) represents a specific electron transfer complex and therefore is a relevant template when considering interactions in mammalian CYP–POR systems, has permitted the construction of mammalian models that are in harmony with findings from previous investigations. These CYP450BM-3-derived models therefore provide a context and platform upon which all prior efforts may be assessed, rationalised and integrated. They also provide a platform upon which further hypotheses regarding this important biomolecular interaction and the ensuing electron transfer may be constructed and subsequently experimentally investigated.

51

Conflict of interest No conflict of interest declared. Acknowledgements A.Z., N.B.-K., and L.C. wish to thank the National Research Foundation (NRF) of South Africa for financial support through bursaries. J.B. wishes to thank the South African Research Chair Initiative (SARChI) for a research chair. M.G. wishes to thank the National Bioinformatics Network (NBN) of South Africa for funding.Author contributions are as follows: A.Z. and J.B. designed research; A.Z., M.G., L.C., N.B.-K., M.K. and P.M. performed research; A.Z., M.G., L.C. and N.B-K. analyzed data; A.Z., M.G., L.C., N.B.-K., M.K. and J.B. wrote the paper. Appendix A. Supplementary data Supplementary data associated with this article can be found, in the online version, at doi:10.1016/j.compbiolchem.2009.10.006. References Bernhardt, R., Kraft, R., Otto, A., Ruckpaul, K., 1988. Electrostatic interactions between cytochrome P-450 LM2 and NADPH-cytochrome P-450 reductase. Biomed. Biochim. Acta 47, 581–592. Bridges, A., Gruenke, L., Chang, Y.T., Vakser, I.A., Loew, G., Waskell, L., 1998. Identification of the binding site on cytochrome P450 2B4 for cytochrome b5 and cytochrome P450 reductase. J. Biol. Chem. 273, 17036–17049. Chan, S.L., Purisima, E.O., 1998. Molecular surface generation using marching tetrahedra. J. Comput. Chem. 19, 1268–1277. Costantini, S., Colonna, G., Facchiano, A.M., 2008. ESBRI: a web server for evaluating salt bridges in proteins. Bioinformation 3, 137–138. Cvrk, T., Strobel, H.W., 2001. Role of LYS271 and LYS279 residues in the interaction of cytochrome P4501A1 with NADPH-cytochrome P450 reductase. Arch. Biochem. Biophys. 385, 290–300. Davydov, D.R., Darovsky, B.V., Dedinsky, I.R., Kanaeva, I.P., Bachmanova, G.I., Blinov, V.M., Archakov, A.I., 1992. Cytochrome C (Fe2+ ) as a competitive inhibitor of NADPH-dependent reduction of cytochrome P450 LM2: locating protein–protein interaction sites in microsomal electron carriers. Arch. Biochem. Biophys. 297, 304–313. Davydov, D.R., Kariakin, A.A., Petushkova, N.A., Peterson, J.A., 2000. Association of cytochromes P450 with their reductases: opposite sign of the electrostatic interactions in P450BM-3 as compared with the microsomal 2B4 system. Biochemistry 39, 6489–6497. Giardina, C.R., Dougherty, E.R., 1988. Morphological Methods in Image and Signal Processing. Prentice-Hall, Inc. Goodman, J.E., O’Rourke, J., Indyk, P., 2004. Handbook of Discrete and Computational Geometry. Chapman & Hall. Hasemann, C.A., Kurumbail, R.G., Boddupalli, S.S., Peterson, J.A., Deisenhofer, J., 1995. Structure and function of cytochromes P450: a comparative analysis of three crystal structures. Structure 3, 41–62. Horovitz, A., Serrano, L., Avron, B., Bycroft, M., Fersht, A.R., 1990. Strength and cooperativity of contributions of surface salt bridges to protein stability. J. Mol. Biol. 216, 1031–1044. Jackson, R.M., Gabb, H.A., Sternberg, M.J., 1998. Rapid refinement of protein interfaces incorporating solvation: application to the docking problem. J. Mol. Biol. 276, 265–285. Jeanmougin, F., Thompson, J.D., Gouy, M., Higgins, D.G., Gibson, T.J., 1998. Multiple sequence alignment with Clustal X. Trends Biochem. Sci. 23, 403–405. Juvonen, R.O., Iwasaki, M., Negishi, M., 1992. Roles of residues 129 and 209 in the alteration by cytochrome b5 of hydroxylase activities in mouse 2A P450S. Biochemistry 31, 11519–11523. Karyakin, A., Motiejunas, D., Wade, R.C., Jung, C., 2007. FTIR studies of the redox partner interaction in cytochrome P450: the Pdx-P450cam couple. Biochim. Biophys. Acta 1770, 420–431. Kelley, R.W., Reed, J.R., Backes, W.L., 2005. Effects of ionic strength on the functional interactions between CYP2B4 and CYP1A2. Biochemistry 44, 2632–2641. Kuznetsov, V.Y., Poulos, T.L., Sevrioukova, I.F., 2006. Putidaredoxin-to-cytochrome P450cam electron transfer: differences between the two reductive steps required for catalysis. Biochemistry 45, 11934–11944. Lee, D.T., Wong, C.K., 1977. Worst-case analysis for region and partial region searches in multidimensional binary search trees and balanced quad trees. Acta Inform. 9, 23–29. Lindahl, E., Azuara, C., Koehl, P., Delarue, M., 2006. NOMAD-Ref: visualization, deformation and refinement of macromolecular structures based on all-atom normal mode analysis. Nucleic Acids Res. 34, W52–56. Lorensen, W.E., Cline, H.E., 1987. Marching cubes: a high resolution 3D surface construction algorithm. Comput. Graph. 21, 163–169.

52

A. Zawaira et al. / Computational Biology and Chemistry 34 (2010) 42–52

Lu, H., Lu, L., Skolnick, J., 2003. Development of unified statistical potentials describing protein–protein interactions. Biophys. J. 84, 1895–1901. Makhatadze, G.I., Loladze, V.V., Ermolenko, D.N., Chen, X., Thomas, S.T., 2003. Contribution of surface salt bridges to protein stability: guidelines for protein engineering. J. Mol. Biol. 327, 1135–1148. Mayuzumi, H., Sambongi, C., Hiroya, K., Shimizu, T., Tateishi, T., Hatano, M., 1993. Effect of mutations of ionic amino acids of cytochrome P450 1A2 on catalytic activities toward 7-ethoxycoumarin and methanol. Biochemistry 32, 5622–5628. Nelson, D.R., Strobel, H.W., 1987. Evolution of cytochrome P-450 proteins. Mol. Biol. Evol. 4, 572–593. Ofran, Y., Rost, B., 2003. Predicted protein–protein interaction sites from local sequence information. FEBS Lett. 544, 236–239. Poulos, T.L., Finzel, B.C., Gunsalus, I.C., Wagner, G.C., Kraut, J., 1985. The 2.6-A crystal structure of Pseudomonas putida cytochrome P-450. J. Biol. Chem. 260, 16122–16130. Poulos, T.L., Finzel, B.C., Howard, A.J., 1987. High-resolution crystal structure of cytochrome P450cam. J. Mol. Biol. 195, 687–700. Pulim, V., Berger, B., Bienkowska, J., 2008. Optimal contact map alignment of protein–protein interfaces. Bioinformatics 24, 2324–2328. Sansen, S., Yano, J.K., Reynald, R.L., Schoch, G.A., Griffin, K.J., Stout, C.D., Johnson, E.F., 2007. Adaptations for the oxidation of polycyclic aromatic hydrocarbons exhibited by the structure of human P450 1A2. J. Biol. Chem. 282, 14348–14355. Schoch, G.A., Yano, J.K., Wester, M.R., Griffin, K.J., Stout, C.D., Johnson, E.F., 2004. Structure of human microsomal cytochrome P450 2C8. Evidence for a peripheral fatty acid binding site. J. Biol. Chem. 279, 9497–9503. Schueler-Furman, O., Wang, C., Baker, D., 2005. Progress in protein–protein docking: atomic resolution predictions in the CAPRI experiment using RosettaDock with an improved treatment of side-chain flexibility. Proteins 60, 187–194. Serrano, L., Horovitz, A., Avron, B., Bycroft, M., Fersht, A.R., 1990. Estimating the contribution of engineered surface electrostatic interactions to protein stability by using double-mutant cycles. Biochemistry 29, 9343–9352. Sevrioukova, I.F., Hazzard, J.T., Tollin, G., Poulos, T.L., 1999a. The FMN to heme electron transfer in cytochrome P450BM-3. Effect of chemical modification of cysteines engineered at the FMN-heme domain interaction site. J. Biol. Chem. 274, 36097–36106. Sevrioukova, I.F., Li, H., Zhang, H., Peterson, J.A., Poulos, T.L., 1999b. Structure of a cytochrome P450-redox partner electron-transfer complex. Proc. Natl. Acad. Sci. U.S.A. 96, 1863–1868. Shen, A.L., Kasper, C.B., 1995. Role of acidic residues in the interaction of NADPHcytochrome P450 oxidoreductase with cytochrome P450 and cytochrome c. J. Biol. Chem. 270, 27475–27480. Shimizu, T., Tateishi, T., Hatano, M., Fujii-Kuriyama, Y., 1991. Probing the role of lysines and arginines in the catalytic function of cytochrome P450d by sitedirected mutagenesis. Interaction with NADPH-cytochrome P450 reductase. J. Biol. Chem. 266, 3372–3375. Stayton, P.S., Poulos, T.L., Sligar, S.G., 1989. Putidaredoxin competitively inhibits cytochrome b5-cytochrome P-450cam association: a proposed molecular model for a cytochrome P-450cam electron-transfer complex. Biochemistry 28, 8201–8205.

Stonehuerner, J., Williams, J.B., Millett, F., 1979. Interaction between cytochrome c and cytochrome b5. Biochemistry 18, 5422–5427. Thompson, J.D., Gibson, T.J., Plewniak, F., Jeanmougin, F., Higgins, D.G., 1997. The CLUSTAL X windows interface: flexible strategies for multiple sequence alignment aided by quality analysis tools. Nucleic Acids Res. 25, 4876–4882. Thompson, J.D., Higgins, D.G., Gibson, T.J., 1994. CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice. Nucleic Acids Res. 22, 4673–4680. Tina, K.G., Bhadra, R., Srinivasan, N., 2007. PIC: protein interactions calculator. Nucleic Acids Res. 35, W473–476. Totrov, M., Abagyan, R., 1996. The contour-buildup algorithm to calculate the analytical molecular surface. J. Struct. Biol. 116, 138–143. Tovchigrechko, A., Vakser, I.A., 2005. Development and testing of an automated approach to protein docking. Proteins 60, 296–301. Tovchigrechko, A., Vakser, I.A., 2006. GRAMM-X public web server for protein–protein docking. Nucleic Acids Res. 34, W310–314. Waldburger, C.D., Schildbach, J.F., Sauer, R.T., 1995. Are buried salt bridges important for protein stability and conformational specificity? Nat. Struct. Biol. 2, 122–128. Wallace, A.C., Laskowski, R.A., Thornton, J.M., 1995. LIGPLOT: a program to generate schematic diagrams of protein–ligand interactions. Protein Eng. 8, 127–134. Wang, C., Schueler-Furman, O., Baker, D., 2005. Improved side-chain modeling for protein–protein docking. Protein Sci. 14, 1328–1339. Wang, M., Roberts, D.L., Paschke, R., Shea, T.M., Masters, B.S., Kim, J.J., 1997. Three-dimensional structure of NADPH-cytochrome P450 reductase: prototype for FMN- and FAD-containing enzymes. Proc. Natl. Acad. Sci. U.S.A. 94, 8411–8416. Wang, Q., Canutescu, A.A., Dunbrack Jr., R.L., 2008. SCWRL and MolIDE: computer programs for side-chain conformation prediction and homology modeling. Nat. Protoc. 3, 1832–1847. Wei, L., Locuson, C.W., Tracy, T.S., 2007. Polymorphic variants of CYP2C9: mechanisms involved in reduced catalytic activity. Mol. Pharmacol. 72, 1280–1288. Wester, M.R., Yano, J.K., Schoch, G.A., Yang, C., Griffin, K.J., Stout, C.D., Johnson, E.F., 2004. The structure of human cytochrome P450 2C9 complexed with flurbiprofen at 2.0-A resolution. J. Biol. Chem. 279, 35630–35637. Yano, J.K., Wester, M.R., Schoch, G.A., Griffin, K.J., Stout, C.D., Johnson, E.F., 2004. The structure of human microsomal cytochrome P450 3A4 determined by X-ray crystallography to 2.05-A resolution. J. Biol. Chem. 279, 38091–38094. Ye, Y., Godzik, A., 2004. FATCAT: a web server for flexible structure comparison and structure similarity searching. Nucleic Acids Res. 32, W582–585. Zawaira, A., Hitchcock, G., 2009. A Primer for Mathematics Competitions. Oxford University Press, Oxford, UK. Zawaira, A., Matimba, A., Masimirembwa, C., 2008. Prediction of sites under adaptive evolution in cytochrome P450 sequences and their relationship to substrate recognition sites. Pharmacogenet. Genomics 18, 467–476. Zhao, Q., Modi, S., Smith, G., Paine, M., McDonagh, P.D., Wolf, C.R., Tew, D., Lian, L.Y., Roberts, G.C., Driessen, H.P., 1999. Crystal structure of the FMN-binding domain of human cytochrome P450 reductase at 1.93 A resolution. Protein Sci. 8, 298–306.