Modeling Human Metabolism of Benzene Following Occupational and Environmental Exposures

We used natural spline (NS) models to investigate nonlinear relationships between levels of benzene metabolites (E,E-muconic acid, S-phenylmercapturic acid, phenol, hydroquinone, and catechol) and benzene exposure among 386 exposed and control workers in Tianjin, China. After adjusting for background levels (estimated from the 60 control subjects with the lowest benzene exposures), expected mean trends of all metabolite levels increased with benzene air concentrations from 0.03 to 88.9 ppm. Molar fractions for phenol, hydroquinone, and E,E-muconic acid changed continuously with increasing air concentrations, suggesting that competing CYP-mediated metabolic pathways favored E,E-muconic acid and hydroquinone below 20 ppm and favored phenol above 20 ppm. Mean trends of dose-specific levels (μmol/L/ppm benzene) of E,E-muconic acid, phenol, hydroquinone, and catechol all decreased with increasing benzene exposure, with an overall 9-fold reduction of total metabolites. Surprisingly, about 90% of the reductions in dose-specific levels occurred below about 3 ppm for each major metabolite. Using generalized linear models with NS–smoothing functions (GLM + NS models), we detected significant effects upon metabolite levels of gender, age, and smoking status. Metabolite levels were about 20% higher in females and decreased between 1% and 2% per year of life. In addition, levels of hydroquinone and catechol were greater in smoking subjects. Overall, our results indicate that benzene metabolism is highly nonlinear with increasing benzene exposure above 0.03 ppm, and that current human toxicokinetic models do not accurately predict benzene metabolism below 3 ppm. Our results also suggest that GLM + NS models are ideal for evaluating nonlinear relationships between environmental exposures and levels of human biomarkers. (Cancer Epidemiol Biomarkers Prev 2006;15(11):2246–52)

Benzene is an important industrial chemical that is also ubiquitous in the environment due to emissions from gasoline and combustion of hydrocarbons and tobacco (1, 2). Occupational exposure to benzene can cause blood disorders, including aplastic anemia, myelodisplastic syndrome, and acute myelogenous leukemia (3, 4). Significant decreases in the numbers of WBC and platelets have recently been reported in workers exposed to <1 ppm benzene (5). These toxic effects are thought to arise from metabolism of benzene, which proceeds along several lines, as illustrated in Fig. 1. Of the various metabolites, 1,4-benzoquinone and the muconaldehydes are regarded as the most toxic species. However, the mechanism by which benzene causes toxicity and the shape of the exposure-response relationship are not well understood (6-8).

We recently reported dose-specific urine concentrations of the major urinary metabolites of benzene (i.e., phenol, catechol, hydroquinone, and E,E-muconic acid) and a minor metabolite [S-phenylmercapturic acid (SPMA)] in 250 benzene-exposed and 139 control workers from Tianjin, China (9). After grouping subjects according to their benzene exposures (30 subjects per group), median metabolite levels increased nonlinearly with increasing median benzene concentrations between 0.03 and 20 ppm, whereas median dose-specific levels of total metabolites (μmol/L/ppm benzene) decreased about 10-fold.

We sought a parsimonious statistical model with which to elaborate on our previous grouped analyses (9) and to determine effects of significant covariates, such as gender, age, and smoking status, on the levels of benzene metabolites. Given the nonlinear relationships involved, we selected NS as basis functions for these models because they use standard (least-squares or maximum-likelihood) methods for estimating variables and for conducting formal tests; they can be used to represent predictors in final models; and they can easily be added to generalized linear models (GLM) for considering covariate effects (10-13). Although GLM + NS models have been used in time-series studies of health effects associated with community air pollution (14), we could find no reports of their applications to characterize exposure-biomarker relationships.

Exposed and control subjects, from two shoe-making factories and three clothes-manufacturing factories, respectively, in Tianjin, China, were recruited with informed consent as described previously (5, 9, 15). Exposed and control subjects were frequency matched by gender. After excluding three control subjects, who had missing values of at least one metabolite, the samples included 250 exposed subjects and 136 control subjects. Table 1 shows summary statistics for the gender, age, and smoking status of participants. Demographic data were obtained by questionnaires at the time of recruitment.

Methods for sampling air and urine were also previously described (5, 9, 15). Briefly, personal full-shift air measurements were matched with post-shift urine samples from exposed and control workers. Of the 386 subjects in this analysis, 139 had repeated measurements of air and urine, making a total of 617 matched air/urine samples. Among subjects with repeated measurements, the median number of paired air and urine samples was three (range, 2-4).

This study was approved by the Institutional Review Boards of the University of North Carolina, the University of California, Berkeley, the U.S. National Cancer Institute, and the Chinese Academy of Preventive Medicine.

The methods of measuring analytes in air and urine were described previously (9, 15). Briefly, benzene and toluene were measured in air using passive personal monitors (Organic Vapor Monitors, 3M, St. Paul, MN) followed by solvent desorption and gas chromatography (15). Urinary benzene was determined by gas chromatography-mass spectrometry using head-space solid-phase microextraction according to Waidyanatha et al. (16). Urinary phenol, catechol, hydroquinone, E,E-muconic acid, and SPMA were measured as trimethylsilyl ether derivatives by gas chromatography-mass spectrometry according to Waidyanatha et al. (17). Quantitation of all urinary analytes was based on peak areas relative to the corresponding isotopically labeled internal standards.

All air samples from control subjects were below the nominal limits of detection of 0.2 ppm for benzene and 0.3 ppm for toluene. In addition, some air measurements from exposed subjects were below the limits of detection (n = 70 for benzene and 67 for toluene) or were missing (n = 23). Air concentrations for these samples were predicted from the simple linear regression of levels of urinary benzene or toluene on the corresponding air levels (in log scale) as described previously for benzene (9). The minor metabolite SPMA was not detected in 30 urine specimens; a value of the limit of detection divided by the square root of two = 0.591 nmol/L was imputed to these samples (18).

For subjects with multiple measurements, the estimated geometric mean air and urine concentrations were used in all statistical analyses.

Relationships between levels of the urinary metabolites and the corresponding air concentrations of benzene were examined using NS models with 6 knots. Because our analyses were done with the (natural) log-transformed air and urine levels, knots represent joints of (logged) air levels of benzene, showing different polynomial trends. They were assigned using equally spaced quantiles of the observations (10). We found that 6-knot models were optimal for these analyses after evaluating preliminary models with 3 to 7 knots. To reduce dimensionality, insignificant knots (P > 0.10) were removed by stepwise elimination (19). Each NS model had the form,

where E[ln(Y_m,j)|ln(X_j)] is the conditional mean of ln(Y_m,j) representing the logged level of the mth metabolite (μmol/L) in the jth subject exposed to benzene at level ln(X_j) (ppm), and ξ_i is the location of the ith knot (in log-scale of benzene exposure). The function

equals

\([\mathrm{ln}(X_{j}){-}{\xi}_{i}]^{3}\)

for positive values and equals zero otherwise. Because metabolite levels from the 60 subjects with the lowest benzene exposures were used to estimate background levels of these metabolites (described later), Model A was applied to the remaining 326 subjects for each benzene metabolite.

Effects of covariates on metabolite levels were determined using GLM + NS models having the form:

where i^* indicates the ith knot in the final NS model for the mth metabolite (m = 1, …, 5, representing E,E-muconic acid, SPMA, phenol, catechol, and hydroquinone, respectively), C_kj is the value of the kth covariate (k = 1, …, K) in the jth subject, and the remaining terms were the same as for Model A. The following covariates were evaluated: gender (0, female; 1, male), age (centered around the estimated mean of 30.2 years, n = 326), smoking status (0, nonsmoker; 1, smoker), body mass index (centered around the estimated mean of 22.5 kg/m², n = 325), co-exposure to toluene [0, low exposure relative to the median concentration of 3.29 ppm (n = 326); 1, high exposure], antibiotics used within 30 days (0, no; 1, yes), and current alcohol consumption status (0, no; 1, yes). Main effects and two-way interactions were evaluated using Proc GLMSELECT of SAS, with backward selection based upon the smallest values of AICc (20), while retaining gender, age, body mass index, and smoking status in all models. With these main effects in the model, no other covariate effects or interactions were retained in final models.

All statistical analyses were done using SAS software for Windows v. 9.12 (SAS Institute, Cary, NC).

Let Y_m,j|Xj be the conditional mean value of Y_m,j, representing the mth metabolite level in the jth subject, given exposure level X_j under Model A. The molar fraction of the mth metabolite, derived from benzene exposure of the jth subject, was estimated as

⁠, where Y_m,b is the background level of the mth metabolite, and the denominator term represents “total” metabolites from benzene. We assigned values to Y_m,b using the median levels of Y_m,j observed in the 60 control subjects with the lowest benzene exposures (E,E-muconic acid, 1.03 μmol/L; SPMA, 0.002 μmol/L; phenol, 54.4 μmol/L; catechol, 11.7 μmol/L; hydroquinone, 6.43 μmol/L; ref. 9). Dose-specific production (μmol/L/ppm benzene) of the mth metabolite in the jth subject was estimated as Y_m,j|Xj−Y_m,b/X_j where X_j is that subject's benzene exposure. Negative values of (Y_m,j|Xj−Y_m,b) in the above computations were replaced by zeros. Because the proportions of negative values increased rapidly with decreasing exposure levels below 0.03 ppm, molar fractions and dose-specific metabolite levels were only evaluated for subjects exposed to benzene at or above 0.03 ppm (n = 267). The following percentages of negative values were observed between 0.03 and 88.9 ppm: E,E-muconic acid, 0.29% and hydroquinone, 14.6%.

Uncertainties in the model predictions of dose-specific metabolite levels were evaluated via bootstrap resampling with 500 iterations (implemented with the SAS macro, %boot). The pool of all observed benzene exposures (each representing a different subject, n = 386) was sampled, with replacement, to select a reference group (the 60 lowest observations) and an exposed group (the 326 remaining observations). Data from the exposed group (n = 326 observations) were then used to construct NS models for the various metabolites, as described above for the original data set.

After removal of nonsignificant terms from Model A, the following final versions of Model A were selected for the five metabolites (values of Y are in μmol/L, whereas those of X are in ppm):

where ξ₁ = ln(0.004 ppm), ξ₂ = ln(0.040 ppm), ξ₃ = ln(0.513 ppm), ξ₄ = ln(1.05 ppm), ξ₅ = ln(2.37 ppm), and ξ₆ = ln(13.8 ppm). These models are shown in Fig. 2A to E along with the corresponding 95% confidence intervals and the individual observations for the 326 subjects.

Effects of covariates, determined under Model B after adjustment for benzene exposure, are summarized in Table 2. Age and/or gender were important explanatory variables for all five metabolites, with males and older subjects typically having lower metabolite levels. Smokers had significantly higher levels of catechol and hydroquinone, whereas lean subjects (lower body mass index values) had significantly higher levels of catechol. Co-exposure to toluene was not significant in any of the models. Likewise, alcohol consumption was not a significant predictor of any benzene metabolite, either in all subjects (n = 326, including 95 current drinkers) or in male subjects (n = 118, including 88 current drinkers).

The background-adjusted urine concentration of the mth benzene metabolite in the jth subject [i.e., (Y_m,j|Xj−Y_m,b) increased monotonically with benzene exposure from 0.03 to 88.9 ppm, as shown in Fig. 3A. The corresponding molar fractions [i.e.,

] are plotted versus benzene exposure in Fig. 3B. Molar fractions for E,E-muconic acid and hydroquinone increased 3-fold (4-12%) and 4-fold (3-11%), respectively, with benzene exposure from 0.03 to 20 ppm, and then decreased with benzene exposure above 20 ppm. The molar fraction of phenol showed the opposite behavior, with reduction from about 88% to 68%, for benzene exposures between 0.03 and 20 ppm, followed by increases above 20 ppm. Molar fractions for catechol and SPMA remained fairly constant, at about 6% and <1%, respectively, over the whole range of benzene levels, with some upwards curvature above 20 ppm.

As shown in Fig. 4, dose-specific metabolite levels [values of (Y_m,j|Xj−Y_m,b)/X_j], decreased with increasing benzene exposure for all but the minor product SPMA. At air concentrations above 20 ppm, shifts in dose-specific metabolism are apparent, with decreasing values for E,E-muconic acid and hydroquinone and increasing values for phenol, catechol, and SPMA. Uncertainties in model predictions of dose-specific metabolite levels are indicated in Fig. 4 by spaghetti plots from the 500 realizations of bootstrap resampling, along with interquartile ranges, and 95% confidence intervals. Although the 95% confidence intervals tended to be large for phenol, hydroquinone, and catechol at low benzene exposures (<0.1 ppm), interquartile ranges were very modest for all metabolites over the entire range of predicted benzene exposures (0.03-88.9 ppm). In addition, the median values from bootstrap analyses were very close to predictions from the NS models applied to the original 386 subjects in our study.

This study of 386 workers in Tianjin, China represents the most extensive set of measurements reported to date for paired air and urine samples from benzene-exposed workers and matched controls. We previously published the empirical relationships between urinary metabolite levels and benzene exposure, based upon groups of these same subjects who had been aggregated by their exposure levels (9). Those preliminary analyses defined crude shapes of the exposure-biomarker relationships but did not permit expected metabolite levels to be predicted at given air concentrations of benzene, nor did they allow effects of age, gender, and other covariates to be estimated, after adjusting for benzene exposure. In the current study, we found that GLM + NS models were ideal for characterizing the continuous relationships between metabolite levels and benzene exposures (Fig. 2) and for testing effects of demographic factors (Table 2). Moreover, using GLM + NS models, we avoided problems that have plagued generalized additive models that apply backfitting algorithms (14, 21, 22).

Although various spline regression models have been applied to investigate covariates in time-to-health effects or survival analyses (11-14, 23-35), we are unaware of any such applications involving human metabolism or exposure-biomarker relationships. Because GLM + NS models are flexible and robust, they offer great promise as tools for investigating human metabolites and other biomarkers of environmental toxicants.

We used GLM + NS models to test for effects upon metabolite levels of several demographic factors (Table 2). Of the covariates considered, only age and/or gender consistently showed significant effects among the five metabolites, after adjustment for benzene exposure. Male subjects had exposure-adjusted metabolite levels that were about 20% lower than those of females. For example, males had 100 × (1 − e^−0.209) = 18.9% less hydroquinone than females at a given exposure level. Likewise, levels of four metabolites (E,E-muconic acid, SPMA, phenol, and hydroquinone) diminished with age, at rates between 1.1% and 1.9% per year of life (Table 2). This range is the same as that reported for albumin adducts of benzene oxide and 1,4-benzoquinone in a different sample of Chinese workers exposed to benzene (36). Although we expected that co-exposure to toluene (a competitive inhibitor of benzene for CYP metabolism) would affect metabolite levels, this was not the case, possibly due to the relatively low toluene concentrations in our study (median, 3.29 ppm; 10th-90th percentile range, 0.012-21.9 ppm; n = 326).

After adjusting for benzene exposure, smoking subjects had about 40% higher levels of hydroquinone and catechol than nonsmokers (Table 2). Because significant smoking effects were not observed for E,E-muconic acid, SPMA, and phenol, we attribute this result to the uptake of hydroquinone and catechol per se from cigarette smoke (37-40). To quantify the contributions of hydroquinone and catechol derived per cigarette, we regressed the logged levels of hydroquinone and catechol on self-reported smoking frequencies in 131 control subjects, who provided this information. This resulted in the following relationships: ln(hydroquinone, μmol/L) = 1.79 + 0.021 (cigarettes per day; P < 0.01) and ln(catechol, μmol/L) = 2.41 + 0.011 (cigarettes per day; P = 0.19). Based upon these models, smoking 20 cigarettes would result in 52% more hydroquinone [i.e., 100 × (1 − e^{(0.021 × 20)})%] and 20% more catechol than observed in nonsmoking control subjects.

Although background-adjusted levels of all metabolites increased monotonically with benzene exposures up to about 30 ppm (Fig. 3A), the molar fractions for phenol, hydroquinone, and E,E-muconic acid changed continuously with increasing air concentrations (Fig. 3B), whereas those for catechol and SPMA remained relatively constant. This indicates that the competing CYP-mediated pathways (Fig. 1) were sensitive to the air levels of benzene inhaled by these subjects. Below 20 ppm, molar fractions of hydroquinone and E,E-muconic acid increased with exposure, whereas those of phenol decreased with exposure; above 20 ppm, the opposite behavior was observed. Because production of hydroquinone and E,E-muconic acid was preferred to that of phenol (below 20 ppm), we infer that phenol and oxepin were either higher-affinity substrates than benzene for the particular CYP enzymes or were more accessible to these enzymes. This conjecture is supported by studies showing that K_Ms for CYP-mediated metabolism of phenol and oxepin were smaller than those of benzene in tissues from humans and/or animals (41-44). Above 20 ppm, the second CYP oxidation steps, leading to hydroquinone and E,E-muconic acid, seem to have become increasingly saturated; this led to the buildup of phenol and, to lesser extents, of catechol and SPMA (other products of a single CYP-oxidation step; ref. 45).

Using mean trends from the NS models (Fig. 4), we investigated dose-specific levels (μmol/L/ppm) of the five benzene metabolites and their sum (total metabolites). The spaghetti plots were dense along the observed mean trends, and interquartile ranges were relatively small. Wide 95% confidence bands were observed below about 0.1 ppm, due to the large relative errors from background subtraction in this region, particularly for phenol, hydroquinone, and catechol, which have important dietary and endogenous sources (39, 46, 47). However, given the narrow interquartile ranges, our conclusions regarding the mean trends should be reasonable. Overall, expected median values of dose-specific levels of total metabolites decreased about 9-fold between 0.03 and 88.9 ppm of benzene. For benzene exposures below 20 ppm, the decreasing trends were more pronounced for phenol and catechol (7- to 11-fold) than for hydroquinone and E,E-muconic acid (2- to 3-fold), reflecting the apparent preference for metabolism of the latter metabolites at low exposures (described above). At benzene exposures above 20 ppm, the trends accelerate downwards for E,E-muconic acid and hydroquinone and turn upwards for phenol and catechol, consistent with results from previous investigations of workers heavily exposed to benzene (16, 17, 48). The mean trend for the minor product SPMA increased over the observed range of exposures and, therefore, displayed completely different behavior than those of the major metabolites.

Interestingly, about 90% of the reductions in dose-specific metabolite levels occurred below about 3 ppm for each major product. This behavior was unexpected, given current toxicokinetic models that indicate that saturable benzene metabolism should not be observed below about 3 to 10 ppm in humans (49, 50). Thus, our results suggest that current toxicokinetic models for benzene are not accurate for air concentrations below 3 ppm.

It is also important to point out that human health risks associated with benzene exposure are based upon linear extrapolation from epidemiology studies involving workers exposed, on average, to air concentrations of tens to hundreds of ppm (3, 4). Our results indicate that persons exposed to air concentrations <0.1 ppm metabolize benzene about nine times more efficiently than such heavily exposed workers (see Fig. 4F). Because the toxic effects of benzene are thought to result from metabolism, this suggests that the health risks associated with low and very low benzene exposures can be considerably greater than those currently predicted from occupational studies.