What are the minimal folding seeds in proteins? Experimental and theoretical assessment of secondary structure propensities of small peptide fragments.

Certain peptide sequences, some of them as short as amino acid triplets, are significantly overpopulated in specific secondary structure motifs in folded protein structures. For example, 74% of the EAM triplet is found in α-helices, and only 3% occurs in the extended parts of proteins (typically ß-sheets). In contrast, other triplets (such as VIV and IYI) appear almost exclusively in extended parts (79% and 69%, respectively). In order to determine whether such preferences are structurally encoded in a particular peptide fragment or appear only at the level of a complex protein structure, NMR, VCD, and ECD experiments were carried out on selected tripeptides: EAM (denoted as pro-'α-helical' in proteins), KAM(α), ALA(α), DIC(α), EKF(α), IYI(pro-ß-sheet or more generally, pro-extended), and VIV(ß), and the reference α-helical CATWEAMEKCK undecapeptide. The experimental data were in very good agreement with extensive quantum mechanical conformational sampling. Altogether, we clearly showed that the pro-helical vs. pro-extended propensities start to emerge already at the level of tripeptides and can be fully developed at longer sequences. We postulate that certain short peptide sequences can be considered minimal "folding seeds". Admittedly, the inherent secondary structure propensity can be overruled by the large intramolecular interaction energies within the folded and compact protein structures. Still, the correlation of experimental and computational data presented herein suggests that the secondary structure propensity should be considered as one of the key factors that may lead to understanding the underlying physico-chemical principles of protein structure and folding from the first principles.

Exhaustive Mapping of the Conformational Space of Natural Dipeptides by the DFT-D3//COSMO-RS Method.

We extensively mapped energy landscapes and conformations of 22 (including three His protonation states) proteinogenic α-amino acids in trans configuration and the corresponding 484 (222) dipeptides. To mimic the environment in a protein chain, the N- and C-termini of the studied systems were capped with acetyl and N-methylamide groups, respectively. We systematically varied the main chain dihedral angles (ϕ, ψ) by 40° steps and all side chain angles by 90° or 120° steps. We optimized the molecular geometries with the GFN2-xTB semiempirical (SQM) method and performed single point density functional theory calculations at the BP86-D3/DGauss-DZVP//COSMO-RS level in water, 1-octanol, N,N-dimethylformamide, and n-hexane. For each restrained (nonequilibrium) structure, we also calculated energy gradients (in water) and natural atomic charges. The exhaustive and unprecedented QM-based sampling enabled us to construct Ramachandran plots of quantum mechanical (QM(BP86-D3)//COSMO-RS) energies calculated on SQM structures, for all 506 (484 dipeptides and 22 amino acids) studied systems. We showed how the character of an amino acid side chain influences the conformational space of single amino acids and dipeptides. With clustering techniques, we were able to identify unique minima of amino acids and dipeptides (i.e., minima on the GFN2-xTB potential energy surfaces) and analyze the distribution of their BP86-D3//COSMO-RS conformational energies in all four solvents. We also derived an empirical formula for the number of unique minima based on the overall number of rotatable bonds within each peptide. The final peptide conformer data set (PeptideCs) comprises over 400 million structures, all of them annotated with QM(BP86-D3)//COSMO-RS energies. Thanks to its completeness and unbiased nature, the PeptideCs can serve, inter alia, as a data set for the validation of new methods for predicting the energy landscapes of protein structures. This data set may also prove to be useful in the development and reparameterization of biomolecular force fields. The data set is deposited at Figshare (10.25452/figshare.plus.19607172) and can be accessed using a simple web interface at http://peptidecs.uochb.cas.cz.

Understanding desaturation/hydroxylation activity of castor stearoyl Δ9-Desaturase through rational mutagenesis.

A recently proposed reaction mechanism of soluble Δ9 desaturase (Δ9D) allowed us to identify auxiliary residues His203, Asp101, Thr206 and Cys222 localized near the di-iron active site that are supposedly involved in the proton transfer (PT) to and from the active site. The PT, along with the electron transfer (ET), seems to be crucial for efficient desaturation. Thus, perturbing the major PT chains is expected to impair the native reaction and (potentially) amplify minor reaction channels, such as the substrate hydroxylation. To verify this hypothesis, we mutated the four residues mentioned above into their counterparts present in a soluble methane monooxygenase (sMMO), and determined the reaction products of mutants. We found that the mutations significantly promote residual monohydroxylation activities on stearoyl-CoA, often at the expense of native desaturation activity. The favored hydroxylation positions are C9, followed by C10 and C11. Reactions with unsaturated substrate, oleoyl-CoA, yield erythro-9,10-diol, cis-9,10-epoxide and a mixture of allylic alcohols. Additionally, using 9- and 11-hydroxystearoyl-CoA, we showed that the desaturation reaction can proceed only with the hydroxyl group at position C11, whereas the hydroxylation reaction is possible in both cases, i.e. with hydroxyl at position C9 or C11. Despite the fact that the overall outcome of hydroxylation is rather modest and that it is mostly the desaturation/hydroxylation ratio that is affected, our results broaden understanding of the origin of chemo- and stereoselectivity of the Δ9D and provide further insight into the catalytic action of the NHFe2 enzymes.

Evidence for H-bonding interactions to the µ-η2:η2-peroxide of oxy-tyrosinase that activate its coupled binuclear copper site.

The factors that control the diverse reactivity of the µ-η2:η2-peroxo dicopper(II) oxy-intermediates in the coupled binuclear copper proteins remain elusive. Here, spectroscopic and computational methods reveal H-bonding interactions between active-site waters and the µ-η2:η2-peroxide of oxy-tyrosinase, and define their effects on the Cu(II)2O2 electronic structure and O2 activation.

Predicting Effects of Site-Directed Mutagenesis on Enzyme Kinetics by QM/MM and QM Calculations: A Case of Glutamate Carboxypeptidase II.

Quantum and molecular mechanics (QM/MM) and QM-only (cluster model) modeling techniques represent the two workhorses in mechanistic understanding of enzyme catalysis. One of the stringent tests for QM/MM and/or QM approaches is to provide quantitative answers to real-world biochemical questions, such as the effect of single-point mutations on enzyme kinetics. This translates into predicting the relative activation energies to 1-2 kcal·mol-1 accuracy; such predictions can be used for the rational design of novel enzyme variants with desired/improved characteristics. Herein, we employ glutamate carboxypeptidase II (GCPII), a dizinc metallopeptidase, also known as the prostate specific membrane antigen, as a model system. The structure and activity of this major cancer antigen have been thoroughly studied, both experimentally and computationally, which makes it an ideal model system for method development. Its reaction mechanism is quite well understood: the reaction coordinate comprises a "tetrahedral intermediate" and two transition states and experimental activation Gibbs free energy of ∼17.5 kcal·mol-1 can be inferred for the known kcat ≈ 1 s-1. We correlate experimental kinetic data (including the E424H variant, newly characterized in this work) for various GCPII mutants (kcat = 8.6 × 10-5 s-1 to 2.7 s-1) with the energy profiles calculated by QM/MM and QM-only (cluster model) approaches. We show that the near-quantitative agreement between the experimental values and the calculated activation energies (ΔH⧧) can be obtained and recommend the combination of the two protocols: QM/MM optimized structures and cluster model (QM) energetics. The trend in relative activation energies is mostly independent of the QM method (DFT functional) used. Last but not least, a satisfactory correlation between experimental and theoretical data allows us to provide qualitative and fairly simple explanations of the observed kinetic effects which are thus based on a rigorous footing.

Mapping Conformational Space of All 8000 Tripeptides by Quantum Chemical Methods: What Strain Is Affordable within Folded Protein Chains?

To gain more insight into the physicochemical aspects of a protein structure from the first principles, conformational space of all 8000 "capped" tripeptides (i.e., N-Ac-X1X2X3-NH-CH3, where Xi is one of the 20 natural amino acids) was investigated computationally. An enormous dataset (denoted P-CONF_1.6M and containing close to 1 600 000 conformers in total) has been obtained by employing a composite protocol combining density functional theory, semiempirical quantum mechanics (SQM), and state-of-the-art solvation methods with 1000 K molecular dynamics (MD) used to generate initial structures (200 snapshots for each tripeptide). This allowed us to present the first rigorous QM-based glimpse at the vast conformational space spanned by small protein fragments. The same computational procedure was repeated for tripeptide fragments taken from the SCOPe database of three-dimensional protein folds, by restraining them to their geometry in a protein. Such complementary data allowed us to compare the distribution of conformational strain energies of unrestrained tripeptidic fragments "in solvent" with those in existing protein chains. Besides providing a rigorous (ab initio) proof of a few well-known concepts and hypotheses concerning protein structures, such as the distribution of (φ, ψ) angles in Ramachandran plots, we have made several observations that came as a certain surprise: (1) distribution of conformational energies does not significantly differ between the "unbiased/unrestrained" conformers obtained from MD sampling in solvent and the biased conformers, i.e., those of a given tripeptide obtained from protein structures; (2) conformational (strain) energy window up to ∼20 to 25 kcal·mol-1 is readily available to tripeptide fragments within the context of a protein chain; (3) overpopulation in certain regions of Ramachandran plot was observed for the unbiased conformers. Last but not least, the massive dataset of accurate (DFT-D3//COSMO-RS) conformational (free) energies of ∼1.6 M peptide conformers, P-CONF_1.6M, obtained throughout this work may serve as excellent dataset for calibrating and benchmarking of popular force fields.

Dinucleoside Polyphosphates as RNA Building Blocks with Pairing Ability in Transcription Initiation.

Dinucleoside polyphosphates (NpnNs) were discovered 50 years ago in all cells. They are often called alarmones, even though the molecular target of the alarm has not yet been identified. Recently, we showed that they serve as noncanonical initiating nucleotides (NCINs) and fulfill the role of 5' RNA caps in Escherichia coli. Here, we present molecular insight into their ability to be used as NCINs by T7 RNA polymerase in the initiation phase of transcription. In general, we observed NpnNs to be equally good substrates as canonical nucleotides for T7 RNA polymerase. Surprisingly, the incorporation of ApnGs boosts the production of RNA 10-fold. This behavior is due to the pairing ability of both purine moieties with the -1 and +1 positions of the antisense DNA strand. Molecular dynamic simulations revealed noncanonical pairing of adenosine with the thymine of the DNA.

Proton-Electron Transfer to the Active Site Is Essential for the Reaction Mechanism of Soluble Δ9-Desaturase.

A full understanding of the catalytic action of non-heme iron (NHFe) and non-heme diiron (NHFe2) enzymes is still beyond the grasp of contemporary computational and experimental techniques. Many of these enzymes exhibit fascinating chemo-, regio-, and stereoselectivity, in spite of employing highly reactive intermediates which are necessary for activations of most stable chemical bonds. Herein, we study in detail one intriguing representative of the NHFe2 family of enzymes: soluble Δ9 desaturase (Δ9D), which desaturates rather than performing the thermodynamically favorable hydroxylation of substrate. Its catalytic mechanism has been explored in great detail by using QM(DFT)/MM and multireference wave function methods. Starting from the spectroscopically observed 1,2-µ-peroxo diferric P intermediate, the proton-electron uptake by P is the favored mechanism for catalytic activation, since it allows a significant reduction of the barrier of the initial (and rate-determining) H-atom abstraction from the stearoyl substrate as compared to the "proton-only activated" pathway. Also, we ruled out that a Q-like intermediate (high-valent diamond-core bis-µ-oxo-[FeIV]2 unit) is involved in the reaction mechanism. Our mechanistic picture is consistent with the experimental data available for Δ9D and satisfies fairly stringent conditions required by Nature: the chemo-, stereo-, and regioselectivity of the desaturation of stearic acid. Finally, the mechanisms evaluated are placed into a broader context of NHFe2 chemistry, provided by an amino acid sequence analysis through the families of the NHFe2 enzymes. Our study thus represents an important contribution toward understanding the catalytic action of the NHFe2 enzymes and may inspire further work in NHFe(2) biomimetic chemistry.

Interplay between Conformational Strain and Intramolecular Interaction in Protein Structures: Which of Them Is Evolutionarily Conserved?

By computing strain energies of peptide fragments within protein structures and their intramolecular interaction energies, we attempt to reveal general biophysical trends behind the secondary structure formation in the context of protein evolution. Our "protein basis set" consisted of 1143 representatives of different folds obtained from curated SCOPe database, and for each member of the set, the strain and intramolecular energy was calculated on the "rolling tripeptide" basis, employing the DFT-D3/COSMO-RS method for the former and the QM-calibrated force field method (MM) for the latter. The calculated data, strain and interactions, were correlated with the conservation of amino acid residues in secondary structure elements and also with the level of the residue burial within the protein three-dimensional structure. It allowed us to formulate several observations concerning fundamental differences between two main secondary structure motifs: α-helices and ß-strands. We have shown that a strong interaction is one of the determining characteristics of the ß-sheet formation, at least at the level of tripeptides (and likely penta- or heptapeptides, too), and that the ß-strand is a prevailing secondary structure in the strongly-interacting regions of the protein folds conserved by evolution. On the other hand, low strain was neither proven to be an important physicochemical property conserved by evolution nor does it correlate with the propensity for the α-helix and ß-strand. Finally, it has been demonstrated that the strong interaction has a certain level of connection with residue burial; however, we demonstrate that these two characteristics should be rather regarded as two complementary factors. These findings represent an important contribution to understanding protein folding from first principles, which is a complementary approach to ongoing efforts to solve the protein folding problem by knowledge-based approaches and machine-learning.

Dinucleoside polyphosphates act as 5'-RNA caps in bacteria.

It has been more than 50 years since the discovery of dinucleoside polyphosphates (NpnNs) and yet their roles and mechanisms of action remain unclear. Here, we show that both methylated and non-methylated NpnNs serve as RNA caps in Escherichia coli. NpnNs are excellent substrates for T7 and E. coli RNA polymerases (RNAPs) and efficiently initiate transcription. We demonstrate, that the E. coli enzymes RNA 5'-pyrophosphohydrolase (RppH) and bis(5'-nucleosyl)-tetraphosphatase (ApaH) are able to remove the NpnN-caps from RNA. ApaH is able to cleave all NpnN-caps, while RppH is unable to cleave the methylated forms suggesting that the methylation adds an additional layer to RNA stability regulation. Our work introduces a different perspective on the chemical structure of RNA in prokaryotes and on the role of RNA caps. We bring evidence that small molecules, such as NpnNs are incorporated into RNA and may thus influence the cellular metabolism and RNA turnover.

Factors Stabilizing ß-Sheets in Protein Structures from a Quantum-Chemical Perspective.

Protein folds are determined by the interplay between various (de)stabilizing forces, which can be broadly divided into a local strain of the protein chain and intramolecular interactions. In contrast to the α-helix, the ß-sheet secondary protein structure is significantly stabilized by long-range interactions between the individual ß-strands. It has been observed that quite diverse amino acid sequences can form a very similar small ß-sheet fold, such as in the three-ß-strand WW domain. Employing "calibrated" quantum-chemical methods, we show herein on two sequentially diverse examples of the WW domain that the internal strain energy is higher in the ß-strands and lower in the loops, while the interaction energy has an opposite trend. Low strain energy computed for peptide sequences in the loop 1 correlates with its postulated early formation in the folding process. The relatively high strain energy within the ß-strands (up to 8 kcal mol-1 per amino acid residue) is compensated by even higher intramolecular interaction energy (up to 15 kcal mol-1 per residue). It is shown in a quantitative way that the most conserved residues across the structural family of WW domains have the highest contributions to the intramolecular interaction energy. On the other hand, the residues in the regions with the lowest strain are not conserved. We conclude that the internal interaction energy is the physical quantity tuned by evolution to define the ß-sheet protein fold.

Preparation of (Pentafluorosulfanyl)benzenes by Direct Fluorination of Diaryldisulfides: Synthetic Approach and Mechanistic Aspects.

Direct fluorination of ortho-, meta- and para-substituted aromatic thiols and disulfides using elemental fluorine afforded substituted (pentafluorosulfanyl)benzenes. This work thus represents the first study of the scope and limitation of direct fluorination for the synthesis of new SF5 -containing building blocks. Fluorinations in batch and flow modes were compared. A comprehensive computational study was carried out employing density functional and wave function methods to elucidate the reaction mechanism of the transformation of ArSF3 into ArSF5 . Eliminating various nonradical pathways, it has been shown that the reaction proceeds by a radical mechanism, initiated by the attack of the F. on the ArSF3 moiety, propagated via an almost barrierless F2 +ArSF4 . →ArSF5 +F. step and terminated by the ArSF4 . +F. →ArSF5 . Most of the calculated data are in very good agreement with experimental observations concerning the ortho-substituent effect on the reaction rates and yields.

Toward Ab Initio Protein Folding: Inherent Secondary Structure Propensity of Short Peptides from the Bioinformatics and Quantum-Chemical Perspective.

By combining bioinformatics with quantum-chemical calculations, we attempt to address quantitatively some of the physical principles underlying protein folding. The former allowed us to identify tripeptide sequences in existing protein three-dimensional structures with a strong preference for either helical or extended structure. The selected representatives of pro-helical and pro-extended sequences were converted into "isolated" tripeptides-capped at N- and C-termini-and these were subjected to an extensive conformational sampling and geometry optimization (typically thousands to tens of thousands of conformers for each tripeptide). For each conformer, the QM(DFT-D3)/COSMO-RS free-energy value was then calculated, Gconf(solv). The Δ Gconf(solv) is expected to provide an objective, unbiased, and quantitatively accurate measure of the conformational preference of the particular tripeptide sequence. It has been shown that irrespective of the helical vs extended preferences of the selected tripeptide sequences in context of the protein, most of the low-energy conformers of isolated tripeptides prefer the R-helical structure. Nevertheless, pro-helical tripeptides show slightly stronger helix preference than their pro-extended counterparts. Furthermore, when the sampling is repeated in the presence of a partner tripeptide to mimic the situation in a ß-sheet, pro-extended tripeptides (exemplified by the VIV) show a larger free-energy benefit than pro-helical tripeptides (exemplified by the EAM). This effect is even more pronounced in a hydrophobic solvent, which mimics the less polar parts of a protein. This is in line with our bioinformatic results showing that the majority of pro-extended tripeptides are hydrophobic. The preference for a specific secondary structure by the studied tripeptides is thus governed by the plasticity to adopt to its environment. In addition, we show that most of the "naturally occurring" conformations of tripeptide sequences, i.e., those found in existing three-dimensional protein structures, are within ∼10 kcal·mol-1 from their global minima. In summary, our "ab initio" data suggest that complex protein structures may start to emerge already at the level of their small oligopeptidic units, which is in line with a hierarchical nature of protein folding.

Breaking Benzene Aromaticity-Computational Insights into the Mechanism of the Tungsten-Containing Benzoyl-CoA Reductase.

Aromatic compounds are environmental pollutants with toxic and carcinogenic properties. Despite the stability of aromatic rings, bacteria are able to degrade the aromatic compounds into simple metabolites and use them as growth substrates under oxic or even under anoxic conditions. In anaerobic microorganisms, most monocyclic aromatic growth substrates are converted to the central intermediate benzoyl-coenzyme A, which is enzymatically reduced to cyclohexa-1,5-dienoyl-CoA. The strictly anaerobic bacterium Geobacter metallireducens uses the class II benzoyl-CoA reductase complex for this reaction. The catalytic BamB subunit of this complex harbors an active site tungsten-bis-pyranopterin cofactor with the metal being coordinated by five protein/cofactor-derived sulfur atoms and a sixth, so far unknown, ligand. Although BamB has been biochemically and structurally characterized, its mechanism still remains elusive. Here we use continuum electrostatic and QM/MM calculations to model benzoyl-CoA reduction by BamB. We aim to elucidate the identity of the sixth ligand of the active-site tungsten ion together with the interplay of the electron and proton transfer events during the aromatic ring reduction. On the basis of our calculations, we propose that benzoyl-CoA reduction is initiated by a hydrogen atom transfer from a W(IV) species with an aqua ligand, yielding W(V)-[OH-] and a substrate radical intermediate. In the next step, a proton-assisted second electron transfer takes place with a conserved active-site histidine serving as the second proton donor. Interestingly, our calculations suggest that the electron for the second reduction step is taken from the pyranopterin cofactors rather than from the tungsten ion. The resulting cationic radical, which is distributed over both pyranopterins, is stabilized by conserved anionic amino acid residues. The stepwise mechanism of the reduction shows similarities to the Birch reduction known from organic chemistry. However, the strict coupling of protons and electrons allows the reaction to proceed under milder conditions.

Mechanisms and Specificity of Phenazine Biosynthesis Protein PhzF.

Phenazines are bacterial virulence and survival factors with important roles in infectious disease. PhzF catalyzes a key reaction in their biosynthesis by isomerizing (2 S,3 S)-2,3-dihydro-3-hydroxy anthranilate (DHHA) in two steps, a [1,5]-hydrogen shift followed by tautomerization to an aminoketone. While the [1,5]-hydrogen shift requires the conserved glutamate E45, suggesting acid/base catalysis, it also shows hallmarks of a sigmatropic rearrangement, namely the suprafacial migration of a non-acidic proton. To discriminate these mechanistic alternatives, we employed enzyme kinetic measurements and computational methods. Quantum mechanics/molecular mechanics (QM/MM) calculations revealed that the activation barrier of a proton shuttle mechanism involving E45 is significantly lower than that of a sigmatropic [1,5]-hydrogen shift. QM/MM also predicted a large kinetic isotope effect, which was indeed observed with deuterated substrate. For the tautomerization, QM/MM calculations suggested involvement of E45 and an active site water molecule, explaining the observed stereochemistry. Because these findings imply that PhzF can act only on a limited substrate spectrum, we also investigated the turnover of DHHA derivatives, of which only O-methyl and O-ethyl DHHA were converted. Together, these data reveal how PhzF orchestrates a water-free with a water-dependent step. Its unique mechanism, specificity and essential role in phenazine biosynthesis may offer opportunities for inhibitor development.

Computational Biochemistry-Enzyme Mechanisms Explored.

Understanding enzyme mechanisms is a major task to achieve in order to comprehend how living cells work. Recent advances in biomolecular research provide huge amount of data on enzyme kinetics and structure. The analysis of diverse experimental results and their combination into an overall picture is, however, often challenging. Microscopic details of the enzymatic processes are often anticipated based on several hints from macroscopic experimental data. Computational biochemistry aims at creation of a computational model of an enzyme in order to explain microscopic details of the catalytic process and reproduce or predict macroscopic experimental findings. Results of such computations are in part complementary to experimental data and provide an explanation of a biochemical process at the microscopic level. In order to evaluate the mechanism of an enzyme, a structural model is constructed which can be analyzed by several theoretical approaches. Several simulation methods can and should be combined to get a reliable picture of the process of interest. Furthermore, abstract models of biological systems can be constructed combining computational and experimental data. In this review, we discuss structural computational models of enzymatic systems. We first discuss various models to simulate enzyme catalysis. Furthermore, we review various approaches how to characterize the enzyme mechanism both qualitatively and quantitatively using different modeling approaches.

PyCPR - a python-based implementation of the Conjugate Peak Refinement (CPR) algorithm for finding transition state structures.

Conjugate peak refinement (CPR) is a powerful and robust method to search transition states on a molecular potential energy surface. Nevertheless, the method was to the best of our knowledge so far only implemented in CHARMM. In this paper, we present PyCPR, a new Python-based implementation of the CPR algorithm within the pDynamo framework. We provide a detailed description of the theory underlying our implementation and discuss the different parts of the implementation. The method is applied to two different problems. First, we illustrate the method by analyzing the gauche to anti-periplanar transition of butane using a semiempirical QM method. Second, we reanalyze the mechanism of a glycyl-radical enzyme, namely of 4-hydroxyphenylacetate decarboxylase (HPD) using QM/MM calculations. In the end, we suggest a strategy how to use our implementation of the CPR algorithm. The integration of PyCPR into the framework pDynamo allows the combination of CPR with the large variety of methods implemented in pDynamo. PyCPR can be used in combination with quantum mechanical and molecular mechanical methods (and hybrid methods) implemented directly in pDynamo, but also in combination with external programs such as ORCA using pDynamo as interface. PyCPR is distributed as free, open source software and can be downloaded from http://www.bisb.uni-bayreuth.de/index.php?page=downloads . Graphical Abstract PyCPR is a search tool for finding saddle points on the potential energy landscape of a molecular system.

Ferrous and ferric state of cytochromes P450 in intact Escherichia coli cells: a possible role of cytochrome P450-flavodoxin interactions.

OBJECTIVES: Cytochromes P450 (CYPs) are heme enzymes oxygenating a broad range of substrates. Their activity is dependent on the presence of a suitable electron donor (eukaryotic NADPH:CYP oxidoreductase or cytochrome b5). The Escherichia naturally contain no CYPs and no NADPH:CYP oxidoreductase, however it was reported that some CYPs heterologously expressed in E. coli may exist in the ferrous form. A small bacterial flavoprotein, flavodoxin is considered to be responsible for reduction some of these CYPs. METHODS: The reduction state of several human CYPs expressed in the intact living E. coli cells was examined. In addition, molecular dynamics and steered molecular dynamics simulations were performed to predict and compare affinity of flavodoxin toward selected CYPs. RESULTS: We determined the reduction state of five human CYPs heterologously expressed in E. coli. The computationally predicted stabilities of CYP-flavodoxin complexes correlate with the percentage of reduced CYPs in bacterial cells. The mean electron transfer distance within optimized complexes was also related to the percentage of reduced CYPs. CONCLUSION: Depending on the resting state, the CYPs heterologously expressed in E. coli could be divided into two groups; CYP2C8, 2C9, 3A4 are in E. coli present mainly in the oxidized form; while CYP1A1, 1A2, 2A6, 2A13, 2B6, 2D6 are found predominantly in the reduced form. We found a significant correlation between the stability of CYP-flavodoxin complexes and the percentage of reduced CYPs in bacteria. Hence, the naturally expressed flavodoxin is probably responsible for reduction of a larger group of human CYPs in bacterial cells.