Evaluation of the biodegradation of Alaska North Slope oil in microcosms using the biodegradation model BIOB

We present the details of a numerical model, BIOB that is capable of simulating the biodegradation of oil entrapped in the sediment. The model uses Monod kinetics to simulate the growth of bacteria in the presence of nutrients and the subsequent consumption of hydrocarbons. The model was used to simulate experimental results of Exxon Valdez oil biodegradation in laboratory columns (Venosa et al., 2010). In that study, samples were collected from three different islands: Eleanor Island (EL107), Knight Island (KN114A), and Smith Island (SM006B), and placed in laboratory microcosms for a duration of 168 days to investigate oil bioremediation through natural attenuation and nutrient amendment. The kinetic parameters of the BIOB model were estimated by fitting to the experimental data using a parameter estimation tool based on Genetic Algorithms (GA). The parameter values of EL107 and KN114A were similar whereas those of SM006B were different from the two other sites; in particular biomass growth at SM006B was four times slower than at the other two islands. Grain size analysis from each site revealed that the specific surface area per unit mass of sediment was considerably lower at SM006B, which suggest that the surface area of sediments is a key control parameter for microbial growth in sediments. Comparison of the BIOB results with exponential decay curves fitted to the data indicated that BIOB provided better fit for KN114A and SM006B in nutrient amended treatments, and for EL107 and KN114A in natural attenuation. In particular, BIOB was able to capture the initial slow biodegradation due to the lag phase in microbial growth. Sensitivity analyses revealed that oil biodegradation at all three locations were sensitive to nutrient concentration whereas SM006B was sensitive to initial biomass concentration due to its slow growth rate. Analyses were also performed to compare the half-lives of individual compounds with that of the overall polycyclic aromatic hydrocarbons (PAHs).

Introduction

The Exxon Valdez oil spill (EVOS) of 1989 in Prince William Sound (PWS) resulted in the contamination of approximately 2000 km of shoreline within the Gulf of Alaska (Bragg et al., 1994). The more recent Deepwater Horizon oil spill contaminated over 1600 km of shoreline within the Gulf of Mexico (Barron, 2012; Lubchenco et al., 2012). Removal of oil spilled in the environment requires approaches that are both effective and environmentally safe, and one of the techniques is bioremediation, which relies on augmenting the natural biodegradation rate of oil (Atlas, 1995; Vogel, 1996; Xu and Obbard, 2004; Cunliffe and Kertesz, 2006; Karamalidis et al., 2010).

Bioremediation can be studied by conducting laboratory scale experiments (Boufadel et al., 1999; Sabaté et al., 2004; Chen et al., 2008; Tian et al., 2008; Beolchini et al., 2010; Lors et al., 2012) or field scale experiments. Boufadel et al. (1999) evaluated the nitrate concentration required for maximum biodegradation of n-heptadecane in sand columns with monitored oxygen consumption and carbon dioxide production. They reported that the maximum biodegradation occurs at a concentration of 2.5 mg nitrate-N/L. They also speculated nitrogen recycling by biomass when the influent nitrate concentration was zero. This occurred by lysing (break apart) of cells to provide nutrients that other cells could use for growth and biodegradation. Wrenn et al. (2006) studied the effects of nutrient source and supply in crude oil biodegradation. They observed that the extent and rate of oil biodegradation remain unchanged regardless of whether the nutrient supply was intermittent or continuous. Desai et al. (2008) investigated the biodegradation kinetics of individual components and PAHs and their studies indicated that the prediction of natural or enhanced biodegradation of PAHs cannot be based on single compound kinetics as this assumption would overestimate the rate of disappearance, and does not account for the inhibition due to competitive interactions between the compounds.

Venosa et al. (2010) conducted laboratory experiments to study the biodegradation of 19 year old lingering weathered oil in PWS due to EVOS at three different islands. They observed that the most weathered oil was the most biodegradable. The experimental results obtained from laboratory scale experiments could be used to gain a better understanding of the kinetic processes involved in biodegradation by the use of numerical models.

Numerical models have been routinely used by researchers to predict the results. Essaid et al. (2003) used the US Geological Survey (USGS) solute transport and biodegradation code BIOMOC coupled with the USGS universal inverse modeling code UCODE to quantify the BTEX dissolution and biodegradation at a site located in Bemidji, MN. Vilcáez et al. (2013) developed a numerical model for the biodegradation of oil droplets as opposed to dissolved oil, and it was applied to estimate the time scale of biodegradation of Deepwater Horizon oil spill in Gulf of Mexico. The model revealed that small oil droplets biodegraded faster due to their larger surface area per unit mass. Herold et al. (2011) modeled the enhanced bioremediation of groundwater contaminated by acenaphthene, methylbenzofurans, and dimethylbenzofurans. The calibrated model was used to explore the feasibility and efficiency of remediation scenarios. It was observed from their simulations that even with several simplifications made in conceptualization, they were able to demonstrate the ability of their model to detect key processes needed for an effective remediation scheme. Geng et al. (2013) developed a mathematical model to simulate the biodegradation of residual hydrocarbon in a variably-saturated sand column. They estimated the biodegradation kinetic parameters by fitting the model to experimental data of oxygen, CO₂ and residual mass of heptadecane obtained from two columns. They were also able to predict accurately the biodegradation for three other columns using the same parameters.

The goal of this study is to present a numerical model, BIOB which is capable of simulating the biodegradation of oil entrapped in the sediments, and to estimate its parameters using the experimental results from Venosa et al. (2010). The estimated parameters will be used to predict biodegradation under different environmental conditions. A parameter estimation tool based on Genetic Algorithm (GA) is employed to automatically conducting the fitting. Finally, the sensitivity of the model prediction to the kinetic parameters of the BIOB model is evaluated.

Biodegradation of Exxon Valdez Oil in Microcosms

Venosa et al. (2010) conducted laboratory experiments to test the biodegradability of 19 year old lingering oil due to the 1989 EVOS in Alaska, USA. They assessed whether or not oil weathered more than 70% could undergo further biodegradation. Oil weathering could occur as a result of its physical (e.g., evaporation, dissolution, washout), chemical (e.g., photooxidation), or microbial (e.g., biodegradation) processes. Oil weathering by biodegradation can be distinguished from the other processes by using biomarkers (e.g., hopane) to normalize the measured oil concentration. Biomarkers are compounds that biodegrade very slowly within the time frame of interest and have been used to evaluate the effectiveness of biodegradation (Bragg et al., 1994; Venosa et al., 1996). In Venosa et al. (2010), the PAH concentrations were normalized using hopane. The degree of weathering is calculated by comparing the constitution of the weathered oil with the fresh oil. Samples were collected from three different sites of Eleanor Island (EL107), Knight Island (KN114A), and Smith Island (SM006B) whose respective mass weathering indices (MWI) were 30, 76, and 60%. These samples of oil-contaminated sediment were collected from each site by digging up to 10–40 cm below the surface. The volume of each sample was about 50 L and the collected samples were used for microcosm experiments at the University of Cincinnati. Seawater samples from each island were also collected for natural attenuation experiments.

The collected samples were then used for experiments to study (a) natural attenuation, the biodegradation due to the nutrients naturally present in seawater/sediment and (b) nutrient-amended treatment, where nutrients were added to the samples and the concentrations of the nutrients were maintained by adding 10 mg/L of KNO₃ to ensure no nutrient limitation occurs (Venosa et al., 1996; Boufadel et al., 1999; Du et al., 1999). Each microcosm contained 800 mL volume of mixed sediment weighing approximately 1.5 kg with glass beads of diameter 2 and 3 mm at the top and bottom of the sediment to prevent sediment fines from exiting the unit. A peristaltic pump was used to withdraw seawater from the bottom of the microcosms at a flow rate of 2 mL/min to an overhead reservoir. The overhead reservoir worked by a siphon mechanism and the seawater is released to the natural attenuation microcosms every 4 h. Samples were collected at 0, 14, 28, 56, 112, and 168 days and the temperature during the experimental study was similar to the conditions in PWS which was 15°C.

The sediments from each microcosm were analyzed for hydrocarbons using gas chromatography/mass spectrometry (GC/MS). The Total Kjeldahl Nitrogen (TKN) was measured using the Hach method (Hach et al., 1985). The results obtained from these experimental analyses indicated that there was significant biodegradation in the nutrient-amended samples, even when the MWI was high. The experimental results were fitted by Venosa et al. (2010) using a first-order decay model. Their model showed that the highest biodegradation rate occurred in the most weathered oil (KN114A). The experimental results also showed that natural attenuation (non-nutrient amended microcosms) resulted in significant oil biodegradation although the biodegradation rates were lower compared to the nutrient-amended microcosms. This shows that supplying nutrients to the contaminated sediment increased the biodegradation potential in the microcosms and thereby this could be used to significantly increase the biodegradation rates of PAHs at the contaminated sites.

Methods

Numerical Model

The numerical model, BIOB, uses Monod kinetics to simulate the biodegradation of hydrocarbon entrapped between the sediments. It also simulates the growth of bacterial biomass in the presence of nutrients. The decay of hydrocarbon can be mathematically expressed as follows (Geng et al., 2014):

$\begin{array}{cc}\frac{dS}{dt}=-\frac{\mu }{Y_X}X& (1)\end{array}$

Where S is the concentration of the polycyclic aromatic hydrocarbons (PAHs) (mg S/kg sediment), X is the concentration of the biomass (mg X/kg sediment), Y_x is the biomass yield coefficient for growth on the hydrocarbon (mg X/mg S) and μ is the growth rate of the biomass (day⁻¹) given by Geng et al. (2014):

$\begin{array}{cc}\mu ={\mu }_{\max }\left(1-\frac{X}{X_{\max }}\right)\frac{S}{K_S+S}\frac{N}{K_N+N}& (2)\end{array}$

Where μ_max is the maximum growth rate (day⁻¹), X_max is the maximum allowable microbial concentration (mg X/kg sediment), K_s is the half saturation concentration of PAH (mg S/kg sediment), N is the nitrogen-based nutrient concentration (nitrate+nitrite+ammonia) (mg-N/L of pore water), and K_N is the half-saturation concentration for nitrogen consumption (mg-N/L of pore water). The term $\left(1-\frac{X}{X_{\max }}\right)$ was introduced by Geng et al. (2014) to account for the decrease in biomass accumulation when the microbial concentration approaches its maximum value. Geng et al. (2014) also included an expression for hopane removal due to physical processes in Equation 1. However, it was noted in Venosa et al. (2010) study that hopane concentration remained constant throughout the experiment and hence this term was made equal to zero for the purpose of this study. This suggests that the hydrocarbon decay in the microcosms occurred due to the biodegradation process and not due to abiotic processes such as leaching or washout. BIOB also includes equations to solve for the consumption of nutrients, O₂ and the production of CO₂ levels. Typically, O₂ consumption and CO₂ production could be used to evaluate the microbial growth. Experimentally, the microbial concentration can be estimated using MPN (most probable number) (Alexander, 1965) or qPCR (quantitative polymerase chain reaction) which is used to quantify a specific targeted DNA (Heid et al., 1996; Smith and Osborn, 2009). In the Venosa et al. (2010) study, the nutrients, O₂ and CO₂ levels were not measured in the experiments and they were always present in excess, and hence these equations were not considered in the current model. The complete set of equations used in the BIOB model is presented in Supplementary Material (Appendix A1).

The biomass growth can be expressed mathematically as follows (Geng et al., 2014):

$\begin{array}{cc}\frac{dX}{dt}=\left(\mu -k_d\right)X& (3)\end{array}$

Where k_d is the endogenous biomass decay rate (day⁻¹).

The model can be solved for the concentration of biomass and PAH at different times by numerically solving the Equations (1) through (3) using an ordinary differential equation (ODE) solver. BIOB uses a Runge–Kutta Felhberg adaptive time-stepping scheme to solve the equations (Chapra and Canale, 1998).

Inverse Problem—Parameter Estimation

The Monod kinetic model requires several parameters that define the system. Therefore, we need to estimate these model parameters before it were modeled using BIOB. Since no microbial concentration was presented in the study, the maximum microbial concentration (X_max) and the initial microbial concentration (X) were also estimated. Therefore, a total of seven parameters were estimated using a parameter estimation tool based on a mathematical tool known as GA (Torlapati, 2013). It should be noted that the GA is only used for estimating the different parameters that are being used in BIOB (Equations 1–3). Table 1 shows the seven parameters that were estimated along with their higher and lower bounds. The general procedure of a GA is described below.

TABLE 1

Table 1. Parameters estimated by the model using individual datasets.

The six key steps involved in a traditional GA are: encoding, population generation, selection, crossover, mutation, and termination (Holland, 1975). The GA starts with a randomly-generated initial set of solutions (also known as chromosomes) between the higher and lower bounds set by the user and this is called the initial population. A solution is a randomly generated parameter between the higher and lower bound set by the user that will used to calculate the objective function. These randomly generated solutions or chromosomes are used to calculate the fitness of the population using the objective function. The objective function value is assigned as the fitness for each chromosome and it is used to assess its ability to survive the current generation. For a minimization problem, a lower value of fitness is desirable. Based on this fitness value, two parents are selected using a selection process. The selected parents undergo a crossover, where the genetic information is exchanged between the parents using a crossover function. Since the genetic information is transferred to the subsequent generation of children, it is always preferable to choose individuals with better fitness in the selection process. It is also possible that an offspring generated from the crossover of the parents could undergo a mutation operation governed by a mutation probability. The fitness of the offspring is calculated and is combined with the entire population. The individuals with poor fitness are removed from the population (death) at the end of the generation. There are several strategies available for the discarding bad solutions, and for implementing the process of encoding, selection, crossover, and mutation. Further details about the GA are available in Torlapati (2013). The specific methods used in this study are discussed in section Results.

Parameter Sensitivity

The sensitivity of the model output to the parameter estimates were evaluated using the following approach. The covariance matrix, V_x, of the parameters was evaluated according to the equation (Bard, 1974):

$\begin{array}{cc}{V}_x={\sigma }^2{H}^{-1}& (4)\end{array}$

Where H is the Hessian matrix whose terms are the second derivative of the objective function with respect to the parameters. The variance of errors, σ², can be estimated by Bard (1974):

$\begin{array}{cc}{\sigma }^2=\frac{F}{n-p}& (5)\end{array}$

Where F is the value of the objective function at the optimum, n is the number of observations and p is number of estimated parameters.

Results

Bioavailability of oil to the degrading microbial communities is an important phenomenon since the biodegradation occurs at the oil-water interface for hydrocarbons with low solubility (Johnsen et al., 2005). Since experimental data was used to calibrate the parameters for the model, the rate of biodegradation for oil with low bioavailability will be lower compared to a location with higher bioavailability. The experimental dataset from the nutrient amended experiments was used to estimate the kinetic parameters required for BIOB. For the inverse problem, an initial population size of 48 was used. This means that 48 random solutions were generated for each parameter between the higher and lower bound given in Table 1. This is to ensure that the estimated parameters do not exceed values observed in the literature. The initial population size determines the solution space in which the GA searches for the parameters. A low initial population might not provide a wide range of solutions and we might be stuck in a local minimum whereas a higher population increases the computational time. It was observed from our simulations that the results did not improve beyond the population size of 48 and hence a population size of 48 was used for all parameter estimation simulations. The number of generations used in our parameter estimation model was 300. The initial hopane-normalized PAH concentrations were (Venosa et al., 2010) 72.65, 44.68, and 58.68 mg/mg hopane for EL107, KN114A, and SM006B, respectively. The concentration of nitrogen (N) in Equation 2 was set to 10 mg/L for the complete duration of 168 days since Venosa et al. (2010) study ensured that there was sufficient nutrient concentration available for nutrient-amended experiments. The input data for the parameter estimation model was experimental data obtained from the study. Since the experiments were performed in triplicates, the average value was obtained at each sample measurement and was provided as an input to the GA. The “fitness” of each solution was then evaluated by calculating the error (difference) between the model concentration and the observed (input) experimental concentration that was provided. This error was normalized by the observed concentration at the corresponding time, and the normalized error was squared, added for all the experimental data, and assigned as fitness for that population. The objective of the GA was to minimize this weighted least squares (WLS) by genetic recombination of the solutions generated in the initial population. This was continued for 300 generations and the solution with the best fitness at the end of the simulation was used for the final parameters. All the simulations were performed on Intel Core i7 desktop computer with 8 cores and a total clock speed of 3.40 GHz. The time step used for all the simulations was 3 h.

Initially, the parameter estimation tool was run with the nutrient amended experimental data for all three islands provided as simultaneous input. The estimated parameters over-predicted the biodegradation for the SM006B dataset. Therefore, the parameters were estimated for each island dataset separately, and the fit improved greatly. The parameters estimated for individual island datasets are presented in Table 1. The estimated parameters that were obtained for each island dataset were compared against each other and it was observed that the EL107 and KN114A datasets had similar parameters. Hence, the parameter estimation tool was run again with the experimental data from the both the islands provided as input to obtain the common solution for these two datasets. The estimated parameters are provided in Table 2 along with the parameters estimated for SM006B. Note that the best fit parameters presented in Table 2 also include the standard deviation for each parameter and the details of these calculations are presented in section Sensitivity Analysis. Subsequently, the best-fit parameters obtained for the nutrient amended experiments (Table 2) were used to estimate the nutrient concentration in the natural attenuation experiments. The natural attenuation microcosms were supplied with seawater collected from each site every 4 h to provide reaeration of the seawater and gentle mixing. The purpose of the parameter estimation tool was to estimate the naturally present nutrients (N) in the seawater that were utilized for biodegradation. The concentrations of nitrogen for natural attenuation estimated by the parameter estimation tool were 0.56, 0.79, and 2.21 mg/L for EL107, KN114A, and SM006B, respectively. The concentration of nutrients in SM006B location was about three to four times higher than the other two locations.

TABLE 2

Table 2. Parameters estimated in the model and their best-fit value.

The time-varying hopane-normalized PAH concentrations for three islands (both nutrient amended experiments and natural attenuation) of EL107, KN114A, and SM006B with the parameters obtained from the GA for the nutrient amended treatments are shown in Figures 1–3, respectively. Details about the performance of parameter estimation tool with increasing generations are presented in Supplementary Material (Appendix A2). Comparisons in these Figures were made against the experimental data collected at different days and the exponential fit. The comparisons were also made against the individually estimated parameters for EL107 (Figure 1A) and KN114A (Figure 2A) locations. It can be observed from these Figures that the results from the model simulations predicted the experimental results well. It can also be observed from the Figures that the results from the combined data parameter estimates and the individual island data parameter estimates are similar.

FIGURE 1

Figure 1. Comparison of different model results and the experimental results for Eleanor Island (EL107) (A) Nutrient amended experiments (B) Natural attenuation.

FIGURE 2

Figure 2. Comparison of model results and the experimental results for Knight Island (KN114A) (A) Nutrient amended experiments (B) Natural attenuation.

FIGURE 3

Figure 3. Comparison of model results and the experimental results for Smith Island (SM006B) (A) Nutrient amended experiments (B) Natural attenuation.

It can be observed from Table 2 that the biomass growth rate (μ_max) for SM006B location is four times lower compared to the other two islands. However, the values of biomass decay rate (k_d) and half saturation constant for nitrogen consumption (K_N) are similar for all three locations. This suggests that the biodegradation kinetics was different at the SM006B location. The differences in biodegradation kinetics could be due to the fact that the interfacial area between the oil and water may have been a limiting factor in SM006B (Geng et al., 2013). The composition of the oil could also be different as the presence of heavier PAH compounds could reduce the overall biodegradation rate. To check for the possible limitation of interfacial area, the specific surface area per unit mass was evaluated for each of the islands by using the grain size distribution obtained from the literature. The average diameter of the sediment was 4.2 mm (Li and Boufadel, 2010), 7 mm (Xia et al., 2010), and 9.4 mm (Xia and Boufadel, 2011) for Eleanor Island (EL056C), Knight Island (KN114A), and Smith Island (SM006B), respectively. The grain size distribution of Eleanor Island presented in this study is not from a same beach as was used in the Venosa et al. (2010) study. The value was used as a representative sample for this study as a qualitative analysis. The grain size distribution can be used to calculate the specific surface area per mass of sediment (A_S) as follows (Geng et al., 2013):

$\begin{array}{cc}{A}_S=\frac{6}{\phi d_{avg}{\rho }_{sediment}}& (6)\end{array}$

Where ϕ is the shape factor, d_avg is the average grain size, and ρ_sediment is the density of sediment.

It can be observed from Equation 6 that the specific surface area per unit mass is inversely proportional to average grain size. It was assumed that all the particles are spherical for simplicity and the value of shape factor was taken equal to 1.0. The specific surface area calculated using the above equation for EL107, KN114A, and SM006B was 0.86, 0.52, and 0.39 m²/kg, respectively. This indicates that the specific surface area per mass of sediment was lower in the SM006B location than the other two locations and could be one of the reasons for low interfacial area between oil and water and hence low biodegradation rates. This is reflected in the low values of growth rate (μ_max) and maximum biomass accumulation values (as observed in Table 2) for the SM006B location.

In order to compare the goodness of fit for the exponential fit and BIOB, the objective function (WLS) was evaluated. The comparison of objective function values for BIOB and exponential fit are presented in Table 3. A lower value of the objective function indicates that the model results are closer to the experimental data. Therefore, for nutrient amended experiments, BIOB provided a better fit than the exponential fit for KN114A and SM006B datasets and the exponential fit is marginally better for EL107 dataset. For the natural attenuation experiments, BIOB provided a better fit for EL107 and KN114A whereas the exponential fit provided a marginally better fit for SM006B. BIOB's prediction for SM006B in the nutrient amended experiments is also more intuitive as there is a slight lag before the actual biodegradation begins whereas the exponential fit indicates an immediate biodegradation. This is not accurate because the biomass requires some time to grow before they are able to biodegrade the PAHs and this observation is consistent with the experimental data and hence a better fit than the exponential model as shown by the objective function values.

TABLE 3

Table 3. Comparison of the objective function values between the BIOB and Exponential models.

Sensitivity Analysis

The sensitivity analyses were performed only for the nutrient amended experiments as the only parameter that was estimated in natural attenuation experiments was nutrient concentration and we present the sensitivity to nutrient concentration in a section Sensitivity to Nutrient Concentration. Based on the Equation 5, the variance of errors (σ²) at the optimum was 0.038 for EL107 and KN114A, and 0.0017 for SM006B for the nutrient amended experiments.

The sensitivity of the parameter estimates was checked by varying each parameter between −20 and +50% of its original value while the other parameters were kept the same and the fitness was re-evaluated with the deviated parameters. Figures 4A,B show the results from analyses for EL107 and KN114A, and SM006B, respectively. It can be observed from Figure 4A (EL107 and KN114A) that the model is not sensitive to parameters K_N and the initial biomass concentration (X). The model is sensitive to yield coefficient (Y_X), biomass decay rate (k_d), half saturation constant (K_S), maximum allowable biomass concentration (X_max) and biomass growth rate (μ_max). In Figure 4B (SM006B), the model is not sensitive to the initial biomass concentration (X), slightly sensitive to half saturation constant of nitrogen (K_N), moderately sensitive to the parameter X_max, and Y_X, and highly sensitive to the parameters K_S, μ_max, and k_d. The difference in sensitivities for the two systems as shown in Figures 4A,B indicates that the biodegradation kinetics of SM006B is indeed different from the remaining datasets. The sensitive parameters are different for each system and one needs to pay attention to these differences when employing a field scale bioremediation process. Researchers who use parameter estimation tools to estimate bioremediation parameters should pay attention to the highly sensitive parameters as small variations can change the fitness values significantly. Finally, it can be observed that the minimum value of the fitness function was obtained at 0% deviation and this shows the robustness of the search algorithm.

FIGURE 4

Figure 4. Sensitivity analysis for kinetic parameter estimates for (A) EL107 and KN114A (B) SM006B. Note that the minimum value of the variance of errors is at 0% deviation (i.e., estimated parameters).

To quantify the findings from Figures 4A,B, the covariance matrix, V_x, of the parameters was evaluated according to the Equation 4. To estimate the values of the Hessian matrix (H in Equation 5), a quadratic function can be fitted to the objective function vs. deviation from the best fit parameter (as shown in Figures 4A,B) and the second derivative of this fitted quadratic function was calculated. For the sake of simplicity, only the diagonal values of the Hessian matrix were evaluated for this study. The standard deviation values evaluated from the covariance matrix for each parameter are presented in Table 2 along with the best fit parameters. The value of standard deviation is large for some parameters like initial biomass concentration and half saturation constant for nitrogen consumption. This is due to low curvature of the objective function for these parameters and this indicates the low sensitivity to these parameters. Also the standard deviation values were quite small for the rest of the parameters and there was no overlap between the two systems except for K_N and k_d values. This further validates the hypothesis that the biodegradation kinetics of SM006B was indeed different from the other two islands.

Sensitivity to Nutrient Concentration

Nutrient concentration is an important parameter in biodegradation, as it has been observed by several researchers that the nitrate concentration should be about 2–10 mg/L for near maximum biodegradation rate (Atlas and Bartha, 1972; Atlas, 1981; Boufadel et al., 1999; Du et al., 1999; Wrenn et al., 2006). The nutrient concentration for the base case scenario was 10 mg/L. Simulations were performed with different nutrient concentrations between 0.1 and 10 mg/L for EL107 and KN114A locations whereas for the SM006B location, the nutrient concentration was varied between 1 and 10 mg/L. The results are presented in Figures 5A–C for EL107, KN114A, and SM006B datasets, respectively. It can be observed from the Figures that for all datasets that, as the nutrient concentration decreases, the biodegradation rate continues to decrease. Table 4 reports the predicted PAH concentrations as a percentage of the initial concentration at 168 days for various nutrient concentrations for all sites. As the nutrient concentration increases, the PAH concentration at 168 days becomes smaller. For EL107 and KN117, the concentration was affected by small increase in the nutrient concentrations for values less than 1.0 mg-N/L. Such was not the case for SM006B, where nutrient concentrations less than 1.0 mg-N/L had no effect on oil biodegradation. For this reason, the lowest nutrient value for SM006B in Table 4 was 1.0 mg-N/L, and such a value resulted in 99.2% of the initial PAH concentration remaining at 168 days. This might provide a partial explanation of the slow biodegradation of oil in the beaches of Smith Island, including SM006B.

FIGURE 5

Figure 5. Sensitivity of the numerical model to nutrient concentration for (A) EL107, (B) KN114A, (C) SM006B. Note that SM006B location has slower biomass growth rate and hence it is more sensitive to nutrient concentration.

TABLE 4

Table 4. PAH concentrations as percentage of the initial concentrations for different nutrient concentrations at t = 168 days.

Sensitivity to Initial Biomass Concentration

Different scenarios were chosen in addition to the base case scenario to study the sensitivity of the system to the initial biomass concentration. The initial biomass concentrations for the base case were 1.5E-3, 4.34E-03, and 2.46E-03 mg X/kg of sediments for EL107, KN114A, and SM006B, respectively. The order of magnitude of the initial biomass concentration was increased 10-fold and decreased 100, 10,000-fold for EL107, KN114A, and SM006B datasets. The results from these simulations are presented in Figures 6A–C for EL107, KN114A, and SM006B locations, respectively. It can be observed from the Figure 6A (EL107) and Figure 6B (KN114A) that the model is not very sensitive to the initial biomass concentration after the first few days. This may be due to a high nutrient concentration and a relatively high growth rate of the bacteria in EL107 and KN114A locations. This allows the bacteria to grow rapidly and reach the maximum allowable concentration in the system. As a result, the oil biodegradation is not adversely impacted after the initial phase in EL107 and KN114A locations. The final hydrocarbon concentration after the 168 day experiment was similar for all initial biomass concentrations for both these locations. However, for the SM006B location (Figure 6C), the initial lag before biodegradation continued to increase as the initial biomass concentration was decreased. This is due to the fact that biomass growth rate is slower and hence the biomass takes longer to grow. There was a decrease in the PAH biodegradation by 14 and 20% when the biomass concentration was decreased by 100 and 10,000-fold. There was no difference in the final PAH concentration when the biomass concentration was increased 10-fold. Therefore, a low initial biomass concentration does not affect the overall biodegradation rate in EL107 and KN114A case due to its faster biomass growth rate whereas the initial biomass concentration seems to play a significant role in the SM006B case when it was decreased. This suggests that the biodegradation kinetics was not limited by the initial biomass concentration as the biodegradation happens almost immediately in the case of EL107 and KN114A (Figure 3A). This further supports the hypothesis that the interfacial area between oil and water plays a significant role in bioremediation than initial biomass concentration.

FIGURE 6

Figure 6. Sensitivity of the numerical model to biomass concentration for (A) EL107, (B) KN114A, (C) SM006B. Note that X represents the initial biomass concentration. * in the legend indicates the base case scenario.

Discussion

To compare the biodegradation rate of the polycyclic aromatic compounds (PAHs) with the biodegradation rate of individual compounds, we identified the compound with first-order decay rate similar to that of the PAHs for each location as reported by Venosa et al. (2010). This compound was designated as the representative compound which means that the rate of biodegradation of PAH is representative of the biodegradation rate of the identified compound. The first-order biodegradation rates were 0.017, 0.0294, and 0.0099 day⁻¹ for nutrient amended experiments and 0.0072, 0.0138, and 0.0043 day⁻¹ for natural attenuation experiments of EL107, KN114A, and SM006B, respectively. Half-life (t_1/2) was evaluated from the first order decay rates (k) as follows:

$\begin{array}{cc}{t}_{1/2}=\frac{\ln 2}{k}& (7)\end{array}$

Table 5 shows the names of the representative compounds and their respective half-life values for nutrient amended experiments. The first-order decay rates for the hopane-normalized individual compounds are available in the supplementary information of Venosa et al. (2010) and Equation 7 was used to evaluate the half-life values. These values are reported in Table 5 for nutrient-amended experiments. In addition to the representative compound, the slowest and fastest biodegradable compounds are also identified for each island. Table 5 also shows the fastest and slowest biodegradable compounds for each island along with their half-life values. Table 6 provides the names and half-life values of the representative compounds as well the fastest and slowest compound for the natural attenuation experiments.

TABLE 5

Table 5. Comparison of biodegradation rates for nutrient-amended experiments.

TABLE 6

Table 6. Comparison of biodegradation rates for natural attenuation experiments.

The parameters in Tables 5, 6 are used to plot the temporal variation of the concentration of hydrocarbon for each island. It was assumed that the respective compounds have the same concentration as the initial PAH concentration at that island to allow for easy comparison. The results from the simulations are shown in Figures 7A–C for EL107, KN114A, and SM006B, respectively, for the nutrients amended experiments and Figures 8A–C show simulation results for natural attenuation experiments for EL107, KN114A, and SM006B, respectively. It can be observed from the Figures that the biodegradation rates vary vastly between different compounds.

FIGURE 7

Figure 7. Comparison of the temporal variation for different compounds in nutrient-amended experiments for (A) EL107, (B) KN114A, and (C) SM006B.

FIGURE 8

Figure 8. Comparison of the temporal variation for different compounds in natural attenuation experiments for (A) EL107, (B) KN114A, and (C) SM006B.

We have also compared the number of compounds with biodegradation rate higher and lower than the representative compound for that island. For the nutrient-amended experiments, there were 10 compounds biodegrading at higher rate than the representative compound for all islands whereas there were 17, 11, and 13 compounds biodegrading at lower rate than the representative compound for EL107, KN114A, and SM006B Islands, respectively. Similarly, for the natural attenuation experiments there were 11, 9, and 9 compounds biodegrading at higher rate than the representative compound for EL107, KN114A, and SM006B, respectively, whereas there were 16, 12, and 15 compounds biodegrading at lower rate than the representative compound for EL107, KN114A, and SM006B Islands, respectively.

The weighted rates at each site were evaluated by multiplying the initial concentration of the compound with its rate and this summation was divided by the concentration of the PAH at the site. These weighted rates were used to compute the half-life using Equation 7 and comparison between the half-life computed using exponential curve fitted data and weighted half-lives are shown in Table 7. It can be observed from the Table that the weighted half-lives for each site are lower than the half-lives computed from the exponential curve fitted data. It should also be noted that BIOB is capable of estimating the biodegradation of individual hydrocarbons where information about these compounds such as the concentrations of the hydrocarbon and the microbes capable of biodegrading that hydrocarbon are available. In Equations 1 and 2, the term S will represent the concentration of the hydrocarbon instead of the PAH concentration if we use BIOB to estimate the biodegradation of individual compounds.

TABLE 7

Table 7. Comparison of the half-lives between the weighted and fitted decay rates.

Figures 9A, B report the concentration of total extractable hydrocarbon (TEH) as function of time in the microcosms of EL107 and KN114A, respectively. Simple exponential decays were fitted to illustrate the decreasing behavior. The Natural Attenuation data for EL107 (Figure 9) show a lag in the decrease of TEH with time, similar to that observed for the hopane-normalized PAH concentration for EL107 (Figure 2B), which suggests that the lag is due to biotic factors, namely the initial phase of microbial growth (or microbial acclimation period). For KN114A, the decrease is more or less immediate for both data sets (Figure 10), but a sudden drop in the TEH concentration occurs at 120 days for all microcosms. We have considered and ruled out any possibility for this decrease, which occurred for both Natural attenuation and nutrient amended microcosms. We believe these data points do not represent correct results, but they don't invalidate the overall trend of decrease with time.

FIGURE 9

Figure 9. Decrease of TEH (total extractable hydrocarbon) as function of time in (A) the EL107 and (B) KN114A microcosms along with fitted exponential decay curves. Diamond symbols = Natural attenuation, and square symbols = Nutrient amended.

FIGURE 10

Figure 10. Hopane concentration normalized by the concentration of TEH (total extractable hydrocarbon) for (A) the EL107 and (B) the KN114A microcosms. The overall increase with time reflects biodegradation.

Figure 10 reports the concentration of hopane normalized by the total mass of TEH for EL107 and KN114A. The increase with time reflects the biodegradation of compounds less resistant to biodegradation than hopane. For KN114A, the decrease between 120 and 180 days reflect the low TEH values at 112 days (Figure 9B). The sudden decrease cannot be due to the biodegradation of hopane, because although hopane does biodegrade, its time scale of biodegradation is longer than a few weeks, more on the order of years (Atlas and Bragg, 2007).

Figures 11A, B report the PAH concentration normalized by the concentration of TEH as function of time for (A) EL107 and (B) KN114A. Fitted straight lines were obtained to determine if there is a correlation between the biodegradation of PAH and that of TEH. For EL107, a low R² was noted along with a small slope, which suggests that the biodegradation of PAH is not directly related to the biodegradation of TEH. We think this has to do also with the variability of the TEH at 112 days. However, the high R² for KN114A suggests a strong linear relation between the biodegradation of PAH and that of TEH. At face value, Figure 11 suggests that tracking the biodegradation of PAHs is not sufficient to infer the biodegradation of the rest of the compounds (i.e., TEH).

FIGURE 11

Figure 11. PAH concentration normalized by the concentration of total extractable hydrocarbon (TEH) as function of time for (A) EL107 and (B) KN114A. Fitted straight lines were obtained to determine if there is a correlation between the biodegradation of PAH and that of TEH.

It should also be noted that the seawater that was renewed to the microcosms was recirculated to the reservoir that supplies the seawater to the microcosms. Therefore, in case of oil leeching, it would probably be replenished to the microcosms when the seawater is renewed. The leaching of oil is also not mentioned in Venosa et al. (2010) and the hopane concentration remains fairly constant throughout the length of the experiment which means the removal of oil due to physical processes was limited. Therefore, we did not account for leaching in the model.

Conclusions

In this study, the biodegradation model BIOB was used to simulate the experimental data from nutrient amended and the natural attenuation treatments from Venosa et al. (2010) study. The kinetic parameters required for the BIOB were estimated using a parameter estimation tool based on GA. It was observed from the estimated parameters that the SM006B location had relatively slower biodegradation rate compared to the other two locations. Grain size analysis from each site revealed that the specific surface area per unit mass of sediment was considerably lower at SM006B, which suggests that the surface area of sediments is a key control parameter for microbial growth in sediments. Subsequently, the best-fit parameters obtained from the nutrient amended experiments were used to estimate the nutrient concentration for the natural attenuation experiments. The estimated parameters were able to predict the experimental results well for all the datasets. In comparison with the exponential decay model, the results from BIOB showed better goodness of fit for KN114A and SM006B datasets while the fit for EL107 dataset for nutrient amended experiments was relatively similar for both exponential decay and BIOB. For natural attenuation treatment, the exponential fit provided a better fit for the SM006B dataset whereas BIOB provided a better fit for EL107 and KN114A datasets. In particular, BIOB was able to capture the initial slow biodegradation due to the lag phase in microbial growth. In addition, the exponential fit does not account for the limitation of nutrients and has no predictive capabilities whereas the parameters estimated for BIOB are still valid at different nutrient concentrations. Therefore, BIOB can be used to predict the bioremediation scenarios at the site.

It was observed from our sensitivity analyses that the initial biomass concentration and the half saturation coefficient of nitrate were not sensitive parameters for EL107 and KN114A locations whereas the biomass growth rate, biomass decay rate, yield coefficient, half-saturation constant for hopane-normalized PAH, maximum allowable biomass were the most sensitive parameters. For the SM006B location, the system was highly sensitive to biomass decay rate, biomass growth rate, half saturation constant for hopane-normalized PAH, yield coefficient and slightly sensitive to maximum allowable biomass concentration and half-saturation constant for nitrogen, and insensitive to small changes in initial biomass concentration. In addition, when the sensitivity analyses were performed on the nutrient concentration, SM006B location was more sensitive to small changes than EL107 and KN114A locations due to the slow biomass growth rate. Even though the nutrient concentration was a sensitive parameter, there was no significant difference in concentration trends for the nutrient concentrations of 5 and 10 mg of N/L in EL107 and KN114A locations. Also, the initial biomass concentration had much lower effects on the final PAH concentration after 168 days in EL107 and KN114A locations than in the SM006B location. The high nutrient concentration allowed for a rapid biomass growth and thus, the PAH biodegradation was not substantially affected after the initial few days in the EL107 and KN114A location whereas in the SM006B location, when the initial biomass concentration was decreased by several orders of magnitude, the biomass required several days before the biomass concentration was sufficient enough to cause hydrocarbon decay. Therefore, it can be concluded that the initial biomass concentration affects oil biodegradation significantly during the initial phase of the biodegradation process. Comparisons were also made to study the biodegradation decay rates of individual compounds against the decay rate of the PAH. It was observed that the weighted average of the individual decay rates was greater than the decay rate of the PAH. In conclusion, we determine that the Monod kinetic model, BIOB is a useful tool that can be used to predict the extent of biodegradation in laboratory scale experiments and in marine environments such as beaches where the oil is entrapped within sediments with no significant transport or dissolution.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

This work was supported by Exxon Valdez Oil Spill Trustee Council under Project No. 11100836. However, it doesn't reflect the views of the Council, and no official endorsement should be inferred. Input from Dr. Yves Robert Personna is greatly appreciated.

Supplementary Material

The Supplementary Material for this article can be found online at: http://www.frontiersin.org/journal/10.3389/fmicb.2014.00212/abstract

References

Alexander, M. (1965). “Most-probable-number method for microbial populations,” in Methods of Soil Analysis. Part 2. Chemical and Microbiological Properties, ed A. G. Norman (Madison, WI: American Society of Agronomy, Soil Science Society of America), 1467–1472.

Atlas, R., and Bragg, J. (2007). “Assessing the long-term weathering of petroleum on shorelines: uses of conserved components for calibrating loss and bioremediation potential,” in Arctic and Marine Oilspill Program Technical Seminar (Canada), 263–289.

Atlas, R. M. (1981). Microbial degradation of petroleum hydrocarbons: an environmental perspective. Microbiol. Rev. 45, 180–209.

Pubmed Abstract | Pubmed Full Text

Atlas, R. M. (1995). Petroleum biodegradation and oil spill bioremediation. Mar. Pollut. Bull. 31, 178–182. doi: 10.1016/0025-326X(95)00113-2

CrossRef Full Text

Atlas, R. M., and Bartha, R. (1972). Degradation and mineralization of petroleum in sea water: limitation by nitrogen and phosphorous. Biotechnol. Bioeng. 14, 309–318. doi: 10.1002/bit.260140304

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Bard, Y. (1974). Nonlinear Parameter Estimation. New York, NY: Academic Press.

Barron, M. G. (2012). Ecological impacts of the deepwater horizon oil spill: implications for immunotoxicity. Toxicol. Pathol. 40, 315–320. doi: 10.1177/0192623311428474

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Beolchini, F., Rocchetti, L., Regoli, F., and Dell'Anno, A. (2010). Bioremediation of marine sediments contaminated by hydrocarbons: experimental analysis and kinetic modeling. J. Hazard. Mater. 182, 403–407. doi: 10.1016/j.jhazmat.2010.06.047

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Boufadel, M., Reeser, P., Suidan, M., Wrenn, B., Cheng, J., Du, X., et al. (1999). Optimal nitrate concentration for the biodegradation of n-heptadecane in a variably-saturated sand column. Environ. Technol. 20, 191–199. doi: 10.1080/09593332008616808

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Bragg, J. R., Prince, R. C., Harner, E. J., and Atlas, R. M. (1994). Effectiveness of bioremediation for the Exxon Valdez oil spill. Nature 368, 413–418. doi: 10.1038/368413a0

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Chapra, S. C., and Canale, R. P. (1998). Numerical Methods for Engineers: With Programming and Software Applications. Boston, MA: WCB; McGraw-Hill.

Chen, J., Wong, M., Wong, Y., and Tam, N. F. (2008). Multi-factors on biodegradation kinetics of polycyclic aromatic hydrocarbons (PAHs) by Sphingomonas sp. a bacterial strain isolated from mangrove sediment. Mar. Pollut. Bull. 57, 695–702. doi: 10.1016/j.marpolbul.2008.03.013

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Cunliffe, M., and Kertesz, M. A. (2006). Effect of Sphingobium yanoikuyae B1 inoculation on bacterial community dynamics and polycyclic aromatic hydrocarbon degradation in aged and freshly PAH-contaminated soils. Environ. Pollut. 144, 228–237. doi: 10.1016/j.envpol.2005.12.026

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Desai, A. M., Autenrieth, R. L., Dimitriou-Christidis, P., and McDonald, T. J. (2008). Biodegradation kinetics of select polycyclic aromatic hydrocarbon (PAH) mixtures by Sphingomonas paucimobilis EPA505. Biodegradation 19, 223–233. doi: 10.1007/s10532-007-9129-3

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Du, X., Reeser, P., Suidan, M. T., Huang, T., Moteleb, M., Boufadel, M. C., et al. (1999). “Optimum nitrogen concentration supporting maximum crude oil biodegradation in microcosms,” in International Oil Spill Conference (Seattle, WA: American Petroleum Institute), 485–488. doi: 10.7901/2169-3358-1999-1-485

CrossRef Full Text

Essaid, H. I., and Bekins, B. A. (1997). BIOMOC, a Multispecies Solute-Transport Model With Biodegradation. Menlo Park, CA: US Department of the Interior, US Geological Survey.

Essaid, H. I., Bekins, B. A., Godsy, E. M., Warren, E., Baedecker, M. J., and Cozzarelli, I. M. (1995). Simulation of aerobic and anaerobic biodegradation processes at a crude oil spill site. Water Resour. Res. 31, 3309–3327. doi: 10.1029/95WR02567

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Essaid, H. I., Cozzarelli, I. M., Eganhouse, R. P., Herkelrath, W. N., Bekins, B. A., and Delin, G. N. (2003). Inverse modeling of BTEX dissolution and biodegradation at the Bemidji, MN crude-oil spill site. J. Contam. Hydrol. 67, 269–299. doi: 10.1016/S0169-7722(03)00034-2

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Geng, X., Boufadel, M. C., Personna, Y., Lee, K., and Tsao, D. (2014). BioB: a mathematical model for the biodegradation of low solubility hydrocarbons. Mar. Pollut. Bull. doi: 10.1016/j.marpolbul.2014.04.007. [Epub ahead of print].

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Geng, X., Boufadel, M. C., and Wrenn, B. (2013). Mathematical modeling of the biodegradation of residual hydrocarbon in a variably-saturated sand column. Biodegradation 24, 153–163. doi: 10.1007/s10532-012-9566-5

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Hach, C. C., Brayton, S. V., and Kopelove, A. B. (1985). A powerful Kjeldahl nitrogen method using peroxymonosulfuric acid. J. Agric. Food Chem. 33, 1117–1123. doi: 10.1021/jf00066a025

CrossRef Full Text

Heid, C. A., Stevens, J., Livak, K. J., and Williams, P. M. (1996). Real time quantitative PCR. Genome Res. 6, 986–994. doi: 10.1101/gr.6.10.986

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Herold, M., Greskowiak, J., Ptak, T., and Prommer, H. (2011). Modelling of an enhanced PAH attenuation experiment and associated biogeochemical changes at a former gasworks site in southern Germany. J. Contam. Hydrol. 119, 99–112. doi: 10.1016/j.jconhyd.2010.09.012

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Holland, J. H. (1975). Adaptation in Natural and Artificial Systems: An Introductory Analysis With Applications to Biology, Control, and Artificial Intelligence. Ann Arbor, MI: U Michigan Press.

Johnsen, A. R., Wick, L. Y., and Harms, H. (2005). Principles of microbial PAH-degradation in soil. Environ. Pollut. 133, 71–84. doi: 10.1016/j.envpol.2004.04.015

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Karamalidis, A., Evangelou, A., Karabika, E., Koukkou, A., Drainas, C., and Voudrias, E. (2010). Laboratory scale bioremediation of petroleum-contaminated soil by indigenous microorganisms and added Pseudomonas aeruginosa strain Spet. Bioresour. Technol. 101, 6545–6552. doi: 10.1016/j.biortech.2010.03.055

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Köpke, B., Wilms, R., Engelen, B., Cypionka, H., and Sass, H. (2005). Microbial diversity in coastal subsurface sediments: a cultivation approach using various electron acceptors and substrate gradients. Appl. Environ. Microbiol. 71, 7819–7830. doi: 10.1128/AEM.71.12.7819-7830.2005

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Li, H., and Boufadel, M. C. (2010). Long-term persistence of oil from the Exxon Valdez spill in two-layer beaches. Nat. Geosci. 3, 96–99. doi: 10.1038/ngeo749

CrossRef Full Text

Lors, C., Damidot, D., Ponge, J.-F., and Périé, F. (2012). Comparison of a bioremediation process of PAHs in a PAH-contaminated soil at field and laboratory scales. Environ. Pollut. 165, 11–17. doi: 10.1016/j.envpol.2012.02.004

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Lubchenco, J., McNutt, M. K., Dreyfus, G., Murawski, S. A., Kennedy, D. M., Anastas, P. T., et al. (2012). Science in support of the deepwater horizon response. Proc. Natl. Acad. Sci. 109, 20212–20221. doi: 10.1073/pnas.1204729109

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Nicol, J. P., Wise, W. R., Molz, F. J., and Benefield, L. D. (1994). Modeling biodegradation of residual petroleum in a saturated porous column. Water Resour. Res. 30, 3313–3325. doi: 10.1029/94WR01879

CrossRef Full Text

Sabaté, J., Vinas, M., and Solanas, A. (2004). Laboratory-scale bioremediation experiments on hydrocarbon-contaminated soils. Int. Biodeterior. Biodegrad. 54, 19–25. doi: 10.1016/j.ibiod.2003.12.002

CrossRef Full Text

Smith, C. J., and Osborn, A. M. (2009). Advantages and limitations of quantitative PCR (Q-PCR)-based approaches in microbial ecology. FEMS Microbiol. Ecol. 67, 6–20. doi: 10.1111/j.1574-6941.2008.00629.x

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Tian, Y., Liu, H. J., Zheng, T. L., Kwon, K. K., Kim, S. J., and Yan, C. L. (2008). PAHs contamination and bacterial communities in mangrove surface sediments of the Jiulong River Estuary, China. Mar. Pollut. Bull. 57, 707–715. doi: 10.1016/j.marpolbul.2008.03.011

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Torlapati, J. (2013). Development and Application of One Dimensional Multi-component Reactive Transport Models. Alabama: Auburn University.

Venosa, A. D., Campo, P., and Suidan, M. T. (2010). Biodegradability of lingering crude oil 19 years after the Exxon Valdez oil spill. Environ. Sci. Technol. 44, 7613–7621. doi: 10.1021/es101042h

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Venosa, A. D., Suidan, M. T., Wrenn, B. A., Strohmeier, K. L., Haines, J. R., Eberhart, B. L., et al. (1996). Bioremediation of an experimental oil spill on the shoreline of Delaware Bay. Environ. Sci. Technol. 30, 1764–1775. doi: 10.1021/es950754r

CrossRef Full Text

Vilcáez, J., Li, L., and Hubbard, S. S. (2013). A new model for the biodegradation kinetics of oil droplets: application to the Deepwater Horizon oil spill in the Gulf of Mexico. Geochem. Trans. 14, 4. doi: 10.1186/1467-4866-14-4

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Vogel, T. M. (1996). Bioaugmentation as a soil bioremediation approach. Curr. Opin. Biotechnol. 7, 311–316. doi: 10.1016/S0958-1669(96)80036-X

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Wrenn, B., Sarnecki, K., Kohar, E., Lee, K., and Venosa, A. (2006). Effects of nutrient source and supply on crude oil biodegradation in continuous-flow beach microcosms. J. Environ. Eng. 132, 75–84. doi: 10.1061/(ASCE)0733-9372(2006)132:1(75)

CrossRef Full Text

Xia, Y., and Boufadel, M. C. (2011). Beach geomorphic factors for the persistence of subsurface oil from the Exxon Valdez spill in Alaska. Environ. Monit. Assess. 183, 5–21. doi: 10.1007/s10661-011-1902-4

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Xia, Y., Li, H., Boufadel, M. C., and Sharifi, Y. (2010). Hydrodynamic factors affecting the persistence of the Exxon Valdez oil in a shallow bedrock beach. Water Resour. Res. 46, W10528. doi: 10.1029/2010WR009179

CrossRef Full Text

Xu, R., and Obbard, J. P. (2004). Biodegradation of polycyclic aromatic hydrocarbons in oil-contaminated beach sediments treated with nutrient amendments. J. Environ. Qual. 33, 861–867. doi: 10.2134/jeq2004.0861

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text