w EPA
   United States
   Environmental Protection
   Agency
velopmeril
  Watershed Classification Framework
      for the State of West Virginia
   L
  J. , *.

                                         / i

-------
                                                          EPA/600/R-03/141
                                                                   July 2004


Watershed Classification Framework for the

                   State of West Virginia



                WV R-EMAP Final Report

                            Naomi E. Detenbeck
                    US EPA Mid-Continent Ecology Division
            National Health and Environmental Effects Research Laboratory
                             Duluth, MN 55804

                              Leslie A. Jagger1
                              Stacey L. Stark1
                             OAO Corporation
                             Duluth, MN 55804

                             Matthew A. Starry2
                        Computer Sciences Corporation
                             Duluth, MN 55804

                   'EPA FAIR Contract, 2EPA FAIR II Contract


                                CR82872001

                              Project Officer

                             Naomi Detenbeck
                        Mid-Continent Ecology Division
            National Health and Environmental Effects Research Laboratory
                             Duluth, MN 55804
                        Mid-Continent Ecology Division
            National Health and Environmental Effects Research Laboratory
                      Office of Research and Development
                      U.S. Environmental Protection Agency

                              Recyded/Recydable
                              Printed with vegetable-based ink on
                              paper that contains a minimum of
                              50% post-consumer fiber content
                              processed chl ori ne free

-------
                                        Notice

The information in this document has been funded wholly by the U.S. Environmental Protection
Agency (US EPA), including support to U.S. Geological Survey EROS Data Center (IAG DW-
14938973) and West Virginia Department of Natural Resources (CR-82872001); and contract
support to OAO Corporation under US EPA FAIR Contract 68-W5-0065, Delivery Order #24, and
Computer Science (CSC) Corporation under US EPA FAIR II Contract 68-W01/W02-032, Task
Order #024. This document has been prepared at the US EPA National Health and Environmental
Effects Research Laboratory, Mid-Continent Ecology Division, in Duluth, Minnesota, with support
from EPA FAIR II Contract 68-W01 /W02-032 to CSC.  It has been subjected to review by the US
EPA National Health and Environmental Effects Research Laboratory and approved for
publication. Approval does not signify that the contents reflect the views of the Agency, nor does
mention of trade names or commercial products constitute endorsement or recommendation for
use.
 Preferred citation:
Detenbeck, N.E., L.A. Jagger, S.L. Stark, and M.A. Starry. 2004. Watershed classification framework for
the state of West Virgnia: WV R-EMAP Final Report. EPA/600/R-03/141. U.S. Environmental
Protection Agency, Office of Research and Development, National Health and Environmental
Effects Research Laboratory, Mid-Continent Ecology Division, Duluth, MN.

-------
                            TABLE OF CONTENTS

LIST OF ACRONYMS	iv

LIST OF FIGURES	  v

LIST OF TABLES	   vii

LIST OF APPENDICES	vui

INTRODUCTION 	  I

PURPOSE OF CLASSIFICATION AND INTENDED USE OF DATA  	  2
      Purpose of watershed classification	  2
      Hydrology-based classification of watersheds	  4
      Relationship between HUCs and watersheds	  5
      Use of watershed classification as a sampling and assessment framework	  6

METHODS	  9
      USGS watershed delineation	  9
      Delineation and coding of hydrologic units for the state of WV	  9
      Watershed characterization	  14
      Creation of sample frame and sample design for 12-digit HUCs  	  17
      Response  threshold derivation	  20

CLASSIFICATION RESULTS	  30
      Hydrology thresholds	  30
      Distribution of land-use/land-cover variables across 12-digit subwatersheds in
             West Virginia 	  33

FUTURE UPDATES AND CLASSIFICATION REFINEMENTS	  45

CONCLUSIONS	  50

ACKNOWLEDGMENTS 	  51

REFERENCES  	  53

APPENDIX A: METADATA for MAPS AND DATABASES, INCLUDING DATA
      QUALITY INFORMATION	A-l

APPENDIX B: AMLS	B-l

APPENDIX C: DATABASE TABLES	C-l
                                                                              111

-------
LIST OF ACRONYMS




DEM        Digital elevation model




DEP         Department of Environmental Protection




DNR        Department of Natural Resources




EDNA       Elevation Dataset with National Applications




EMAP       Environmental Monitoring and Assessment Program




FGDC       Federal Geographic Data Committee




GIS         Geographic Information System




MAIA        Mid-Atlantic Integrated Assessment




NHD        National Hydrography Dataset




NWBD      National Watershed Boundary Database




NWI        National Wetlands Inventory




PRISM       Parameter-elevation Regressions on Independent Slopes Model




R-EMAP     Regional Environmental Monitoring and Assessment Program




US EPA      U.S. Environmental Protection Agency




USGS        U.S. Geological Survey
i\

-------
LIST OF FIGURES

Figure 1.   Use of watershed classification	  3

Figure 2.   Definition of a subwatershed region	  7

Figure 3.   Generic example of Pfafstetter coding for an 8-digit HUC	  8

Figure 4.   USGS gaging stations and associated watershed boundaries	  10

Figure 5.   Map of 8-digit hydrologic cataloging units (HUCs), 10-digit watersheds, and 12-
           digit subwatershed units produced for West Virginia	  11

Figure 6.   Automated derivation of main channel slope	  18

Figure 7.   Major drainage basins associated with the state of West Virginia and surrounding
           states, as defined by 2- and 4-digit HUCs	  19

Figure 8.   Example of unequal weighting of interbasin HUCs	  21

Figure 9.   Map of ecoregions overlapping with West Virginia watersheds	  22

Figure 10.  Map of land-use covering West Virginia watersheds	  27

Figure 11.  Map of palustrine and lacustrine classes from National Wetlands Inventory for the
           state of West Virginia	  29

Figure 12.  Predicted 2-year flood normalized to watershed area as a function of
           watershed storage	  31

Figure 13.  Results of Classification and Regression Tree (CART) analysis	  32

Figure 14.  Map of WV 12-digit HUCs coded by fraction watershed storage	  34

Figure 15.  Map of WV 12-digit HUCs coded by estimated fraction impervious surface area ....  35

Figure 16.  Map of WV 12-digit HUCs coded by estimated fraction surface mining area	36

Figure 17.  Map of WV 12-digit HUCs coded by fraction agricultural area	  37

Figure 18.  Map of WV 12-digit HUCs coded by fraction forested area	  38

Figure 19.  Map of Omernik Level III ecoregions within West Virginia with land-use summaries.  40

Figure 20.  Map of WV 12-digit subwatershed units by land-use and watershed storage cksses.. .  41

-------
Figure 21. Watershed cksses associated with 12-digit HUCs in sample population for Central
          Appakchkn Pkteau and Central Ridge and Valley ecoregions in West Virginia.  ...   43

Figure 22. Watershed classes associated with 12-digit HUCs in sample population for Western
          Allegheny Plateau ecoregion in West Virginia	 44

Figure 23. Effect of exponent term for storage on shape of relationship between watershed
          storage and normalized peak flow  	 48

Figure 24. Effect of peak flow breakpoint, exponent term, and prediction error on  correct
          classification rate for high peak flow ckss or misckssification rate for low peak flow
          ckss	 49
 VI

-------
LIST OF TABLES

Table 1.   Guidelines for hydrologic unit subdivision for National Watershed Boundary
          Database, according to interagency protocol (FGDC 2002)	  13

Table 2.   Geographic information system databases used to characterize WV watersheds	  15

Table 3.   Definition of sample frame for R-EMAP wadeable stream sampling in
          WV, 2000-2001	  20

Table 4.   Potential watershed characteristics related to peak flows  	  24
                                                                                       VII

-------
LIST OF APPENDICES (CD-ROM, back pocket)




APPENDIX A: INVENTORY OF MAPS AND DATABASES	 A-l




APPENDIX B: AMLS	B-l




APPENDIX C: DATABASE TABLES  	 C-l
VIM

-------
INTRODUCTION





       The Environmental Monitoring and Assessment Program (EMAP) is a research program to




develop the tools necessary to monitor and assess the status and trends of national ecological




resources. These tools include probability-based survey designs and indicators of biological




condition.  Regional EMAP (R-EMAP) is designed to evaluate how these tools can be applied at




local and regional scales to meet the management needs of the states and regions. One of the




emerging issues in monitoring programs for the states and tribes is the need to develop designs to




meet multiple objectives outlined in different sections of the Clean Water Act. To meet the




requirements of Section 305(b), probabilistic survey designs are needed to produce estimates of




regional condition with a known level of confidence.  In addition, monitoring programs must




identify impaired waters and associated causes of impairment to meet the listing requirements of




Section 303 (d) of the Clean Water Act (US EPA 2002). Ideally, monitoring programs will provide




information to allow states to more efficiently plan the next round of monitoring, in order to




identify as many impaired waters of the state as possible.  One approach which can accomplish both




of these needs is to incorporate risk-based categories into sampling designs, either as strata in a




random-stratified design, or as  categories to which different probability-weights are assigned. This




approach allows managers to summarize not just information  about regional condition, but also




information about the risk of impairment associated with  different classes of aquatic resources.





       A R-EMAP project was designed for the state of West Virginia to accomplish multiple




objectives, including the evaluation of a risk-based random-stratified  sampling process.  An




additional goal of the WV R-EMAP project was to refine the fish index of biotic integrity (IBI)




developed through EPA's Mid-Atlantic Integrated Assessment (MALA) project to produce separate




indices for different thermal classes of streams and specific to the state of West Virginia. The




purpose of this report is to describe the methods involved in producing a preliminary watershed




classification and incorporating the risk-based watershed classification into a monitoring design for

-------
 the state of West Vkginia. Validation and refinement of the watershed dassification scheme, based




 on analysis of pending monitoring data, will be addressed in a subsequent report.










 PURPOSE OF CLASSIFICATION AND INTENDED USE OF DATA




 Purpose of watershed classification




       Historically, classification systems have been used for inventory purposes, for stratifying




 landscapes prior to characterizing reference or regional condition, and for facilitating




 communication with natural resource managers and the public  (Omernik  and Gallant 1988;




 Heiskary and Wilson 1990; Herlihy et al. 2000; Pan et al. 2000; Waite et al. 2000). More recently, it




 has become necessary to develop classification systems that can explain differences in vulnerability




 of aquatic ecosystems to stressors, aid in regionalizing water quality criteria, and predict probability




 and causes of impairment of aquatic systems (US EPA 2002).




       Most sources of stream impairment contained in state 303(d) listings are related to nonpoint




 source pollution (US EPA 2001, 2003). To more efficiently deal with water quality management




 issues, an integrated approach to small watershed assessment, diagnosis, and restoration planning is




 needed (Figure 1). Current guidance from US EPA supports the  development of a consolidated




 assessment and listing approach to enable joint application of Clean Water Act Sections 305(b) and




 303(d) (US EPA 2003).  Monitoring strategies developed for 305(b) regional assessments also must




 inform the 303(d) listing process.  In the 2001-2003 WV Regional Environmental Monitoring and




 Assessment Project (R-EMAP), use of a new watershed classification technique specific for West




Vkginia should allow the state to improve its ability to discern anthropogenic causes of impakment




 from non-anthropogenic sources of variability in reference condition (Cincotta 2000; Brazner et al.




2003) by producing a series of stressor-response relationships specific to different watershed classes.

-------
                      threshold identification
                    strategies: monitoring,
                    modeling, restoration
                           i
                           I
                            cost-benefit
                            analyses
                 watershed prioritization
                                                 classification
monitoring framework
                                              class condition and cumulative
                                              condition assessment
                            screenm
                                            threshold confirmation
                                            [translators]
       Figure 1. Use of classification in sequence of monitoring, assessment of condition,
       and prioritization of watersheds for further monitoring, modeling, and restoration.

In addition, use of the watershed classification scheme within a monitoring strategy will allow

estimation of the condition of classes of watersheds and prediction of ecological risk for

unmonitored systems.

       Gradients of land-use can be studied in combination with hydrology moderating factors to

establish relationships of changes in runoff, and associated water quality and biological responses,

with land-use activities (Verry 1986; Detenbeck et al. 2000).  Since multiple non-point stressors

potentially are affecting biological condition, it is appropriate to use watersheds as the sampling

and/or assessment unit.  In order to establish an appropriate stratified random sampling regime to

diagnose causes of impairment, watersheds within the state of West Virginia must be properly

classified. The classification developed and described in this document is based on identifying the

levels of hydrology-moderating factors such as watershed storage and watershed development at

which rapid degradation occurs. Watershed storage capacity is operationally defined as the fraction

-------
 of watershed area present in lake + wetland area, and can be used to predict the magnitude of floods




 (Jennings et al 1993). Peak flows are generally associated with a large fraction of nutrient and clean




 sediment yields from watersheds, so that watershed storage can also be used to predict attenuation




 of nonpoint-source loadings. Wetlands and lakes increase retention time within the watershed,




 allowing processes such as nutrient uptake, denitrification, and sedimentation in wetlands to




 improve water quality. Watershed storage can also be a good predictor of baseflow water quality




 (Johnston et al. 1990; Detenbeck et al. 2000; Detenbeck et al. 2003a,b).




 Hydrology-based classification of watersheds




       Watershed properties such as depressional storage volume and surface-water




 hydrogeomorphic types can be used to classify lakes or streams according to relative risk of impact.




 Our approach to predicting stream sensitivity to nonpoint source pollution is based on the nonlinear




 response of hydrologic regimes and associated loadings of non-point source pollutants to watershed




 properties (Richards 1990). Selected hydrologic thresholds are related to 1) variation in levels of




 watershed storage (either natural or anthropogenic) and 2) land-use activities affecting runoff and




 the hydrologic regime. In common usage, the term "threshold" refers to "the point at which a




 physiological or psychological effect begins to be produced" or "the point at which something




 starts." We define a hydrologic thresholds a breakpoint or inflection point in a nonlinear relationship




 between a watershed property and hydrologic response variable such as peak flows. For example,




peak flows per unit area often increase exponentially once the storage capacity of a watershed has




been exceeded. We have identified these thresholds based on a space-for-time substitution, or




comparison across watersheds, as inadequate historical records are available for analysis of time




series of individual watersheds.




       The  United States Geological Survey (USGS) has defined a series of empirical nonlinear

-------
equations relating watershed properties such as watershed area, channel slope, watershed storage,




and land-use (percent forested, percent urbanization, or percent impervious surface area) to peak




flows of given recurrence intervals (Q2, Q5,... Q100; Jennings et al. 1993). Peak flows increase




exponentially as watershed storage decreases below a given threshold (Detenbeck et al. 2000).




These thresholds form the basis for a watershed classification system.




Relationship between HUCs and watersheds




       The first step in developing a classification strategy is to define  the population of interest.




Theoretically, an infinite number of watersheds could be defined, starting at any random point on a




stream network (Detenbeck et al. 2003a). Benefits of defining a fixed set of watershed units include:




a) the ability to tie assessment units to management plans and actions, b) ease in mapping




classification units and classes, and c) the significance of fixed management units for stakeholders




who form watershed management organizations to facilitate local environmental protection.




       The USGS is coordinating with other Federal agencies and the  states to develop a National




Watershed Boundary Database (NWDB), which will sequentially divide existing 8-digit Hydrologic




Units (HUCs) into 10-digit HUCs, and 10-digit HUCs into 12-digit HUCs (Legleiter 2001). Use of a




seamless national boundary database will facilitate communication among neighboring states, and




planning and management of basins that cross interstate boundaries. In addition, creation of




seamless nationwide GIS kyers will allow development of a series of elevation, hydrography, and




derived coverages that are internally consistent with one another (Franken et al. 2001).




       We chose  12-digit HUCs as the base map for our classification  framework.  Watersheds




associated with 12-digit HUCs span the size range of wadeable streams in West Virginia (400 -




40,000 ha), correspond to the  size range of most management units in the state, and fit within a




hierarchical scheme of larger management units (i.e., 10- and 8-digit HUCs). A total of 36 8-digit

-------
HLJC units drain into the state of West Virginia; these 8-digit HUCs contain a total of 242 10-digit




HUCs (USGS watershed units) and 1427 12-digit HUCs (USGS subwatershed units).  Of these 12-




digit HUCs, 883 fall predominandy within the borders of West Virginia.




       Not all HUCs are equivalent to watersheds. HUCs can be defined as either basin or




interbasin.  Noncoastal basin HUCs are generally equivalent with watersheds at the 12-digit scale,




while interbasin HUCs  divide the main channel of the next larger size class of HUCs into segments.




Coastal HUCs can contain multiple parallel watersheds draining to a large water body (Great Lake or




ocean), but these are not an issue in West Virginia.  Interbasins may include a middle segment of a




larger river that encompasses many smaller side tributaries. Interbasin HUCs can be aggregated with




upstream units to define full watersheds (Figure 2). Throughout this report, derivation of watershed




attributes for 12-digit HUCs is based on the boundary and characteristics of the full aggregated




watershed regions.




       Use of Pfafstetter codes facilitates the discrimination between basin and interbasin HUCs,




and identification of all upstream or downstream HUCs.  Pfafstetter codes are not the same as HUC




codes, but are produced by the USGS during the process of creating elevational derivatives such as




watershed boundaries.  Pfafstetter codes are hierarchical, with a 2-digit segment added at each new




level of subdivision. The last digit of the Pfafstetter code increases consecutively from mouth to




headwater sub-basins, with even numbers associated with tributaries (basin HUCs) and odd




numbers associated with HUCs along the mainstem (interbasin HUCs) (Verdin and Verdin 1999;




Figure 3).





Use of watershed classification as a sampling and assessment framework




       Once a fixed set of watersheds is identified, delineated, characterized, and classified, this set




can be used to define a sample population and survey design. At least two options exist for

-------
                                Subwatershed Region
HUC 050701010103
Region with National
Hydrography Data
and Shaded National
Elevation Data
                                                                  Legend

                                                                  |    | 050701010103 Region
                                                                  -/v— NHD stream reaches
Map Projection Alters Equal Area Conic
 Figure 2. Definition of a subwatershed region for 12-digit HUC 050701010103 through
 aggregation of 05 additional upstream HUCs.

-------
Figure 3. Generic example of Pfafstetter coding for an 8-digit HUC (red outline). The base of the
main channel is coded as "1" while the headwater portion of the main channel is coded as "9".
Each of the four largest tributaries (basin HUCs) is assigned an even number between 2 and 8, in
ascending sequence from the mouth towards the headwaters.  Likewise, each of the segments
draining to the mainstem between the tributaries (interbasin HUCs) is assigned an odd number
between 3 and 7, in ascending sequence from the mouth to the headwaters. Coding of subdivisions
within one of these 9 HUCs would proceed in a similar fashion, with the addition of a second digit
to the code.

-------
 incorporating watershed classes into a survey design framework: 1) classes can be used as strata in a




stratified random design, or 2) classes can be assigned unequal probability weights to influence the




distribution of selected watershed units across classes (Detenbeck et al. 2003a). Both options allow




a rare but significant watershed class (e.g., one with high associated probability of impact) to be




sampled at a higher frequency than would be possible if a simple random survey design were




implemented.  Monitoring data can then be assessed for classes of watersheds as well as for an entire




region, and a series of appropriate management strategies devised.










METHODS




USGS watershed delineation




       To define hydrologic thresholds, gaging stations with long-term flow records and derived




peak flow statistics were identified from previous studies in the state of West Virginia (Frye and




Runner 1970; Runner 1980; Figure 4). Watersheds associated with USGS gaging stations were




delineated onscreen in Arclnfo using digital raster graphic (DRG) backdrops (1:24,000) and




National Hydrography Dataset (NHD) coverages.




Delineation and coding of hydrologic units for the state ofWV




       In contrast to manual delineations for USGS gaged watersheds, hydrologic units in West




Virginia were delineated through a two or three stage process as part of the development of the




Elevation Dataset with National Applications (EDNA, Franken et al. 2001;




http://edcntsl2.cr.usgs.gov/ned-h/). The EDNA process produces HUC boundaries suitable for




inclusion in the National Watershed Boundary Dataset (Figure 5). Delineation tasks were




accomplished through an interagency agreement with the USGS EROS Data Center and




collaborative efforts of US EPA Region III staff.  Stage I of the EDNA process previously had been

-------
                                                                              USGS Gaging Station
                                                                              Watersheds
Map Projection AJbar* Equal ATM Conic
                                                                                • Gagng Slatnn Ports
                                                                                  Gagng Station Walenhedi
                                                                               Level III Ecorrgions
                                                                               Name
                                                                                  Wteslem Allegheny PbKau
                                                                                  Central Appaladnam
                                                                                  Cenlial Appalachun Ridgn and Va
 Figure 4. Location of USGS gaging stations and associated watershed boundaries used for
 derivation of hydrologic thresholds in West Virginia.
 10

-------
                                                                         WV8,10, &12-Digit
                                                                         HUCs with Contributing
                                                                         8-Digit HUCs
                                                                           |   J 8-Digit HUCs (EDNA)
                                                                               8-Digit HUCs (contributing)
                                                                               10-DigitHUCs
                                                                               12-Digit HUCs
                                                                           Processed HUCs
                                                                           Data Type
                                                                               USGS1:250K
                                                                               EDNA Stage II
                                                                               EDNA Stage III
Map Projection. Alters Equal Area Conic
 Figure 5.  Map of 8-digit hydrologic catologing units (HUCs), 10-digit watersheds, and 12-digit
 subwatershed units produced for West Virginia.  Stage III 12-digit HUCs are shown only for three
 of the 13 8-digit HUCs completed as of June 2003, i.e., those three 8-digit HUCs requiring
 additional processing beyond Stage II to correct flow directionality.
                                                                                               11

-------
accomplished through a memorandum of understanding between the USGS and US Weather





Service. Stage I consists of an automated watershed delineation process using digital elevation




models (DEMs) with GIS to produce a nationwide coverage of catchments at the scale of 1-2 mil2.




Catchments produced at this stage reflect the predicted movement of surface water over surface




topography; elements such as karst features are not included. Stage I produces several products,




including catchment boundaries, synthetic streamlines (predicted location of streams and river




courses), and elevational derivative raster coverages such as slope, aspect, flow direction, and flow




accumulation.




       During Stage II, Stage I catchments are aggregated to the scale of 10- and 12-digit HUCs,




based on a set of national guidelines (FGDC 2002; Table 1). At this point, GIS coverages are




evaluated for internal consistency and potential errors are flagged.  Indicators of potential errors




either in digital elevation models or in NHD include divergence of synthetic streamlines from




mapped NHD streams and inconsistencies between HUC boundaries and other sources of mapped




watershed boundaries. During Stage III, a series of Arc View tools are used to correct small regions




of the original digital elevation models near the flagged errors to represent the actual flow of surface




water.  For example, a highway overpass could be detected as a barrier to flow in an original DEM,




when in fact, water can freely flow underneath the underpass.  Creating a small "notch" in the




original DEM at that neighborhood allows the stream course to be correctly predicted. After the




hydrologically-corrected DEMs have been created, boundaries for 10- and 12-digit HUCs are




regenerated, as are raster coverages for other hydrologic derivatives.
12

-------
 Table 1. Guidelines for hydrologic unit subdivision for National Watershed Boundary
 Database, according to interagency protocol (FGDC 2002).
Unit digits
8-digit
10-digit
12-digit
Name Size range (acres)
cataloging unit
watershed 40j000 tQ 250>m
subwatershed -, Q QQQ ^ QQQ
Subdivisions/unit
5-15
5-15

       For the West Virginia R-EMAP project, we used Stage II HUCs to develop the watershed

classification framework.  Currently, Stage III HUCs are only partially completed for the state, and

were only used where Stage II HUCs required substantive correction, e.g., to correct direction of

flow (Figure 5).


       The assignment of Pfafstetter codes during Stages I and II of the process allows ready

identification of upstream and downstream 10- and 12-digit HUCs within 8-digit HUCs, and thus

aggregation of interbasin HUCs to define watersheds using REGION processes in Arclnfo.

Definition of Arclnfo regions is desirable considering the nested nature of watersheds. Regions are

model areal geographic features described by one or more polygons.  Regions are efficient ways to

model nested or overlapping data because of the concept of shared geometry. For example, one

polygon feature can belong to multiple regions without the topology being duplicated. Two files

store the region-arc relationship and the region-polygon relationship "behind the scenes." Regions

are stored as subclasses within a polygon coverage and each region subclass has its own separate

attribute table.


       The bulk of the work  was completed using individual Arc Macro Language (AML) scripts

developed for each component of the characterization process. The scripts were used to create and

combine region subclasses into an integrated coverage for further geoprocessing using the


                                                                                         13

-------
REGIONQUERY command in Arclnfo (ESRI, Redknds CA). The REGIONQUERY command





mimics traditional oveday commands on region subclasses in an integrated coverage.  So, for




example, if a single coverage had a 12-digit watershed subclass and a wetland subclass, the




REGIONQUERY command could be used to create a new subclass computed from the geometric




intersection of the two subclasses (watersheds and wetlands). The region tables keep track of which




wetland polygons belong to each watershed polygon while only physically representing the polygons




once. One table could then be created that summarized wetland type by area for each watershed




region.





       Due to the size of input datasets and software/hardware limitations, 8-digit HUCs




containing 12-digit watersheds were processed individually. This presented an issue when a 12-digit




watershed required additional upstream 8-digit HUC(s) to be considered a true watershed. A




method was needed to identify diese "interbasin" watersheds and combine the output for their




required counterparts. Each interbasin watershed was given a special identifier code in the




watershed region attribute table.  A table was made listing all possible interbasin codes and the 12-




digit watersheds needed to complete the full watershed.  This table was imported into a Microsoft




Access database. The individual output tables for 12-digit watersheds within an 8-digit HUC were




combined into one table for each characterization process and imported into the same database.  A




series of queries were performed to combine the individual output values to determine the total




output values for interbasin watersheds for each characterization process.





Watershed characterisation






       Watershed regions associated with all 12-digit  HUCs in the state of West Virginia and long-




term gaging stations were characterized for attributes expected to be related to generation of peak




flows (Fne and Runner 1970; Runner 1980): watershed area, main  channel length, main channel




14

-------
slope, watershed storage, mean January minimum temperature, snowpack, and percent forest cover.




In addition, watersheds were characterized for other major land-use and land-cover characteristics




related to generation of runoff and nonpoint source pollution: percent mined land, percent




impervious surface area, and percent agriculture. GIS databases used to calculate watershed




characteristics are described in Table 2.
Table 2. Geographic information system databases used to characterize West Virginia watersheds. EDNA = Elevation Dataset with
National Applications, NHD = National Hydrography Database, NWI = National Wetlands Inventory, OWI = Ohio Wetlands Inventory,
NYWI = New York Wetlands Inventory, PRISM = Parameter-elevation Regressions on Independent Slopes Model.
Database
EDNA



NHD

NLCD-
mining
modified


NWI
OWI
PRISM
Layers and attributes
10-digit watershed boundaries
12-digit watershed boundaries
Pfafstetter codes
elevation grids

stream reaches
stream level
land-use, with updates for surface
mining


palustrine wetland polygons
lacustrine polygons
wetland and open water types by
grid cell
monthly snowfall
annual precipitation
minimum January temperature
Derived variables
sample units
watershed boundaries
watershed areas
upstream connections
main channel slope
main channel slope
upstream connections
percent land-cover



percent wetland + lake area
percent wetland + lake area
snowfall
mean annual precipitation
mean minimum January
temperature
Last update Scale Source
2003 1:24,000 USGS EROS Data
,..„ ... Center
(30 m gnd)



1999 1:100,000 USGS

2001 1:24,000 NLCD - EPA
(30 m grid) surface mining grid -
Tennessee Valley
Authority
Thru:
US EPA Region 3,
Wheeling, WV
1981 to 1:24,000 US FWS, WV DEP
present composite coverage
based on 1:24,000 Ohio Department of
1985-1987 Natural Resources
TM scenes
1998, based 1:250,000 Climate Data Source,
on 1961-1990 Corvallis, OR
climate data,
                                                                                          15

-------
       The PRISM (Parameter-elevation Regressions on Independent Slopes Model; Climate





Source, Corvallis, OR) GIS coverage was used to generate climatic summaries across watersheds.




PRISM is a continuous grid (raster coverage) of climatic variables generated from weather station




data and regional nonlinear regressions that take into account regional precipitation averages as well




as orographic effects associated with altitude, and proximity/influence of large surface water bodies




(Daly et all 994).





       We used National Wetlands  Inventory (NWI) coverages, supplemented by the Ohio




Wetlands Inventory coverages to define watershed storage. Total areal coverage of palustrine and




lacustrine polygons was added and divided by watershed area. For the Ohio Wetlands Inventory




coverage, polygons coded as 34 to 38 (wetland and open water cksses) were included.  For the New




York Wetlands Inventory, cover classes 1-4 and 9 were included (all classified or unclassified




wetland polygons).





       Fraction impervious surface  area was derived as a function of knd-use intensity as:





fimperv = (0.55*flwintrs)+(0.90*£hintres)+fcomind,





       where fimperv =  fraction impervious area in watershed,





              flwintrs = fraction low intensity residential area in watershed,





              fhintres = fraction high intensity residential area in watershed, and





              fcomind = fraction commercial, industrial, and transportation knd-use area in





              watershed.






based on median of estimated constructed materials associated with different land-use cksses in





NLCD (http://www.epa.gov/mrlc/dennitions.html). These weightings are higher than those
16

-------
reported for northern Virginia (NVPDC 1980), because the latter included only effective impervious




surface area (unconnected rooftops were not included), but are consistent with classification




guidelines for NLCD.






       Main channel slope was defined first by selecting main channel reaches from the NHD




coverage based on the minimum value for the reach LEVEL attribute  within a watershed.  Main




channel reaches were buffered by 30 meters, and this buffer was used to clip a subset of elevation




grids from the DEM coverages. Frequency analyses were used to estimate the 10th and 85th




percentile of elevation values from the stream elevation buffer files; combined with main channel




length, these could be used to calculate main channel slope (Figure 6).





Creation of sample frame and sample design for 12-digit HUCs





       Although we developed the classification framework for the entire state of West Virginia, in




practice, we restricted membership in the sample frame for the 2001-03 WV R-EMAP project.




Sample frames were defined based on drainage basin, ecoregion, and 12-digit HUC watershed region




size (Table 3).  The Potomac River basin was excluded from study because of species differences




related to  biogeography, and because the fish IBIs developed for adjacent regions in Maryland




should be adequate for the state of WV to apply in that basin (Figure 7). In practice, we also




excluded watersheds with ill-defined drainage networks related to karst topography (e.g., no single




outlet, underground streams). Selection of 12-digit HUCs from within the sample frame each year




was based on a random-stratified design. Strata were based on watershed classes defined by




hydrologic thresholds and knd-use intensity.  After defining 12-digit HUCs as single independent




watersheds (basin HUCs) or groups of hydrologically adjacent (dependent) interbasin units,
                                                                                          17

-------
                                HUC 060701010103
                                Region with NHD and
                                Shaded NED
/




c

"Ajf-WI. 1
0 Ot5 01






N
A
                         MUC 050701010103
                         Region wtm NMD nd
                         ShxtedNED
                                              0.8
                                              0.6
                                              0.4
                                              0.2
                                                      Derivation of main channel slope
                                                           85%
                                                    10%
   472 - 436 = 36 m
/ change in elevation
                    36 m elev/ 29.928 km
                    main channel length
                                                0
                                                420  430  440  450  460  470 480 490  500  510
                                                          Elevation (meters above sea level)
Figure 6. Automated derivation of main channel slope. The main channel is identified based on
finding stream reaches with the minimum "LEVEL" attribute in NHD (top left), then is buffered
and intersected with a digital elevation model (top right). A frequency analysis of output elevation
grid cells yields estimates of the 10th and 85lh percentiles. Combined with main channel length,
these can be used to calculate main channel slope.

-------
                                                                           Major Drainage Basins
Map Projection: Alters Equal Area Conic
                                                                            |    | 2-Digit HUCs
                                                                            |    | 4-Digit HUCs
HUC 02
02
05
HUC 04
0207
0208
0501
0502
0503
0504
0505
0507
0509
REGION NAME
Mid Atlantic Region
Ohio Region
SUBREGION NAME
Potomac
Lower Chesapeake
Allegheny
Monongahela
Upper Ohio
Muskingum
Kanawha
Big Sandy-Guyandotte
Middle Ohio
 Figure 7.  Major drainage basins associated with the state of West Virginia and surrounding states,
 as defined by 2- and 4-digit HUCs.
                                                                                                 19

-------
unequal probability-weighting was used to limit the chance of selecting subwatersheds directly

upstream or downstream of one another (Figure 8).


 Table 3. Definition of sample frame for R-EMAP wadeable stream sampling in
 West Virginia, 2000-2001.	

 Year   Ecoregion(s)                   Drainage basin              Watershed size (ha)

 2000   Central Appalachian Pkteau      inclusive                     400-40000

 2000   Central Ridge and Valley         all excluding Potomac River   400-40000
         Ecoregion                      basin

 2001    Western Allegheny Pkteau	inclusive	400 - 40000	


Response threshold derivation

       We carried out analyses and applied sampling designs separately by ecoregion for studies in

2000 and 2001 (Table 3, Figure 9).  To define hydrologic thresholds related to natural watershed

attributes, we identified gaging stations with long-term flow records and derived peak flow statistics

from previous studies in the state of West Virginia (Frye and Runner 1970; Runner 1980; Figure 4).

We delineated watersheds for each gaging station and characterized these watersheds based on

attributes recorded in Frye and Runner (1970) and Runner (1980).  We updated watershed storage

estimates from Frye and Runner's values using GIS analysis of National Wetlands Inventory

coverages and calculation of fraction watershed area covered by palustrine and kcustrine systems.

In earlier analyses  conducted before GIS was widely avaikble, Frye and Runner calcukted watershed

storage by overkying a point grid on USGS topographic quadrangles and counting intersections

with wetlands, lakes, and ponds. With a rare aquatic resource such as wetknds occurring in West

Virginia or other regions with well-developed drainage systems, the accuracy of such methods is

limited by the temporal accuracy of maps, the density of grid points used, and number of grid
20

-------
                                         Kilometers
                           25  50
                                       100
A
                                                                          Random-Stratified
                                                                          Selection Process
                                                                          Weighting Procedure
Map Projection Alters Equal Area Conic
                                                                            Legend
                                                                            |    | Independent 12-digitHUCs
                                                                            |    | Dependent 12-digit HUCs
                                                                            i     Central Appalachians
                                                                            Inclusion Probabilities
                                                                               | 0.013
                                                                                 0.018
                                                                                 0026
                                                                                 0.053
 Figure 8. Example of unequal weighting of interbasin HUCs for one watershed class in the Central
 Appalachian Plateau ecoregion of West Virginia.
                                                                                                21

-------
Map Pro)*cton Albert Equal Aiea Cone
                                                                                 US EPA
                                                                                 Level III Ecoregions
                                                                                  Level III Ecoregions
                                                                                  Name
                                                                                     Wo3»m Mogfwny Fbt.au
                                                                                     C«n0.il AppatichKuij
                                                                                     Cones! AppaLxrhoo Rrioo, and Vifay^
  Figure 9.  Map of ecoregions overlapping with West Virginia watersheds.
22

-------
points counted per unit area. For West Virginia, percent coverage had only been recorded to a




single significant unit (as 0 or 1).





       We applied two different methods to identify hydrologic thresholds, i.e., regions in plots of




area-normalized peak discharge (e.g., Q2/watershed area) showing nonlinear discrete shifts along a




gradient of watershed attributes. We obtained values for 2-year peak flows from Frye and Runner




(1970) and Runner (1980). We chose a recurrence interval of two years as an endpoint because this




is the frequency of flooding associated with development of channel morphology and sediment




delivery (Rosgen 1996).  Our first method, applied in 2000, involved screening potential parameters




of interest identified by Frye and Runner (1970; see Table 4) using Linear regression analysis of log-




transformed Q2/area, with observations weighted by period of record to account for differences in




level of uncertainty (Tasker and Stedinger 1989). We chose the best number and combination of




predictive variables using Mallow's Cp  statistic (SAS 1990). We identified thresholds through visual




graphical analysis, by plotting (nontransformed) Q2/area as a function of single variables or




combinations of variables and observing where nonlinearities in relationships were observed.




       In 2001, we applied a second technique, based on classification and regression tree (CART)




analysis of potential variables in SYSTAT (Wilkinson 1999). Nonlinear regression analysis followed




by visual analysis of graphic bi-plots was less successful in 2001 because of collinearities in




independent variables.  CART is a nonparametric technique that sequentially bisects observations




into subpopulations based on a single response variable, and identifies dependent variables




associated with breaks in the distribution.  In this way, CART is able to identify nonlinearities in




relationships, as well as interactions among predictor variables, based on a minimum set of





assumptions (Wilkinson 1999).
                                                                                         23

-------
 Table 4. Potential watershed characteristics related to peak flows, Frye and Runner (1970).
          Code
          A
          S
          St
          F
          L
          E
          S,
Units
mil2
feet -mil"1
          Si
          T,
Definition
watershed area
channel slope
storage (area in lakes and ponds +1)
forest area
main channel length
mean basin elevation
snowfall index
average annual precipitation - 20
precipitation intensity/24 hours, 2-year     inches
recurrence interval
soil infiltration index
average minimum January temperature     deg F
miles
feet above sea level
inches
inches
       We identified potential thresholds of response related to land-use change from the literature

because insufficient data were available from gaged sites with significant periods of record.  We

based potential land-use thresholds of response on previous analyses of data from the Maryland

Biological Stream Survey (MBSS, Boward et al. 1999). For the MBSS data set, species loss in fish

communities had been found at levels of impervious surface area above 2% of watershed area.

Nonlinear shifts in water quality were associated with levels of agricultural land-use above 25%.

Other studies have  shown a change in peak snowmelt associated with loss of mature forest cover

greater than 40-60% (Verry  1986). This approach is similar to trying to bracket expected response

diresholds in bioassays for suspected toxicants.

       Less information is available from the literature to detect thresholds of impact rekted to

mining activity in the Appalachian region.   Scott (1984) identified percent disturbed area as the best
24

-------
predictor of changes in stream hydrology in a subset of mined watersheds of West Virginia,




presumably associated with changes in water quality (Scott 1984).  However, Scott did not have long




hydrologic time series for these watersheds, but instead based his conclusions on modeled results;




his analysis showed an apparent threshold at 30% disturbed area.  We chose to use a conservative




estimate of 10% for the mining threshold, similar to that shown as a threshold for percent




impervious surface area in most literature reviews (Schueler 1994), with the knowledge that the




magnitude of mining thresholds would have to be refined in the future.




       We defined land-use coverage for watersheds associated with 12-digit HUCs by intersecting




watershed boundaries with a mining-modified NLCD coverage (Table 2, Figure 10). We defined




percent watershed storage by intersecting watershed boundaries with wetland inventory coverages




(Figure 11).
                                                                                          25

-------
                                                                                                Mining-Modified NLCD
                                                                                                Land-Use/Land-Cover
                                                                                                 Mining-Modified NLCD
                                                                                                 Grouped Land-Use Codes
                                                                                                       Mining
                                                                                                       Agriculture/Orchards
                                                                                                       Forest
                                                                                                       Grassland
                                                                                                       Urban/Residential
                                                                                                       Open Water/Wetlands
                                                                                                       Barren
    Map Projection: Albers Equal Area Conic
Figure 10. Map of land-use covering West Virginia watersheds, based on National Landcover Characterization Database, updated for mining category.
                                                                                                                                             27

-------
                       0   25  50
Map Projection: Abets Equal Area Conic
                                                                     Lakes and Wetlands
                                                                     in National Wetland
                                                                     Inventory (NWI)
                                                                       Legend

                                                                          12-digit HUCs
                                                                       NWI SYSTEM
                                                                          Lacustrine
                                                                          Palustrine
 Figure 11.  Map of palustrine and lacustrine classes from National Wetlands Inventory for the state
 of West Virginia.
                                                                                        29

-------
  CLASSIFICATION RESULTS





  Hydrology thresholds





       In year 2000, the nonlinear equation we derived to describe 2-year peak flows (Q2) for the




Central Appalachian Plateau and Central Ridge and Valley ecoregions, was:




       Q2 = Aa Pp T Snsn  St5'  (p < 0.05)




       where Q2 — 2-year peak flow




         A — watershed area




         P = average annual precipitation




         T = average minimum January temperature




         Sn = snowfall




         St = watershed storage.




         Percent forest was not included in the equation for 2-year peak flow, probably because of




a lack of significant variation in percent forest cover in these two ecoregions. When gaging station




data were combined for the full state, percent forest cover was included as a significant predictor




variable in peak flow equations for 5- and 10-year recurrence intervals, but not for 2-year recurrence




intervals. Visual graphical analysis showed that percent forest alone showed no threshold effect on




peak flows normalized to watershed area, but the product of percent forest cover, snowfall, and




minimum January temperature did show a noisy threshold with a lot of scatter.  The sharpest




threshold for Q2/A identified based on visual graphical analysis was for percent watershed storage




(Figure 12), with a threshold identified at ~0.5% storage (fraction storage = 0.005).




       We identified variables and associated threshold levels for Q2/A for the Western Allegheny




Plateau in 2001  through CART analysis. At each step, the CART procedure sequentially divided




observations into two categories based on the magnitude or category of the predictor variable (in








30

-------
           0    0.02  0.04  0.06   0.08   0.1   0.12   0.14  0.16

                   Fraction storage (lake + wetland/wshd)


 Figure 12. Predicted 2-year flood normalized to watershed area (cfs/mi!2) as a function of
 watershed storage (lake + wetland area/watershed area) for USGS gaging station watersheds in the
 Central Appalachian Plateau and Central Ridge and Valley ecoregions of West Virginia.


this case fraction storage or main channel length) that best separated the response variable,

Q2/watershed area, into low and high magnitude subpopulations. Each box in the output is

basically a frequency plot, with each observation represented by a colored dot, and each color

representing the final subclass assignment. Separation of the dataset based on combination of two

variables accomplished a 67% reduction in total variance  (Figure 13).  Flashy hydrologic regimes

(significantly larger Q2/A) for the Western Allegheny Plateau ecoregion occurred in low storage

watersheds (< 0.03% watershed storage)  or relatively short drainage basins (main channel length <

13.5 km).
                                                                                        31

-------
         Q2/watershed area
                                                                I
                                                      %storage<0.031
                          Main channel length (km)<13.5
Figure 13. Results of Classification and Regression Tree (CART) analysis of 2-year flood
normalized to watershed area (cfs/mil2) as a function of watershed attributes for USGS gaging
station watersheds in Western Allegheny Plateau ecoregion of West Virginia (percent reduction in
error = 67 %). At each step, the CART procedure sequentially divided observations into two
categories based on the magnitude or category of a predictor variable (in this case fraction storage
or main channel length) that best separated the response variable, Q2/watershed area, into low and
high magnitude subpopulations. Each box represents a frequency plot, with each observation
represented by a colored dot, and each color representing the final subclass assignment.

-------
  Distribution of land-use/ land-cover variables across 12-digit subwatersheds in West Virginia




       Overall, the level of watershed storage in West Virginia is very low. For all watersheds




associated with 12-digit HUCs in the 8-digit HUCs characterized (those within or overlapping with




WV state borders), the top 20th percentile of HUCs had a range of 0.8 to 18.6% watershed storage




(Figure 14). Most 12-digit HUCs in the top 20th percentile were actually located near the border or




even outside of WV state boundaries, in the headwaters of the Upper and Middle Ohio River basins.




Nonetheless storage thresholds associated with a marked change in hydrologic regime were low




enough so that both high and low storage classes could be identified within each set of ecoregions




(Figure 14).





       The range of percent impervious surface area in watersheds  associated with the 12-digit




HUCs characterized is also relatively low (0-18.4 %), with the top 20th percentile covering a range of




only 1.2 - 18.4 %. The spatial distribution of the top 20th percentile of percent impervious surface




area is similar to that of the top 20th percentile of watershed storage, paralleling development along




road or river corridors, and probably the creation of reservoirs (Figure 15). As for urban




development, the highest concentration of agriculture within West Virginia drainage basins (top 20th




percentile = 48-83 %) occurs outside of the state boundaries in the Upper Ohio or Potomac




drainages (Figure 16).  In contrast, the highest concentration of mining activity (top 20th percentile =




1.4-16.4% watershed area) occurs within the Central Appalachian Plateau, in the southwestern





portion of the state (Figure 17).  Forest  cover across the state is relatively high and constant; most




instate HUCs are in the top 40th percentile of HUCs characterized (79 - 100 % forest cover; Figure





18.)
                                                                                             33

-------
                                                                         Fraction
                                                                         Watershed Storage
                                                                         By 12-DigitHUC
                                                                           12-Digrt HUCs
                                                                           FractionStorage
                                                                               0.0000-0.0010
                                                                               00011 -00023
                                                                               00024 - 0.0042
                                                                           HH 0 0043 - 0 0082
                                                                           •• 00083-0.1860
klip Protection Albett tqujl Aiea Come
 Figure 14.  Map of 12-digit HUCs within 8-digit HUCs overlapping with West Virginia, coded by
 fraction watershed storage ([lake + wetland area/watershed area] x 100).  Classes mapped are not
 equal interval classes, but represent population quintiles.

-------
                                                                          Fraction
                                                                          Impervious Surface
                                                                          By12-DigitHUC
                                                                           12-Digit HUCs
                                                                           Fraction Impervious Surface
                                                                                0.0000 - 0.0008

                                                                                0.0009 - 0.0028

                                                                                0.0029 - 0.0053

                                                                                0.0054-0.0121

                                                                                0.0122-0.1837
Map Projection. Abers Equal Area Conic
 Figure 15. Map of 12-digit HUCs within 8-digit HUCs overlapping with West Virginia, coded by
 estimated fraction impervious surface area in associated watersheds. Classes mapped are not equal
 interval classes, but represent population quintiles.
                                                                                               35

-------
                                                                          Fraction Mining
                                                                          By12-DigitHUC
                                                                            12-Digit HUCs
                                                                            Fraction Mining
                                                                                 0.0000
                                                                                 00001-0.0012
                                                                                 00013-0.0064
                                                                            IB} 0.0065-0.0138
                                                                            Ij^B 0.0139-0.1644
Map Pto)*ctKm *Jb«f* Equal Aiea Cone
  Figure 16. Map of 12-digit HUCs within 8-digit HUCs overlapping with West Virginia, coded by
  estimated fraction surface mining area in associated watersheds.  Classes mapped are not equal
  interval classes, but represent population quintiles.
36

-------
                                                                           Fraction Agriculture
                                                                           By 12-DigitHUC
                                                                            12-Digit HUCs
                                                                            Fraction Agriculture
                                                                                 0.0002 - 0.0822
                                                                                 0.0823-0.1910
                                                                                 0.1911 -0.3140
                                                                               .  0.3141-0.4814
                                                                            ^^1 0.4815-0.8315
Map Projection: AJbers Equal Area Conic
 Figure 17.  Map of 12-digit HUCs within those 8-digit HUCs overlapping with West Virginia,
 coded by fraction agricultural area in associated watersheds.  Classes mapped are not equal interval
 classes but represent population quintiles.
                                                                                                 37

-------
                                                                          Fraction Forest
                                                                          By 12-DigitHUC
                                                                            12-Digrt HUCs
                                                                            Fraction Forest
                                                                                0.1400 - 0.4788
                                                                                0 4789 - 0.6605
                                                                                0 6606 - 0 7870
                                                                            IB 07871 -08937
                                                                            ^H 0 8938 - 0.9989
M«p Pioftctlon Albef* Equal Ai«« Come
 Figure 18.  Map of 12-digit HUCs within those 8-digit HUCs overlapping with West Virginia,
 coded by fraction forested area in associated watersheds. Classes mapped are not equal interval
 classes but represent population quintiles.
38

-------
       Because of the uneven distribution of watershed storage, watershed morphology, and land-




use across the state of West Virginia, the distribution of 12-digit HUCs is not even across watershed




classes; in fact some combinations of watershed classes are not even represented (Figures 19-20).




For example, in the Western Allegheny Plateau ecoregion, only one 12-digit HUC watershed




occurred in the high urbanization/low storage class, and no 12-digit HUC watersheds occurred with




combinations of moderate or high agriculture and low storage. Thus, not all potential watershed




classes can be compared within each ecoregion (Figures 21-22).
                                                                                           39

-------
Mdp Pto|*t.t)on Alb*ft Equal Ai»« Cone
                                                                               Aggregated Land-Use
                                                                               by Ecoregion
                                                                                 [    | HUC Bomdwy
                                                                                         Agriculture/Orchardi
                                                                                    BH Forest
                                                                                         Grasslands
                                                                                         Urban/fteskfartfal
                                                                                    H Open Water/Wetlands
                                                                                         Barren
                                                                                 Level III Ecoregions
                                                                                 Name
                                                                                      Western Allegheny Plateau
                                                                                      Central Appalachians
                                                                                      Central App. Ridges & Valeys
  Figure 19. Map of Omemik Level III ecoregions within state of West Virginia, and associated
  land-use summaries at ecoregion level for region within WV watersheds.
41)

-------
        Western
        Allegheny
        Plateau
         Central Appalachian
         Ridges & Valleys
                                                                                                                        12-DigitHUC
                                                                                                                        Watershed Classes
Watershed Classes
     • Low Impact/Low Storage
     j Low Impact/High Storage
      High Urbanization/Low Storage
      High Urbanization/High Storage
      Moderate Agriculture/Low Storage
      Moderate Agriculture/High Storage
      High Agriculture/Low Storage
      High Agriculture/High Storage
          Ecoregions
             Ecoregions
                                                                                                                          Level
                                                                                                                                Level
       Map Projection: Albers Equal Area Conic
Figure 20. Map of 12-digit subwatershed units within West Virginia's 8-digit HUCs, classified by land-use and watershed storage classes. Storage classes are based on empirical hydrologic thresholds for
watershed storage (CA, CARV ecoregions) or combinations of storage and main channel length (WAP ecoregion).
                                                     41

-------
Map Projection: Abere Equal Area Conic
                                                                          Random-Stratified
                                                                          Selection of HUCs
                                                                          for Sampling
                                                                           Legend
                                                                           |    | HUCs selected for sampling
                                                                           Thresholds
                                                                                Reference, low storage
                                                                           HI Reference, high storage
                                                                                High mining, low storage
                                                                                High mining, high storage
                                                                           BH High agriculture, low storage
                                                                                High agriculture, high storage
                                                                                High impervious, low storage
                                                                           i	| High impervious, high storage
                                                                           Level III Ecoregions
                                                                           Name
                                                                           |    | Central Appalachians
 Figure 21. Watershed classes associated with 12-digit HUCs in sample population for Central
 Appalachian Plateau and Central Ridge and Valley ecoregions in West Virginia.  Sample population
 included 12-digit subwatersheds with associated watersheds in size range of wadeable streams and
 outside of Potomac River basin.
                                                                                                43

-------
                                                                              Watershed Classes in
                                                                              Western Allegheny Plateau
    Map Projection. A!b« r_ Equal Area Conic
                                                                               Legend
                                                                                    Sample Ponts
                                                                               Watershed Classes
                                                                                    Low Impact/Low Storage
                                                                                    Low Impact/High Storage
                                                                                  | High Urbanization/Low Storage
                                                                                    High Urbanization/High Storage
                                                                                    Moderate Agriculture/Low Storage
                                                                                    Moderate Agncuttvjre/High Storage
                                                                                  | High Agriculture/Low Storage
                                                                                    High Agriculture/High Storage
                                                                               Level III Ecoreglons
                                                                               Name
                                                                                   1 Western Allegheny Plateau
  Figure 22. Watershed classes associated with 12-digit HUCs in sample population for Western
  Allegheny Plateau ecoregion in West Virginia. Sample population included 12-digit HUCs with
  associated watersheds in size range of wadeable streams.
44

-------
  FUTURE UPDATES AND CLASSIFICATION REFINEMENTS




       Potential future improvements to the watershed classification system for West Virginia




include the following:





       • Examination and documentation of differences in reference condition among low and




        high storage reference ("low impact") classes based on WV R-EMAP data, EMAP data




        from the MidAtlantic Integrated Assessment (MAIA), and macroinvertebrate and water




        quality data collected as part of the WV Department of Environmental Protection (DEP)




        rotating basin monitoring program





       • Examination and documentation of differences in condition between reference and




        "impacted" watershed classes within storage categories based on WV R-EMAP data,




        EMAP MAIA data, and macroinvertebrate and water quality data collected as part of the




        WV DEP rotating basin monitoring program




       • Examination and adjustment of thresholds of impact associated with watershed storage




        and knd-use intensity based on the analyses described above.




       • Potential refinement of watershed impact classes based on  empirical analysis of existing




        data with the addition of other watershed attributes such as development within stream




        and river corridor buffers, road density, and density of stream or river road crossings.




        Land-use impacts in a mountainous region might be better reflected by estimation of valley





        development, as compared to land-use intensity on a whole-watershed basis.




       • Improved estimation of percent impervious surface area in WV watersheds. The method





        used for estimation in this project was indirect.  However, the Chesapeake Bay program is





        currently mapping impervious surface area directly based on satellite imagery from Landsat




        7, IKONOS, and MODIS platforms (http://chesapeake.towson.edu/data/), and








                                                                                      45

-------
        additional studies are underway to calibrate percent impervious surface area estimation




        equations based on actual measurements from aerial photographs (David Jennings, US




        EPA, Reston, VA, personal communication).




      • Improved estimation of percent mined area (surficial disturbance). Current estimates were




        based on interpretation of TM + SPOT satellite imagery by  Tennessee Valley Authority




        under contract to Canaan Valley Institute for permitted mining areas (Hope Childers, US




        EPA Region III). The state of West Vkginia is currently constructing coverages of mining




        activity that might eventually allow better definition of area! extent and different sources of




        mining impacts (e.g., valley fill vs mountain top mining,




        http://www.nrac.wvu.edu/hollowfill/, http://gis.dep.state.wv.us/data/omr.html);




        however, current coverages are based on permit boundaries which do not always reflect




        mining activity.




       • Slight adjustments in boundaries and characterizations from Stage II to Stage III  10- and




        12-digit HUCs as Stage III EDNA process is completed for the state of West Vkginia.




       There are several ways in which derivation of thresholds could be improved, depending on




the objectives of a monitoring project. CART analysis objectively maximizes the separation of two




subpopulations of a dependent variable (e.g., Q2/A) as a function of the best predictive variable,




e.g., fraction watershed storage. Thresholds derived from graphical analysis could be verified




through application of a statistical procedure such as CART  (Wilkinson 1999). Reanalysis of USGS




data for the Central Appalachian ecoregions using CART yielded an apparent  threshold of fraction




storage = 0.001, or 0.1%, very close to the threshold of 0.5% derived through visual examination).
46

-------
       Alternatively, thresholds could be chosen by optimizing the chance of correctly identifying




watersheds with high normalized peak flows, i.e., those with greater potential for bank erosion and




transport of nonpoint source pollutants, while minimizing the chance of misidentifying watersheds




expected to have low normalized peak flows.  If watershed storage is used as the indicator of




normalized peak flows, then the probability of the two types of misclassification error will depend




on the strength of the relationship between normalized peak flows and watershed storage, the error




associated with the relationship between the two variables, and the underlying distribution of




watersheds across the range of watershed storage. The absolute value of the exponent term in the




equation relating storage and peak flows is a measure of the "strength" of the relationship.  An




exponent of zero would indicate that there is no relationship between storage and normalized peak




flows, while a highly negative value indicates a rapid rate of increase in peak flows as storage is




reduced (Figure 23).




       Several examples are  shown below to illustrate the combined effect of the magnitude of the




exponent (-2 versus -0.01), the magnitude of a desired breakpoint in Q/A relative to its underlying




distribution (median versus 95th percentile), and the effect of differing levels of error of prediction




(standard error of estimate = 0 to 100%; Figures 24a-h).  In these simulations, an even distribution




of watersheds across storage  cksses was used. In all cases there is a tradeoff between maximizing




the proportion of correct predictions for high normalized peak flows (Figures 24a,c,e,g) and




minimizing the proportion of incorrect predictions for the low normalized peak flow class (Figures





24b,d,f,h). As we reduce the chosen threshold for percent watershed storage, our ability to predict





watersheds with high normalized peak flows improves, but the misclassification error rate for the





low peak flow ckss also increases (e.g., compare Figure 24a with Figure 24b). Our chances for
                                                                                         47

-------
           Effect of exponent term for storage
                                                              — STA-5
                                                              — STA-2
                                                              — STA-1
                                                                 STA-0.5
                                                              — STA-0.2
                                                              — STA-0.1
                                                              — STA-0.01
                                                              — STA0
                       20     40     60     80    100    120
                            % watershed storage
        Figure 23.  Effect of magnitude of exponent term for storage on shape of
        relationship between watershed storage and normalized peak flow (Q/area =
        ST).

making correct high flow predictions and minimizing errors in predicting low flow watersheds arc

generally better the stronger the relationship is between storage and peak flows (e.g., compare Figure

24a with Figure 24c and Figure 24b with Figure 24d). The greater the error in the relationship

between storage and normalized peak flows (higher SEE), the worse our predictions will be.

(Compare different color lines within Figure 24a or 24b.) Flowever, prediction error appears to

have little influence if we wish to distinguish only those watersheds with the most extreme peak

flows (> 95lh percentile of Q/A; Figures 24e and f).  In this example, our ability to avoid making

misclassificanon errors  for low peak-flow watersheds is relatively small if our breakpoint for Q/A is

high (Figures 24f and h). If at the same time, the strength of the relationship is weak (e.g., exponent

= -O.U1), however, we also have little chance of correctly identifying the high peak-flow watersheds

(l-igurc 24g).
48

-------
a)
c)
  e)
 g)
Fraction LO storage dass with Q/A > median
. s s e s - s
Eft
I c<
1 in
4 «
| 1
i»
1"
| 0.4
S 0.2
I o
S.
Effect of standard error of estimate on
correct prediction of Q/A > threshold
in LO storage class, exponent = -2
K^V^
- SEE = 0%
- SEE = 1%
- SEE=10%
SEE=50%
- SEE=7S%
- SEE=100%
> 20 40 SO 80 100 120
Chosen threshold, % watershed storage
ect of standard error of estimate on
>rrect prediction of Q/A > threshold
LO storage class, exponent = -0.01
)hf^~~^ 	 ^ — • —
SEE=60%
- SEE=7S%
- SEE=100%
0 20 40 60 80 100 120
Chosen threshold, % watershed storage

| Eff
£ CC
< in
a «
1 i
1«
r
§04
I"
o
fo
| Eff
1 «
» in
i«
1 '
§0.8
•D
|«
£0.4
2 0.2
1
ect of standard error of estimate on
irrect prediction of Q/A > threshold
LO storage classes, exponent = -2
V 	
- SEE = 0%
- SEE = 1%
- SEE=10%
SEE=SO%
- SEE=75%
- SEE=100%
9 20 40 60 80 100 120
Chosen mreshold.% watershed storage
ect of standard error of estimate on
rrect prediction of Q/A > threshold
LO storage class, exponent = -0.01

SEE=SO%
- SEE=76%
- SEE=100%
1 ° 	
£ 0 20 40 60 SO 100 120
Chosen threshold, % watershed storage
                                               b)
                                               d)
                                               f)
Effect of standard error of estimate on
c incorrect inclusion of Q/A > threshold
I in HI storage class, exponent = -2
§ 0.7
f 0.6
S 0.6
I**
|0.3
V 0.2
V
- SEE = 0%
- SEE = 1%
- SEE=10%
SEE=60%
- SEE=76%
- SEE=100%
£ 0 20 40 60 60 100 120
Chosen threshold, % watershed storage
Eff
c in<
TO .
1 ln
f 0.8
5 °-7
|0.6
S 0.6
CD
« °A
i 0.2
ts o
C
ect of standard error of estimate on
correct inclusion of Q/A > threshold
HI storage class, exponent = -0.01
__^A
SEE=60%
- SEE=75%
- SEE=100%
9 20 40 60 80 100 120
Chosen threshold, % watershed storage

S Eff
8 in(
1 0.07
£ 0.06
1 0.05
•§ 0.04
| 0.03
£ 0.02
* 0.01
ect of standard error of estimate on
correct inclusion of Q/A > threshold
n HI storage class, exponent = -2
I
- SEE = 0%
- SEE = 1%
- SEE=10%
SEE=SO%
- SEE=7S%
- SEE=100%
2 0 20 40 60 80 100 120
Chosen threshold, % watershed storage

. Efl
g in
? in
< 0.07
1 0.06
S 0.04
f 0.03
a 0.02
S 0.01
I °
ect of standard error of estimate on
correct inclusion of Q/A > threshold
HI storage class, exponent = -0.01
/V\A
SEE=50%
- SEE=76%
- SEE=100%
0 20 40 60 80 100 120
Chosen threshold, % watershed storage
Figure 24. Effect of peak flow breakpoint (Q/A > median [a-d] versus Q/A >95th percentile [e-h]),
magnitude of exponent term in relationship, Q/area = STa (a = -2 [a,b,e,fj versus a = -0.01[c,d,g,h])
and prediction error on correct classification rate for high peak flow class [a,c,e,g] or
misclassification rate for low peak flow class [b,d,f,h].
                                                                                            49

-------
CONCLUSIONS




       In some regions such as West Vkginia, the watershed classification framework for




monitoring will be limited initially by available regional data and information on relationships




between land-use and hydrologic or biological response. However, the process is iterative and




classification schemes can be validated and improved upon based on initial monitoring results.  In




previous approaches to watershed classification, we have developed thresholds based on land-use




attributes calculated as a fraction of watershed area. Land-use variables described as a fraction of




watershed area might be less useful in montane regions where development is concentrated in




valleys; in these cases it might be more appropriate to express land-use as a fraction of stream buffer
zone area.
       Watershed storage was one to two orders of magnitude lower in West Vkginia than in Great




Lake regions for which we have derived hydrologic thresholds for storage.  However, the hydrologic




threshold for watershed storage was also an order of magnitude lower for WV as compared to the




upper Midwest, suggesting that storage is still a critical component in determining hydrologic




regime. Both removal and addition of watershed storage elements is affected by human activities.




Thus, it is possible that not all combinations of land-use and watershed storage classes will occur




within a region, as was found for the Western Allegheny Plateau ecoregion.  If examination of all




potential interactions between hydrologic and land-use classes is necessary or desired, we will need




to include data from adjacent states with a wider gradient of watershed conditions.




        The use of nonlinear regression analysis followed by visual analysis of bi-plots to identify




variables of interest for hydrologic thresholds is less objective and will be of limited use in cases of




collinearity of independent variables, or where classes have alternative definitions, e.g.,  Condition B




 < X OR Condition C < Y, as was observed for the Western Allegheny Plateau ecoregion. In these
 50

-------
cases, CART analysis should be a more powerful approach to define thresholds.  It is possible,




however, that graphical analysis will be more useful when interactive terms (Condition B < X and




Condition C < Y) are important in denning thresholds.  In the current analysis, visual graphical




analysis was used because of the form of the equation relating peak flows and watershed variables;




no inflection point can be defined in an exponential curve because the second derivative is constant.




Graphical analysis can be followed by CART analysis to confirm the placement of the threshold in




an unbiased fashion.









  ACKNOWLEDGMENTS




       We gratefully acknowledge contributions by Federal and contract staff of the USGS EROS




Data Center (Sue Greenlee, Kris Verdin, Jay Kost, Dean Tyler) under interagency agreement DW-




14938973 in developing Stage II and Stage III tools used for catchment aggregation and creation of




hydrologicaUy-corrected digital elevation models, without which this project would not have been




possible.  We also gratefully acknowledge the contributions of Sharon Batterman (US EPA MED-




Duluth, IAG Project Officer) and GIS staff of US EPA Region III, in particular, Don Evans and




Carmen Constantine, in helping to create Stage II and Stage III 10- and 12-digit HUC boundaries




for the state of West Virginia; and of Wendy Bkke-Coleman, US EPA Office of Environmental




Information, for facilitating this arrangement. We thank Hope Childers, contract staff at US EPA




Region III, Wheeling, WV office for providing the mining-modified NLCD land-use coverage for




West Virginia. We thank US EPA staff at Region III- Philadelphia, PA and Fort Meade, MD (Tom




DeMoss, Tom Pheiffer) and Wheeling, WV offices (Maggie Passmore, Frank Borsuk); state of WV




DNR and DEP staff (Dan Cincotta, WV DNR; Pat Campbell, Dave Montali, John Wilts, WV




DEP); and staff of the Canaan Valley Institute (Tom DeMoss, Ron Preston) for consultation on all
                                                                                       51

-------
aspects of this project We thank Debra Taylor at the US EPA MED-Duluth for help in data entry




and USGS gaging station watershed delineations. We gratefully acknowledge the GIS support




provided by OAO Corporation under US EPA FAIR Contract 68-W5-0065, Delivery Order #24;




and Computer Sciences Corporation under US EPA FAIR II Contract 68-W01/W02-032, Task




Order #024 in watershed characterization and map production. We also thank Peg Pelletier and




Jonathan Smith for providing helpful review comments on an earlier draft of this report.

-------
                                       REFERENCES

Boward, D.M., P.P. Kazyak, S.A. Stranko, M.K. Kurd, and T.P. Prochaska. 1999. From the
       Mountains to the Sea: The State of Maryland's Freshwater Streams. EPA/903/R-99/023.
       Maryland Department of Natural Resources, Monitoring and Nontidal Assessment
       Division, Annapolis, MD.

Brazner, J.C., D.K. Tanner, N.E. Detenbeck, S.L. Batterman, S.K. Stark, L.A. Jagger, and V.M.
       Snarski. 2003. Landscape influences on fish assembkge structure and function in western
       Lake Superior tributaries. Submitted to Environmental Management.

Cincotta, D. 2000. A Small Watershed Characterization, Classification, and Assessment for West
       Virginia Utilizing EMAP Design and Tools. R-EMAP proposal submitted to US EPA
       Mid-Continent Ecology Division, Duluth, MN and US EPA Region III, Philadelphia, PA.

Daly, C., R.P. Neilson, and D.L. Phillips. 1994. A statistical-topographic model for mapping
       climatological precipitation over mountainous terrain. Journal of Applied Meteorology,
       33:140-158.

Detenbeck N.E., D. Cincotta, J.M. Denver, S.K. Greenlee, and A.R. Olsen. 2003a. Watershed-
       based survey designs. Submitted to Environmental Monitoring and Assessment (accepted with
       minor modifications).

Detenbeck N.E., C.M. Elonen, D.L. Taylor, L.E. Anderson, T.M. Jicha, and S.L. Batterman.
       2003b. Effects of hydrogeomorphic region, watershed storage and mature forest on
       baseflow and snowmelt stream water quality in second-order Lake Superior Basin
       tributaries. Freshwater Biology 48(5):911-27.

Detenbeck, N.E., S.L. Batterman, V.J. Brady, J.C. Brazner, V.M. Snarski, D.L. Taylor, J.A.
       Thompson, and J.W. Arthur. 2000. A test of watershed classification systems for ecological
       risk assessment. Environmental Toxicology and Chemistry 19(4):1174-1181.

FGDC. 2002. Federal Standards for Delineation of Hydrologic Unit Boundaries.  Federal
       Geographic Data Committee. http://www.ftw.nrcs.usda.gov/huc_data.html

Franken, S.K., D.J.Tyler, and K.L. Verdin. 2001. Development of a National Seamless Database
       of Topography and Hydrologic Derivatives, Paper 730 in Proceedings, 2001 ESRI
       International User's  Conference, San Diego, CA.

Frye, P.M. and G.S.  Runner. 1970. A Proposed Streamflow Data Program for West Virginia, U.S.
       Dept. of the Interior, Geological Survey Water Resources Division, Open-file report,
       Charleston, WV.

Heiskary, S.A.  and C.B. Wilson 1990.  Minnesota lake water quality assessment report. Minnesota
       Pollution Control Agency, St. Paul, MN.
                                                                                      53

-------
 Herlihy A.T., D.P. Larsen, S.G. Paulsen, N.S. Urquhart, and BJ. Rosenbaum. 2000. Designing a
        spatially balanced, randomized site selection process for regional stream surveys: The
        EMAP Mid-Atlantic pilot study. "Environmental Monitoring and Assessment, 63:95-113.

 Jennings, M.E., W.O. Thomas, Jr., and H.C. Riggs. 1993. Nationwide Summary of U.S. Geological
        Survey Regional Regression Equations for Estimating Magnitude and Frequency of Floods
        for Ungaged Sites, 1993. U.S Geological Survey, Reston,VA. Water Resources
        Investigations Report WRI 94-4002.

 Johnston, C.A., N.E. Detenbeck, and GJ. NiemL 1990. Cumulative impacts of wetland loss on
        stream water quality and quantity. Biogeochemistry 10:105-141.

 Legleiter, KJ. 2001. Interagency Development of National Watershed and Subwatershed
        Hydrologic Units, Paper 492 in Proceedings, 2001 ESRI International User's Conference,
        San Diego, CA.

 NVPDC (Northern Virginia Planning District Commission). 1980. In Schueler, R.T. 1984.  The
        importance of imperviousness. Watershed Protection Techniques 1(3):100-111.

 Omemik, J.M. and A.L. Gallant. 1988. Ecoregions of the upper Midwest states. EPA/600/3-
        88/037, U. S. Environmental Protection Agency, Environmental Research Laboratory,
        Corvallis, OR.

 Pan, Y., R.J. Stevenson, B.H. Hill, and A.T. Herlihy. 2000. Ecoregions and benthic diatom
        assembkges in Mid-Atlantic Highlands streams, USA. Journal of the North American
        Benthological Society, 19:518-540.

 Richards, R.P. 1990. Measures of flow variability and a new flow-based classification of great lakes
        tributaries. J. Great Lakes Res. 16:53-70.

 Rosgen, D. 1996. Applied river morphology. Wildland Hydrology, Pagosa Springs, CO.

 Runner, G.S. 1980. Hydrologic Data for Runoff Studies on Small Drainage Areas, West Virginia
        Department of Highways, Research Project 16, Open-File Report 80-560, Charleston,
        WV.

 SAS. 1990. SAS/STAT Users Guide. Version 6, Fourth edition. SAS Institute, Gary, NC.

 Schueler, T.R. 1994. The  importance of imperviousness. Watershed Protection Techniques 1 (3):100-111.

 Scott, A.G. 1984.  Analysis of characteristics of simulated flows from small surface-mined and
        undisturbed Appalachian watersheds in the Tug Fork basin of Kentucky, Virginia, and
        West Virginia. U.S. Geological Survey, Reston, VA. Water Resources Investigations
        Report WRI 84-4151.
54

-------
Tasker, G.D., and J.R. Stedinger. 1989. An operational GLS model for hydrologic regression.
      Journal of Hydrology 111:361-375.

US EPA. 2001. The National Costs of the Total Maximum Daily Load Program (Draft Report).
       U.S. Environmental Protection Agency, Office of Water, Washington, DC EPA/841/D-
       01/003.

US EPA. 2002.  Consolidated Assessment and Listing Methodology: Toward a Compendium of
       Best Practices. First Edition. U.S. Environmental Protection Agency, Office of Wetlands,
       Oceans, and Watersheds. Washington, DC
       (http://www.epa.gov/owow/monitoring/calm.html).

US EPA. 2003.  National 303(d) List Fact Sheet.
       (http: / /oaspub.epa.gov/waters/national_rept. control).

Verdin, K.L., and J.P. Verdin. 1999. A topological system for delineation and codification of the
       earth's river basins. /. Hydrol. 218:1-12.

Verry, E.S. 1986. Forest harvesting and water: The Lake States experience. Water Resources Bulletin
       22:1039-1047.

Waite, I.R., A.T. Herlihy, D.P. Larsen, and D.J. Klemm. 2000.  Comparing strengths of geographic
       and nongeographic classifications of stream benthic macroinvertebrates in the
       Mid-Atlantic Ffcghlands, USA. Journal of the North American Onthological Society, 19:429^-41.

Wilkinson, L. 1999. Systat 9.0 Statistics I. SPSS Inc. Chicago, IL.
                                                                                      55

-------
SEPA
     United States
     Environmental Protection
     Agency
Please make all necessary changes on the below label,
detach or copy, and return to the address in the upper
left-hand comer.

If you do not wish to receive these reports CHECK HEREfJ ;
detach, or copy this cover, and return to the address in the
upper left-hand corner.
PRESORTED STANDARD
 POSTAGE & FEES PAID
         EPA
   PERMIT No. G-35



-------