Figure 1: Workflow of the docking study with the 9 successive stages (cf. 9 boxes). Swiss PDB viewer, VEGA ZZ and FlexX under Sybyl were used at stages 3, 4 and 5. At stages 6 and 7 a combination of Autodock 4 and docking program Yeti 8 with Biograf R were used in a more systematic way. Abbr. “Ms” for metabolite(s)
Figure 2: The four options (a through d) of possible molecular geometries and the atom types required to formulate the input for docking. While options a and b have an intact aromatic ring with sp2 hybridization and aromatic bond types, options c and d operate with a carbon atom with sp3 hybridization in tetrahedral bond geometry. Option d circumvents the necessity of replacing a (undefined) coordination with (parameterized) covalent bonding, and proofed superior results during docking
Figure 3: The substrate (S) with its aromatic ring approaches the catalytic site with the heme group to form a tetrahedral transition state (2). A ternary coordination complex is formed and embraces the heme iron IV and the activated oxygen moiety which in turn is bound to an sp3 hybrid carbon atom (step 2). Two reaction pathways are possible (steps 3a and 3b) to yield the hydroxylated metabolite [39]. The mechanism is known as “NIH shift”, i.e. a hydrogen atom is replaced by deuterium (D) to reflect the atom shift. The docked ligand poses adopted a sp3 hybrid geometry of the attacked carbon atom in order to simulate the reaction transition state into the product (step 4)
Figure 4: Geometrical comparison between final poses of 7-hydroxylated PPC (left) and the observed case of hydroxylated camphor complexes (PDB codes: 3CPP, 1NOO) [59,42]. Color code: O red, N blue, C grey, Fe orange, S yellow, H white [43]. The views show each ligand in the center and the heme group leftmost on top of the ligands
Figure 5: Geometrical comparison between final poses of 2’-hydroxylated PPC and the observed case of hydroxylated flurbiprofen (white carbon atoms) complexes (PDB codes: 1R9O) [38]. Color code: O red, N blue, C grey, F light blue, Fe orange, S yellow, H omitted for better viewing. The two heavy atoms which make the distinction between warfarin and phenprocoumon are colored in magenta (O=) and green (-CH3). The oxygen species attached to heme iron is presented as a sphere with full atom radius. The proximal cysteine (Cys435) is displayed, too (leftmost yellow stick). The observed distances: Fe - hydroxyl O (2.4 Å), Iron Fe-C atom of SoM (4.9 Å), hydroxyl O–C atom of SoM (3.0 Å) [101]. The required coordination distances had the same lengths as found in the crystal structures: Fe-O between 1.6 and 2,4 Å or O—C2’ between 3.0 and 3.1 Å, respectively [38,42]
Figure 6: Minima (MS-n) and transition state structures (TS-n) located on the potential surface for the catalytic reaction (Figure 7)
Figure 7: Potential energy surface profile for the reaction. Structures of minima and transition states are indicated in accordance with Figure 6

Year

Description and key words

References

2006

Overcoming the receptor's rigidity of its backbone by assessing its structural flexibility during molecular dynamics (MD) studies with explicit water. But MD with explicit solvent simulation has two mayor flaws:
(i) the limited computational resources and
(ii) the addition of even tiniest calculation errors (oversimplifications by force fields, parameters and electronic and entropic effects) into relevant model artifacts during large time periods of simulations during tens of nanoseconds.

[21]

2007

Combining softdocking from receptor ensembles and MD creates a new approach called Fleksy.

[22]

2007

To circumvent the missing conformational flexibility during rigid body docking the so-called soft-docking approach can be applied. It reduces the impact of the repulsion term in the force field equations between ligand and binding site.

[22]

2008

The application of grid-based docking approaches was reviewed in a seminal work.

[23]

2009

Defining a few binding-relevent amino acids at the binding site to simulate conformational changes of their side chains. For example Autodock 4.2.

[17]

2010

Use of Linear Intereaction Energie (LIE) with Merck Molecular Force Field and ligand-based virtual screening fot identification of CYP1A2 ligands by structure-based

[24]

2011

The proposition of an in-silico method that predict the SoMs as “metabolophoros”

[15]

2013

Docking with GOLD was used to describe the impact of CYP1 inhibitors on regioselectivity for hydroxylation reactions.

[25]

2015

Grid docking by Glide software was combined with 3D QSAR techniques to elucidate the CYP2C19 inhibition activity of certain alkaloids.

[26]

2015

A special case is the simulation of covalently bound ligands to form receptor complexes.

[27]

2017

Using many different crystal structures of the same protein (related PDB entries) for docking to describe protein flexibility. This is known as docking by multiple receptor structures, i.e. a receptor ensemble is used as target to reflect the protein flexibility, e.g.
Glide ensemble docking.

[28]

Table 1: Listing of published concepts for ligand - receptor docking and examples applied to CYP biotransformation during the last decade

PPC2´

Total E
[Kcal / mol]

Electro-static IA

VdW IA

H-bonding

Metal IA

Total E minus Mtl IA

a

-9

-1

-8

-0.1

0

-9

c

-9.5

-1.5

-8

-0.1

0

-9.5

d

-20

-0.01

-7

0

-13

-7

b

-19

-4

-1

-0.1

-13

-6

PPC4’

 

 

 

 

 

 

a

-7

0.2

-7

0

0

-7

c

-7

-0.1

-7

0

0

-7

d

-30

-7

-6

0

-16

-14

b

-29

-8

-6

-0.3

-15

-14

PPC6

 

 

 

 

 

 

a

-3

1

-3

-0.002

0

-3

c

-3

-0.1

-3

-0.002

0

-3

d

-24

-5

-4

-0.004

-16

-8

b

-41

-18

-6

-2

-27

-14

PPC7

 

 

 

 

 

 

a

-0.7

1

-1

-0.003

0

-0.7

c

-7

-0.002

-7

-0.001

0

-7

d

-37

-10

-2

0

-26

-11

b

-34

-14

6

-0.01

-26

-8

PPC10

 

 

 

 

 

 

a and c

-21

0

-21

0

0

-21

d

-24

0

-24

-0.005

0

-24

b

-24

0

-24

-0.004

0

-24

FLP

 

 

 

 

 

 

a

-31

-18

-9

-5

0

-31

c

-31

-18

-9

-5

0

-31

d

-48

-25

-7

-3

-14

-34

b

-41

-24

-3

-2

-12

-29

Table 2: Listing of docking solutions at the catalytic site. 1st column refers to the four options (a through d) of possible molecular geometries (Figure 2). The hydroxylation site was accessed through the oxygen-iron (III)-heme-Cys435 octaeder. In the case of circle 3 as the missing structural datum, all 4 remaining possibilities (not hydroxylated sites) were evaluated (9R10S, 9R10S, 9S10R, 9S10S). “No” implies necessary spatial changes which reflects induced fit, e.g. side chain rotations on Val113, Leu365 and Thr301 made the 2’ position accessible to hydroxylation [49]

 

Fe(II)P

 

Triplet

Singlet

Quintent

Erelat

0.00

12.88

17.14

Fe-Np

2.01

1.99

1.97

Fe(II)

0.738

0.737

0.971

Np

-0.655

-0.621

-0.710

 

Fe(III)P

 

Quartet

Doblet

Sextet

Erelat

321.46

326.63

333.62

Fe-Np

1.98

1.97

1.98

Fe(III)

1.028

1.026

1.013

Np

-0.703

-0.706

-0.701

 

Fe(IV)P

 

Quintent

Singlet

Triplet

Erelat

721.87

737.29

737.75

Fe-Np

1.96

1.99

1.99

Fe(IV)

1.065

0.857

1.064

Np

-0.696

-0.615

-0.695

Table 3: Relative energies respect to the ground state (Erelat) in kcal/mol optimized geometrical parameters (Fe-Np) in Å and Mulliken atomic charges of Fe and NP atoms corresponding to the lowest-lying spin state of Fe(II)P, Fe(III)P and Fe(IV)P at B3LYP/6-31G level plus the Land2DZ basis set was used for the Fe atom (only 10 inner core electrons were described with Land2DZ pseudopotential)

Structures

ΔG

ΔG#

n

MS-5

-157.7336

 

 

TS-3

 

-157.7042

337.15 i

MS-6

-157.7080

 

 

TS-7b

 

-157.6596

2097.50 i

TS-7a

 

-157.6781

211.94 i

MS-7a

-157.7206

 

 

MS-7b

-157.7483

 

 

Table 4: Thermodynamic properties (ΔG and ΔG#) in kcal/mol and imaginary frequencies (n) in cm-1 of the structures in the reaction path calculated at 298.15 K

Structures

Fe-Np

Fe-O

NpFeO

TS-3

2.03

1.92

105.0

MS-6

2.01

1.77

103.6

TS-7a

2.02

1.97

87.7

TS-7b

2.00

4.18

97.8

MS-7a

2.01

2.03

88.8

MS-7b

2.01

3.54

74.5

Table 5: Selected geometrical parameters of the minima and transition states optimized with B3LYP/6-31G plus the LanlDZ (bond lengths in Å and angles in degrees)