Enhancing CYP4B1 catalytic performance through an orthologous sequence-guided engineering approach

We describe a sequence-guided approach to identify and engineer key amino acids in orthologous CYP4B1 enzymes for selective ω-hydroxylation of saturated fatty acids and sequential hydroxylation of trans -stilbene to trans -4,4′- and trans -4,3′-dihydroxystilbene. Orthologous enzymes generally share high sequence identities and catalyze the same reactions, but often exhibit subtle yet distinct activity and selectivity profiles. Understanding the sequence determinants that underlie this observation is not only essential for elucidation of evolutionary adaptations, but can also be used for rational biocatalyst design. Within this study, sequence comparison of a set of orthologous CYP4B1 enzymes identified four amino acids as key activity and selectivity drivers for catalysis. These were systematically engineered by introducing potentially beneficial single amino acid exchanges into less active/selective CYP4B1 orthologs, and vice versa. The exchanges shifted activity and selectivity trends for ω-hydroxylation of lauric acid by the mutated orthologs in a manner consistent with those observed for the corresponding wildtype enzymes across multiple orthologous backgrounds. Notably, combinatorial mutagenesis of two amino acids resulted not only in additive, but also in orthologous- and substrate-dependent adverse influences on activity and selectivity, which likely reflects the enzyme function optimized by nature in the context of the individual orthologs evolutionary background. For conversion of trans -stilbene, two beneficial amino acid exchanges between two CYP4B1 orthologs resulted in 95% conversion of trans -stilbene and 85% dihydroxystilbenes. In summary, we demonstrate that systematic introduction of only one or two sequence-guided mutations can be sufficient to substantially alter activity and selectivity profiles of orthologous enzymes for a given substrate, which highlights the usefulness of this approach for rational design of biocatalysts.

Cytochrome P450 monooxygenases (P450s) are heme thiolate enzymes widespread in all domains of life (Lamb et al., 2019; Nelson, 2018). They are best known for their ability to hydroxylate non-activated C-H bonds in organic compounds through reductive cleavage of molecular oxygen. Other reaction types include deamination, epoxidation of C=C double bonds, and O- and N-dealkylation (Bernhardt and Urlacher, 2014; Guengerich and Munro, 2013). As one of the most versatile enzyme superfamilies, the substrate scope ranges from small molecules, such as n -alkanes, over aromatic compounds to more sophisticated substrate groups, such as steroids and lignans. Their ability to catalyze reactions in a regio-, chemo-, or stereoselective manner, together with their versatility, makes P450s powerful biocatalysts (Urlacher and Girhard, 2019).

Based on sequence identities, P450s are classified into families and subfamilies in a systematic nomenclature (Nelson, 2018; Nelson et al., 1993): The term ‘CYP’ is followed by a number indicating the families that share >40% sequence identity; the subsequent letter indicates the subfamily sharing >55% sequence identity. Within a subfamily, individual enzymes are numbered consecutively, e.g., CYP102A1 is the systematic nomenclature for the well-investigated P450 BM3 enzyme from Bacillus megaterium (Fansher et al., 2024). One exception of these general nomenclature rules occurs in mammals; here, the same enzyme identifier (e.g., CYP4B1) is used for different species even though the respective P450s slightly differ in their amino acid sequence. These P450s are orthologous enzymes that share a common ancestor and therefore generally exhibit high sequence identities of up to 99% across species, depending on the degree of evolutionary conservation of the respective genes.

Rational design of biocatalysts aiming to enhance activity or selectivity is usually carried out by mimicking evolution through introduction of structure-guided mutations at amino acid positions assumed to interact with the respective substrate of interest (Chubb et al., 2026; Xu and Du, 2018). This approach usually requires introduction of a larger variety of amino acids at multiple positions resulting in a considerable body of work, especially if single mutations are combined, and is also highly dependent on the quality of the homology model if high resolution crystal structures are not available. Generally, in the first step of this approach a homology model is generated based on available crystal structures using methods such as SwissModel or AlphaFold (Abramson et al., 2024; Waterhouse et al., 2018). Usually, the properties of the obtained model depend strongly on the resolution and sequence similarity to that of the template structure. Furthermore, not only the orientation of the bound substrate in the active site needs to be suitable for the desired reaction, but also the accessibility of a ligand into the binding pocket plays a crucial role (Kokkonen et al., 2019). Therefore, several prediction tools can be applied to analyze these pathways, including CAVER, CaverDock or MOLEonline, which provide insights to static tunnels and substrate channelling of a given crystal structure or homology model, respectively (Chovancova et al., 2012; Pravda et al., 2018; Filipovič et al., 2020). Next, for interactions with substrates underlying dynamic processes in enzymes, different conformations of the enzymes and motion of the substrates need to be considered. However, one main limitation of homology modelling is that static enzyme structures are predicted (Abramson et al., 2024). If a structure prediction is based on crystal structures in a closed conformation, also the predicted homology model is accordingly closed. Moreover, predictions based on open conformations can also end up in a closed state of a homology model as exemplified shown by Abramson et al. on the ubiquitin ligase cereblon (Abramson et al., 2024). In consequence, the analysis of dynamic properties and engineering of the involved positions can be challenging if based solely on homology models.

In contrast to structure-guided analysis, the identification of amino acids that underlie the distinct activity and selectivity profiles of orthologues enzymes should be possible by a rather easy and straightforward sequence-guided approach: The activity and selectivity profiles of a set of orthologous wildtype enzymes against the substrate of interest is determined, accompanied by analysis of a multiple sequence alignment of the chosen orthologs (if applicable in combination with additional orthologs). Due to their generally high evolutionary sequence conservation, often only few or even single amino acid positions can be elucidated that are responsible to substantially alter catalytic efficiency, substrate scope, or product distribution of an individual ortholog, respectively. Thus, identifying these key amino acids is highly valuable for rational biocatalyst design. Enzyme families that exhibit conserved catalytic mechanisms, but functional diversity across species are particularly suitable for such comparative analyses, and one prominent example is the P450 enzyme family.

CYP4B1 belongs to the mammalian P450 family 4 (CYP4) acting mainly as ω-hydroxylases in the metabolism of fatty acids and eicosanoids. For such compounds containing saturated alkyl chains (e.g., n -alkanes, n -alkanols or saturated fatty acids), a multitude of P450s have been described that exclusively attack either the subterminal positions (e.g., CYP102A1 acting at ω-1, ω-2, ω-3 (Jang et al., 2016; Padayachee et al., 2025)), or positions close to the primary alcohol or carboxylic group (e.g., CYP152A1 preferring the β-position, or CYP152A2 preferring the α-position (Girhard et al., 2007; Matsunaga et al., 1999)). In contrast, several members of the mammalian CYP4 family share the ability to hydroxylate the terminal ω-position, which is chemically more challenging (Eschenfeldt et al., 2003; Ortiz de Montellano, 2008; Scott, 2017). This ability has been attributed to an autocatalytic post-translational modification, whereby a glutamic acid in the protein backbone forms a covalent ester bond with the 5-methylpropionate of the heme cofactor, leading to a more rigid structure of the heme itself. Best investigated in this regard is the CYP4B1 ortholog from rabbit (rCYP4B1): While the tertiary structure and major elements of secondary structure of the CYP4B1 catalytic domain resemble that of other mammalian microsomal P450s, p.E310 was found responsible for covalent heme binding (Hsu et al., 2017), and this position is highly conserved within the underlying sequence of almost every cyp4b1 gene described in the literature.

Strikingly, the situation is completely different in humans: Previous studies revealed that the human wildtype ortholog (hCYP4B1) is catalytically inactive, which was attributed from structural point of view to one amino acid exchange in the meander region where a serine is present instead of an otherwise highly conserved proline (p.S427) (Röder et al., 2023; Zheng et al., 1998).

However, our own research revealed that the picture of CYP4B1 emerging from the conclusions of previous studies was incomplete: First, when we reinvestigated oxidation of saturated fatty acids and fatty alcohols with C9–C16 carbon chain length, we found that upon conversion of lauric acid 1 by wildtype rCYP4B1 the ω-position is not exclusively hydroxylated, but instead seven hydroxylated products are formed, and only 40% of all products is the ω-hydroxylated one (Thesseling et al., 2020).

Second, we investigated the activity profiles of 12 CYP4B1 orthologs from non-human primate origin to precisely pinpoint if, and at which time points in evolution these orthologs were rendered inactive (Hüsken et al., 2025). Sequence-guided amino acid analysis of the chosen CYP4B1 orthologs revealed four additional evolutionary changes, each leading to loss of catalytic function of ape and human orthologs: p.V71G specific for Denisovans, p.R106C, p.R244H, and an exon deletion found exclusively in gorilla CYP4B1. Systematic functional analyses proved that these changes rendered the respective orthologs catalytically inactivate, but we also demonstrated that reversion by mutagenesis restored their catalytic activity indicating that their structural framework remained intact. Despite this evolutionary aspect, we noticed that the catalytically active CYP4B1 orthologs showed distinct activity and selectivity profiles for conversion of lauric acid 1 (Figure 1A), and that some orthologs exhibit considerably higher selectivity for terminal ω-hydroxylation (1a) in comparison to rCYP4B1 (Figure 1C). This suggested adaptations of the orthologs in response to an evolutionary requirement for higher selectivity.

Third, in a follow-up study, we described a systematic substrate screening with two representative catalytically active CYP4B1 orthologs from gray mouse lemur (mlCYP4B1) and green monkey (chCYP4B1), and reported their product and activity profiles in comparison to that of rCYP4B1 (Röder et al., 2025). In summary it turned out that in addition to differences in regioselectivity toward lauric acid 1, the orthologs also displayed distinct substrate profiles: Whereas rCYP4B1 generally favored smaller molecules (shorter chain-lengths fatty acids and cyclic terpenes, such as carvone), the primate orthologs showed higher acceptance of longer chain-lengths fatty acids and sterically more demanding substrates. For instance, mlCYP4B1 and chCYP4B1 were more efficient than rCYP4B1 in catalyzing hydroxylations of trans -stilbene 2 (∼80% conversion versus 27%), and producing trans -4,4′-dihydroxystilbene 2b and trans -4,3′-dihydroxystilbene 2c in two consecutive hydroxylations via the mono-hydroxylated intermediate trans -4-hydroxystilbene 2a (Figure 1B).

In conclusion, our studies on CYP4B1 suggest a more complex picture of the metabolic roles of individual CYP4B1 orthologs, yet the structural and mechanistic details underlying this hypothesis are not well understood and further investigation is required to close this knowledge gap. Within this study we describe a sequence-guided approach to identify key amino acids that influence the activity and selectivity of orthologous CYP4B1 enzymes. Through systematic engineering of these positions by mutagenesis, we were able to link the observed variations in activity and selectivity profiles to individual sequence variations of four representative CYP4B1 orthologs. Our experimental results are underpinned by computational analyses elucidating mechanistic details of tunnel formation, substrate guidance and orientation in the active site of individual orthologs.

Despite ω-hydroxylation 1a of lauric acid 1 (Figure 1A) we have investigated the sequential hydroxylation of trans -stilbene 2 to trans -4,4′-dihydroxystilbene 2b and trans -4,3′-dihydroxystilbene 2c. Recent studies by other research groups have reported a pronounced potential of 2b and 2c in cancer treatment that has also been described for resveratrol (trans -3,5,4′-trihydroxystilbene 2d) (Kimura et al., 2020; Savio et al., 2022; Savio et al., 2016). Since dihydroxystilbenes (DHS) 2b and 2c occur only rarely in nature, while 2 is an abundant and cheap precursor, producing 2b and 2c with biocatalysts, such as CYP4B1, presents a viable option (Figure 1B).

Plasmid construction for cyp4b1 expression in E. coli

Primate cyp4b1 -genes (chCYP4B1, Chlorocebus sabaeus, XM007978860.2; mlCYP4B1, Microcebus murinus, XM_012776101.1; rhCYP4B1, Rhinopithecus bieti, XM_017881706.1; human codon-optimized synthetic genes obtained from Geneart AG, Regensburg, Germany) were amplified via PCR and ligated into pET22b(+) as described previously (Röder et al., 2025). Recombinant expression and modification of rabbit cyp4b1 (NM_001082103, isoform 1) was achieved using construct 4B1#7 (N-terminal truncation of residues 3–19; (Cheesman et al., 2003)) in the pET22b(+) vector, as described in (Roellecke et al., 2017).

The Q5 Site-Directed Mutagenesis Kit was used according to the manufacturer’s manual (catalog no. E0554S, NEB, Ipswich, MA, United States); the respective primers are listed in Supplementary Table S1.

Expression of recombinant enzymes in E. coli

Recombinant expression of cyp4b1 in strain Escherichia coli OverExpress C43(DE3) [F − ompT hsdS B (r B − m B −) gal dcm (DE3)] (Lucigen Corporation, Middleton, WI, United States) and purification via affinity chromatography were carried out as described previously (Röder et al., 2025). SDS-gels of all CYP4B1 variants investigated within this study are provided in Supplementary Figure S1.

The P450 concentration of each sample was determined by recording CO-difference spectra, using an extinction coefficient of ε 450–490 nm = 91 mM -1 cm -1 (Guengerich et al., 2009; Omura and Sato, 1964). CO-difference spectra of all CYP4B1 variants investigated within this study are provided in Supplementary Figure S2.

Expression and purification of flavodoxin (YkuN) from Bacillus subtilis and flavodoxinreductase (Fdr) from E. coli JM109 were carried out as described previously (Girhard et al., 2010).

For expression of bovine cytochrome b 5 (bCYB5), the gene carrying plasmid pET17 b (kindly provided by Prof. Rita Bernhardt as described in (Neunzig et al., 2014)) was used in E. coli OverExpress C43(DE3) following an expression protocol described by Mulrooney and Waskell (2000). Purification via IMAC was carried out following the protocol described for bovine CPR (Zöllner et al., 2008).

Glucose dehydrogenase (GDH) from Bacillus megaterium (gdhIV; GenBank D10626) was expressed in E. coli BL21(DE3) from pET22b(+) and partially purified by a combination of salt precipitation and heat denaturation of endogenous E. coli proteins as described elsewhere (Kranz-Finger et al., 2018).

Reconstitution of CYP4B1 activity

Substrate conversions were carried out in vitro in 50 mM potassium phosphate buffer, pH 7.5, in a total volume of 100 μL. All reaction mixtures contained 0.25 μM CYP4B1, 1 μM Fdr, 10 μM YkuN, 0.5 μM bCYB5, 1.000 U mL −1 catalase, 25 U mL −1 GDH, 20 mM glucose, 25 μg mL −1 1,2-didodecanoylsn-glycero-3-phosphocholine (DLPC), 200 μM substrate (from a 10 mM stock solution dissolved in dimethyl sulfoxide), and 200 μM NADP +. Samples were incubated at 30 °C with 600 rpm shaking in a thermomixer C device (Eppendorf, Hamburg, Germany) for 90 min in case of 1, or 180 min in case of 2 and 2a.

GC/MS and LC/MS analysis

Conversions and product distributions of samples containing 1 were analyzed via GC/MS (GC/MS-QP2010 plus, Shimadzu, Germany), equipped with FS-Supreme-5ms column as described in our previous study (Röder et al., 2025). Briefly, samples were acidified with 2 µL of 37% HCl, followed by the addition of the internal standard undecanoic acid (100 µM). The substrates and products were extracted twice in 500 µL of tert-butyl methyl ether (MTBE), dried over MgSO 4, and MTBE was evaporated. The residues were dissolved in 30 µL of N,O-bis(trimethylsilyl)trifluoroacetamide with 1% trimethylchlorosilane, incubated at 80 °C for 30 min, and analyzed by GC/MS.

Samples containing 2 or 2a were analyzed via LC/ESI-MS (LCMS-2020, Shimadzu, Germany), equipped with Chromolith Performance RP-18e column, as described in our previous study (Röder et al., 2025). Briefly, the internal standard flavone (100 µM) was added. The substrates and products were extracted in 300 µL ethyl acetate. The organic solvent was evaporated, and the residues were resolved in 30 µL acetonitrile, and then subjected to LC/ESI-MS analysis. Product yields were calculated based on calibrations with authentic standards (see Supplementary Data).

Computational methods

Homology models of chCYP4B1, rhCYP4B1 and mlCYP4B1 were generated by SwissModel (Waterhouse et al., 2018) on the basis of the crystal structure of rCYP4B1 (5t6q.pdb), and mutations inserted using the software Swiss Pdb-Viewer (Windows version 4.1.0). For substrate channel predictions, the Pymol plugin CAVER 3.03 was applied using a minimum probe radius of 1.0 Å, a shell depth of 4 Å, a shell radius of 3 Å, and a clustering threshold of 3.5 Å (Chovancova et al., 2012). The heme iron was selected as starting point with a maximum distance of 3 Å in the desired radius of 5 Å. Further analysis of the substrate channels was carried out with the software CAVER Analyst 2.0 beta using the same settings (Jurcik et al., 2018).

For in silico docking, the P450 models were prepared with AutoDock Tools (Windows version 1.5.7) by addition of polar hydrogens and Kollmann charges to the amino acids. The heme was prepared by the addition of Gasteiger charges, and the charge of the heme iron was manually set to +0.400e to correspond to Fe(II). This was counterbalanced by modifying the partial charges of the nitrogen atoms that were set to −0.348e. The ligand molecules 1, 2, or 2a were energetically optimized employing the build-in MM + force field of HYPERCHEM version 8.0.10 (Hypercube Inc.,Gainsville, FL, United States). AutoDock atom types and partial charges were computed with an in-house PERL script to be consistent with our previous reports (Khatri et al., 2015; Röder et al., 2025). In silico docking was conducted using AutoDock Vina (Windows version 1.2.7) applying default settings, in which the protein was kept rigid and all rotatable bonds of the ligands were treated as flexible.

Initial determination of activity and selectivity driving positions in CYP4B1 orthologs

Key amino acid positions hypothesized to influence activity and/or selectivity of CYP4B1 were identified by a straightforward multiple sequence alignment of seven orthologous amino acid sequences sharing high sequence identities ranging from 85% up to 99% to each other (Supplementary Figure S3). The multiple sequence alignment revealed several positions with sequence differences, out of which 14 displayed systematic variations of individual orthologs from the otherwise highly conserved evolutionary consensus sequence (Supplementary Figure S4); e.g., every ortholog had uniformly a conserved phenylalanine at position 113/118 (p.F113/118), with the exception of two orthologs containing a leucin at the respective position.

Initial investigation of these 14 positions (Supplementary Table S2) revealed four key positions showing consistent trends toward the previously observed properties of the respective wildtype enzymes on activity and/or selectivity profiles for hydroxylation of lauric acid 1 across all tested orthologs and were therefore chosen for detailed analysis. These positions are namely, p.113/118, p.138/143, p.221/226, and p.483/488 (Figure 2A). According to the crystal structure of rCYP4B1 three out of these four positions are located within substrate recognition sites 1 (p.113/118), 2 (p.221/226), and 6 (p.483/488), whereas p.138/143 is located at the end of the C helix at the transition to the C-D loop (Figure 2B; Supplementary Figure S4).

Position 113/118

p.113/118 consists of an evolutionary conserved phenylalanine, whereas in rhCYP4B1 (and piCYP4B1; not investigated within this study) leucine is present instead. In our previous study, we have observed that these two CYP4B1 orthologs exhibited the lowest activities toward 1 and 4-ipomeanol as compared to the other active orthologs from primate origin and rCYP4B1 (Hüsken et al., 2025). As hypothesized, when mutation p.F113L was introduced into rCYP4B1, a drop in conversion of 1 from 65% to 14% was observed, while the ω-hydroxylation selectivity increased from 40% to 85%. This initial result indicated that p.113/118 influences both, activity and selectivity.

Position 138/143

p.138/143 is located outside the active site and displays another exchange exclusive for rhCYP4B1 (and piCYP4B1), where an alanine is present instead of a conserved proline at the end of the C helix. The initial mutation p.A143P in rhCYP4B1 resulted in a distinct increase in conversion of 1 from 41% to 86%, while selectivity for ω-hydroxylation remained unchanged.

Position 221/226

p.221/226 consists of an evolutionary conserved leucine, whereas in rCYP4B1 an isoleucine, and in mlCYP4B1 a methionine is found, respectively. Initial introduction of mutations p.I221L or p.I221M into rCYP4B1 resulted in both cases in increased selectivity for ω-hydroxylation of 1, but also a slight decrease in conversion.

Position 483/488

Finally, at the beginning of the substrate channel in SRS-6, a conserved p.L488 was identified, except for rCYP4B1 containing p.P483 instead. The initial mutation p.P483L in rCYP4B1 led to an increase in both, increase in selectivity for 1a from 40% to 65%, and increase in conversion of 1 from 65% to 87%.

Systematic single mutagenesis of key positions

For detailed investigation of the individual influences of p.113/118, p.138/143, p.221/226 and p.483/488 systematic single mutated variants of the four investigated CYP4B1 orthologs were generated and their activity and selectivity profiles toward lauric acid 1 (Figure 3, left panel) and toward trans -stilbene 2 (Figure 3, right panel) were determined.

Position 113/118

p.F113 located in the active site at SRS-1 is strongly conserved in CYP4B1 orthologs along evolution. The presence of p.L118 in rhCYP4B1 was hypothesized as one reason for the observed low activity toward 1 and 4-ipomeanol in our previous study (Hüsken et al., 2025). This assumption is supported by the observation that mutations p.F113L in rCYP4B1, and p.F118L in chCYP4B1 and mlCYP4B1 resulted in a 3- to 11-times decreased conversion of 1, albeit accompanied by an increase in selectivity for ω-hydroxylation (Figure 3A, left panel). The same behavior of single-mutated variants was also observed for conversion of 2 (Figure 3A, right panel). Unexpectedly, the mutation p.L118F in rhCYP4B1 did not increase activity toward 1 and 2 as hypothesized, but instead also resulted in a decrease of conversion (from 41% to 6% for 1, from 57% to 36% for 2). This indicated that potential interactions of p.113/118 with additional amino acids occur that cannot be evaluated solely by sequence comparison (see chapter 3.3.1).

Position 138/143

The highly conserved p.P138/143 is located at the end of the C helix and initiates the transition to the C-D loop in rCYP4B1. Sequence alignments of this position with CYP102A1 and CYP3A4 reveals its high conservation as a proline even across non-orthologous P450 enzymes (Sevrioukova et al., 1999; Yano et al., 2004). The mutations p.P138A in rCYP4B1, as well as p.P143A in chCYP41 and mlCYP4B1 resulted in a decreased conversion of 1, while preserving ω-hydroxylation selectivity (Figure 3B, left panel). More importantly, and as hypothesized, the mutation p.A143P in rhCYP4B1 led to increased conversions of 1 (from 41% to 86%) and 2 (from 57% to 82%).

Position 221/226

p.221/226 mainly impacts the ω-hydroxylation selectivity for 1. Whereas the wildtype CYP4B1 orthologs with higher selectivity contain either p.L226 or p.M226, the less selective wildtype rCYP4B1 contains p.I221 instead. As expected, the mutations p.L226I and p.M226I in chCYP4B1 and mlCYP4B1 resulted in a decreased selectivity toward the ω-position of 1 (Figure 3C, left panel). Conversely, introduction of p.I221M or p.I221L in rCYP4B1 increased selectivity for ω-hydroxylation, and it can generally be stated that in case of 1, the selectivity with methionine was higher than with leucine (with exception of chCYP4B1; Supplementary Table S3).

Furthermore, p.221/226 also influenced the ability to form 2b and 2c from 2. This reaction is assumed to occur via two consecutive hydroxylation steps with 2a as intermediate. The second hydroxylation of 2 occurs at its opposite 4′- or 3′-site, but a direct re-orientation in the active site is highly unlikely due to space-limitations. Instead, 2a has to leave and re-enter the active site in its opposite orientation. Consequently, the ability of CYP4B1 to produce 2b and 2c depends on its acceptance of both compounds, 2 and 2a. The calculated conversion values reflect the initial substrate depletion of 2 with selective formation of 2a. The depletion of 2 is either not, or only slightly changed by mutations at p.221/226 - with an exception for rhCYP4B1, where p.L226I led to an increase, but p.L226M to a decrease in conversion of ∼35% (Figure 3C, right panel). The activity with 2a is reflected by the calculated values for formation of 2b and 2c. These values indicate the quantity of 2b and 2c in relation to all reaction products (=100%), whereby the remainder to 100% is almost exclusively the reaction intermediate 2a (value not included to keep clarity). From these values it can be concluded that among all orthologs, those with p.L221/226 promote the second hydroxylation step with highest activity, followed by p.I221/226 (Figure 3C, right panel), while orthologs with p.M221/226 produce mostly only 2a. This hypothesis is supported by the observation that the mutations p.M226I and p.M226L in mlCYP4B1, increased the yields of 2b and 2c from 33% to 59% and 83%, respectively.

Position 483/488

p.483/488 is located in the 4-1 beta sheet near the entrance of the CYP4B1 substrate channel, and consists of a conserved leucine. In rCYP4B1, however, p.P483 is present instead. The mutation p.P483L in rCYP4B1 resulted in a higher selectivity for ω-hydroxylation of 1 (65% compared to 40%). The reverse mutation p. L488P in the other orthologs did either not or only slightly decrease their selectivity (Figure 3D, left panel; Supplementary Table S3), but led to an increased conversion of 1 by rhCYP4B1 (82% compared to 41%).

In case of 2, mutation p.L488P in chCYP4B1 and mlCYP4B1 resulted in an increased formation of 2b and 2c (Figure 3D, right panel). This indicates substrate dependent effects for this position that are different compared to those of 1.

Subsequent computational analysis of key positions

Position 113/118

The influence of p.113/118 became only visible upon subsequent computational analysis of the interactions with other amino acid residues (Figure 4). In rCYP4B1, p.F113 is part of a ‘network’ consisting of methionine p.M217 and three aromatic residues (p.F113, p.W116 p.F309), which share Met-sulphur-π (Met-aromatic) and π-π interactions with each other (Figure 4A). This network connects the helices F, B′ and I, and stabilizes the active site. The Met-aromatic interaction between one methionine in a ‘bridging motif’ with two aromatic residues provides additional stabilization energy to the active site and is found in numerous oxidoreductases (Waters, 2004; Weber and Warren, 2018; Weber and Warren, 2019). With the exception of rhCYP4B1 (and piCYP4B1), the involved residues are highly conserved in CYP4B1 orthologs. Importantly, in rhCYP4B1, the p.F118L substitution is not the only alteration of that motif, but the methionine is also substituted by alanine (p.M222A), which results in different interactions in that region (Figure 4B). Since no interaction between p.L118 and p.A222 is predicted in that case, this leads to the formation of an additional channel between these two residues that is not possible within rCYP4B1 where the aforementioned network (p.F113, p.M217) forms a barrier over the active site. Furthermore, stabilizing π-π interactions with p.W121 (B′ helix) are not possible with p.L118 in place of p.F113; instead, p.L118 shares weaker van der Waals forces with p.V128 of the B′-C loop. These amino acids substitutions in the active site are most probably one reason for the lower conversions by rhCYP4B1 observed in our experimental data, and furthermore explain the activity diminishing effects of the p.F113L mutation in rCYP4B1, as well as p.F118L mutations in chCYP4B1 and mlCYP4B1 (Figure 3A). In the in silico docking with 2 and 2a stabilizing edge-to-face T-form π-π interactions between the aromatic group and p.F113/118 are hypothesized. This is also not possible with p.L113/118, and seems to be another reason for the low activity, especially toward aromatic compounds (Supplementary Figure S5).

To prove this hypothesis by experiment, we produced double-mutated variants of exemplary chosen rhCYP4B1 and rCYP4B1. Indeed, the additional mutation of p.217/222 in both directions (p.F113L, p.M217A in rCYP4B1; p.L118F, p.A222M in rhCYP4B1) resulted in a recovery of activity toward 1 (Figure 5, left panel), and even higher activity for sequential hydroxylation of 2 (Figure 5, right panel). In contrast, the single mutation of p.217/222 resulted in a complete loss of activity in rhCYP4B1, and reduction of activity in rCYP4B1 with both substrates (drop from 65% to 16% for 1 with unchanged selectivity, and from 89% to 46% for 2).

Position 138/143

Computational analysis of p.138/143 revealed no direct interaction with any substrate, albeit our experimental data demonstrated its influence on CYP4B1 activity (Figure 3B). The presence of p.A143 in rhCYP4B1 instead of the conserved p.P143 is predicted to prolong the C helix as compared to that of rCYP4B1 (Figure 3B; middle panel; rhCYP4B1 in orange). Together with the helices D and H, the C helix is located on the proximal side of the heme, where the electron transfer from the FMN of flavodoxin occurs (Supplementary Figure S6). Structural changes in this region can influence the electron transfer between the P450 and its redox partner proteins, as shown e.g., for CYP102A1, where mutations of p.C104 in the C helix led to decreased FMN-to-heme electron transfer rates by 10- to 20-fold (Sevrioukova et al., 1999; Su et al., 2020). Therefore, this is suggested as most obvious cause of the declined activity observed for the p.P138/143A variants.

Positions 221/226 and 483/488

In silico docking suggests that p.221/226 is involved in the stabilization of 2 and 2a in the active site via van der Waals forces by sharing hydrophobic interactions with their aromatic group located further away from the heme (Supplementary Figure S5). For further analysis of p.221/226 and p.483/488 we applied CAVER and CAVER Analyst 2.0 (Chovancova et al., 2012; Jurcik et al., 2018) to obtain insights about the access of the substrate into the active site, even though the P450 homology models were generated based on the substrate bound crystal structure of rCYP4B1 (5t6q.pdb) that likely reflects the closed conformation. In silico prediction of the individual substrate channels and their contour revealed that all wildtype CYP4B1 orthologs displayed at least two different possibilities for substrates to enter the active site (following the prediction obtained with the set parameters; Figures 6A–D), but in every observed model the most probable substrate channel (with highest priority in the prediction) involves the amino acids p.221/226, p.483/488, as well as p.113/118.

The first distinct region is located around 10–13 Å distance to the heme. Here, the channel describes a curve in the direction to the active site. Within this region p.221/226 has the largest influence on the contour of the channel (Figure 6E): The radius of CYP4B1 orthologs with the smaller amino acid residues (p.I221 in rCYP4B1 and p.L226 in chCYP4B1 and rhCYP4B1) is 2.1–2.3 Å, whereas that of mlCYP4B1 with the largest residue p.M226 is clearly narrowed down to 1.8 Å. This observation would explain the experimentally observed behavior, where p.M221/226 showed the highest selectivity for ω-hydroxylation of 1 (Figure 3C; Supplementary Table S3), because there is less space available for movements of the substrate.

For the most probable channel, the path of substrates from the surface to the heme iron was predicted as follows: The substrate enters the channel by passing p.483/488 (β 4-1 sheet). At p.221/226 (F helix) it is directed along a curved path that leads to the active site, where the channel becomes perpendicular to the heme flanked by p.113/118 (B′ helix) on one side. By comparison of the channel contours between the wildtype orthologs, it became clear that they share a bottleneck region at around 17 Å distance to the heme iron. At this distance, p.483/488 has the biggest impact on the channel’s architecture (Figure 6E): Orthologs bearing p.L488 have a more confined bottleneck radius (BR) of 1.3–1.4 Å, whereas rCYP4B1 with p.P483 exhibits the widest BR of 1.7 Å, suggesting that the leucine causes a narrowing of the channel at this position. Conversely, the mutation p.L488P is hypothesized to widen the bottleneck which leads to an increase in catalytic activity. However, our experimental data revealed that this mutation only contributed to a higher activity of rhCYP4B1 toward 1 and 2, while the activities of chCYP4B1 and mlCYP4B1 remained unchanged (Figure 3D), and the ‘reverse mutation’ p.P483L in rCYP4B1 also increased conversion of 1 in addition to a higher ω-selectivity, even though leucine should narrow the entrance of the channel.

Upon closer inspection of the surrounding amino acid residues and by comparing the computational data with the multiple sequence alignment, it transpired that the majority of amino acids in this region are conserved across all CYP4B1 orthologs. However, despite the above investigated positions, we noticed another highly variable amino acid position that was not included in our investigation so far, since it appeared to be completely non-conserved in the multiple sequence alignment. This residue is located at p.222/227 as p.D222 in rCYP4B1, p.V227 in chCYP4B1, p.A227 in rhCYP4B1, and p.E227 in mlCYP4B1. p.222/227 forms the bottleneck together with p.483/488, p.221/226 and the conserved residues p.Q377/382 and p.Q218/223. Interactions between these residues therefore may explain the experimentally observed apparently contradictory effects of solely mutating p.483/488 (Figure 3D). Since alternative substrate channels are also predicted that do not pass via p.483/488, it is also possible that mutations at that position do not, or only barely affect the activity toward the respective substrate if an alternative channel is used.

Combinatorial mutagenesis of key positions

Systematic double-mutagenesis of rCYP4B1

We next aimed to analyze whether the positive and negative effects observed with the single mutated orthologs are additive when combined in double mutated variants. rCYP4B1 was chosen as exemplary ortholog for generation of the 9 possible double-mutated combinations, and their activity and selectivity profiles against 1 (Figure 7A) and 2 (Figure 7B) were determined.

We found that the combined mutations in rCYP4B1 generally resulted in additive effects with regard to the conversion of 1; e.g., the activity declining effect of both single mutations p.F113L and p.P138A resulted in only 2% conversion of 1 by the respective double-mutated variant. Notably, the activity loss by the single-mutated variants at p.113, p.138 and p.221 could in all cases partially be recovered by additional introduction of the activity increasing mutation p.P483L (Figure 7A). Regarding ω-hydroxylation selectivity, it can generally be stated that additive effects of double-mutated variants are not as obvious as compared to the respective single-mutated variants. Albeit almost every double-mutated variant had (up to five-fold) higher ω-hydroxylation selectivity when compared to the wildtype rCYP4B1, the selectivity increasing effects of the single mutations were either not or only partially additive in double-mutated variants.

The best overall performance was observed with the p.I221L, p.P483L variant that reached four-fold higher selectivity without loss of activity. Notably, the activity and selectivity profile of this variant is similar to that of the wildtype chCYP4B1 ortholog, which has the same residues at the respective positions (p.L226, p.L488).

Regarding conversion of 2 the already high activity of the wildtype rCYP4B1 ortholog (89% conversion) was preserved in most double-mutated variants, with exception of those carrying the p.F113L mutation (Figure 7B). However, none of the double-mutated variants showed an additive effect with regard to formation of 2b and 2c; instead decreased formation was observed as compared to both, wildtype and respective single-mutated variants (with exception of the p.I221L, p.P483L variant). This indicated that the intermediate 2a was less accepted, which became evident in reactions containing 2a as substrate (Supplementary Figure S7A).

Combinatorial mutagenesis for biocatalytic production of dihydroxystilbenes

Finally, we aimed to increase production of 2b and 2c from 2 by combinatorial introduction of the assumed beneficial mutations at p.221/226 and p.483/488 into the most active CYP4B1 ortholog. Wildtype chCYP4B1 produced 40% and mlCYP4B1 33% 2b and 2c from 2, whereas wildtype rCYP4B1 and rhCYP4B1 produced only 12% and 7%, respectively (Supplementary Figure S7B). When 2a was used as substrate under the same reaction conditions, highest conversion was observed with mlCYP4B1 (81%), followed by chCYP4B1 (54%), whereas rCYP4B1 and rhCYP4B1 converted only 30% and 3% of 2a, respectively (Supplementary Figure S7B).

As expected, the introduction of p.M226L and p.L488P into mlCYP4B1 resulted in almost complete conversion of 2 with formation of 161 µM 2b and 2c (≡ 85% conversion of 2a); this corresponds to a calculated substrate turnover frequency (TOF) of 7.8 µM substrate µM P450 −1 min -1 (Figure 8A). Computational analysis revealed that the bottleneck radius is decreased in each mutated mlCYP4B1 variant compared to the wildtype enzyme (Figures 8C–F). The mutation p.M226L creates more space at the region in the active center where the hydroxyl group of 2a is located during the second hydroxylation reaction (Figure 8D). This can also be seen on the contour of the channel around 11 Å length in both mlCYP4B1 variants with p.M226L (Figure 8B). The mutation p.L488P influences the region around 17 Å distance to the heme by widening the channel (Figure 8E). This results in a higher conversion of 2, as well as higher yield of 2b and 2c, that is, however, accompanied by a decrease of selectivity toward the 4′ position (reflected in formation of 2b): For wildtype mlCYP4B1 a 2b/2c ratio of 5.8 was calculated, but all mutated variants have a lower ratio (p.M226L = 1.9, p.L488P = 1.3, double mutated variant = 1.9). This decrease in selectivity (reflected by increases in 2c) is probably caused by the extended space in the active site allowing slight movements of 2a.

Four key activity and selectivity drivers in CYP4B1 were identified by a sequence-guided approach, namely, p.113/118, p.138/143, p.221/226, and p.483/488. Albeit these positions are generally conserved, they differ in one or two orthologs. As hypothesized, systematic single mutagenesis of the respective residues in a specific ortholog resulted in most cases in activity and/or selectivity changes consistent with that of the wildtype ortholog from which the amino acid exchange was extracted. Supported by subsequent computational analysis, we found that three of these positions (p.113/118, p.221/226, and p.483/488), which were exclusively found by the sequence-guided approach, are assumingly involved in the tunnel formation in CYP4B1. Therefore, these positions influence substrate channelling towards the heme, as well as orientation of the substrates in the active site.

Through combinatorial mutagenesis of two amino acids in rCYP4B1 (p.I221L and p.P483L) the activity and selectivity of wildtype chCYP4B1 toward 1 could be mirrored; this suggests that these two positions play a crucial role in the evolutionary divergent activity with regard to fatty acid metabolism. However, combinatorial mutagenesis resulted not only in additive, but also in orthologous- and substrate-dependent adverse effects on activity and selectivity, which likely reflects the enzyme function optimized by nature in the context of the individual orthologs evolutionary background.

The conclusions drawn from the data of this approach allowed us to improve the biocatalytic performance of mlCYP4B1 for sequential double hydroxylation of trans -stilbene 2 to achieve almost complete conversion and 85% yield of trans -dihydroxystilbenes 2b and 2c by introduction of only two beneficial mutations elucidated from the sequence of the other orthologs (p.M226L from chCYP4B1 and rhCYP4B1, as well as p.L488P from rCYP4B1).