In total, seven strains of E. coli N4830-1 (F⁻suo thi-1 thr-1 leuB6 lacY1 fhuA21 supE44 rfbD1 mcrA1 his ilv galK8 Δ(hemF-exp) Δ(bio uvrB) [λ ΔBam N⁺cI857 Δ(CroattR)]) harboring pEX-CYT plasmid were kindly provided by Dr. Naheed Kaderbhai. The strains carry different copy numbers of cytb₅ gene ranging from 0 to 6 which are labeled as N0 to N6 (E. coli pλ-0cyt to E. coli pλ-6cyt). The vector pEX-CYT also contains an ampicillin-resistant gene as a selectable marker (Supplementary Figure S1). All the strains were streaked onto LB agar (Formedium, UK) prior to growing overnight at 37°C on 75 μg/mL ampicillin-supplemented LB agar for selection of the desired colonies.

A single colony was inoculated into 5 mL of LB broth supplemented with 75 μg/mL ampicillin in a 50 mL sterile conical tube and incubated overnight at 37°C with 200 rpm shaking. The inoculum was then diluted with the medium to a final OD_600nm of 0.1 prior to transferring as 200 μL aliquots to a sterile 100 well-plate (10 replicates). Bacterial growth was then monitored with a Bioscreen C analyzer (Oy Growth Curves Ab Ltd, Finland) using the following settings: 5 min preheating, continuous medium shaking, OD_600nm measurement with 10 min intervals, and incubation at 37°C for 24 h.

The cells were cultured for the Cyt b₅ production following the procedure described by Kaderbhai et al. (1992). Briefly, 1 mL overnight culture was inoculated into a 250 mL Erlenmeyer flask containing 25 mL of LB medium plus 75 μg/mL ampicillin. The cells were incubated at 150 rpm and 30°C until reaching the mid-log phase (~2.5 h), then the incubation temperature was immediately increased to 38.5°C to activate the heat-sensitive promoter, λP_L. The culture was induced and incubated for 8 h to allow for the protein production, while a control set was kept uninduced and continuously incubated at 30°C. The final biomass for each sample was determined by measuring the OD_600nm afterward. To optimize the production condition, the bacteria were also cultured as 3 mL aliquots in 15 mL tubes for various lengths of time ranging from 8 to 14 h.

The normalized biomass (~1 mL) of each sample were collected by centrifugation (5,000 g) for 10 min at room temperature and the supernatant was then kept separately to examine any protein leakage. The pellets were then prepared for analysis using sodium dodecyl sulfate-polyacrylamide gel electrophoresis (SDS–PAGE) following the protocol described previously (Kaderbhai et al., 1992). Subsequently, 15 μL of each sample was loaded into sample wells of a 15% polyacrylamide gel and 5 μL of PageRuler plus Prestained Protein Ladder (Thermo Fisher Scientific Inc., USA) was used as a protein marker. Electrophoresis was performed at 180 V for 50 min. The protein bands were subsequently subjected to western blot analysis, and transferred to a polyvinylidene fluoride (PVDF) membrane by applying a constant current at 400 mA for 1 h. The antibodies for the Cyt b₅ target protein and the loading control, β-actin were purchased from Invitrogen (Thermo Fisher Scientific Inc., USA) and the detection reagents from Pierce ECL Western Blotting Substrate (Thermo Fisher Scientific Inc., USA). The target protein bands were detected chemiluminescently using ImageQuant LS 4000 (GE Healthcare UK Ltd, UK).

Overall, the growth behavior of all bacterial strains under 37°C were similar throughout the incubation period, with all strains reaching the stationary phase at around 10 h (Figure 1). At the end of the incubation period (24 h), a decrease in OD_600nm has been observed that may indicate the start of death phase. The bacterial strain N1 showed the lowest OD_600nm value when they reached the stationary phase, while strain N6 displayed as the highest OD_600nm. However, these are insignificantly different as illustrated by the error bars showing the SD of each strain (calculated from 10 biological replicates, n = 10).

All the bacterial strains reached an average final OD_600nm of 3.5. The pH of the media for all samples were specifically monitored to ensure that any changes detected in the FT-IR metabolic fingerprint of the cells were because of Cyt b₅ production and not a result of any environmental changes. The pH values displayed similar changes for all the examined strains and changed from ~6.8 to around 8.2 (average values) after induction (Supplementary Table S1). This alkalinization of the medium is expected due to the cells being cultured in a rich medium where amino acids are catabolized (deamination process), to be used as a carbon source followed by the release of the generated ammonium into the surrounding medium. However, the harvested bacterial cells were visible as pale pink color in some strains (N3–N6) as shown in Supplementary Figure S2.

The FT-IR spectra of induced E. coli N4830-1 strain after pre-processing displayed the commonly known vibrational regions of bacterial cells (Table 1), i.e., the region between 1,200 and 900 cm⁻¹ absorption bands of nucleic acids (PO2-) and bacterial cell wall component such as polysaccharides (C-O-C) (Zarnowiec et al., 2015). The 1,700–1,500 cm⁻¹ region is particularly important as it is usually the most intense region and provides information on the amide I (C=O stretching) and II (C–N stretching, N–H bending) vibrational bands that are the vital chemical characteristic of peptide bonds in proteins (Zarnowiec et al., 2015). The second derivative of the FT-IR spectral data illustrated the underlying peaks as shown in Supplementary Figure S3 that can be assigned to the secondary structure of proteins, i.e., the 1,631 cm⁻¹ region belongs to the β-sheet structure of amide I and 1,658 cm⁻¹ represents the α-helix structure (Holloway and Mantsch, 1989; Kong and Yu, 2007). The vibration in amide II region at about 1,567 cm⁻¹ also relates to the β-sheet structure of the protein (Brian et al., 2014). Amide bands of proteins were also observed in the uninduced set; however, its second derivative spectra (Supplementary Figure S3A) showed similar intensities among all strains, while there was variation in the induced set (Supplementary Figure S3B). This band intensity variation among the 7 strains in Cyt b₅ production condition (38.5°C), especially at 1,658 cm⁻¹ (Supplementary Figure S3C), indicates a difference in the level of the protein content and their structure under the examined condition.

PCA scores plot of FT-IR data from various incubation time (8–14 h) still showed the separation between the control and induced samples, which agreed with the previous results (Figure 2A). The PC–DFA scores plot, using 25 principal components (PCs) accounting for 99.79% of the TEV, in Figure 2B further illustrated the discrimination of 3 main groups, i.e., control samples, induced samples with low (N0–2) and high (N3–6) cytb₅ copy numbers. The control and induced sample groups were also analyzed separately by PC–DFA, using 25 PCs accounting for 99.82 and 99.79% of TEV, respectively. The PC–DFA scores plot of the control set (Figure 2C), only displayed the separation of the samples according to the incubation time, which is because of the biochemical changes that are generally linked to the cells going through the various growth phases. Unlike the control samples, the PC–DFA scores plot (Figure 2D) of the induced samples displayed a gradient of discrimination among all 7 strains from N0 to N6 according to DF1 axis. Moreover, the 8-h incubation set also discriminated from the other incubation times according to DF2 axis, indicating the similarity within the groups of 9–14 h induction. Therefore, 9 h was selected as the optimum-inducing time for all the following experiments to study the metabolic changes resulting from the protein production, as it requires less incubation time and will also be more cost-effective for potentially future large-scale studies.

The production of Cyt b₅ was confirmed in the relevant strains using SDS–PAGE technique. The results showed that there were strong visible bands (13 kDa) in the thermally induced cultures, as expected (Figure 3B), while these were not detected in the control samples (Figure 3A). The detected protein bands were further identified by the western blot analysis (Supplementary Figure S4), which illustrated the western blot result of bacterial samples (N0–N6) with a protein ladder (m) showing bands of Cyt b₅ at the expected position in all 6 producing strains (N1–N6). A non- Cyt b₅ producing strain such as N0 showed no band in that region.

All the bacterial strains illustrated comparable growth profiles at 37°C incubation temperature (Figure 1) and did not show any adverse effects from carrying various copies of the cytb₅ gene. It was noteworthy that the different number of cytb₅ genes seems not to affect their overall growth behavior, and the same was also observed under Cyt b₅ producing condition (Supplementary Table S1). However, the harvested cells illustrated different shades of color where strains N3–N6 showed a visible pale pink color while the other strains did not (Supplementary Figure S2). This is perhaps not only due to the level of Cyt b₅, as was also reported in a previous study (Kaderbhai et al., 1992), but also to do with the limited hemin availability and its incorporation into the protein structure. Therefore, the shades of color from the pellets were not directly indicative of the level of Cyt b₅ produced. Subsequently, a confirmation by alternative methods for relative protein quantification was required. Hence, densitometric data from SDS–PAGE and colorimetric method using the Soret peak was employed, which are discussed in the following sections.

To test the reproducibility of the process, these bacteria were cultured three times in different batches and on three different days. For all batches, the PCA scores plots showed two main clusters according to the conditions, control, and induced groups (Supplementary Figure S5). The PC–DFA scores plots also displayed a clear trend for the induced samples according to the gene copy number from N0–N6, which was confirmed and illustrated in the findings of all the three separate batch cultures, highlighting the reproducibility of the whole process (Supplementary Figure S6). These multivariate analyses also indicated a clear change in the metabolic fingerprints of the control vs. induced cells.

To interpret the chemometrics data further, the loadings plots (Supplementary Figure S7) were investigated. The PC–DFA loadings plot of all samples (Supplementary Figure S7A) including the control and induced sets showed important vibrations that discriminated the samples, these included vibrations in the amide I and II regions (i.e., 1,569 and 1,626 cm⁻¹), lipids (1,744 cm⁻¹), and nucleic acids (1,082 cm⁻¹), which highlighted significant metabolic changes due to the protein production process. While the PC–DFA scores plot of control samples (Figure 2C) displayed a separation based on the incubation time (8–14 h) according to DF1 axis, its loadings plot (Supplementary Figure S7B) highlighted the vibration in amide region, i.e., 1,628 cm⁻¹ and lipids vibration at 1,749 cm⁻¹ (C=O stretching) to be significantly changing. The latter vibration was only detected as significant in the control set (Supplementary Figure S7B) compared with that of the induced group (Supplementary Figure S7C), suggesting these differences to be mainly as a result of the biochemical changes in the cells because of the various incubation lengths and growth phases on bacterial membrane. In addition, the λP_L promoter that is used for induction of this recombinant protein uses a temperature shift from 37 to 38.5°C and it is known that temperature shifts alter membrane fluidity in bacterial cell walls (Mejía et al., 1995) and this is likely to be the cause as to why there are changes in these lipid vibrations. This, we shall explore further when we analyze these bacteria using lipidomics. By contrast, the PC–DFA loadings plot of induced samples (Supplementary Figure S7C) demonstrated that the most significant changes in DF1 axis is accounted for by the variation in amide vibrations relating to the protein secondary structures, such as the vibrational bands at 1,626 and 1,569 cm⁻¹ (β-sheet) and 1,658 cm⁻¹ (α-helix) (Holloway and Mantsch, 1989; Kong and Yu, 2007; Brian et al., 2014). These findings suggested that the 7 bacterial strains examined in this study are mainly discriminated based on their Cyt b₅ content, which is directly correlated to the their cytb₅ gene copy number. To test these findings further, the peak height ratios of the significant amide vibrations, 1,631 cm⁻¹ (β-sheet) and 1,658 cm⁻¹ (α-helix) were calculated and found to be showing an increasing trend from bacterial strains N0–N6, as expected (Supplementary Figure S8). Since the majority structures of the Cyt b₅ are β-sheets (Holloway and Mantsch, 1989), this result also indicated the different level of Cyt b₅ production among the 7 strains.

The first approach was to quantify the Cyt b₅ via western blot analysis, using the β-actin band for normalization. However, the Cyt b₅ was expressed at very high levels, and resulted in saturation of the bands. Several attempts were made to correct this issue, including dilution of the protein extracts and reducing the volume of samples loaded into the gel; however, this also reduced the control bands and compromised the reproducibility of this approach (Supplementary Figure S4). Therefore, to quantify the Cyt b₅ more reliably, the total protein bands in SDS–PAGE was used along with the random ordering of samples to avoid any bias resulting from samples' positions on the gel. With 7 biological replicates, this technique provided an increasing trend of Cyt b₅ content similarly to that from the colorimetric method. As expected, the absorbance at Soret band increased with the concentration of supplemented hemin (0–100 μM) in Cyt b₅-producing strains before reaching the saturation level. This indicated the lack of heme incorporation into the produced protein providing the Soret absorption, especially, in strains with higher copy numbers of the gene. The findings of densitometric analysis from SDS–PAGE and of the Soret absorption using spectroscopy for the relative Cyt b₅ quantification are consistent as suggested in a combination plot (Figure 6A) with R² value of 0.978. Furthermore, a plot of DF1 scores (Figure 6B) from FT-IR fingerprinting data of induced samples also provided a consistent result to those from the wet-lab protein quantitative techniques.