# A framework for developing ecosystem-specific nutrient criteria: Integ rating biological thresholds with predictive modeling

Patricia A. Soranno 1 

### Department of Fisheries and Wildlife, Michigan State University, East La nsing, Michigan 48824



## Kendra Spence Cheruvelil, 2 

## R. Jan Stevenson, and Scott L. Rollins 3 

### Department of Zoology, Michigan State University, East Lansing, Michiga n 48824



## Sarah W. Holden and Sylvia Heaton Water Bureau, Michigan Department of Environmental Quality, Lansing, Mi chigan 48909



## Eric Torng Department of Computer Science and Engineering, Michigan State Universi ty, East Lansing, Michigan 48824

Abstract
We present a novel ecosystem-specific framework for developing nutrient criteria from biological thresholds and predictive modeling (BTPM) and an application of this framework to lakes i n Michigan, U.S. The four main components for the BTPM framework are: (1) to predict each ecosystem’s ‘ex pected nutrient concentration’ in the absence of human effects using a predictive model, (2) to identify importa nt biological thresholds along a nutrient gradient (i.e., biological [BIO] benchmarks), (3) to determine each ecos ystem’s current nutrient concentration, and (4) to use the above information to derive a nutrient criterion for each ecosyste m using the BTPM algorithm. The BTPM framework is extremely flexible in that it can be applied to any aquatic eco system type or nutrient and the four components can be implemented in a variety of ways. Our BTPM framework has two additional features: it recognizes that prior to human disturbance, ecosystems varied in their natural nutri ent concentrations, and it incorporates risk into the decision-making process. In the simplest scheme, a nutrient criterio n is set at a BIO benchmark greater than the expected nutrient concentration. However, to protect ecosystems more co nservatively, a criterion is set at current lake nutrient concentrations if current is less than the BIO benchmark. In our a pplication of the BTPM framework, we developed total phosphorus (TP) criteria for a diverse set of 374 lakes in M I. The expected lake TP concentrations in the absence of human effects ranged from 3 m g L

2 1 to 24 m g L

2 1 , suggesting that a single criterion approach would not be appropriate. We found two predominant benchmarks in the biological data a long the TP gradient, one for zooplankton metrics at 8 m g L

2 1 , and one for phytoplankton metrics at 18 m g L

2 1 . We present the sequence of analyses and decisions that could be used to apply this approach in a management context using Michigan lakes as an example.



### Although the primary role that phosphorus plays in lake eutrophication has been known for decades, setting management goals for phosphorus has been challenging. In the U.S., the measurement of nutrient levels in water bodies serves as the main basis for assessing progress towards meeting the goals of the Clean Water Act to protect the designated uses of water bodies (USEPA 2000; Reckhow et al. 2005). In the European Union, the Water Framework Directive adopted in 2000 mandates that ‘good ecological status’ be achieved in all waterbodies in European Union countries by 2015, but provides little guidance on methods for achieving this goal (European Union 2000; Søndergaard et al. 2005). In the U.S., states and tribes consistently identify excessive levels of nutrients as a major cause of why as much as half of the surface waters surveyed in the country do not meet water qualitystandards (USEPA 2000). The U.S. Environmental Protection Agency has mandated that states either adopt the ambient nutrient criteria developed by that agency, or develop scientifically defensible numerical criteria for 1 Corresponding author (soranno@msu.edu).

2 Present address: Lyman Briggs College, Michigan State University, East Lansing, Michigan 48824.

3 Present address: Environmental Studies, University of California, Santa Cruz, California 95064.
Acknowledgments Funding for this research was provided by a Michigan Department of Environmental Quality (MI-DEQ) grant to R.J.S. and P.A.S., an MSU–Michigan Agricultural Experiment Station Disciplinary grant to P.A.S., and an EPA grant to P.A.S., K.S.C., R.J.S., and others. We thank James Breck for providing some of the GIS base layers, Michael Belligan, Justin Booth, and David Lusch for GIS processing of the landscape databases, and Katherine Webster and two anonymous reviewers for helpful comments on earlier drafts. Jason Smith and Steve Holden provided valuable database management help. We also thank Carrie Scheele, Lesley Knoll, Aaron Jubar, Sherry Martin, and Mary Bremigan for access to their databases for the biological analyses, and Chelsea Stephens for zooplankton processing. Craig Stow and Jim Bence provided statistical guidance, and Howard Wandell provided help with the historical STORET water chemistry database. Finally, we thank the MI-DEQ Nutrient Criteria Technical committee and the MI-DEQ Nutrient Criteria Implementation Committee for many critical and insightful discussions during the development of this approach, and Candice Bauer from USEPA for suggestions and critical evaluation throughout the development of this approach. Limnol. Oceanogr., 53(2), 2008, 773–787 

### E 2008, by the American Society of Limnology and Oceanography, Inc.



### 773

### nutrients (phosphorus and nitrogen) that will protect the designated uses of water bodies (USEPA 2000). In fact, without nutrient criteria, a water body typically must reach nuisance algal conditions before an action is taken, which is often subjective at best. Despite this need, few states have established numerical nutrient criteria on a state-wide scale, and as of this writing, less than ten states have adopted numeric nutrient criteria for the majority of their lakes (C. Bauer, pers. comm.). Several approaches for setting nutrient criteria have been proposed, and an important assumption in each method is that a given nutrient criterion acts as an indicator of whether designated uses are being met (Reckhow et al. 2005). If we consider the ‘aquatic life’ designated use as measured by biological integrity, then the different approaches to setting nutrient criteria can be classified as either implicit or explicit in this regard. In an implicit approach, biological integrity is assumed to be protected at minimal human disturbance levels defined by some human disturbance gradient and associated nutrient value, but the biological integrity is not measured. In an explicit approach, changes in specific biological responses are used as a surrogate for the designated use more explicitly by examining how they change along a nutrient gradient either through analytical approaches, expert judgment, or some combination of the two (Stevenson et al. 2004; Reckhow et al. 2005). Surprisingly, biological response gradients have only sometimes been used for setting nutrient criteria (King and Richardson 2003; Heiskary and Wilson 2005; Reckhow et al. 2005). This situation may occur because many of the relationships between lake biology and phosphorus arelinear which provide little support for establishing a nutrient criterion. However, although not always explicitlyquantified, there are many examples of nonlinear relationships between biological responses and nutrients in lakes (Watson et al. 1997; Jeppesen et al. 2000; Downing et al.2001; Vadeboncoeur et al. 2003) and in wetlands (Stevenson et al. 2002; King and Richardson 2003; Lougheed et al. 2007). In cases where thresholds can be clearly quantified, a nutrient criterion can be justifiably set at or below that threshold value (Stevenson et al. 2004). Although there has been a wealth of research on eutrophication, its causes, and its associated effects on aquatic foodweb dynamics (Smith 1998), there are still few integrated approaches available for setting nutrient criteria in aquatic ecosystems that incorporate these insights into a framework that can be easily applied at the state (or national) level by agency managers and biologists. Such a framework must obviously strike a balance between model and analytical complexity and implementationpracticality in a management context. Management agencies are charged with managing often thousands of aquaticecosystems within social, political, and logistical constraints such as limited funding, limited data availability and limited personnel for implementation. Our goals were to develop a framework that can work within suchconstraints by integrating biological thresholds and predictive modeling (BTPM) of nutrient expected conditions, and to present a practical implementation of the framework that can be readily adopted by management agencies. In our approach, the biological thresholds can be viewed as an indicator of whether designated uses are being met, which are integrated into the overall criteria development process, and which in different applications can take different form.

Overview of the BTPM framework— Our BTPM framework integrates nutrient modeling to predict ecosystemspecific expected condition (using any model form), biological thresholds (i.e., nonlinear biological responsesalong a nutrient gradient), and current nutrient concentrations into an algorithm to derive ecosystem-specificnutrient criteria. Our framework has four main components (Fig. 1). The framework is flexible in that it can be applied to any aquatic ecosystem type or nutrient, it canincorporate any approach to predict an individual ecosystem’s expected condition or to identify thresholds, and these data can be integrated into an algorithm to determine nutrient criteria in a variety of different ways. Our framework has two additional important features. It recognizes that prior to human disturbance, water bodies differed in their natural nutrient concentration as a result of their hydrogeomorphic (HGM) features, which should be reflected in the establishment of nutrient criteria, and it incorporates risk into the decision-making process tochoose biological thresholds. Application of this framework requires data on as many ecosystems as possible (to obtain a wide range in HGM features). Ideally, these ecosystems would be sampled weekly or monthly during an important index period (such as summer) and over the course of multiple years. However, in Michigan and many other states, there are few lakes with monitoring data from multiple time periods within or across summers. Our framework can accommodate data from single time periods or multiple time periods. Similarly, if lacking the necessarydata on relevant biological responses of interest, informaFig. 1. The four main components of an ecosystem-specific framework for developing nutrient criteria. 

### 774 Soranno et al.

### tion from the literature can be used, preferably in combination with sampled data. Below, we present an application of our BTPM framework.

Application of the BTPM framework to Michigan lakes— We applied our BTPM framework to a diverse set of 374 lakes in Michigan, in consultation with the Michigan Department of Environmental Quality (MI-DEQ), the statemanagement agency responsible for setting nutrients standards for Michigan’s aquatic resources. We present results for lake total phosphorus (TP) criteria. A key advantage of our Michigan lake application is that these individual lake criteria can be derived efficiently given standard and widelyavailable limnological monitoring data sets including variables such as: lake morphometry, lake chemistry (total nutrients), water color, and geographic information system (GIS) data that are becoming widely available at the state (and national) level such as geology, ecoregion, and land use/ cover (LULC). Next we describe the main steps that we used in this application (Fig. 2). There were two steps to predict the lake-specific expected TP (also called lake reference condition). First, we built a hydrogeomorphic-land use model (HGM-LU) that models lake TP as a function of commonly used predictors such as HGM features and human LULC (Norvell et al. 1979; Ha  ̊kanson 1996). The predictive model for this step can take a variety of different forms; we used a multiple regression approach of cross-sectional data on several hundreds of lakes. Because we did not have phosphorus loading estimates for these lakes, we used variables that are related to phosphorus loading. In the second step, we used the HGM-LU model to predict expected lake TP (EXP 0 

### ) inthe absence of human LULC effects through backcalculation of each lake by setting the human LULC coefficients to zero (Dodds and Oakes 2004; Baker et al. 2005). EXP 0 

### is the expected condition for each lake defined as the least disturbed condition given the state of today’s landscape (Stoddard et al. 2006). We then calculated expected condition with allowance (EXP A 

### ) by adding an ‘allowance’ to EXP 0 

### . The allowance can serve any numberof different purposes such as representing model uncertainty in the estimation of EXP 0 

### , or some low or minimal level of allowable human disturbance to the lake (Fig. 2). The allowance should be chosen based on the goals of the nutrient criteria, the definition of expected condition, ecological justification, or some estimate of prediction uncertainty, all of which are likely to vary from state to state. We chose to represent model uncertainty. To identify biological thresholds, we analyzed lake biological responses (such as zooplankton, phytoplankton and macrophyte metrics) along a TP gradient to identify critical thresholds where major changes in lake biology occur (King and Richardson 2003; Lougheed et al. 2007), which we define as BIO benchmarks (Fig. 2). These benchmarks were used to define BIO zones, which are defined as the TP values between two BIO benchmarks (as well as between the lowest benchmark and the lowest TP concentration, and the highest benchmark and the highest TP concentration in the dataset). We determined current TP by using existing sampling data from state databases. The fourth component was to integrate each lake’s expected TP, current TP, and the BIO benchmarks to derive a nutrient criterion for each lake using the BTPM algorithm (Fig. 2). The BTPM algorithm is made up of four rules for setting individual lake criteria (Fig. 2): (1) If modeled expected condition (EXP A 

### ) is greater than current TP, then set criterion at current TP. This rule is particularly important because it ensures that there is no further degradation in lakes that have lower TP than is predicted from the model. (2) If EXP A 

### is less than current TP, but is in the same BIO zone as current TP (but not the highest BIO zone), then set criterion at current TP. (3) If EXP A 

### is less than current TP, but in the highest TP BIO zone, then set criterion at EXP A 

### . (4) If EXP A 

### is in a lower TP BIO zone than current TP, then set criterion at the higher TP BIO benchmark for the zone EXP A 

### is in. The integration of the BIO zones with the EXP A 

### recognizes that model predictions have some uncertainty associated with them. Our application of the BTPM framework to Michigan lakes is based on six key assumptions that we evaluated and developed during our criteria development process in a work-group setting with researchers, managers, and biologists. First, we assumed that phosphorus is the main ‘stressor’ to lakes. Second, we assumed that there is important natural variation in TP among lakes becauseof HGM features such as catchment and lake morphometry, geology, and regional scale factors such as climate, which can be modeled. Third, we assumed that benchmarks should be established to sustain desired levels of one or more biological attributes of lakes which are related to a lake’s designated use. Fourth, we assumed that importantlake biological responses should include integrative measures of lake biology from both the pelagic and littoral zones and lake water clarity, which is related to lakebiology through phytoplankton biomass. Fifth, we assumed that human disturbance of a lake can be reasonably approximated as the proportion of human LULC in the lake catchment. Although LULC does not represent all of the possible human effects on lakes, it is commonly used as an ‘indicator’ of human effects and is the most widely available dataset for this purpose. Finally, because lakes in the U.S. are managed at the state level, we chose the state as an appropriate spatial scale to build models.



### Methods

Databases— Water chemistry data: For water chemistry, we used data for 374 natural lakes ( $ 0.20 km 2 

### ) with publicaccess sampled by the Michigan Department of Environmental Quality (MI-DEQ) from 1975–1982 to match the time period for our LULC data (1978–1985). We obtained the data from the Storage and Retrieval database (STORET). Natural lakes are those with or without a water control structure on them, but that had an open-water basin prior to human modification. Data were collected from each lake on a single date during the summerstratified season (July–September; Table 1) for TP, alkalinity, and true water color (Table 1).HGM-LULC data: We created a GIS-based lake HGMLULC database for all 374 lakes (Table 2). Lake morEcosystem-specific nutrient criteria 775Fig. 2. The four main components in our application of the BTPM framework fo r Michigan lakes. LULC is land use/cover, morph is morphometry, and geol is geology. For the fourth component, the letters represent all of the eight possible relationships among lake expected condition (EXP A ), current condition (CUR), and the two biological response (BIO) benchma rks for the 374 Michigan lakes. The

‘ % of lakes’ is the number of lakes in the Michigan lake dataset that fall withi n each of the eight possible relationships for the BTPM algorithm. 

### 776 Soranno et al.

### phometry was measured as: mean depth, maximum depth, shoreline development factor (the ratio of shoreline perimeter divided by the circumference of a circle of the same area), and lake basin slope [(surface area) 1/2 

### /mean depth] (Table 2). Catchment areas were delineated as cumulative catchments that include all upstream stream and lake catchments. Stream length included the length of all upstream streams directly connected to lakes. Bedrock geology data was grouped into five categories: carbonate, clastic, hard rock, salt, and iron (Geologic Survey Division, MI-DEQ). Surficial geology data was grouped into five categories: dune sand, glacial till, lacustrine, moraine, and outwash (Michigan Natural Features Inventory and Michigan Department of Natural Resources [MDNR]). We obtained LULC data from the Michigan Resource Information Service (MDNR), where urban, agriculture,forest, upland field, and open water classes were determined from aerial photographs taken between 1978 and 1985 at a resolution of 0.025 km 2



### .

We quantified the proportions of LULC, bedrock geology, and surficial geology in the 500-m buffer around each lake, which we use as an indicator of the lake’s catchment. Although these 500-m buffers are not true catchments, analyses on these Michigan lakes (Soranno, unpubl. data) have shown that LULC from the 500-m buffer is highly correlated ( p , 0.001) to LULC of the lake’s catchment ( % urban, r 5 0.63, % agriculture, r 5 0.83), and these data are readily available for all lakes in a region with a GIS coverage. Biological data: We combined biological response datasets from six recent studies of lakes sampled from 1998–2004 (Table 3). To maximize sample sizes, we combined data across studies when lakes were sampled using the same methods. Data were based on single visits to a lake with the exception of the study by Bremigan et al. (2005), from which we had data for six lakes and two years of contrasting TP levels. The biological responses examined included: phytoplankton biomass and clarity metrics, phytoplankton community and toxin metrics, micro- and macrozooplankton community metrics, and macrophyte cover metrics. (For a complete list of metrics, see Tables 6 and 7 and see citations in Table 3 for further details on methods and variables.)

Step 1: Predict lake-specific expected TP— We built the HGM-LU model that predicts lake epilimnetic TP from HGM variables, limnological variables that influence TP (water color and alkalinity; Chow-Fraser and Duthie 1983; Ha

 ̊kanson et al. 2005), ecoregion, and LULC (Tables 1 and 2). To build the model, we used stepwise multiple regression (both forward and backward) and included all predictor variables that were significantly univariately correlated to TP (correlation coefficient, p # 0.10). We transformed variables as needed to meet assumptions of normality and homogeneity of variances. We chose the best model based on the Akaike information criterion (AIC) because it takes into account the number of parameters in the model and penalizes models with more parameters; models with AIC values of 10 units lower than other models were assumed to have the most support (Burnham and Anderson 2002). We took several additional steps to minimize the potential multicollinearity problems associated with using multiple linear regression. First, we conducted both forward and backward stepwise regression and used a conservative significance level ( p # 0.05) to retain predictor variables inthe models. Second, we quantified the level of multiTable 1. Median, minimum, and maximum values for the lake chemistry variables. n 5 374 lakes for all variables except alkalinity ( n 5 372).

Lake chemistry Median Minimum Maximum
Total phosphorus ( m g L

2 1 ) 13.0 1.0 155 Alkalinity (mg L

2 1 CaCO 3 ) 105 1.0 225 Water color (Pt Co) 9.0 1.0 99 Table 2. Median, minimum, and maximum values for the hydrogeomorphic (HGM) and land use/cover (LULC) features. n

5 374 for all variables. Stream length is the length of upstream streams connected to the given lake. Precipitation is the 30-yr annual average. Drainage area is the catchment area plus the lake area. Shoreline development factor, all ratios, and lake basin slope are unitless. Land use/cover percentages, geology percentages and precipitation were calculated on the 500-m buffer surrounding each lake.

HGM features and LULC Median Minimum Maximum
Lake and catchment morphometry

Lake area (km 2 ) 0.61 0.20 75.8 Mean depth (m) 4.5 0.8 42.5 Maximum depth (m) 12.8 3.0 86.9 Shoreline development factor 1.8 1.0 6.3 Lake basin slope 194 46.4 1,825 Catchment area (km 2 ) 7.1 0.21 1,948 Catchment area : lake area 10.0 0.7 2,655 Drainage area (km 2 ) 8.4 0.44 1,949 Drainage area : lake area 11.0 1.7 2,656 Stream length (km) 1.4 0 2,805 Climate Precipitation (mm) 835 706 1,008 Surficial geology

% dune 0 0 53 % outwash 30 0 100 % moraine 0 0 100 % exposed bedrock 0 0 100 % peat and muck 0 0 71 % lacustrine 0 0 100 % glacial till 0 0 100 Bedrock geology

% carbonate 0 0 100 % clastic 100 0 100 % hard rock 0 0 100 % salt 0 0 100 % iron 0 0 100 LULC

% agriculture 6 0 81 % urban 13 0 87 % forest 50 1 98 % upland vegetation 6 0 54 % wetland 6 0 49 % open water 1 0 17 

### Ecosystem-specific nutrient criteria 777

### collinearity in the predictor variables in the final model by calculating variance inflation factors and using a cutoff of 30 (Hair et al. 1995), above which predictor variables are assumed to suffer from high multicollinearity. For the best model, we calculated partial R 2 

### ’s for each predictor variable to assess the worth of collecting each piece of information for practical relevance and application of the model. For the best HGM-LU model, we also tested whether including an ecoregion framework explained additional variation in lake TP. Because ecoregions are also correlated to a variety of HGM and LULC variables (Wickham et al. 2005), we examined the residuals of the best HGM-LU model and performed an ANOVA with ecoregion as the predictor. We used Ecological Drainage Units (Higgins et al. 2005) as our ecoregion because it has been shown to partition variation in TP for Michigan lakes (Cheruvelil et al. in press). Using the best HGM-LU model, we then calculated EXP 0 

### by setting the model coefficients for agriculture and urban LULC to zero. We calculated EXP A 

### by adding the allowance to EXP 0 

### , which we estimated by calculating the upper 75th confidence interval from the HGM-LU model to account for model uncertainty.

Step 2: Identify biological thresholds along TP gradient— We examined BIO response variables that are commonly sampled by agencies, are related to attributes that the public value (e.g., algal biomass influences water clarity) and are good indicators of ecological integrity in lakes (e.g., zooplankton are important component of lake pelagic foodwebs). We quantified the relationship between each BIO response and TP in two ways; we used linear regression on transformed variables (where necessary) to quantify linear relationships and regression tree analysis (RTA) on untransformed variables to quantify nonlinear relationships using Systat 11.0 (Systat Software, Inc.). RTA is a recursive data partitioning algorithm that initially splits the response variable into two subsets based on the value of the predictor variable (TP) that maximizes the reduction in total residual sum of squares from the parent group to the two daughter groups (Breiman et al. 1984). Each node is then split using the same algorithm. At each split, we examined the ‘proportional reduction in error’ (PRE), which is a goodness of fit statistic and similar to an r 2 

### value. We considered that there was scientifically supportable evidence for a shift in the BIO response if the PRE was

$ 0.30 (moderate-strong) and if the PRE was greater than the r 2 

### from the linear model. For these moderate-strong relationships, we quantified uncertainty by defining the ‘TP benchmark risk’ using a bootstrap simulation (Qian et al. 2003). We resampled the data (with replacement) 1,000 times and calculated the single strongest threshold for each simulated dataset. The ‘TP benchmark risk’ was defined as the cumulative probability of a threshold based on the relative frequency of each threshold in the distribution (King and Richardson 2003), in other words, the proportion of bootstrapped thresholds that are equal to or lower than a given TP threshold. The bootstrapped RTA was conducted in the R statistical environment using the rpart and bootstrap libraries (R Development Core Team, 2007). We then plotted cumulative probability curves for each metric, which described the cumulative risk of a change in the response variable (BIO) with changing TP values (Qian et al. 2003).To inform our decision of choosing the BIO benchmarks, we calculated BIO response factors by first calculating the median value of each BIO response variable above and below a BIO threshold, and then dividing the larger median by the smaller median. We quantified the BIO response factor at three risk levels: 1 % , 10 % , and 50 %

. We assumed that there was an ecologically important change in the BIO response when the medians above and below a threshold differed by a factor of two (although any factor can be chosen). Because the BIO response can either increase or decrease with increasing TP dependent on the relationship (e.g., Secchi depth decreases and algal biomass increases with increasing TP), we calculate the factor by dividing the larger median by the smaller regardless of whether the BIO response increased or decreased. We use these BIO response factors to help identify the point at Table 3. The biological datasets used for the biological response threshold analy sis. n refers to the number of lakes in each dataset, years refers to the year that the lakes were sampled in, and TP range is the mi nimum and maximum TP value in the dataset.

Biological datasets Response variables n Years TP range ( m g L

2 1 )
(a) Martin and Soranno 2006 Chl a 64 2003 2.6–34.0 Secchi depth 63 (b) Cheruvelil and Soranno. In press. Macrophyte cover 53 2000–2001 2002–2003 4.4–65.6 4.8–64.0 (c) Knoll 2004 Chl a 64 Secchi depth 64 Extinction coefficient 64 Phytoplankton 56 Microzooplankton 46 (d) Scheele unpubl. data Macrozooplankton 15 2002–2003 5.3–19.3 (e) Bremigan et al. 2005 Chl a 13 1998–1999 11.8–32.3 Secchi depth 13 Macrozooplankton 13 (f) Jubar 2004 Secchi depth 15 2004 10.0–24.6 Macrozooplankton 15 

### 778 Soranno et al.

### which an important change in the BIO response has occurred, which we define as a BIO benchmark and is one of the inputs into the BTPM algorithm.

Step 3: Determine current TP in each lake— We used the existing sampling data from the 374 Michigan lakes that were sampled during one time period in the summer ( see Databases for details ).

Step 4: Derive lake-specific nutrient criteria using the BTPM algorithm—We integrated the results from step 1 to step 3 into the BTPM algorithm to develop lake-specific nutrient criteria (see Fig. 2 for the graphical representation of this algorithm in step 4). The inputs to the algorithm are the measured TP concentration for each lake (current), the expected lake TP concentration with allowance (EXP A 

### ), and the BIO benchmarks. We define ‘BIO zones’ as TP values between any two adjacent BIO benchmarks, and between a BIO benchmark and either the maximum or minimum TP value in the dataset. The BTPM algorithm is made up of the four rules stated earlier that defines all the possible relationships among these variables ( see Application of the BTPM framework to Michigan lakes ).



### Results

Step 1: Predict lake-specific expected TP— The lakes used to build the HGM-LU model ranged in TP from 1 m g L

2 1 

### to 155 m g L

2 1 

### , with a relatively low median of 13 m g L

2 1 

### (Table 1). The lakes also ranged widely in alkalinity and water color (Table 1). Of the many HGM variables univariately related to lake TP (Table 4), the strongest included common predictors of TP such as lake depth, LULC, catchment morphometry, precipitation, and water color. The best HGM-LU model included variables known to be correlated with P loading to lakes (Table 5, Fig. 3). Although significant, outwash geology only adds a small amount to the overall R 2 

### , therefore, if a management agency were to use this model, they could choose whether to retain this variable based on cost and effort of including it. To assess the importance of ecoregions for predicting reference TP, we performed an ANOVA with the residuals from the best model as the response and Ecological Drainage Units as the predictor. Because ecoregion was not significant ( p . 0.05), it was not incorporated into the final model. Using the best HGM-LU model, we estimated EXP 0 

### and EXP A 

### for each lake. The upper TP range for current lake condition exceeded the upper range for both EXP 0 

### and EXP A 

### as expected (Fig. 4). In the absence of human disturbance (EXP 0 

### ), TP in our study lakes ranged from 3 m g L

2 1 

### to 24 m g L

2 1 

### , demonstrating the wide range in TP expected condition. When incorporating an allowance (EXP A 

### ), TP expected condition ranged from 4 m g L

2 1 

### to 33 m g L

2 1



### . Step 2: Identify biological thresholds along TP gradient— There were 36 BIO response relationships with TP as the predictor. Fourteen of them had a PRE $ 0.20, and seven of those had a PRE $ 0.30 (moderate–strong) and a higher fit than the r 2 

### from linear regression models (Tables 6 and 7). We grouped these BIO responses into two ecological classes: 1) phytoplankton and water clarity and 2) zooplankton. Figure 5 shows the raw data for eachmoderate-strong BIO response relationship and the location of the RTA threshold(s). This figure also depicts the cumulative frequency plot of the 1,000 bootstrapped thresholds. The interquartile ranges of the bootstrapped thresholds ranged from 5 m g L

2 1 

### to 37 m g L

2 1 

### TP for phytoplankton and clarity metrics, and from 6 m g L

2 1 

### to 13 m g L

2 1 

### TP for zooplankton metrics (Table 7). Because the thresholds calculated using RTA were high relative to the range of bootstrapped thresholds (Fig. 5 and Table 7), we suggest caution in using RTA values for benchmark Table 4. Univariate correlations between hydrogeomorphic (HGM), land use/cover (LULC), and water chemistry predictor variables and TP. Units and descriptions are the same as for Table 2. Correlations were performed on transformed variables when necessary. Correlation coefficients that are significant at alpha # 0.10 are in bold and were included in the predictive modeling step. n 5 374 for all variables, except for alkalinity ( n

5 372).

Predictor variables r p Lake and catchment morphometry

Lake area (km 2 )

2 0.050 0.336 Mean depth (m)

2 0.292 , 0.001 Maximum depth (m)

2 0.301 , 0.001 Shoreline development factor 0.044 0.394 Lake basin slope 0.199 , 0.001 Catchment area (km 2 ) 0.107 0.038 Catchment area : lake area 0.160 0.002 Drainage area (km 2 ) 0.096 0.063 Drainage area : lake area 0.093 0.073 Stream length (m) 0.101 0.051 Climate Precipitation (mm) 0.218 , 0.001 Surficial geology

% dune 0.095 0.066 % outwash

2 0.097 0.062 % moraine 0.107 0.039 % exposed bedrock

2 0.032 0.531 % peat and muck 0.053 0.310 % lacustrine

2 0.062 0.233 % glacial till 0.052 0.313 Bedrock geology

% carbonate

2 0.104 0.044 % clastic 0.171 0.001 % hard rock

2 0.077 0.137 % salt 0.013 0.799 % iron

2 0.083 0.109 LULC

% agriculture 0.287 , 0.001 % urban 0.148 0.004 % forest

2 0.318 , 0.001 % upland vegetation 0.093 0.072 % wetland 0.007 0.894 % open water

2 0.058 0.265 Water chemistry

Alkalinity 0.084 0.131 Water color 0.410 , 0.001 

### Ecosystem-specific nutrient criteria 779

### identification without consideration of the uncertainty or risk associated with those values. We defined risk as the proportion of bootstrapped thresholds that were equal to or lower than a given TP threshold, which we incorporated into the next step that identifies the BIO benchmarks. To choose the BIO benchmarks for the BTPM algorithm, we used our estimates of risk in the seven moderately-strong BIO threshold relationships at three levels (1 % , 10 % , and 50 % ) and the BIO response factors (Table 8). For phytoplankton and water clarity metrics, medians for both chlorophyll a (Chl a ) and phytoplankton biomass differed by more than a factor of two (i.e., biologically important) at the 1 % risk level, which was exacerbated at higher risk values. For extinction coefficient, the medians differed by a factor of two at the 10 % risk level. In contrast, medians for Secchi depth never differed by a factor of two. This suggests that although the thresholds derived for Secchi depth werestatistically significant, they were not biologically meaningful. In addition, the PRE for Secchi depth was the lowest of all of the moderate-strong relationships (0.32, Table 7), and in fact, was very similar to the r 2 

### for the linear regression Table 5. Multiple regression results for the best hydrogeomorphic-land use model (HGM-LU) for predicting lake TP ( R 2 5

2 0.39, p , 0.001, F 5, 368 5 46.6). SE refers to the parameter estimate standard errors and partial R 2 ’s provide the contribution of each variable to the overall model.

Predictor variables Parameter estimate SE Partial R 2

Constant 1.80 0.138 Mean depth

2 0.24 0.050 0.05 Surficial geology

( % outwash)

2 0.18 0.065 0.02 Water color 0.36 0.037 0.18 % agriculture 1.10 0.146 0.13 % urban 1.24 0.193 0.10

Fig. 3. Predicted versus observed (current) TP for the besthydrogeomorphic-land use model. See Table 5 for model description. Fig. 4. Frequency histograms of lake TP for (A) current TP, (B) the expected TP with allowance (EXP A ), and (C) the expected TP without allowance (EXP 0 ). The x -axis ends at 75 m g L

2 1 , thus two lakes are not shown at 118 m g L

2 1 and 155 m g L

2 1 for (A). 

### 780 Soranno et al.

### (0.29). Because the nonlinear relationship was not well supported by the data, thresholds for Secchi depth did not contribute towards setting criteria in Michigan lakes. For the zooplankton metrics, the difference in the medians above and below the benchmarks was biologically important for proportional biomass of Daphnia and zooplankton biomass at the 1 % risk level. Based on these analyses, we identified two BIO benchmarks, one associated with phytoplankton metrics, and the other associated with zooplankton metrics. BIO benchmarks must be chosen with care to minimize the risk of losing biological integrity in an ecosystem by setting a benchmark either too high or too low. Although choosing the lowest possible BIO benchmark may seem protective and conservative, this may not actually be the case once the benchmarks are incorporated into the BTPMframework. For example, consider what happens in a twobenchmark algorithm with 20 m g L

2 1 

### TP as the upper BIO benchmark and either 10 m g L

2 1 

### or 12 m g L

2 1 

### TP as the lower benchmark. If we use 12 m g L

2 1 

### as the lower BIO benchmark, lakes with EXP A 

### of 11 m g L

2 1 

### will be managed to maintain a benchmark of 12 m g L

2 1 

### . Although choosing 10 m g L

2 1 

### as the lower BIO benchmark seems as if it would be more protective, doing so means that lakes with EXP A 

### of 11 m g L

2 1 

### will be managed for the upper BIO benchmark of 20 m g L

2 1 

### . As a result, choosing excessively low or high benchmarks may not be prudent. Therefore, we chose intermediate benchmark values, using risk levels from either 1 % or 10 % (depending on the BIO response factor) for thresholds that were associated with an important ecological change. For phytoplankton metrics, we chose 18 m g L

2 1 

### (range from 16.0 m g L

2 1 

### to 20.9 m g L

2 1 

### ); for zooplankton metrics, we chose 8 m g L

2 1 

### (range from 7.5 to m g L

2 1 

### 9.1 m g L

2 1 

### ).

Step 3: Determine current TP in each lake— We show the histogram of the current TP data in Fig. 4. Table 6. Results of the regression tree analysis (RTA) of the relationship between TP and the untransformed biological responses (BIO) where the fit statistic, the proportional reduction in variance (RT A-PRE), is , 0.30; and r 2 and p values for the linear regression between TP and the BIO responses (transformed where necessary). For the da taset description, letter designation, and citation, see Table 3. n is the number of lakes in each analysis. For relationships that had more tha n one threshold, the values in bold are the thresholds that occur first, with the greatest reduction in deviance. Tox in is microcystin toxin.

Dataset Biological responses TP range ( m g L

2 1 ) n TP threshold(s) ( m g L

2 1 ) RTA-PRE Linear regression

r 2 p c Phytoplankton community and toxins Microcystis biomass ( m g L

2 1 ) 4.8–64.0 56 34.6 0.11 0.06 0.07 Toxin concentration per Chl a ( m g m g 2 1 ) 4.8–64.0 64 8.1, 9.8 0.21 0.00 0.91 Environmental concentration of toxin ( m g L

2 1 ) 4.8–64.0 64 34.6 0.25 0.17 0.001 Toxin per microcystis biomass ( m g m g 2 1 ) 4.8–64.0 56 34.6 0.17 0.06 0.08 Cyanobacteria biomass ( m g L

2 1 ) 4.8–21.0 51 — 0 0.03 0.26 Microzooplankton Ciliate biomass ( m g L

2 1 ) 4.8–64.0 54 — 0 0.04 0.15 Rotifer biomass ( m g L

2 1 ) 4.8–44.0 61 20.7 0.16 0.02 0.23 Rotifer species number ( # ) 4.8–21.0 57 — 0 0.02 0.26 d, e, f Macrozooplankton community Bosmina spp. mean length (mm) 5.3–25.6 31 19.3 0.07 0.01 0.68 Ceriodaphnia spp. mean length (mm) 5.3–18.2 20 12.3 0.23 0.18 0.06 Daphnia spp. mean length (mm) 5.3–32.3 42 11.2 , 12.7 0.25 0.02 0.44 d, e Copepod biomass ( m g L

2 1 ) 5.3–32.3 28 11, 13.4 0.27 0.27 0.004 Daphnia pulicaria total biomass ( m g L

2 1 ) 5.3–25.6 17 11.2 0.16 0.02 0.63 Bosmina spp. total biomass ( m g L

2 1 ) 5.3–25.6 18 14.0 0.25 0.11 0.18 b Macrophyte cover ( % , unless otherwise indicated) Macrophyte edge (perimeter : area, m

2 1 ) 4.4–65.6 53 8.3, 22.7 , 28.4 0.27 0.00 0.90 Lake cover 4.4–65.6 53 8.3, 22.7 , 28.4 0.18 0.00 0.69 Lake dense cover 4.4–65.6 53 8.4 , 18.7 0.14 0.06 0.74 Lake emergent cover 4.4–65.6 53 19.9 0.08 0.03 0.23 Lake floating cover 4.4–65.6 53 — 0 0.01 0.50 Lake submersed cover 4.4–65.6 53 8.3, 22.7 0.14 0.01 0.44 Lake EWM cover 4.4–65.6 53 12.1 0.10 0.07 0.05 Lake dense EWM 1 cover 4.4–65.6 53 28.4 0.14 0.07 0.05 Littoral cover 4.4–65.6 53 19.9 0.10 0.05 0.10 Littoral dense cover 4.4–65.6 53 8.4 , 10.8 0.13 0.00 0.97 Littoral emergent cover 4.4–65.6 53 10.8 0.06 0.03 0.19 Littoral floating cover 4.4–65.6 53 — 0 0.00 0.71 Littoral submersed cover 4.4–65.6 53 19.9 0.11 0.07 0.06 Littoral EWM cover 4.4–65.6 53 12.2 0.08 0.05 0.11 Littoral dense EWM cover 4.4–65.6 53 — 0 0.05 0.11

1 EWM is the non-native Eurasian watermilfoil ( Myriophyllum spicatum L.) 

### Ecosystem-specific nutrient criteria 781

### Step 4: Derive lake-specific nutrient criteria using the BTPM algorithm— We entered the following into the BTPM algorithm: two BIO benchmarks of 8 m g L

2 1 

### and 18 m g L

2 1 

### , expected lake TP with allowance for each lake (EXP A 

### ), and current lake TP (CUR) for each lake. Figure 2 shows the graphical representation of this algorithm showing all possible relationships (A–H) among the inputs. For the 374 Michigan lakes, we calculated how many lakes fell within each category (Fig. 2). For the 23 % of lakes for which there was a BIO benchmark between the CUR and EXP A 

### , criteria were set at those benchmarks (B and D). For the 26 % of lakes that did not exceed a BIO benchmark, criteria were set at CUR (A and C). The only situation in which criteria were set at EXP A 

### was when both CUR and EXP A 

### were greater than the largest BIO benchmark and CUR was greater than EXP A 

### (E), which only occurred in 3 % of lakes in our database. For 48 % of the lakes, CUR TP was less than EXP A 

### TP (F–H). However, recall that we added an ‘allowance’ to each lake’s EXP 0 

### , including those lakes with relatively low human disturbance, and consequently low TP. Therefore, EXP A 

### is likely to be higher than CUR for many minimally disturbed lakes that we have in our dataset. To accommodate this situation and to prevent degradation, one of the rules in the algorithm states that under these conditions, criteria are set at CUR rather than EXP A 

### (F–H).



### Discussion

We have presented a framework to develop ecosystemspecific nutrient criteria in a scientifically-defensible way that is flexible enough to incorporate many different context-specific situations that might occur in different political jurisdictions or geographic regions. We have also presented an application of our framework to set TP criteria in a diverse set of Michigan lakes. We will first discuss implications of the BTBM framework and general topics that must be addressed in any application of this framework using our application to Michigan lakes as a concrete example, and then we discuss some features and issues that arose from our specific application of the BTPM framework to Michigan lakes.

Implications and general features of the BTPM Framework— One critical advantage of our BTPM framework, compared to other approaches for setting nutrient criteria,is that our framework explicitly addresses several interdependent questions that should be answered when setting nutrient criteria: What is the expected condition? What is the effect of nutrients on biological responses? What level of nutrients protects biological integrity? And, what is a reasonable level of protection? Our BTPM framework addresses these questions by: (1) integrating multiple pieces of information for setting nutrient criteria (expected conditions, BIO response thresholds, and current nutrient concentrations) into a single algorithm, (2)allowing nutrient criteria to be developed on an ecosystem-specific basis, instead of a single criterion for allecosystems within a given geographic region, (3) incorporating uncertainty into the benchmark-identification process which allows policy makers to make risk-level decisions based on scientific evidence, and (4) allowing for updating as additional data become available through refinements of the predictive model or the BIO benchmark identification procedure. Although we have presented a fairly detailed description of a framework to set nutrient criteria using Michigan lakes as an example, each application of this framework will likely result in a different set of decisions and justifications based on available data. In addition, there are still many decisions that need to be made at the implementation stage. For example, in our case, lake managers must still decide the lake sampling period for assessment, the management actions for lakes that exceed its nutrient criterion, how ‘close’ a lake has to be to its criterion to be considered attaining or nonattaining (i.e., how much measurement error is acceptable), and how to regulate surface water discharges of phosphorus to ensure that lakes meet the nutrient standards. These decisions will vary region to Table 7. Results of the regression tree analysis (RTA) of the relationship between TP and the untransformed biological responses (BIO) where the fit statistic, the proportional reduction in variance (RT A-PRE), is $ 0.30; and r 2 values for the linear regression between TP and the BIO responses (transformed where necessary); all linear regres sions are significant at p # 0.05, unless noted otherwise. For the dataset description, letter designation, and citation, see Table 3. n is the number of lakes in each analysis. For relationships that had more than one threshold, the values in bold are the thresholds that occur first, with the greatest reduction in deviance. The IQ-range is the interquartile range of the bootstrapped thresholds.

Dataset Biological responses TP range ( m g L

2 1 ) n TP threshold(s) ( m g L

2 1 ) RTA - PRE Linear regression r 2 IQ-range in bootstrapped

thresholds
Phytoplankton and clarity a, c, e Chl a ( m g L

2 1 ) 2.6–64.0 140 34.1 0.59 0.42 16–37 c Phytoplankton biomass ( m g L

2 1 ) 4.8–64.0 56 26.7 0.78 0.44 16–31 c Extinction coefficient (m

2 1 ) 4.8–64.0 64 26.7 0.71 0.47 12–37 a, c, e, f Secchi depth (m) 2.6–64.0 153 5.5, 11.4 , 22.0 0.32 0.29 5–21

Zooplankton d, e, f Cladoceran mean length (mm) 5.3–32.3 43 11.2 , 15.3 0.43 ns 6–11 d, e Proportion Daphnia biomass 5.3–32.3 28 11.7 , 16.5 0.64 0.53 9–12 Zooplankton biomass ( m g L

2 1 ) 5.3–32.3 28 13.0 , 15.6 0.30 0.21 8–13 

### 782 Soranno et al.

### region and are beyond the scope of our study here, but are necessary to address when applying this or any other framework. Another challenging issue when developing nutrient criteria is determining an acceptable allowance (i.e., deriving EXP A 

### ) for the expected nutrient concentrations. In our application, TP criteria were set higher than ‘reference condition’ (EXP 0 

### ) through the incorporation of an allowance term. However, lakes should be managed tomaintain resilience in response to potentially high interannual variability in nonpoint P loading (Carpenter and Cottingham 1997) and resulting high variability in lake TP Fig. 5. (A–D) Biological (BIO) responses for phytoplankton and clarity me trics versus TP for the variables with moderate-strong regression tree analysis (RTA) fits ( $ 0.30). Each data point represents a single lake. The plotted line is the cum ulative frequency of the 1,000 bootstrapped thresholds (right-hand y -axis, proportion). The arrows refer to the thresholds identified from th e individual RTA analysis. Where there is more than one threshold identified, the solid arr ow is the first threshold and dotted arrows are subsequent thresholds. n is the number of lakes in each dataset. (E–G) Same as above except for zoopla nkton metrics. 

### Ecosystem-specific nutrient criteria 783

### and Chl a (Knowlton and Jones 2006). There is evidence from modeling studies of single ecosystems that suggests that we should take a precautionary approach in light ofthe large uncertainties in limnological relationships, stochastic interannual P inputs, and lags in implementing policy (Carpenter et al. 1999). An important limitation to our application of the BTPM framework to Michigan lakes is the lack of quantification of interannual P variability associated with CUR and EXP A 

### , which can influence resulting management recommendations (Stow et al. 1997). Unfortunately, quantification of prediction uncertainty associated with interannual variability requires a large amount of data typically not available for substantial numbers of lakes. Future applications of the BTPM framework should incorporate these important factors where data are available (Knowlton and Jones 2006) or develop novel approaches, such as examining lake variance as an indicator of thresholds (and regime shifts) inecosystems (Brock and Carpenter 2006) or taking a Bayesian approach (Borsuk et al. 2004; Lamon and Stow 2004). Finally, an important component in establishing any standards for management is the quantification and communication of uncertainty (Harris et al. 2003). In the BTPM framework, uncertainty can be quantified for each of the four major steps. However, there are tradeoffs among complete modeling of uncertainty, implementation ease, communication between managers and the public, and data availability. In our application of the BTPMframework for Michigan lakes, we incorporated uncertainty in the first two steps. In our predictive HGM-LU model, we incorporated model error by using the 75th confidence interval for EXP A 

### . For the BIO benchmarks, we quantified uncertainty (defined as risk) using a bootstrap simulation (King and Richardson 2003). Although we did not formally propagate errors through to the final step of assigning each nutrient criterion, we found that the BTPM algorithm is somewhat robust for uncertainties from the previous two steps, particularly the prediction of EXP A 

### . For example, as long as the predicted value of EXP A 

### is in the same BIO zone as its true value, then the nutrient criterion (and consequently, the management action) remains the same. We chose the simplest case in our Michigan lakes example where we treat both the BIO benchmarks and the expected conditions as single point values for demonstration purposes and due to lack of data. To handle more complex methods of quantifying uncertainty including Bayesian analyses, one could treat both the BIO benchmarks and expected conditions as ranges of values with probabilitiesassociated with those ranges. In sum, the BTPM framework can be adapted to incorporate uncertainty at any level of complexity.

Application of the BTPM framework to Michigan lakes— In our application of the BTPM framework to Michigan lakes, we used an HGM-LU model to predict TP in Michigan lakes, and we generated BIO benchmarks bysearching for nonlinear relationships of commonly measured biological responses with TP. Next, we discuss issues that arose during our application of these two specific steps in more detail, pointing out limitations of our application,Table 8.The median of the biological (BIO) responses above and below a given TP threshold and the BIO response factor (the larger median divided by the smallermedian).nis the number of lakes below the threshold (top) and the number of lakes above the threshold (bottom). Units of BIO response metrics are the same as Table 7.Resp. is response, and daph. biom. isDaphniabiomass. Values in bold have factors.2.0, which we define as biologically important. The medians and factors are calculated atthree risk levels (1%, 10%, and 50%).BIOresponsemetric1%risk: 10%risk: 50%risk:TPthresholdMedianbelowMedianaboveBIOresp.factornTPthresholdMedianbelowMedianaboveBIOresp.factornTPthresholdMedianbelowMedianaboveBIO resp.factornPhytoplankton and clarity 110 132 137Chla16.02.82 6.682.43025.83.00 20.206.7834.33.00 51.3017.13Phytoplank- 45 50 52ton biomass16.437.85 121.43.21120.438.10 500.613.1630.438.90 658.816.94Extinction 42 58 60coefficient 11.920.2420.33 1.4 2220.920.2420.783.3633.120.2421.124.74Secchi 10 18 75depth 5.1 4.55 3.40 1.3 143 5.8 4.40 3.38 1.3 135 11.3 3.90 3.00 1.3 78ZooplanktonCladoceran 2 4 7mean length 6.1 1.31 0.78 1.7 41 7.6 1.27 0.78 1.6 39 11.1 1.21 0.76 1.6 36Proportion 5 6 9daph. biom9.10.48 0.0195.62311.10.49 0.0197.42211.60.48 0.00119.519Zooplankton 4 5 13biomass7.539.50 132.43.4248.439.60 138.63.523 12.8 85.80 142.8 1.7 15 

### 784 Soranno et al.

### topics that need further study, and implications for other applications of this framework. For modeling lake TP, we took a hydrogeomorphic perspective that recognizes that lake nutrients are a function of both local and regional hydrogeomorphic features; and we used data readily available to state agencies. Predictive ability would be improved if we had sufficient data to use well-established process-based models that predict lake TP concentrations from loading rates and lake hydrology, but this is rarely the case for large numbers of lakes that need to be managed. For example, MI has 6,595 lakes . 0.04 km 2



### , most of which are public and must be managed by the state, often without any lake-specific data available. Thus, building models from GIS databases that are available statewide has obvious appeal and warrants further research to improve predictive ability. We found that lake HGM predictors captured variation in TP better than the regional predictor ecoregion. This result agrees with previous studies that have shown that ecoregions alone do not capture TP variation well for streams (Dodds and Oakes 2004, Wickham et al. 2005) or lakes (Cheruvelil et al. in press), most likely because ecoregion ignores important local, site-level features such as lake depth, catchment area or finer-scaled patterns in surficial geology. Furthermore, the ability of ecoregion to explain variation in TP is even lower once the effect of LULC is removed (Wickham et al. 2005; Cheruvelil et al. in press). Nevertheless, ecoregions should still be considered when developing a HGM-LU model because ecoregions are useful for a variety of management applications (Stoddard 2004), continue to be used by state and federal agencies, and may become more relevant for larger spatial extents than studied here (Cheruvelil et al. in press). In addition, ecoregion is likely to become important for study areas with more heterogeneous landscapes than Michigan (Heiskary and Wilson 2005) or for other variables besides TP (such as total nitrogen) (Cheruvelil et al. in press). One way to test our HGM-LU model is to compare EXP 0 

### of our lakes with expected TP from presettlement times obtained from paleolimnological reconstruction of TP using diatom communities. In a study of northeastern U.S. lakes, Dixit et al. (1999) found that 47 % of their lakes had TP , 10 m g L

2 1 

### , and 49 % were between 10 m g L

2 1 

### and 30 m g L

2 1 

### . For our study lakes, EXP 0 

### was , 10 m g L

2 1 

### TP for 66 % of our lakes, and the remaining lakes had TP

, 30 m g L

2 1 

### . Although lakes in northeastern U.S. are not identical to our study lakes, they share many similar HGMfeatures, including climate. Paleolimnological reconstruction of historical TP levels for a limited number of lakes within the region under study is a good supplement to our approach and can be integrated into our framework. An important challenge in setting nutrient criteria is theselection of BIO responses to use for identifying benchmarks. The phytoplankton metrics that we used for the Michigan lakes analysis have obvious relevance to lakemanagement because these metrics measure some component of algal biomass which subsequently influences water clarity. Cladoceran mean length has also been shown to be strongly related to Secchi transparency because cladocerans favor algal particles that have high light attenuation (Stemberger and Miller 2003). Similarly, the proportion of daphnid biomass is relevant given the dominance of daphnids in filtering particles, their importance in the zooplankton community, and their critical role in pelagic lake foodwebs as a source of food for planktivorous and juvenile fishes. Depending on data availability, there may be additional, highly relevant biological responses that have threshold responses to changes in TP for other management units. Our results suggest that nonlinear relationships of BIO responses with TP may be more common than previously thought, especially when examined in the context of the risk of nuisance conditions, such as cyanobacteria blooms (Downing et al. 2001). Visual inspection of publishedrelationships suggests the existence of nonlinear relationships, even though they have not been quantified (Jeppesen et al. 2000). For example, several taxonomic groups of phytoplankton exhibit nonlinear relationships to increasing TP, where the thresholds appear to occur between 10 m g L

2 1 

### and 30 m g L

2 1 

### TP (Watson et al. 1997), or in the case of % cyanobacteria biomass at

, 30 ug L

2 1 

### and

, 100 ug L

2 1 

### (Downing et al. 2001). In another study, the proportion of picophytoplankton declined dramatically above a TP value of

, 8 m g L

2 1 

### (Watson et al. 1992).Finally, Vadeboncoeur et al. (2003) examined the proportion of total lake primary production due to periphyton versus pelagic algae, and visual inspection of their data suggests a breakpoint at

, 10 m g L

2 1 

### TP, above which the proportion of total primary production due to benthic production dramatically declines. Having relatively high benthic production could be viewed as a desirable state because it signifies that there is enough light for production to occur on sediment surfaces in the littoral zone, facilitating both periphyton and aquatic plant growth and consequently littoral zone foodwebs. Evidence such as the above suggests that it will be fruitful to reexamine and quantify many existing datasets for nonlinear responses of lake biological and ecosystem dynamics along a TP gradient. These nonlinear responses are important to consider in a management context as they are likely points where there is an abrupt change in an ecosystem, or where small changes in the driver produces a large response, which is conceptually linked to the idea of resilience (Groffman et al. 2006). An assumption that is often implicit in nutrient criteria development is that during years of unusually high phosphorus loads due to climate and anthropogenic variability, lake biology and phosphorus dynamics may respond, but not irreversibly so, as might occur with a regime shift (Carpenter 2003). Thus, during these times, the lake would not shift to a ‘new state’ beyond which it would likely not return, as occurs in regime shifts between clear and turbid states in lakes (Genkai-Kato and Carpenter 2005). If threshold responses are in fact tied to regime shifts, we should manage lakes in such a way to avoid them because reversal is either difficult or impossible over relevant time scales (Carpenter 2003). Unfortunately, phosphorus is just one component that pushes a lake to another regime; factors such as food web structure, lake physical modifications, climate change, Ecosystem-specific nutrient criteria 785

### and non-native species introductions also determine where a lake is in relation to regime shifts. Therefore, all of these factors need to be managed concurrently withphosphorus and provide further support for a precautionary approach to setting nutrient criteria to ensure that the criteria correctly represent the designated uses they are intended to protect.



### References

B AKER , E. A., K. E. W EHRLY , P. W. S EELBACH , L. W ANG , M. J. W ILEY , AND T. S IMON . 2005. A multimetric assessment of stream condition in the Northern Lakes and Forests Ecoregion using spatially explicit statistical modeling and regional normalization. Trans. Am. Fish. Soc. 134: 697–710. B ORSUK , M. E., C. A. S TOW , AND K. H. R ECKHOW . 2004. A Bayesian network of eutrophication models for synthesis, prediction, and uncertainty analysis. Ecol. Model. 173: 219–239. B REIMAN , L., J. H. F RIEDMAN , R. O LSHEN , AND C. J. S TONE . 1984. Classification and regression trees. Wadsworth International Group. B REMIGAN , M. T., S. M. H ANSON , P. A. S ORANNO , K. S. C HERUVELIL , AND R. D. V ALLEY . 2005. Aquatic vegetation, largemouth bass and water quality responses to low-dose fluridone two years post-treatment. J. Aquat. Plant Manage. 43: 65–75. B ROCK , W. A., AND S. R. C ARPENTER . 2006. Variance as a leading indicator of regime shift in ecosystem services. Ecol. Soc. 11: 9 [accessed on October 12, 2006]. Available from http://www. ecologyandsociety.org/vol11/iss2/art9/ B URNHAM , K. P., AND D. R. A NDERSON . 2002. Model selection and multimodel inference, 2nd ed. Springer-Verlag. C ARPENTER , S. R. 2003. Regime shifts in lake ecosystems: Pattern and variation. V. 15. Excellence in Ecology Series. Ecology Institute.

———, AND K. L. C OTTINGHAM . 1997. Resilience and restoration of lakes. Conservation Ecology 1: 2 [accessed on October 12, 2006]. Available from http://www.consecol.org/vol1/iss1/art2. ———, D. L UDWIG , AND W. A. B ROCK . 1999. Management of eutrophication for lakes subject to potentially irreversible change. Ecol. Appl. 9: 751–771. C HERUVELIL , K. S., AND P. A. S ORANNO . In press. Relationships between lake macrophyte cover and landscape features. Aquatic Botany doi:10.1016/j.aquabot.2007.10.005. ———, ———, M. T. B REMIGAN , T. W AGNER , AND S. L. M ARTIN

. In press. Grouping lakes for water quality assessment and monitoring: The roles of regionalization and spatial scale. Environ. Manage. C HOW -F RASER , P., AND D. H. D UTHIE . 1983. An interpretation of phosphorus loadings in dystrophic lakes. Arch. Hydrobiol. 97: 109–121. D IXIT , S. S., J. P. S MOL , D. F. C HARLES , R. M. H UGHES , S. G. P AULSEN , AND G. B. C OLLINS . 1999. Assessing water quality changes in the lakes of the northeastern United States using sediment. Can. J. Fish. Aquat. Sci. 56: 131–152. D ODDS , W. K., AND R. M. O AKES . 2004. A technique forestablishing reference nutrient concentrations across watersheds affected by humans. Limnol. Oceanogr. Methods 2: 333–341. D OWNING , J. A., S. B. W ATSON , AND E. M C C AULEY . 2001. Predicting Cyanobacteria dominance in lakes. Can. J. Fish. Aquat. Sci. 58: 1905–1908. E UROPEAN U NION . 2000. Directive 2000/60/EC of the European Parliament and of the council establishing a framework for the community action in the field of water policy. European Commission, Off. J. Eur. Commun. L327 (2000) 1. G ENKAI -K ATO , M., AND S. R. C ARPENTER . 2005. Eutrophication due to phosphorus recycling in relation to lake morphometry, temperature, and macrophytes. Ecology 86: 210–219. G ROFFMAN , P. M., AND o THERS . 2006. Ecological thresholds: The key to successful environmental management or an important concept with no practical application. Ecosystems 9: 1–13. H AIR , J. F., R. E. A NDERSON , R. L. T ATHAM , AND W. C. B LACK

. 1995. Multivariate data analysis with readings. Prentice Hall. H A

 ̊ KANSON , L. 1996. Predicting important lake habitat variables

from maps using modern modeling tools. Can. J. Fish. Aquat. Sci. 53 (Suppl. 1): 364–382. ———, A. C. B LENCKNER , B. H ELLSTROM , AND S. S. H ELLSTROM

. 2005. The influence of calcium on the chlorophyll–phosphorus relationship and lake Secchi depths. Hydrobiologia 537: 111–123. H ARRIS , G. P., AND o THERS . 2003. The role of models in ecosystem management, p. 300–306. In C. D. Canham, J. J. Cole and W. K. Lauenroth, [eds.], Models in ecosystem science. Princeton University Press. H EISKARY , S. A., AND C. B. W ILSON . 2005. Minnesota lake water quality assessment report: Developing nutrient criteria. 3rd ed. Water Assessment and Environmental Information Section Environmental Analysis & Outcomes Division Watershed Section Regional Division. Minnesota Pollution Control Agency. September. H IGGINS , J. V., M. T. B RYER , M. L. K HOURY , AND T. W. F ITZHUGH . 2005. A freshwater classification approach for biodiversity conservation planning. Conserv. Biol. 19: 432–445. J EPPESEN , E., J. P. J ENSEN , M. S ØNDERGAARD , T. L AURIDSEN , AND F. L ANDKILDEHUS . 2000. Trophic structure, species richness and biodiversity in Danish lakes: Changes along a phosphorus gradient. Freshw. Biol. 45: 201–218. J UBAR , A. K. 2004. Quantifying effects of residential lakeshoredevelopment on littoral fishes and habitat: Toward a framework for lake ecosystem conservation. M.S. thesis. Michigan State University. K ING , R. S., AND C. J. R ICHARDSON . 2003. Integrating bioassessment and ecological risk assessment: An approach todeveloping numerical water-quality criteria. Environ. Manage. 31: 795–809. K NOLL , L. B. 2004. Understanding zebra mussel ( Dreissena polymorpha ) impacts on autotrophic and heterotrophic plankton of inland lake ecosystems. M.S. Thesis. Michigan State University. K NOWLTON , M. F., AND J. R. J ONES . 2006. Natural variability in lakes and reservoirs should be recognized in setting nutrient criteria. Lake and Reserv. Manage. 22: 161–166. L AMON , E. C., AND C. A. S TOW . 2004. Bayesian methods for regional-scale eutrophication models. Water Res. 38: 2764– 2774. L OUGHEED , V. L., C. A. P ARKER , AND R. J. S TEVENSON . 2007. Using nonlinear responses of multiple taxonomic groups to establish criteria indicative of wetland biological condition. Wetlands 27: 96–109. M ARTIN , S. L., AND P. A. S ORANNO . 2006. Defining lake landscapeposition: relationships to hydrologic connectivity and landscape features. Limnol. Oceanogr. 51: 801–814. N ORVELL , W. A., C. R. F RINK , AND D. E. H ILL . 1979. Phosphorus in Connecticut lakes predicted by land use. Proc. Natl. Acad. Sci. 76: 5426–5429. 

### 786 Soranno et al.Q IAN , S. S., R. S. K ING , AND C. J. R ICHARDSON . 2003. Twostatistical methods for the detection of environmental thresholds. Ecol. Model. 166: 87–97. R D EVELOPMENT C ORE T EAM . 2007. A language and environment for statistical computing. R Foundation for Statistical Computing. R ECKHOW , K. H., AND o THERS . 2005. A predictive approach to nutrient criteria. Environ. Sci. Technol. 39: 2913–2919. S MITH , V. H. 1998. Cultural eutrophication of inland, estuarine, and coastal waters, pp. 7–49. In M. L. Pace, and P. M. Groffman [eds.], Successes, limitations, and frontiers in ecosystem science. Springer-Verlag. S ØNDERGAARD , M., E. J EPPESEN , J. P. J ENSEN , AND S. L ILDAL A MSINCK . 2005. Water framework directive: Ecological classification of Danish lakes. J. Appl. Ecol. 42: 616–629. S TEMBERGER , R. S., AND E. K. M ILLER . 2003. Cladoceran body length and Secchi disk transparency in northeastern U.S. lakes. Can. J. Fish. Aquat. Sci. 60: 1477–1486. S TEVENSON , R. J., AND o THERS . 2004. Interpreting results of ecological assessments, p. 85–111. In M. T. Barbour, S. B. Norton, H. R. Preston and K. W. Thornton [eds.], Ecological assessment of aquatic resources: Linking science to decision-making. Society of Environmental Toxicology and Chemistry. ———, Y. P AN , AND P. V AITHIYANATHAN . 2002. Ecological assessment and indicator development. in wetlands: the case of algae in the Everglades, USA. Verh. Int. Ver. Theor. Angew. Limnol. 28: 1248–1252. S TODDARD , J. L. 2004. Use of ecological regions in aquatic assessments of ecological condition. Environ. Manage. 34: S61–S70. ———, D. P. L ARSEN , C. P. H AWKINS , R. K. J OHNSON , AND R. H. N ORRIS . 2006. Setting expectations for the ecological condition of streams: the concept of reference condition. Ecol. Appl. 16: 1267–1276. S TOW , C. A., S. R. C ARPENTER , AND R. C. L ATHROP . 1997. A Bayesian observation error model to predict cyanobacterial biovolume from spring total phosphorus in Lake Mendota, Wisconsin. Can. J. Fish. Aquat. Sci. 54: 464–473. U.S. E NVIRONMENTAL P ROTECTION A GENCY [USEPA]. 2000. Ambient water quality criteria recommendations information supporting the development of state and tribal nutrient criteria: Lakes and reservoirs in nutrient ecoregion VI. Office of Water. Washington, D.C. V ADEBONCOEUR , Y., E. J EPPESEN , M. J. V ANDER Z ANDEN , H.-H. S CHIERUP , K. C HRISTOFFERSEN , AND D. M. L ODGE . 2003. From Greenland to green lakes: Cultural eutrophication and the loss of benthic pathways in lakes. Limnol. Oceanogr. 48: 1408–1418. W ATSON , S. B., E. M C C AULEY , AND J. A. D OWNING . 1997. Patterns in phytoplankton taxonomic composition across temperate lakes of differing nutrient status. Limnol. Oceanogr. 42: 487–495.

———, ———, AND ———. 1992. Sigmoid relationships between phosphorus, algal biomass, and algal community structure. Can. J. Fish. Aquat. Sci. 49: 2605–2610. W ICKHAM , J. D., K. H. R IITTERS , T. G. W ADE , AND K. B. J ONES

. 2005. Evaluating the relative roles of ecological regions and land-cover composition for guiding establishment of nutrient criteria. Landscape Ecol. 20: 791–798.



### Received 23 October 2006 Accepted 23 September 2007

Amended 29 October 2007 Ecosystem-specific nutrient criteria 787