Comparison of Two Alternative Methods for Developing TMDLs to Address Sediment Impairments

While excessive sediment is a leading cause of aquatic life use impairments in free-flowing rivers in Virginia, there is no numeric sediment-water quality criterion. As a result, total maximum daily load (TMDL) sediment loads are often established using a comparable, nonimpaired reference watershed. Selecting a suitable reference watershed can be problematic. This case study compared the reference watershed approach (RWA) which uses the Generalized Watershed Loading Function and the disaggregate method (DM) which uses output from Phase 5.3 of the Chesapeake BayWatershed Model. In this case study, the two methods were used to develop sediment TMDLs for three impaired watersheds in Virginia (Taylor Creek, Turley Creek, and Long Meadow Run). In this case study comparison, the RWA required between 12.8 and 14.7 times greater sediment load reductions (t=year) to reach the TMDL load (Taylor Creek > Long Meadow Run > Turley Creek) when compared to the reductions called for using the DM. While each TMDL development method has inherent limitations, the DM uses output from the Chesapeake Bay Watershed Model to establish TMDL target loads. This means that the application of the DM is restricted to the Chesapeake Bay Watershed. DOI: 10.1061/(ASCE)HE.1943-5584.0001728. This work is made available under the terms of the Creative Commons Attribution 4.0 International license, http://creativecommons.org/licenses/by/4.0/. Author keywords: Benthic impairments; Chesapeake Bay; Reference watershed approach; Sediment; Total maximum daily load (TMDL); Generalized watershed loading function (GWLF).


Introduction
The Clean Water Act (Section 303d) requires that total maximum daily loads (TMDLs) be developed for impaired waters (USEPA 2002). A TMDL defines the "allowable" load of a specific pollutant that the waterbody can assimilate and still meet water quality standards: where TMDL = allowable load (mass=time); WLA = waste load allocation from point sources; LA = load allocation from nonpoint sources and naturally occurring background sources; and MOS = margin of safety (which accounts for uncertainty). Virginia's narrative general water quality standard (Section 9- VAC-25-260-20;SWCB 2007) specifies the need to maintain the biological integrity of a waterbody. Compliance with this standard is based on assessment of benthic macroinvertebrate inventories and habitat evaluation. Benthic macroinvertebrates are small, bottom-dwelling organisms which are large enough to see with the naked eye. The quantity and diversity of these benthic organisms are affected by changes in water quality. Benthic macroinvertebrates are relatively immobile (compared to fish) and are therefore more vulnerable to the effects of sediment and other pollutants introduced into a water body (Lenat et al. 1981;Wood and Armitage 1997;Zweig and Rabeni 2001). The diversity of the benthic community and their relative sensitivity to pollution make them responsive to changes in water quality. Benthic organisms are also typically long-lived, making them an important indicator of past and present water quality conditions despite being a major intermediate constituent of the aquatic food chain (Benham et al. 2009).
In 2006, the Virginia Department of Environmental Quality (VDEQ) adopted a standardized, multimetric, bioassessment macroinvertebrate index to assess the aquatic life use status of wadeable freshwater streams and rivers in noncoastal areas of Virginia. The Virginia stream condition index (VSCI) (Burton and Gerritsen 2003) includes a series of biological metrics that are regionally calibrated to an appropriate reference condition and combines them into a single value that is sensitive to a wide range of stressors. Free flowing streams and rivers with a VSCI score of less than 60 are deemed impaired. If an impairment is identified, a stressor analysis is performed that further assesses the status of the benthic community (e.g., assemblage heterogeneity, quantity) in light of potential stressors (e.g., chemical pollution, excess organic matter, flow/ hydrologic modifications, excess sediment) ). This analysis attempts to establish linkages between potential offending pollutant(s) and the status of a benthic community. Once the offending pollutant(s) is (are) identified, pollutant fate and transport computer simulation watershed models are often used to determine the level of pollutant source reductions required to address the impairment.
Per the 2010 integrated report, approximately 28,571 km, or 34%, of Virginia's rivers and streams were assessed for compliance with the appropriate designated uses. Of the assessed rivers and stream reaches, 19,473 km (68%) were impaired for at least one designated use. Aquatic life use impairments constituted approximately 3,640 km (19%) of the total length of impaired rivers and stream reaches (VDEQ 2011a).
The US Environmental Protection Agency identifies excessive sedimentation as a leading cause of aquatic life use impairments in the United States (USEPA 1999(USEPA , 2002. For example, excess sediment can fill the interstitial spaces in gravel and cobble stream substrate, eliminating benthic macroinvertebrate habitat. Sources of sediment include residential runoff, forestry and mining operations, agricultural practices, construction sites, stream bank erosion, and in-stream disturbances (USEPA 2002). A review of TMDLs completed in Virginia between 1999 and 2007 found that sediment was the primary stressor in 46 of the 76 TMDLs developed to address aquatic life use impairments (VDEQ 2011b).
When developing a TMDL for a case in which sediment has been identified as the primary pollutant, the target TMDL sediment load is often determined using a nonimpaired reference watershed, i.e., the reference watershed approach (Younos et al. 2007). The reference watershed approach (RWA) has historically been used in Virginia to establish TMDL target sediment loads because there is no numeric ambient water quality criterion for sediment (VT-BSE 2003, 2004, 2007, 2015, 2016. With the RWA, a nonimpaired watershed having similar characteristics (e.g., similar physiographic region and land use distribution) to the impaired watershed is modeled and the simulated load is used as the target TMDL sediment load. The difference between the impaired and RWA-modeled sediment loads determines the sediment load reduction required to achieve the TMDL target load. The assumption is that if the sediment load in the impaired watershed is reduced to that of the RWA load (i.e., the TMDL target load), the impaired benthic community will (eventually) recover and the stream will be delisted. When applying the RWA, selecting a suitable reference watershed can be challenging because of the limited availability of nonimpaired watersheds with physiographic features similar to impaired watersheds.
In a review of models used for simulating sediment and nutrients in watersheds and receiving waters, Borah et al. (2006) categorized models having the potential for use in TMDLs. Watershed models like Hydrological Simulation Program-FORTRAN (HSPF, Bicknell et al. 2001) are capable of estimating loads (water, sediment, and chemical) from different parts of a watershed as well as simulating water quality and quantity of a receiving waterbody. Loading models such as the Generalized Watershed Loading Function (GWLF) are used to estimate the amount of water, sediment, and nutrients being delivered at a watershed outlet (Haith 1985;Haith and Shoemaker 1987). GWLF has typically been used when applying the RWA in Virginia (VT-BSE 2003, 2004, 2007. The quantification of the Chesapeake Bay TMDL in 2010 addressed both nutrients (nitrogen and phosphorus) and sediment using Phase 5.3 of the Chesapeake Bay Watershed Model (CBWM), which simulates sediment fate and transport using HSPF. Using the output from the CBWM offers a potentially simpler and potentially more consistent method of calculating TMDL target loads for sedimentimpaired watersheds within the Chesapeake Bay watershed. Termed the "disaggregate method" (Yagow et al. 2012, p. 49), this alternative to the RWA relies on land-use-specific unit-area loads derived from the CBWM coupled with a detailed land-use assessment of an impaired watershed to determine a target TMDL sediment load. The disaggregate method (DM) eliminates the need to select and model a reference watershed. This case study used three watersheds in two different physiographic regions in central and western Virginia to compare the TMDL target loads established when using two different procedures-the RWA and the DM.

Impaired Watersheds
Three Virginia watersheds that lie within the Chesapeake Bay watershed were selected for this study (Fig. 1). Details about each impaired watershed are listed in Table 1. Based on the VSCI impairment threshold value of 60, Taylor and Turley Creeks are slightly impaired, while the impairment on Long Meadow Run is more severe.

Reference Watershed Approach
A list of similar, nonimpaired, candidate reference watersheds and data considered most relevant for comparison with the impaired watersheds was assembled. Watershed characteristics evaluated in selecting the reference watershed included mean elevation, land-use distribution, physiographic region, slope, soil erodibility, and VSCI score ( Table 2). The 2000 Mid-Atlantic Regional Earth Science Applications Center (RESAC) land-use data (Goetz et al. 2000) were used when comparing impaired and candidate reference watershed land-use distributions. Soils data used for selecting the candidate reference watersheds came from the USDA Natural Resources Conservation Service Soil Survey Geographic Database (USDA-NRCS SSURGO) (USDA 2011a).
While no two watersheds are identical, a candidate reference watershed that was in the same eco-region and most closely matched the slope and land-use distribution of each impaired watershed was selected as the corresponding reference watershed. Robinson River was selected as the reference watershed for Taylor Creek, and Upper Opequon Creek and Toms Creek were selected as reference watersheds for Long Meadow Run and Turley Creek, respectively (Table 2 and Fig. 2).
GWLF is a lumped parameter model developed to simulate monthly sediment and nutrient loadings in nongauged watersheds (Haith and Shoemaker 1987). GWLF simulates surface and subsurface flows, sediment yield and sediment delivery, and dissolved and attached nitrogen and phosphorous loads from rural, urban, and mixed-land-use watersheds. It can also simulate septic system loads and accommodate point source discharge data (Evans et al. 2003). GWLF assumes that all the sediment generated within a given year exits the watershed during the same year (no net annual sediment deposition). The GWLF model-year runs from April 1st to March 31st (Borah et al. 2006). GWLF requires that land uses be divided into rural (predominantly pervious areas) and urban (predominantly impervious areas) categories. The sediment available for transport from pervious areas is multiplied by a sediment delivery ratio based on watershed size and a transport capacity based on average daily runoff (Yagow et al. 2004). Sediment loads from impervious areas are calculated using the Sartor and Boyd (1972) equation, an exponential function that describes the buildup and wash-off of sediment from impervious surfaces. GWLF estimates sediment loading for a given watershed by aggregating loads from all landuse areas into a single watershed sediment yield. Flow routing and streambank and channel erosion were not considered in the original GWLF model (Haith 1985). However, to improve the model's efficiency, researchers at Penn State University developed a regression equation that can be used to calculate streambank and channel erosion (Evans et al. 2003). This equation was added to the Arc-View version of GWLF (AVGWLF, Evans et al. 2003). AVGWLF was further modified to consider monthly sediment yield by land use and to allow sediment loads from impervious areas to be included in the total sediment load generated by the watershed. The AVGWLF model uses a daily time step to estimate sediment yield, which is then aggregated into a monthly yield. Because GWLF was originally developed for use in ungauged watersheds and is typically not calibrated, the default values suggested by Haith et al. (1992) were used for this application. However, parameters related to sediment and nutrient transport may be adjusted if necessary (Qi et al. 2017).   For the RWA modeling, watershed boundaries for the impaired and reference watersheds were delineated using the 30-m National Elevation Dataset (NED) obtained from the NRCS-USDA Geospatial Data Gateway (USDA 2011b). Soil data were also obtained from the NRCS-USDA Geospatial Data Gateway. The outlet of each impaired watershed was coincident with the downstream limit of the impaired stream segment as defined by VDEQ. Outlets for each referenced watershed were designated as the closest tributary confluence downstream from the corresponding biological monitoring station. The watershed areas, land-use polygons, and flow lengths were determined using ArcGIS 9.3. Soil erodibility factors (SSURGO K-factor), percent slope, and mean elevations were calculated for each land-use category using the zonal statistics function in ArcGIS. Land-use information for the impaired and reference watersheds was derived from the 2000 RESAC land-use dataset (Goetz et al. 2000). The RESAC land-use classes and their distributions within the impaired and reference watersheds are shown in Table 3.
The Biological Systems Engineering (BSE) department at Virginia Tech uses 12 standardized land-use categories when modeling with GWLF (G. Yagow, personal communication, 2011). The RESAC land-use data categories were aggregated into the land-use categories used by BSE ( Table 4). Several of the BSE land-use categories included two or more of the RESAC land uses, while others split some of the RESAC data into two or more BSE land-use categories (e.g., the RESAC land-use category cropland was split between high-till and low-till). The grouping or splitting of land-use data is based, in part, on data from the 2011 National Land Cover Database (NLCD) (Homer et al. 2015), the Conservation Technology Information Center (CTIC 2016), the agricultural census (USDA 2014), and the best professional judgment of the modeler. Note: SSURGO = soil data from the Soil Survey Geographic Database; and K-factor = soil erodibility factor. VSCI scores are based on an average of the last 3 years of VSCI scores (prior to the 303d listings for the impaired watersheds and prior to this study for the candidate reference watersheds). AVSCI score less than 60 signifies impairment.  Yagow et al. (2002) created two Microsoft Excel spreadsheets ("WATERSHED" and "LANDUSE") to facilitate GWLF parameterization. The spreadsheets contain lookup tables with parameter values based on characteristics such as watershed location, soil type, and land-use categories. Use of these spreadsheets provided consistency when parameterizing GWLF.
Often, reference and impaired watersheds are not the same size. The size of a reference watershed is typically adjusted to match the impaired watershed. This area adjustment maintains the land-use distribution in the reference watershed. The area adjustment for the reference watershed is calculated as where A adj = area-adjusted reference watershed (ha); A imp = area of the impaired watershed (ha); A ref = area of the reference watershed (ha); and A l = area of each land use within the reference watershed (ha). Because the sediment delivery ratio (SDR) and the mean channel depths are both functions of the watershed area, the area adjustment ensures that the SDR and the mean channel depth for the impaired watershed and its companion area-adjusted reference watershed are the same. The TMDL target sediment load for the RWA procedure is the simulated load from the area-adjusted reference watershed. The required sediment reduction is calculated as where S red = required sediment reduction (t=year); S ref = areaadjusted reference watershed sediment load (t=year); and S imp = impaired watershed sediment load (t=year).

Disaggregate Method
The disaggregate method used land-use-specific unit-area loads (UALs) from two CBWM model runs (an existing-condition run and the 2010 Chesapeake Bay TMDL target load run) and finerscale, locally assessed land-use inventories of the impaired watersheds of interest in order to calculate TMDL target sediment loads for the impaired watersheds. The CBWM simulated annual sediment loads for the period 1986-2000 for land-based source categories (i.e., nonpoint source loads) for the 2009 progress and 2010 watershed implementation plan (WIP) simulations for the Chesapeake Bay land-river segments that encompassed the impaired and reference watersheds used in this study (USEPA 2010b). The 2009 progress scenario (USEPA 2010b) included pollution controlling best management practices (BMPs) in place in 2009 based on data reported by Virginia to the Chesapeake Bay Program (CBP). The 2010 WIP scenario included those BMPs that were proposed to be implemented by the state to meet the Chesapeake Bay TMDL. Table 5 shows land-use-specific unit-area loads for Albemarle County land-river segment (A51003_to_JL1_6770_ 6850). The UALs were calculated by dividing the annual sediment load for a given land use by the corresponding area: where UAL l = unit-area load for a given land use/source category (t=ha=year); S l = simulated annual sediment load for the land use/source category (t=year); and A l = area of the land use/source category within a given land-river segment (ha). For consistency between the RWA and the DM, the same RESAC land uses used in the DM were used to parameterize the GWLF model. To accomplish this, it was necessary to associate locally assessed RESAC land-use classes with the 24 sedimentgenerating land-use categories used in the CBWM. The CBWM land-use categories were created from a combination of data, such as NLCD imagery and statistical data from the USDA's agricultural census (USDA 2014). The 24 sediment-generating land-use categories from the Phase 5.3 CBWM were associated with the appropriate corresponding land-use categories from the GWLF modeling in the RWA. Using the percentage of each CBWM land-use category associated with each CBWM model run, the areas of the 12 broader GLWF land uses in the impaired and reference watersheds were redistributed, as shown in Table 6. The land-use distribution for the 2010 WIP run may be different from the 2009 progress run because some pollution control measures simulated in the 2010 WIP run required a land-use change in the CBWM.
Taylor Creek watershed spanned two CBWM land-river segments (A51125_JL1_6770_6850 and A51003_to_JL1_6770_6850; Turley Creek and Long Meadow Run both lay within one land-river segment (A51165_PS2_5560_5100). Existing sediment loads were calculated by multiplying the 2009 progress run UALs with corresponding locally assessed land-use areas, and the TMDL target loads were calculated by multiplying the 2010 WIP UALs with corresponding locally assessed land-use areas. For example, if the Table 5. Land-use-specific unit-area loads for Albemarle County, Virginia  Table 6. UAL for high-till w/o manure land use was 2.4 t=ha=year based on the 2009 progress run, and the total area of high-till w/o manure land use within the watershed was 30 ha, the total existing load generated from that land use is (2.4 times 30) 72 t=year.

Sediment Load Comparison
Long-term average annual loads are typically used to define the TMDL load endpoints in sediment TMDLs. This averaging period is used to account for both wet and dry periods in the hydrologic cycle and to account for variations and other important factors (such as snowmelt) that may affect sediment loading. The annual sediment load reductions (t=year and percent) calculated using the RWA and the GWLF model were compared to the annual sediment load reductions calculated using the DM and the CBWM model for each of the study watersheds. The sediment load reductions estimated by the two methods were considered independent of each other, and a Wilcoxon nonparametric test (Bauer 1972) was used to evaluate the pairwise differences in the annual sediment reductions called for by the two methods (α ¼ 0.05). A significance level of 0.05 indicated that there was only a 5% risk of concluding that a difference existed when there was no difference. For each impaired watershed, target TMDL sediment loads were subtracted from the corresponding existing sediment load to determine the required load reductions both in terms of t=year [Eq. (5)] and percent reduction [Eq. (6)]: where S existing = impaired watershed existing sediment load (t=year); and S target = TMDL target sediment load (t=year).

Results
This study compared two different approaches for developing sediment TMDLs to address localized benthic impairments-the RWA and the DM, which use land-use-specific, unit-area loads derived from CBWM output (Phase 5.3) and locally derived land-use areas to calculate existing and TMDL target loads. Existing and TMDL target sediment loads for three impaired and three reference watersheds were calculated for a 15-year simulation period (1986)(1987)(1988)(1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999)(2000) using the two the different procedures. When compared to the DM, the RWA consistently required greater sediment load reductions in order to achieve target sediment loads (Table 7). For Taylor Creek, the RWA predicted greater existing sediment loads than the DM in all but two years (1993 and 1995). The RWA also predicted greater target sediment loads than the DM for Robinson River (the Taylor Creek reference watershed) in all but three years (1992, 1993, and 1995). These inconsistencies were due to a combination of higher UALs obtained from the CBWM for land-river segment A51003_to_JLI_6770_6850 forest land use and lower UALs generated by GWLF simulations for 1992, 1993, and 1995. The average existing sediment load predicted using the RWA for Taylor Creek was 813.2 t=year (1,267.1-453.9) greater than that predicted using the DM, while the average target load predicted using the RWA was 408.7 t=year (841.8-433.1) greater than that predicted using the DM. For Turley Creek, the RWA predicted greater existing sediment loads than the DM in all but one year (1993); the RWA also predicted greater target sediment loads than the DM in eight of the 15 years (1986-1991, 1996, and 1999). The average existing load predicted using the RWA was 1,391.4 t=year (2,672.9-1,281.5) greater than that predicted using the DM, but the average target load predicted using the RWA (i.e., for Toms Creek) was 420.2 t=year (1,060.1-639.9) less than that predicted for Turley Creek using the DM. For Long Meadow Run, the RWA predicted greater existing sediment loads than the DM in all but one year (1993), while the RWA predicted lesser target sediment loads than the DM in eight of the 15 years (1992-1998 and 2000). These inconsistencies were due to higher UALs obtained from CBWM land-river segment A51165_PS2_5560_5100 for the years in question. The average existing load predicted using the RWA was 4960.8 t=year (8,948.7-3,987.9) greater than that predicted using the DM, but the target load predicted using the RWA (i.e., for Upper Opequon Creek) was 1,445.6 t=year (3,304.3-1,858.7) less than the average target load for Long Meadow Run predicted using the DM. The RWA simulated larger sediment loads for the impaired watersheds than for their respective reference watersheds. For the DM, the sediment loads simulated using the 2009 progress run UALs were greater than those simulated using the 2010 WIP run UALs. Table 8 contains the annual and average sediment load reductions predicted by each method, expressed in magnitude (t) and on a percentage basis. The pairwise and average differences (RWA-DM) for each year of the 15-year simulation period are also presented. For Taylor Creek, the RWA required a greater sediment load reduction (t) and a larger percent reduction than the DM in all 15 years. The average sediment load reduction called for by the RWA was 425.4 t=year, but the reduction called for by the DM was 20.8 t=year-a difference of 404.5 t=year. The average percent reduction required by the RWA was 33.6%, which is approximately seven times greater than the 4.6% reduction required by the DM. For Turley Creek, the RWA required a greater sediment load reduction and a larger percent reduction than the DM in all years. The average sediment load reduction called for by the RWA was 2,033.0 t=year, but the average reduction called for by the DM was 221.4 t=year-a difference of 1,811.6 t=year. The average percent reduction was 76.1% for the RWA and 17.3% for the DM. Again, the average percent reduction called for by the RWA was approximately four times greater than that called for by the DM. For Long Meadow, the RWA required a greater sediment load reduction and a larger percent reduction than the DM in all years. The average sediment load reduction called for by the RWA was 7,090.0 t=year, but the average reduction called for using the DM was 683.5 t=year-a difference of 6,406.5 t=year. The average percent reduction was 79.2% for the RWA and 17.1% for the DM. The average percent reduction called for by the RWA was approximately five times greater than that called for by the DM.
Uncertainty analysis is recommended as the base for margin of safety determination in TMDLs (NRC 2001;Reckhow 2003;Hantush and Chaudhary 2014). However, TMDLs developed for benthic impairments in Virginia using the RWA typically include a 10% explicit margin of safety due to modeling-related uncertainties. The DM, on the other hand, assumes an implicit margin of safety, because conservative parameter values and estimates of BMP efficiencies were used as inputs into the CBWM and the Chesapeake Bay TMDL. Consequently, if the RWA TMDL target loads were reduced by 10%, the RWA would still require greater reductions in sediment loads than the DM. A statistical analysis using the Wilcoxon nonparametric test confirmed that at the α ¼ 0.05 level, the load reductions called for by the RWA were all significantly greater than those called for by the DM.

Discussion
The two alternative TMDL development methods compared in this study resulted in significant differences in estimated sediment loads, which led to significant differences in required sediment load reductions for the three impaired watersheds. The differences in sediment load reductions required by the two TMDL development methods occurred for several reasons. First, when the RWA was used, differences in watershed characteristics such as land-use distribution and slope and soil erodibility between the impaired and reference watersheds were manifested in the load reductions called for in the impaired watersheds. The greater the differences between the impaired and reference watershed characteristics, the greater the required sediment load reduction. However, with the DM, sediment load reductions were based on the simulated implementation of BMPs within the impaired watersheds alone. Second, the GWLF model used with the RWA predicted greater unit-area sediment loads for cropland and forested land uses compared to the DM, but the DM unit-area loads for the pasture/hay land-use categories were greater than in the GWLF. Since the impaired watersheds used in this study were primarily comprised of cropland and forest land uses, the RWA produced larger total watershed sediment loads than the DM.
Differences between the RWA and DM sediment loads can also be attributed to differences between the two sediment simulation models and the degree to which each simulates BMPs. In the RWA, the GWLF model does not explicitly account for any BMPs that may be in place except those that are comprised of a land-use change. In the DM, however, the CBWM incorporates a time series of best management practices based on available data for the different model segments. Both of the CBWM scenarios used in the DM method included BMPs that are designed to reduce sediment loading. The 2009 progress scenario included BMPs in place in 2009 based on data reported by Virginia to the Chesapeake Bay Program. The 2010 WIP scenario included the BMPs that were proposed by the state to meet the Bay TMDL. The simulation of BMPs in the DM reduced the amount of sediment delivered to the stream.
Despite the explainable differences between the two approaches, the consistently high sediment load reductions called for when using the RWA in watersheds with only minor benthic impairments raises concerns about the use of the RWA in establishing TMDLs for these watersheds. VSCI scores indicated that Taylor Creek and Turley Creek were only slightly impaired and that the Long Meadow Run impairment was more severe. However, when the RWA was used, the average annual sediment load reduction required for Turley Creek (76%) was similar to that of Long Meadow Run (79%). A reference watershed is selected on the basis of having a VSCI score of 60 or above and having characteristics similar to the impaired watersheds. Finding suitable reference watersheds with characteristics similar to the impaired watersheds can sometimes be difficult to achieve. Differences in watershed characteristics provide a major source of uncertainty in the RWA by affecting the TMDL target load and, ultimately, the relative required reduction. This challenge facing the RWA is greater in urban areas where most watersheds are impaired.
A major limitation of this study was the application of the GWLF with the RWA. Applying another sediment fate and transport model like Soil and Water Assessment Tool (SWAT) or Hydrologic Simulation Program-FORTRAN may have resulted in significantly different results. As has been discussed in the literature, choice of model, which is often a matter of user preference, can have a significant impact on the TMDL outcome (Borah et al. 2006;Shoemaker et al. 2005). Another limitation of this study was the application of the DM using output from Phase 5.3 of the CBWM. This model was designed to simulate drainages and river reaches with discharges of 2.83 cubic meters per second (cms) or greater. Exceptions were made to the 2.83-cms cutoff rule for some river reaches that had sufficient monitoring data but had discharges less than 2.83 cms. Scaling the CBWM output to develop local TMDLs for, in some cases, smaller watersheds did introduce a potential source of error. However, both the RWA and the DM rely on relative comparisons for TMDL development. As such, absolute errors in simulated sediment loads are less consequential. This case study was intended to illustrate possible alternative methods for developing local sediment TMDLs in the Chesapeake Bay watershed.

Summary and Conclusions
This case study compared two alternative methods for establishing target TMDL loads to address aquatic life use impairments caused by excess sediment. The case study assessed three impaired/ reference watershed pairs, Taylor Creek/Robinson River, Turley Creek/Toms Creek, and Long Meadow Run/Upper Opequon Creek. The results indicated that sediment load reductions called for by using the reference watershed approach, both on a t=year and percentage basis, were greater than those called for by using the disaggregate method. The RWA required between 12.8 and Table 8. Required sediment reduction for the impaired watersheds (magnitude and percentages) predicted by the reference watershed approach and the disaggregate method and the pairwise differences between the two methods 14.7 times greater sediment load reductions (t=year) to reach the TMDL load (Taylor Creek > Long Meadow Run > Turley Creek) compared to the reductions called for using the DM. These differences mean that, compared to the DM, the RWA would require the implementation of more pollution control measures (best management practices) in order to achieve the sediment load reduction required to meet the TMDL.
The differences in sediment load reductions required by the two methods occurred for several reasons. When the RWA is used, selected differences in characteristics between the impaired and reference (e.g., land-use distributions, slope and soil erodibility) can lead to large differences in the simulated loads that can produce large required load reductions for the impaired watershed. The DM, however, predicts reductions based on the simulated implementation of BMPs within a given watershed. The GWLF model used with the RWA predicts greater sediment loads per unit area in cropland and forested areas compared to the DM, which predicts greater sediment loads in pasture/hay land-use categories. Because the impaired watersheds used in this study were comprised primarily of cropland and forest land uses, the RWA yielded greater sediment loads than the DM.
There was little or no evidence to support claims that either method will result in attainment of the water quality standards; therefore, additional research is needed to evaluate how effective TMDLs developed using both approaches are at addressing underlying impairment. Additional research is also needed to evaluate whether the DM is appropriate for developing local sediment TMDLs elsewhere in the Chesapeake Bay Watershed.