AI-Guided Prediction of Aflatoxin B1 Bioactivation Validated by Toxicity and DNA Adduct Formation in VERO E6 Cells

1. Introduction

VERO-E6 cells are derived from the kidney epithelium of the African green monkey (Chlorocebus sabaeus) and are widely used as a mammalian cell model in virology (Ammerman et al., 2008) [1], vaccine production (Barrett et al., 2009) [2], and toxicological studies (Costa et al., 2016) [3]. ,They provide a useful system for studying how cells respond to chemical toxicants because they grow consistently in culture while maintaining intact apoptotic and stress-response pathways (Osada et al., 2014, Sène et al., 2021) [4,5]. As epithelial cells, they serve as a model for the first point of contact for xenobiotics, making them important for studying chemical uptake, cellular injury, and downstream stress responses. This allows investigation of genotoxic mechanisms as well as broader toxicant-induced cellular effects (Gómez-Lechón et al., 2010) [6].

Studies using VERO-E6 cells have already shown how different toxins affect cellular function. For example, microcystin-LR (MCLR) exposure leads to dose- and time-dependent cytotoxicity, mitochondrial dysfunction, oxidative DNA damage, and activation of ERK1/2 signaling pathways (Alverca et al., 2009, Dias et al., 2010, Menezes et al., 2013) [7,8,9]. In a similar way, the mycotoxin ochratoxin A (OTA) has been reported to increase reactive oxygen species (ROS) levels, induce apoptosis, and cause micronucleus formation in VERO cells. These observations suggest that oxidative stress plays a significant role in its toxicity in this system (Costa et al., 2016; Costa et al., 2015) [3,10].

VERO cells have also been used to investigate the combined effects of aflatoxin B1 (AFB1) and OTA, where non-linear dose–response relationships and enhanced genotoxicity were observed (El Golli-Bennour et al., 2010)[11]. Additionally, cellular assays in VERO-E6 cells have been used to detect industrial and environmental pollutants, such as 2,4,6-trichlorophenol, and to assess leachates from drinking water materials, demonstrating that morphological and membrane integrity changes can serve as sensitive indicators of chemical toxicity (Liao et al., 2010) [12]. More recently, VERO-based cytotoxicity datasets have been used to develop naïve Bayesian models for predicting compound toxicity, highlighting their utility in computational toxicology pipelines (Perryman et al., 2018) [13]. These models have been applied to classify structurally diverse small molecules as toxic or non-toxic and to identify chemical features associated with cytotoxic responses. Together, these studies demonstrate that VERO cells are a useful platform for detecting diverse toxicological endpoints.

While the Chlorocebus sabaeus genome has been sequenced, revealing orthologs of major Homo sapiens (human) cytochrome P450 enzymes involved in xenobiotic metabolism, including CYP1A and CYP3A family members, the metabolic capabilities of VERO cells remain incompletely understood. (Warren et al., 2015) [14]. However, many xenobiotics, including mycotoxins, polyaromatic hydrocarbons, heterocyclic aromatic amines, and small hydrophobic molecules, require metabolic activation by cytochrome P450 (CYP) enzymes to generate reactive intermediates capable of forming DNA adducts and inducing mutagenesis. For example, AFB1 is bioactivated primarily by CYP1A2 and CYP3A family enzymes to form the highly reactive exo-8,9-epoxide (Figure 1), which covalently binds DNA and initiates carcinogenic processes (Kensler et al., 2011, Guengerich et al., 2002) [15,16]. While hepatic systems are typically used to study such metabolic activation, limited evidence suggests that non-hepatic cell lines, including VERO cells, may express certain CYP enzymes (Gómez-Lechón et al., 2010; Bader et al., 2019) [6,17]. However, individual cytochrome P450 enzymes from VERO cells have not been purified or systematically characterized at the biochemical level, and their catalytic activity toward well-characterized human CYP substrates remains largely undefined (Guengerich et al., 2008, Zanger & Schwab., 2013) [18,19]. While in vitro studies have demonstrated active cellular metabolism in VERO cells, including glycolytic and mitochondrial pathways (Petiot et al., 2010) [20], direct evidence of functional CYP-mediated xenobiotic metabolism is limited. Consequently, the extent to which VERO cells support CYP-dependent metabolic activation of pro-carcinogens remains unclear.

Recent advances in AI, driven structural biology, and computational toxicology provide new opportunities to address key gaps, particularly with respect to protein structure, substrate binding, and catalytic positioning within CYP enzymes, thereby enabling an understanding of the structural and mechanistic basis of xenobiotic metabolism in VERO cells. Methods such as AlphaFold are now widely used to predict protein structures with high accuracy, and molecular docking can be used to examine how enzymes and ligands may interact (Jumper et al., 2021) [22]. In the context of toxicology, these approaches can help us understand how compounds like AFB1 might interact with different CYP isoforms, including how the ligand is positioned and how close it comes to the catalytic site. At the same time, these methods have limitations. They do not provide information about gene expression, enzyme levels, or whether the enzymes are active in the cell. Because of this, computational results alone cannot confirm whether a predicted CYP enzyme is functionally involved under biological conditions. For this reason, it is better to use both computational and experimental approaches, particularly in systems like VERO cells where detailed biochemical information is not available. While computational analysis can give us an idea of what might happen, experimental data are still needed to confirm what is occurring in the cell.

In this study, we used a combined experimental and computational approach to evaluate the metabolic potential of VERO-E6 cells. We first quantified the transcriptional expression of CYP genes and then used AI-based structural modeling and molecular docking to examine possible interactions between CYP enzymes and AFB1. To assess these predictions, we measured DNA adduct formation and assessed cytotoxic effects after AFB1 exposure. By combining gene expression data with the modeling results, we were able to better understand how VERO cells metabolize AFB1.

2. Results

2.1. Expression of C. sabaeus CYP3A4, CYP3A5, and CYP3A7 RNA Transcripts

We used RT-qPCR to determine which CYP genes were expressed in C. sabaeus (VERO E6) cells (Figure 2), and CYP expression was compared to the housekeeping gene GAPDH, where the Ct values are inversely proportional to mRNA abundance. The Ct value for GAPDH was around 25. CYP3A5 showed a similar Ct (≈ 25), which indicates that it has relatively high basal expression. CYP3A4 and CYP3A7 gave Ct values of about 27 and 36, suggesting lower expression compared to CYP3A5. Statistical comparison of Ct values relative to GAPDH demonstrated significant differences for CYP3A5 (p < 0.01), CYP3A4 (p < 0.05), and CYP3A7 (p < 0.01), with CYP3A7 showing the largest Ct difference, consistent with comparatively lower transcript abundance. In contrast, CYP1B1 and CYP1A2 had much higher Ct values (≈ 39 and 42), which points to extremely low levels of transcription. We also could not detect CYP2D6, CYP1A1, or CYP2C9 under these conditions, likely because their expression was below the detection limit of the assay. Overall, what we see is that VERO E6 cells do express CYP genes, with CYP3A5 and CYP3A4 showing the strongest expression in our dataset.

2.2. Comparison of the CYP Amino Acid Sequences from H. sapiens to Those of C. sabeus

We observed close similarities between the amino acid sequences of the CYP3A enzymes from C. sabaeus (VERO E6) and those of human versions (Supplementary Figure S1), which were all approximately500 amino acids in length.. Among the individual isoforms, CYP3A4 showed the highest similarity, with about 93 % identity and 97 % similarity between C. sabaeus and human sequences. CYP3A5 and CYP3A7 isoforms were also quite similar, with roughly 91 % identity and 95 % similarity. \ Aligned sequences exhibited at most one gap, indicating few insertions or deletions . Some similarity is due to conservative substitutions, where the amino acids changes do not affect charge or hydrophobicity. The catalytic regions were the most conserved, including key motifs and the substrate recognition sites (SRS1–SRS6), which are important for enzyme function. For example, the heme-binding motif (FGxGPRNCIG), which contains the cysteine that coordinates the heme iron, was conserved in all three CYP3A proteins. In C. sabaeus, the “x” position is threonine (T), while in humans it is serine (S), which is still a conservative change. Other motifs, such as PERF and EXXR, were also present in both species, suggesting that the catalytic core is maintained. On the other hand, most of the differences were seen in the N-terminal region and in some surface-exposed loops. The apparent flexibility of the surface-expose loops suggests that these changes have minimal effects. The sequences are still quite similar overall, and the key catalytic features are present in both species. Based on this, CYP3A enzymes in C. sabaeus and humans likely share functional similarities.

2.3. Conservation of the H. sapiens and C. sabaeus Three-Dimensional CYPs Structures

To determine whether C. sabaeus CYP enzymes have the structural features needed for catalysis, we compared AlphaFold-predicted structures with experimentally determined human CYP structures including CYP3A5 (PDB: 5VEU) (Hsu et al., 2018) [23], CYP3A4 (PDB: 9BV6) (Wang et al., 2025) [24] and CYP3A7 (PDB: 8GK3) (Liu et al., 2023) [25] (Figure 3). All C. sabaeus CYP3A isoforms exhibited the typical cytochrome P450 fold, including the heme-binding catalytic core. Overall, the fold and the position of the heme were similar to those seen in human enzymes, suggesting that the basic catalytic structure is preserved.

Notable differences were observed in substrate access regions, particularly the B–C and F–G loops, which appeared slightly more flexible and extended in C. sabaeus enzymes. Since these loops are known to control substrate entry into the active site, this added flexibility may facilitate ligand entry and accommodation without affecting the core catalytic structure. A comparison of human and C. sabaeus CYP3A isoforms is summarized in Table 1, which highlights a conserved catalytic core and minor differences in loop flexibility and active-site openness. Another difference was observed in the N-terminal region between human and C. sabaeus CYP3A isoforms. This region mainly helps keep the enzyme associated with the endoplasmic reticulum membrane. In many human CYP crystal structures, this part is often missing or only partially resolved because it is flexible and difficult to crystallize. In contrast, the AlphaFold models include the full-length protein, so the N-terminal region looks more extended and flexible. Because of this, its orientation appears different from that seen in human structures. That said, this region is mainly involved in anchoring the protein to the membrane and is not close to the catalytic site. So, any differences here are unlikely to have a major effect on the heme-binding pocket or enzyme activity. From this, it seems that C. sabaeus CYP3A enzymes still retain the key structural features needed for substrate binding and catalysis.

2.4. Multiple Docking Orientations Reveal Accessible Binding Landscapes

We used molecular docking to look at how AFB1 might interact with CYP enzymes. In the simulations, AFB1 did not stay in a single position within the CYP3A active site (Figure 4). Instead, it adopted several different binding orientations. We saw differences.

In how deeply the ligand sat in the pocket, how it was oriented with respect to the heme, and where it was positioned along the access channel. Overall, this suggests that the CYP binding pocket is quite flexible. Across docking runs, AFB1 was observed in both the central catalytic region and the nearby access channel. In some cases, such as CYP3A5 Mode 1 and CYP3A4 Mode 1 (Figure 4), the difuran moiety of AFB1 was oriented toward the heme group. In other poses, including CYP3A4 Mode 2 and CYP3A7 Mode 3 (Figure 4), the coumarin ring was directed toward the catalytic center. These differences indicate that AFB1 can reposition within the active site, allowing distinct parts of the molecule to approach the heme iron during docking.

Analysis of binding energies indicated that AFB1 interacts well with all CYP3A isoforms (Table 2). Among the isoforms, CYP3A4 generally had the lowest binding energies, suggesting that it accommodates AFB1 more tightly. In comparison, CYP3A5 and CYP3A7 also had favorable interactions, but to a lesser extent. These differences are probably due to minor changes in active-site shape, residue composition, and the flexibility of the access channel in each isoform.

A more detailed look at individual docking modes is summarized in Table 3. In the productive poses, AFB1 stabilization within the active site involved several residues located along the catalytic cavity and access channel. For CYP3A4, we kept seeing the same set of residues come up as Phe108, Phe213, Phe215, Phe241, and Phe304. These are mostly aromatic and hydrophobic, and they seem to form a kind of pocket that holds the planar rings of AFB1 in place, probably through hydrophobic and π–π interactions. Leu482 is also part of this region and adds to the hydrophobic surface. Residues like Arg105 and Arg212 are nearby as well and seem to help with how the ligand sits in the pocket. CYP3A5 did not look the same. Thr309 is close to the catalytic cavity and might be involved in positioning the ligand, possibly through polar interactions or hydrogen bonding. Arg105 and Phe213 are still there, but the interactions here seem less dominated by hydrophobic contacts and more by polar ones. In CYP3A7, the favorable poses were mainly associated with arginine residues, such as Arg105 and Arg372, as well as hydrophobic residues, such as Ala305, Ala370, and Leu482. In this case, AFB1 seems to be stabilized more by a mix of electrostatic interactions and hydrophobic contacts than by the strong aromatic interactions we see in CYP3A4.

Importantly, a high binding energy alone does not necessarily indicate that catalysis will occur. Several poses with favorable docking scores placed AFB1 too far from the heme center (Table 3), making them unlikely to support epoxidation. For catalytic activation, the reactive C8–C9 double bond of the difuran ring needs to be positioned close to the heme iron–oxo species. Only those orientations that bring this region near the heme are likely to lead to the formation of the AFB1-8,9-exo-epoxide, which is responsible for DNA adduct formation. The differences in residues involved in stabilizing these poses likely reflect small variations in active-site structure and flexibility among the CYP3A isoforms. Although the overall fold and heme-binding architecture are conserved, even minor differences in side-chain composition can influence how the ligand is positioned (Ekroos & Sjögren., 2006) [26]. As a result, each CYP isoform appears to use a slightly different combination of aromatic, hydrophobic, and polar residues to orient AFB1 for catalysis. The fact that we observe multiple favorable orientations aligns with what is already known about CYP enzymes, which are quite flexible and allow substrates to move within the active site before settling into a productive position. Overall, this suggests that CYP3A enzymes can bind AFB1 in multiple ways. This flexibility likely helps the ligand adjust its position until the reactive difuran region is aligned with the heme catalytic center.

2.5. Catalytically Favorable Binding Orientations Support Epoxidation Potential

Among the docking results, we selected representative poses with both low binding energies and favorable catalytic geometries for closer analysis (Figure 5). For AFB1 to be metabolically activated, the reactive C8–C9 double bond needs to be positioned close to the heme iron so that oxygen transfer can occur and the AFB1-8,9-epoxide can form. Based on this, we examined the docking poses by considering the ligand’s orientation, its distance from the heme, and its fit within the active site. In the selected poses, AFB1 was generally positioned with its difuran ring facing the heme, bringing the C8–C9 bond close enough to potentially support catalysis.

Among the CYP isoforms, CYP3A4 placed AFB1 closest to the heme iron, which likely makes it the most favorable for epoxidation. CYP3A5 and CYP3A7 exhibited slightly larger Fe–C8/C9 distances, but the ligand was still positioned in a way that could support catalytic activation. This suggests that all three enzymes may be capable of forming epoxide.

The active-site views in Figure 5 indicate that these productive orientations are stabilized by interactions between AFB1 and residues lining the catalytic cavity. Aromatic residues appear to form hydrophobic and π–π interactions with the planar aflatoxin structure, helping to hold the ligand in place above the heme. We also noticed that nearby polar and charged residues seem to affect how the ligand sits in the pocket and limit how much it moves. Because of that, the reactive double bond tends to stay pointed toward the catalytic center.

Taken together, this kind of positioning looks suitable for oxygen insertion. The details of the different modes, including Fe–C8/C9 distances and interacting residues, are listed in Table 3. From this, it is clear that catalytic potential is not determined by binding energy alone. For example, some poses with good docking scores, such as CYP3A4 Mode 2 and CYP3A5 Mode (Table 3), placed the ligand away from the heme and are therefore unlikely to support epoxidation at the C8–C9 position. These poses may reflect alternative ways the ligand can bind within the active site, but the orientation does not appear suitable for oxidation at the main reactive site. In contrast, the poses shown in Figure 5 not only have good binding energies but also place the ligand closer to the heme with supportive interactions, keeping the C8–C9 bond near the catalytic center. These poses are therefore more likely to reflect orientations that can support catalysis. The way CYP3A enzymes position AFB1 here also fits with what is already known about their flexibility and ability to work with different substrates. The minor differences in Fe–C8/C9 distances between the isoforms likely come from differences in active-site shape or how flexible the surrounding loops are, which may influence how efficiently the reaction proceeds. Overall, these results suggest that C. sabaeus CYP enzymes are structurally capable of positioning AFB1 in orientations that support formation of the reactive epoxide intermediate responsible for DNA adduct formation.

2.6. DNA Adduct Formation Confirms AFB1 Bioactivation

The most direct evidence of metabolic activation in VERO E6 cells is the formation of AFB1-derived DNA adducts. We quantified AFB1-derived DNA adducts following short-term exposure. Cells were treated with 200 nmol/L AFB1 for 2 h, 4h, and 16 h, and total DNA adducts, including the AFB1 N7Gua and AFB1 Fapy adducts, were quantified (Figure 6C). After 4 h exposure to 200 nmol/L AFB1, total DNA adduct levels increased approximately fourfold relative to 2 h treatment. The maximum adduct burden reached approximately 1.2 adducts per 1000 kb of DNA at 16 h following a single exposure. Statistical analysis (one-way ANOVA followed by Dunnett’s test, p < 0.001 at 4 h and p < 0.05 at 16 h) confirmed significant increases in adduct formation relative to the 2 h reference time point. These results provide direct biochemical evidence that AFB1 is rapidly bioactivated in VERO-E6 cells, resulting in the formation of canonical genotoxic DNA lesions.

2.7. AFB1 Induces Time- and Dose-Dependent Cytotoxicity

Consistent with the DNA adduct data, AFB1 exposure caused pronounced cytotoxicity in VERO-E6 cells. Dose- and time-dependent reductions in viability were observed, with approximately 50 % lethality by 14 d after prolonged exposure to 200 nmol/L AFB1 (Figure 6A-B). The data indicated that AFB1 exposure not only killed cells but also reduced cell number. This delayed cytotoxic response is similar to that reported in hepatocyte-based models, in which toxicity is mainly linked to the accumulation of reactive metabolites rather than to direct chemical damage (Kensler et al., 2011) [15]. We used a two-tailed Welch’s t-test to compare the AFB1-treated groups with their respective DMSO controls at each time point. At 5 nmol/L, AFB1 caused a small but significant decrease in cell numbers at 6 d (p = 0.029) and again at 14 d (p = 0.027), while no clear change was observed at 10 d (p = 0.469). At the higher dose of 100 nmol/L, a reduction in cell number was seen at later time points, with significant decreases at 10 d (p = 1.46 × 10⁻⁴) and 14 d (p = 4.90 × 10⁻⁵). In contrast, there was no noticeable effect at 6 d (p = 0.806). The strongest effects were seen at 200 nmol/L, where cell number was significantly reduced at all time points—6 d (p = 7.65 × 10⁻⁵), 10 d (p = 5.54 × 10⁻⁴), and 14 d (p = 1.71 × 10⁻⁶). A similar pattern was observed for percent viability. At 5 nmol/L AFB1, we did not see much difference from DMSO at 6 d or 10 d, but by 14 d, there was a small decrease in viability (p = 0.006). In contrast, both 100 nmol/L and 200 nmol/L AFB1 reduced viability at all exposure times (p < 0.05). Overall, the effect clearly depends on both dose and time. One thing we noticed is that the cell number started to drop before we saw a clear decrease in viability, which suggests that AFB1 may first slow down cell growth and then lead to more obvious cell death.

3. Discussion

The present study integrates molecular, structural, computational, and biochemical approaches to evaluate the metabolic competence of VERO-E6 cells for AFB1 bioactivation. From our results, these cells express key xenobiotic-metabolizing genes, including CYP3A7, CYP3A4, and CYP3A5. From the predicted structures using AI-based tools, we observed that these enzymes still maintain the typical cytochrome P450 fold and can accommodate AFB1 in orientations that may support catalytic activation. While AI-based modeling can identify binding orientations that are compatible with the generation of the reactive exo-8,9-epoxide intermediate, it cannot predict the formation of DNA adducts since the epoxide can be catalyzed to less reactive intermediates if additional detoxification enzymes are present. The formation of the reactive intermediate was supported by the experimental detection of AFB1-derived DNA adducts, which provides indirect but definitive evidence that the reactive epoxide intermediate was formed in VERO-E6 cells. This study indicates that combining AI-based characterization of CYP-AFB1 interactions within silico toxicology and experimental validation provides an effective approach for determining xenobiotic metabolism in systems where direct enzymatic characterization is limited. Based on our results, VERO-E6 cells can support CYP-mediated activation of AFB1 and may serve as a useful model for studying xenobiotic metabolism. To our knowledge, this is the first study that combines AI-based structural prediction with experimental detection of AFB1-derived DNA adducts in VERO-E6 cells.

Here, we used AI-based structural modeling and molecular docking to explore CYP-mediated metabolic activation, especially since detailed biochemical data for this system are limited. In VERO-E6 cells, although CYP gene expression can be detected, the lack of purified enzymes and limited functional data make it difficult to directly study enzyme–substrate interactions using conventional methods. Using AI-predicted protein structures, we were able to model the three-dimensional organization of VERO CYP enzymes and examine the active site in more detail. With molecular docking, we looked more closely at how AFB1 might interact with these enzymes—specifically how it sits in the binding pocket, which residues are nearby, and how close it gets to the heme iron. For example, we determined whether the reactive C8–C9 double bond of AFB1 could be positioned within the substrate binding site in a manner that would preferentially favor epoxide formation.. At the same time, there are some limitations with this approach. In most docking methods, the protein is treated as rigid, which does not reflect all dynamic structures (Lee et al., 2025) [27]. We also noticed that scoring functions can sometimes rank poses as favorable just based on energy, even if those poses are not likely to be relevant for catalysis (Desta et al., 2020) [28]. Because of this, binding affinity alone is not enough to predict metabolic activation.

Another limitation is that these computational methods do not give any information how the biological context of the enzyme can affect activity Zanger & Schwab., 2013) [19]. In cells, CYP activity is influenced by factors that are not accounted for by the folded structure predicted by the primary sequence, including-translational modifications, orientation and of the enzyme within the endoplasmic reticulum membrane, its interaction with cytochrome P450 oxidoreductase (POR), and the local membrane environment. From the structural and docking analysis, CYP3A4 is best suitable to bindAFB1 in the most favorable position for epoxidation. However, additional functional studies would be needed to determine which CYP isoform is mainly responsible for this activity in VERO-E6 cells..

When we compared the human and VERO CYP3A enzymes, we observed high sequence and structural similarity, particularly in the catalytic core and heme-binding region. We also saw those important motifs for enzyme activity, such as the conserved cysteine for heme binding and the substrate recognition sites, were present in both species. This points to the enzymes working in an equivalent way. We also noticed some differences in the N-terminal region, which engages in membrane anchoring and may influence how the enzyme is positioned in the endoplasmic reticulum (Denisov et al., 2005) [29]. We also noticed variations in the B–C and F–G loop regions. These regions are quite flexible and seem to be involved in how substrates get into the active site. Because of that, any differences here are more likely to affect how the substrate enters, rather than the catalytic step itself (Yan & Hirao., 2026, Li et al., 2022, Urban et al., 2018)[30,31,32].

The CYP3A family is important for xenobiotic metabolism, and CYP3A4, in particular, is involved in the metabolism of many commonly used drugs (Zanger & Schwab., 2013) [19]. Based on the level of structural similarity we observed, it seems likely that VERO CYP3A enzymes can also interact with different xenobiotics in a way that is comparable to human CYPs. AI-predicted structures and docking can be used to explore how different compounds, such as mycotoxins, environmental contaminants, pharmaceuticals, and industrial chemicals, might interact with metabolic enzymes (Jumper et al., 2021, Stokes et al., 2020) [22,33]. This makes it possible to quickly screen enzyme–substrate interactions, identify binding modes that could be relevant for catalysis, and decide which compounds are worth testing experimentally.

Going forward, this approach could be improved by including molecular dynamics simulations, using better scoring methods, and combining these results with transcriptomic or proteomic data to get a more complete picture of metabolic capacity.. From the RT-qPCR data, we detected expression ofCYP3A5, CYP3A7 and CYP3A4 but not CYP1A2. This underscores the importance of these CYP3A in AFB1 activation but does not preclude that CYP1A2 activity is important at low AFB1 exposure levels (Guengerich et al., 2002, Gallagher et al., 1996)[16,34]. However, we do not know which additional detoxification enzymes, such as glutathione S-transferases, interact with these CYPs to modulate AFB1 adduct formation in VERO-E6 cells, nor the full range of AFB1 metabolites that are generated that may contribute to toxicity.

4. Conclusions

Here, we showed that AFB1 is metabolically activated in VERO-E6 by detecting AFB1-derived DNA adducts. Although VERO-E6 cells do not fully represent hepatic metabolism, they still serve as a useful and responsive system for studying how cells respond to toxicants. Based on the combination of AI computational modeling of enzyme structure and substrate binding, and CYP gene expression, we suggest that bioactivation is mediated by encoded CYP enzymes. Our results indicate that metabolic activation cannot be explained by enzyme expression or binding affinity alone. Instead, the orientation of the substrate within the catalytic site, especially the positioning of the reactive C8–C9 bond relative to the heme iron, plays a critical role. This approach can also be used for other compounds, including environmental toxicants and pharmaceuticals, to better understand how they are metabolized. This study thus provides a way to study xenobiotic metabolism in cell systems where detailed biochemical information is limited.

5. Materials and Methods

5.1. Chemicals and Media

AFB1 was purchased from Sigma and dissolved in dimethyl sulfoxide (DMSO) at a concentration of 10 mg/mL and stored at 4 °C. The Dulbecco’s Modified Eagle Medium (DMEM) medium was purchased from CellTreat. Fetal bovine serum (FBS) was purchased from MidSci.

5.2. Cell Lines and Toxicity Assays

VERO-E6 cells were purchased from the American Type Culture Collection (ATCC). Cells were cultured in DMEM high glucose, 5 % Fetal Bovine Serum, and penicillin/streptomycin and incubated at 37 °C, 5 % CO₂. To assess the toxicity in VERO E6 cells with different toxicants, 50 000 cells were seeded per well in a 6-well plate and allowed to attach overnight. Cells were then treated with AFB1 (DMSO, 5 nmol/L,100 nmol/L, 200 nmol/L) for 6 days, 10 days, or 14 days, with media changes every third day, and fresh AFB1 at the indicated concentrations was re-added with each media change. After exposure, cells were rinsed with PBS, trypsinized, and counted using a hemocytometer to determine the total cell count. Cell lethality was measured using trypan blue (0.4 %, 1:1), a dye that stains non-viable cells. Each experiment was conducted in triplicate to ensure consistent results. To evaluate the effects of AFB1 concentration and exposure time on cell proliferation and viability, statistical analysis was performed at each time point (6 days, 10 days, and 14 days). AFB1-treated groups (5 nmol/L, 100 nmol/L, and 200 nmol/L) were compared with the corresponding DMSO control. Differences between groups were analyzed using Welch’s two-tailed t-test, which is suitable for comparing independent groups when variances may differ. Statistical significance was considered at p < 0.05. (Motulsky et al., 2018)[35]

5.3. Quantitative RT-PCR

Total RNA was isolated from cells using Trizol reagent (1 mL per 5 × 10⁶ cells) following standard procedures (Chomczynski and Mackey., 1995)[36]. Briefly, phase separation was conducted using chloroform, followed by centrifugation at 11 500 rpm for 15 min at 4 °C. The aqueous phase was collected, and RNA was precipitated using 0.5 mL of isopropanol. The resulting pellet was washed with 75 % ethanol, air-dried, and then resuspended in nuclease-free water. To remove any remaining genomic DNA, the samples were treated with DNase I (Thermo Scientific, Lot# 2893971).

Quantitative RT-PCR was then conducted using the Power SYBR Green RNA-to-CT One-Step Kit (Applied Biosystems, Lot# 1707101) on a Step OnePlus Real-Time PCR System (Applied Biosystems). Each reaction contained 300 ng RNA in a final volume of 20 µL and was performed in 96-well plates according to the manufacturer’s instructions. The cycling program consisted of reverse transcription at 48 °C for 30 min, an initial denaturation at 95 °C for 10 min, followed by 50 cycles of 95 °C for 15 s, 56 °C for 15 s, and 60 °C for 1 min. A final extension was performed at 60 °C for 10 min.

During amplification, SYBR Green fluorescence was measured at each PCR cycle using the StepOnePlus Real-Time PCR System. The cycle threshold (Ct) was defined as the PCR cycle number at which the fluorescence signal from the amplified product crossed a fixed threshold above background fluorescence. Ct values were automatically calculated by the instrument software using the same threshold settings across samples. Lower Ct values indicate higher starting mRNA abundance because fewer amplification cycles are required for the fluorescence signal to reach the threshold. GAPDH was used as the housekeeping reference gene for comparison of relative transcript abundance.

cDNA amplification was performed using gene-specific primers. The following forward (F) and reverse (R) primers were used: CYP3A4 (F: 5′-AAACAAAAGCACCGAGTGGA-3′; R: 5′-GGAGAGTGGTGCTAGTTGTG-3′), CYP3A7 (F: 5′-ACTGTCTACAGCCTGTGCCTGG-3′; R: 5′-TGAGCTCTAGAGAGGGCACTT-3′), CYP1A2 (F: 5′GCACAACAAGGGACACAACA-3′; R: 5′-TCCAGTTGCTGTAGCAGGATG-3′), CYP1B1 (F:5′ATCCCAATTCAAGCGCTCCT-3′; R: 5′-TGGTGCTCATGCTGCGG-3′), CYP3A5 (F: 5′-CAGCCTGGTGCTCCTCTATC-3′; R: 5′-ACCATCTTGCGTTCTCCACAT-3′), The resulting cDNA products are in range of 150 bp to 200 bp for Glyceraldehyde-3- phosphate dehydrogenase (GAPDH), CYP3A5, CYP1B1, CYP3A7, CYP3A4 and 270 bp for CYP1A2.

Statistical analysis of RT-qPCR data was performed using paired two-tailed Student’s t-tests comparing Ct values of CYP genes with the endogenous control GAPDH. Results are presented as mean ± standard deviation from three technical replicates, and p < 0.05 was considered statistically significant.

5.4. Detection and Quantification of Aflatoxin B1 (AFB1)–Derived DNA Adducts

The metabolically active form of AFB1, the AFB1 exo-8,9-epoxide, reacts with the N⁷ group of guanine to form the AFB1 N⁷-Gua adducts (Figure 1). The AFB1 N⁷-Gua adduct is highly unstable and converts to trans-AFB1-FAPY-Gua adducts, with the combined levels representing the total DNA adducts (Dohnal et al., 2014)[37]. To assess metabolic activation of AFB1, we quantified total AFB1-associated DNA adducts in exposed cells. Approximately 2 × 10⁷ VERO-E6 cells were cultured in T75 flasks and exposed to 200 nmol/L AFB1 for 2 h, 4 h, or 16 h. After treatment, cells were rinsed twice with phosphate-buffered saline (PBS), detached using trypsin, and washed again with PBS. Genomic DNA was purified using a magnetic bead-based extraction procedure according to the manufacturer’s recommendations (Gen Find V3 mammalian DNA extraction kit). DNA was precipitated with 0.1 mol/L ammonium acetate in 70 % ethanol and resuspended in ethanol for concentration determination. The isolated DNA was then subjected to acid hydrolysis, and AFB1-associated adducts were quantified by liquid chromatography–tandem mass spectrometry as previously reported by (Pawel et al., 2023)[38]. Data represents replicate measurements from independent biological experiments. For time-course studies, total DNA adduct levels were analyzed using one-way analysis of variance (ANOVA) to assess the effect of time. When significance was detected (p < 0.05), Dunnett’s post hoc test was performed to compare each time point with the 2 h reference condition. (Zar et al., 2010, Montgomery et al., 2019)[39,40].

5.5. Computational Modeling and Docking Methodology

Protein structures for CYP3A4, CYP3A5, AND CYP3A7 from C. sabeus were predicted using Alphafold3 to enable structural analysis in the absence of experimentally resolved VERO CYP structures (Jumper et al., 2021)[21]. We used AlphaFold3 to generate structural models of the VERO CYP enzymes, which we then used for computational analysis. Without these structures, it would have been difficult to do the docking or make sense of the interactions.

These models also helped with the docking in AutoDock Vina, since they gave us a better idea of what the proteins look like. As a result, we were able to obtain more consistent predictions for both binding affinity and ligand orientation within the active site. This also allowed us to look more closely at protein–ligand interactions, including identifying key residues involved in binding and checking how close the ligand is to the heme iron (Fe), which is important for catalytic activity in cytochrome P450 enzymes. Details of the individual steps, including protein preparation, ligand processing, docking setup, and interaction analysis, are described in Sections 2.5.1–2.5.6.

5.5.1. Protein Sequence Retrieval and Structural Prediction

We obtained the amino acid sequences of cytochrome P450 enzymes from C. sabaeus from the NCBI protein database. To compare them with human sequences, we ran pairwise alignments between the C. sabaeus and H. sapiens CYP3A isoforms using BLASTP from the NCBI BLAST suite. This allowed us to look at sequence identity, similarity, and whether there were any insertions or deletions (Altschul et al., 1990)[41].

For the structural predictions, we used the default settings and selected the models with the highest confidence scores (pLDDT) for further analysis. The structures were downloaded into CIF format and then converted to PDB format using Python-powered Molecular olet-modeling system (PyMOL). After that, we looked at the structures to make sure the overall folding and key features appeared reasonable.

5.5.2. Heme Transfer and Template Preparation

We retrieved the experimentally determined crystal structures of the corresponding human CYP enzymes from the Protein Data Bank (PDB) (Berman et al., 2000)[42] and used them as a reference for placing the heme group, since their active sites and iron coordination are well characterized. During preparation, we removed other ligands and solvent molecules, while keeping the heme and the key residues around it. The AlphaFold-predicted VERO CYP structures were then aligned with the human CYP structures in PyMOL using the Cα backbone atoms. This helped bring both structures into the same frame, so we could compare conserved features more easily, including the cysteine residue that coordinates the heme iron. Following alignment, the heme cofactor was transferred from the human CYP structure to the Vero CYP model. The Fe–S distance between the heme iron and coordinating cysteine was measured to confirm proper coordination geometry. The resulting heme-bound structures were saved as docking-ready models.

5.5.3. Protein and Ligand Preparation for Docking

We used a structure-based docking approach to look at how Aflatoxin B1 (AFB1) interacts with CYP3A4, CYP3A5, and CYP3A7. During protein preparation, we removed most of the crystallographic water molecules, except for those close to the heme or involved in ligand interactions. The heme group was kept in place. We also fixed duplicate side chains using PyMOL and MGLTools (Rosignoli & Paiardini., 2022)[43]. After that, polar hydrogens were added, Kollman charges were assigned, and the final structures were converted to Protein Data Bank, Quantum and Torsions (PDBQT) format for docking in AutoDock Vina.

The three-dimensional structure of AFB1 was obtained from the PubChem database. Hydrogen atoms were then added, and Gasteiger charges were assigned. The ligand was then converted to PDBQT format, and rotatable bonds were defined for docking in AutoDock Vina.

5.5.4. Docking Site Definition

For docking, we centered the search space on the heme iron in the active site. We used a cubic grid box of 4 nm × 4 nm × 4 nm (40 Å × 40 Å × 40 Å) so that both the catalytic pocket and nearby regions involved in ligand entry and stabilization were included.

This grid size was chosen so that both the active site and nearby regions, such as subpockets and possible access pathways, were included, without making the calculations unnecessarily large. Using this setup, we could check how the ligand sits in the active site, how it interacts with surrounding residues, and whether its orientation is suitable for catalytic activity.

5.5.5. Docking Protocol

Docking simulations were conducted using AutoDock Vina (Eberhardt et al; 2021) [44]. To get better sampling, we ran the docking multiple times using different random seeds. In each run, up to 20 poses were generated within an energy range of 16.73 kJ/mol (4 kcal/mol). The exhaustiveness parameter was increased to 32 (default = 8) so that the program could explore a wider range of ligand conformations and provide more consistent docking results.

5.5.6. Post-Docking Structural Analysis

We looked at the docking poses in PyMOL to see how the ligand sits in the active site, which residues are around it, and how close the reactive C8–C9 atoms are to the heme iron (Rosignoli & Paiardini., 2022) [43]. The C8 and C9 atoms were marked so we could track their positions, and distances from these atoms to the heme iron were measured to determine whether the orientation could support CYP-mediated activation. We also looked at which residues were within about 0.5 nm (5 Å) of the ligand to better understand the binding environment. The distances between the C8–C9 region and the heme iron were measured and, where needed, averaged so we could compare different poses. Dashed lines were mainly used for visualization, which made it easier to see how the ligand sits relative to the catalytic center. To narrow down the most relevant poses, we considered a combination of docking scores, distance measurements, and the surrounding residues. In practice, we focused on poses that showed both a reasonable binding score and a C8–C9 distance of around 0.6 nm (6 Å) or less from the heme iron, since this distance is generally needed for epoxide formation.