2013). MD simulations and MMGBSA or the corresponding Poisson-Boltzmann approach MMPBSA have also been previously used for CYP ligand binding mode and CYP selectivity prediction (Bello et al. 2014, Sato et al. 2017, Juvonen et al. 2020).
Energetics, accessibility and stability metrics using MMGBSA, atomic distances to the heme iron and their STDs, are suggested for the analysis and prediction of CYP ligand binding modes and substrate SOMs. The analysis of binding pose stability and substrate accessibility were crucial additions to the MMGBSA-based binding mode and SOM prediction workflow. This combination of metrics was also applied with good performance for analyzing the isoform selectivity of 7-hydroxylation of novel profluorescent 3-phenylcoumarin-based tool molecules by CYP1 enzymes. Both the stability and accessibility were depicted by the atomic distances of the ligand to the heme iron.
Another approach for stability analysis would be clustering the binding poses from multiple MD simulation trajectories (Sato et al. 2017). This would also be a more automated method for the binding pose comparison. Ligand binding poses can diverge drastically from the starting poses during the simulations (Bello et al.
2014). This was also the case in the simulations here, and thus a partly manual analysis was required for the correct binding mode and SOM prediction.
Accordingly, binding mode clustering would likely enhance the predictions with the suggested metrics. Incorporating the MD-derived metrics with a ligand-based reactivity descriptor would likely further enhance the methodology (Sato et al. 2017) as the intrinsic properties of the ligand also affect SOM selectivity.
Novel profluorescent tool compounds were identified for CYP enzymes, and computational approaches provided insights to the isoform-selectivity of the fluorescence-producing reaction. The most potent CYPs, the CYP1 family and CYP2A6, for catalyzing a specific fluorescence-producing 7-hydroxylation reaction of 3-phenylcoumarin derivatives were identified with a simple comparison of CYP binding sites. All 21 compounds were indeed 7-hydroxylated by at least one of the CYP1 family members and the smallest compounds by CYP2A6, although other isoforms catalyzed the reaction as well. Rigid NIB docking provided clues for the ligand-CYP interactions and binding modes that would facilitate the 7-hydroxylation reaction in CYP isoforms 1A1, 1A2, 1B1, 2A6, 2D6 and 2C19. In MD simulations of 3-phenylcoumarins in complex with the CYP1 family enzymes, the binding modes differed from the initial docking-based hypotheses as the enzymes could relax and water molecules could emerge to the binding sites to mediate ligand-enzyme interactions. The difference between the docking and MD simulations results supports the concept that considering CYP flexibility in ligand binding mode analysis can greatly affect the obtained results (Hritz et al. 2008, Sheng et al. 2014).
In the CYP1 family, water interactions and access channels to the binding site are suggested to have an important role in substrate and SOM selectivity. In the MD simulations here, water molecules mediated ligand-CYP interactions in all CYP1 forms and the access channel compositions differed between the isoforms. Crystal structures of the CYP1 forms all include water molecules inside or at the immediate proximity of the binding site, and a water molecule has been
54
reported to mediate H-bond between CYP1A2 and the co-crystallized ligand (Sansen et al. 2007, Wang et al. 2011, Walsh et al. 2013). In previous MD simulations, water molecules have been shown to readily exchange between the binding site and bulk solvent through access channels across xenobiotic-metabolizing CYP enzymes of varying flexibility (Rydberg et al. 2007, Hendrychova et al. 2012). The binding site water networks of CYP1A2 have also been identified to differ from one ligand to another (Watanabe et al. 2017). Here, the access channel composition varied between the CYP1 enzymes, which also affected the water networks. Access channel composition has been suggested to affect CYP substrate recognition by many prior studies (Urban et al. 2018). A more detailed analysis of the access channel composition and water networks in the CYP1 family with a series of substrates would be beneficial. On one hand, it could provide more information about the effect of the access channel and water network composition to substrate and SOM selectivity. On the other hand, the bound substrate can also have an impact on the access channels and water networks. Such investigation could be useful in shedding light to the structural mechanisms of CYP mediated metabolism and thus aid in future predictions of CYP mediated metabolism. The variation and importance of the hydration sites in the binding site could be especially useful for binding mode and SOM predictions.
The workflows of molecular modelling, docking, and MD simulations applied to CYP enzymes were fairly simple, but the results offer ideas for utilizing the methods in computational cytochrome P450 mediated metabolism.
Compared to docking, MD simulations can provide additional information about the ligand binding pose stability, and insights to the protein flexibility and solvent interactions. The MD simulations of CYP2A6 and 2A5 in complex with small coumarins, and the comparison of 3-phenylcoumarin derivative binding modes in CYP1 enzymes demonstrated how small differences in the highly similar substrates CYP binding sites could affect the catalytic activity and both substrate and SOM selectivity of these enzymes. The MD-based approaches and metrics for CYP mediated metabolism prediction here and in previous studies (Bello et al. 2014, Panneerselvam et al. 2015, Sato et al. 2017, Watanabe et al. 2017, Juvonen et al. 2020) could be applied on a larger set of known CYP isoforms and ligands. It would be beneficial to evaluate the method with a diverse data set, and even an automated workflow could be built. For now, due to their computational demands, the MD/MMGBSA or MMPBSA studies to predict or evaluate CYP selectivity, binding modes and SOMs have been usually performed on a few example cases (Bello et al. 2014, Panneerselvam et al. 2015, Sato et al.
2017, Watanabe et al. 2017, Juvonen et al. 2020). In addition to metabolism prediction, MD simulations have been widely applied to analyze the uniform and differing characteristics of CYP enzymes such as flexibility, channel composition and behaviour in an aqueous solvent or lipid membranes (Rydberg et al. 2007, Hendrychová et al. 2011, Hendrychova et al. 2012, Berka et al. 2013, Jeřábek et al.
2016).
Depending on the available computational resources and the number of hit or lead compounds, the hit-to-lead or lead optimization phase of drug discovery
(Hughes et al. 2011) may be most suitable for the utilization of MD simulations in the optimization of CYP mediated metabolism of the compounds in development. At these stages, the number of compounds have already been filtered with other methods and include only the most potent compound series (Hughes et al. 2011). MD simulations could be especially useful in optimization cases where a drug lead compound is desired to be metabolized at a specific site, or the formation of a certain metabolite is not wanted. As is the case with other computational methods for the prediction of CYP metabolism (Brändén et al.
2014), the application of MD simulations may not be practical in cases where CYP inhibition is the main issue. This is because the mechanisms of inhibition are diverse, including competitive inhibition, binding to an allosteric site, or mechanism-based inactivation (Brändén et al. 2014, Raunio et al. 2015). Issues related to CYP induction should also be assessed with other methods because the induction pathways are not related to ligand-CYP interactions but to other proteins that act as transcription factors (Crivori and Poggesi 2006).
The Panther NIB methodology for rigid docking, VS and rescoring of molecular docking results is well-established and has demonstrated high early enrichment rates in VS applications (Niinivehmas et al. 2015, Kurkinen et al. 2018, 2019, Jokinen et al. 2019). Here, the different aspects of the method protocol were discussed and a practical workflow has been provided to be used by both experts and new users in the field of NIB VS. The provided workflow in the original publication is excellent for demonstrative purposes as it utilizes free or academically free software, although the utilized ligand preparation workflow did not produce the best results amongst the compared software for the step. In the study, it was demonstrated that the VS success is affected by the model creation settings the used protein 3D structure. The effects of the protein conformation to NIB VS have also been demonstrated in studies utilizing protein structures from MD simulations (Virtanen and Pentikäinen 2010) and is a common phenomenon across different docking approaches. The performance of ligand preparation software in reproducing experimental conformers and probing the conformational space have been studied previously (Ebejer et al.
2012, Friedrich et al. 2017). Here, the ligand 3D conformer generation was found to greatly affect the NIB VS performance, and special attention should be paid for the chosen protocol. Considering the excellent performance and speed of NIB docking and rescoring in VS (Niinivehmas et al. 2015, Kurkinen et al. 2018, 2019, Jokinen et al. 2019) and the detailed discussion and workflow provided for the method, the NIB Panther methodology is a great option for structure-based VS.
Thanks to its speed, the NIB Panther VS rigid docking may be used for initial VS of large molecular databases, or as a rescoring method in combination with flexible docking.
Computational prediction of CYP mediated metabolism and VS have wide attention in drug and chemical development. The utilized computational tools for these objectives are diverse, including fast ligand-based methods and protein structure-based methods that, in turn, provide mechanistic insights into ligand-protein interactions. Here, structure-based methods were used for computational prediction and analysis CYP metabolism and VS. MD simulations were used for CYP ligand binding mode and SOM prediction, and to shed light on the isoform-selectivity of novel inhibitors, substrates and profluorescent tool molecules.
Compared to docking, MD simulations provided an important view to the enzyme flexibility, ligand stability and water interactions in ligand-CYP complexes. Metrics of binding free energy, substrate accessibility and binding pose stability were found to be useful in the scope of binding mode and SOM prediction of small CYP2A6 and CYP2A5 ligands. They also showed good performance in the evaluation of reaction isoform-selectivity of 3-phenylcoumarin-based ligands in the CYP1 family. In the following studies, the metrics could be incorporated into a more straightforward binding mode and SOM prediction protocol with, for example, ligand-based reactivity and binding mode clustering. The MD simulations suggested that water interactions and channel composition have an impact on substrate and SOM selectivity in the CYP1 family. In the scope of VS, the priorly developed highly fast structure-based NIB Panther VS methodology was explored. Crucial notes were found and discussed regarding the effects of the used protein structure, options of NIB model generation and ligand 3D conformer generation. A detailed workflow for the NIB method has been provided to be utilized by new users to the protocol. The reviewed methods and protocols are potent tools to be used for the computational prediction and analysis of CYP mediated metabolism and VS in future applications.