The Distribution of Lingering Subsurface Oil from the Exxon Valdez Oil Spill

The study approach consisted of the following steps:

1. Build a preliminary model based on data from Short et al. (2004; 2006). They sampled 124 locations in PWS that were randomly selected from segments that were described during Shoreline Cleanup Assessment Technique (SCAT) surveys as heavily or moderately oiled at any time during the period from 1990 to 1993 and beaches described as heavily oiled during 1989 but that had only light to no oil impact during subsequent years. It is important to note that all of these segments would have a high probability of having lingering oil because of their initial degree of oiling. The data for those segments with subsurface oil were related to available geomorphic variables to generate an ensemble of logistic regression models via Bayesian Model Averaging (Hoeting et al. 1996; 2002; Raftery et al., 1997) to predict the presence of subsurface oil at all locations in PWS and GOA.

2. Use a map of the preliminary model results to randomly select sites for field investigation in 2007 that:

   • Were located in PWS and the GOA, to expand the geographic scope,

   • Included all segments that were documented as having ever being oiled, instead of only heavily oiled areas.

3. Conduct field sampling at 106 new sites in PWS and 32 new sites in the GOA using the same methods as Short et al. (2004) with sites from all locations stratified by predicted preliminary model results. At each site, 48 randomly placed pits were dug to assess presence of subsurface oil.

4. Test the model performance with the new field data.

5. Review field data from pooled data set of all 264 sites surveyed in 2001, 2003, and 2007 to identify the geomorphic and hydrologic factors that appear to be controlling the presence and absence of subsurface oil. Use these data to develop new predictor variables.

6. Refine the model using the new predictor variables, a higher spatial resolution, and the more flexible modeling methodology of boosted classification trees (Ridgeway, 2006; Elith et al., 2008).

7. Validate the refined model through additional field sampling in 2008, focusing on sites predicted to have the most oiling.

8. Use the refined model to generate maps and statistics for different degrees of oiling at different ranges of accuracy.

All analyses were carried out with ArcGIS and the R statistical computing language.

The preliminary model (described in Michel et al., 2010) was used to randomly select 140 sites for investigation, using the methods developed by Short et al. (2004). We visually detected oil in one or more pits on 13 of the 140 (9%) sites investigated in 2007. All of the segments where oil was detected were located in PWS. Of the sites where oil was detected in PWS, two had only surface oil, seven had only subsurface oil, and four had both. The intra-site oil detection frequency, or percentage of pits with any recorded oiling, ranged from 0% to over 30%. Figure 1 depicts subsurface oiling results at the 106 sites surveyed in PWS in 2007.

A systematic review of all collected field data was conducted to summarize geomorphic controls on the persistence of subsurface oil. A number of factors were identified as being associated with increased likelihood of retaining lingering subsurface oiling, including:

• Low topographic slope

• Low exposure to wave action

• Armoring of gravel beaches

• Tombolos/natural breakwaters

• Rubble accumulations

• Limited shallow groundwater flow

• Edge effects (transitions between permeable and impermeable shoreline types)

Similarly, a number of factors and morphologies were identified as being associated with reduced likelihood of retaining lingering subsurface oiling, including:

• Impermeable bedrock

• Platforms with thin sediment veneer

• Wide fine-grained gravel beaches

• Extensively treated bayhead beaches

• Steep exposed gravel beaches

• Low-permeability raised bay-bottom beaches

• Strong shallow groundwater flow

• Sandy tidal flats

• Proximity to a stream outlet

Figure 2 shows the average thickness and depth of burial for the different degrees of oiled sediments, for the 509 oiled pits (out of 12,357 pits dug) from all PWS surveys in 2001, 2003, and 2007. The heavier oiled sediments (HOR and MOR; see definitions in Figure 2) were thicker and nearer to the surface than the lightly oiled sediments (LOR, OF, and SH). Half of the oiled pits contained MOR. The number of pits and their oiling degree varied by tidal elevation: 25% of the pits with subsurface oil were in the upper intertidal zone (2.8–3.8 m above MLLW), 70% were in the middle intertidal zone (1.8–2.8 m above MLLW), and 5% were in the lower intertidal zone (0.8–1.8 m above MLLW).

The refined model included the following predictor variables:

1. Intertidal slope and topographical complexity. A fine-scale topographic/bathymetry digital elevation model was generated and used to calculate intertidal slope and an index of local topographic complexity, as a measure of the influence of slope and natural breakwaters (Hobson, 1972; USGS, 2006; USDA, 1996; NOAA, 2007b).

2. Shoreline type. The shoreline type (based on the Environmental Sensitivity Index [ESI] data) was used to calculate an index of geomorphic complexity based on distance to an interface between a permeable and impermeable shoreline type (e.g., transition from a gravel beach to a rocky shore), as a measure of the edge effect and local sheltering.

3. Substrate permeability. An index of substrate permeability was developed based on the ESI shoreline types.

4. Distance to stream mouths. Distance to stream mouths was calculated, to reflect surface water flow and the tendency for finer-grained sediments and higher sediment reworking (USGS, 2007; USDA, 2002).

5. Shoreline convexity. Indices of shoreline convexity were calculated at scales of 50, 100, and 500 m, to reflect small-scale sheltering from wave energy.

6. Distance from spill source. Overwater distance to the release site was calculated, to reflect a general trend in decreased degree of oiling with distance from the source.

7. Shoreline orientation to spill. An oil-approach angle index was calculated, to reflect the orientation of the shoreline relative to the angle of the approach of floating oil slicks.

8. Oiling history. A quantitative index of oiling history based on SCAT oiling categories (No Oil, Very Light, Light, Moderate, and Heavy) for 1989, 1990, and 1991 was calculated (Neff et al., 1995; Teal, 1991). This index is a circular moving average using a 50 m radius of the back-calculated estimate of shoreline surface oil coverage (m²/m).

9. Wave exposure. Exposure to wave action was represented by fetch and an exposure index (NOAA, 2007a; USACE, 1984; Findlayson, 2005; Rohweder et al., 2008; Hayes, 1996).

The refined model was constructed using an ensemble of tree-based classifiers via Generalized Boosted Modeling (GBM) that can inherently handle non-linearities and interactions among multiple predictor variables (Ridgeway, 2006; Elith et al., 2008). The refined model was constructed with data from this study, Short et al. (2004; 2006), and limited data from Irvine et al. (1999, 2006).

The refined model was also evaluated for accuracy both with the data used to construct it and validated via a limited data set collected in 2008. Figure 3 shows the Receiver Operator Characteristic (ROC) for the preliminary model from both internal and validation data. The Area Under the Curve (AUC) of the refined model ROC curve is 0.99. Though likely over-optimistic, this indicates excellent overall model performance. The AUC of the final model ROC curve via the limited validation data is 0.80 (MOR or > only). Though this is not as good as the internal model validation statistics, the validation data were limited, and this still represents good model performance. The oiling histories (e.g., SCAT oiling descriptors) were the most significant drivers, though all of the predictor variables played some role, as in Figure 4.

Different models were run for the following:

• Any subsurface oil

• Any oil equal to or greater than lightly oiled residues (LOR)

• Any oil equal to or greater than moderately oiled residues (MOR)

• Any oil equal to or greater than heavily oiled residues (HOR)

• Any oil equal to or greater than MOR and greater than 15% cover (representing oil in 1 of 6 pits in the column)

• Any oil equal to or greater than MOR and greater than 30% cover (representing oil in 2 of 6 pits in the column)

• Any oil equal to or greater than MOR and greater than 50% cover (representing oil in 3 of 6 pits in the column)

The results are shown in Tables 1–3 for different degrees of Positive Predictive Value (PPV) – a measure of acceptable error rates. The tables show number of cells (about 10 m in shoreline length due to the rasterization of a linear shoreline), the total length in kilometers, and the number of “sites” defined as a group of cells within 100 m of each other. Figure 5 depicts a map of a subset of model output. In Table 1, under the “known” category, there are 393 cells, 2.85 km, and 67 sites with any kind of subsurface oil. These statistics were derived from only the actual field survey. Based on the model predictions for the entire oil-impacted area (PWS and GOA), using the 90% accuracy cutoff, there are 2,668 cells, 19.36 km, and 167 sites with any subsurface oil. Similar statistics are provided in Table 3 for the three different geographic regions: PWS, the Kenai Peninsula (KEN), and Shelikof Strait (SHL). Using a treatment threshold of any oil equal to or greater than MOR and the 90% accuracy cutoff, there are 492 cells, 3.57 km, and 64 sites. In Table 2, using treatment thresholds of any oil equal to or greater than MOR and greater than 15% frequency of oiled pits, at the 90% accuracy threshold, there are 361 cells, 2.62 km, and 52 sites. Using treatment thresholds of any oil equal to or greater than MOR and greater than 50% frequency of oiled pits, at the 90% accuracy threshold, there are 114 cells, 0.83 km, and 31 sites. The model could be used to predict the number of cells, kilometers of shoreline, and number of sites for different screening criteria and accuracy thresholds.

Figures 6 and 7 show histograms of cumulative length of oiled shoreline per discrete contiguous site. These plots provide a descriptive summary of the estimated along-shore length of the subsurface oil patches meeting the type and relative amount criteria within these discrete sites. As can be seen in Figure 6, about 80% of all the sites with any subsurface oil are greater than 10 m in length. For sites with MOR or greater and 15% coverage, 62% are greater than 10 m in length. The relatively small sizes of the heavier lingering oil sites is not of particular concern in terms of treatability, particularly considering the results of Boufadel et al. (2011) that indicate there are effective treatment methods for remediating small patches of subsurface oil.

We developed a suite of geospatial models to identify potential areas where subsurface oil is still present on the shorelines of PWS and GOA affected by the Exxon Valdez oil spill. It was also important to be able to quantify the relative amount of subsurface oil, so that the models could be used as screening tools to prioritize shorelines for different remediation methods. The models were based on data collected at 314 shoreline segments surveyed between 2001 and 2008. These data allowed identification of a number of geomorphologic and hydrologic factors that have likely contributed to the persistence of subsurface oil within PWS and GOA two decades after the spill.

Because detailed data layers for each of these factors were not available, the model used existing data sets as surrogates to represent these factors, such as distance to a stream mouth or shoreline convexity. While the linkages between the data used and the physical phenomena that drive persistence are not clearly understood in all cases, the performance of these models was remarkably good. The model simultaneously evaluates all identified variables to make a similar assessment in a rigorous, unbiased manner. The model results suggest there are a limited but significant number of as-yet unsurveyed locations in PWS and GOA that are likely to contain lingering subsurface oil. Furthermore, the model results may be used to quantitatively prioritize shoreline for investigation with known uncertainty.

It should be noted that regional differences in physical environment and oiling history, and sparser recent data, mean that extra caution should be used in interpreting model results in the GOA region. With understanding of these limitations, these model results can be used by educated field practitioners in concert with the knowledge of the factors that increase or decrease the likelihood of subsurface oil within a given shoreline segment. Used as such, these models may be capable of predicting the distribution of the lingering oil across the entire spill area with sufficient accuracy and resolution to perform as useful tools for evaluating potential ongoing impacts to users of these shorelines and prioritizing shoreline locations for potential remediation based on various screening criteria, or predicting potential for lingering oil at future spills.

This study was funded by the Exxon Valdez Oil Spill Trustee Council via Contract No. AB133F07CN0097 to RPI through the National Oceanic and Atmospheric Administration (NOAA) with additional funding to the U.S. Geological Survey and NOAA. We appreciate the very helpful coordination of Peter Hagen of NOAA. Mandy Lindeberg, Jeep Rice, and Mark Carls of the NOAA Auke Bay Laboratory (ABL) are acknowledged for their help in all phases of the work. ABL provided the chemical analyses of the samples collected in 2007. Jerry Pella, a retired ABL statistician, provided invaluable assistance in study design. Linos Cotsapas of RPI provided excellent logistical support. The field work was conducted by a dedicated team of scientists: David Betenbaugh, Christine Boring, Ion Cotsapas, Linos Cotsapas, Heidi Dunagan, Joe Holmes, Tom Freeman, Daniel Mann, Mike Tetreau, and Mark White, with industrious pit digging by numerous field assistants. Michel Boufadel of Temple University provided input on the importance of beach hydrology. David and Annette Janka of the Auklet, and their crews, provided safe vessel operations and great knowledge of PWS. John Rogers and John Whittier, captains of the Waters, also provided excellent support for crews working in the challenging Gulf of Alaska environment. Any use of trade names is for descriptive purposes only and does not imply endorsement by the federal government.