THE EFFECTS OF LAND USE DISTURBANCE AT LOCAL AND REGIONAL
SCALES ON WETLAND RELATIONSHIPS WITH LAKE TOTAL PHOSPHORUS
AND WATER COLOR
By
C. Emi Fergus

A THESIS
Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of
MASTER OF SCIENCE
Fisheries and Wildlife
2010

i

ABSTRACT
THE EFFECTS OF LAND USE DISTURBANCE AT LOCAL AND REGIONAL
SCALES ON WETLAND RELATIONSHIPS WITH LAKE TOTAL PHOSPHORUS
AND WATER COLOR
By
C. Emi Fergus
We quantified wetland relationships with lake total phosphorus (TP) and water color
in relation to landscape features at local and regional scales using multilevel mixed
effects models on 1,790 north-temperate lakes. Our two research questions were: 1) do
human land use disturbances at local scales affect wetland relationships with lake TP and
water color? And 2) are wetland relationships with lake TP and water color different
across regions, and if so, what regional characteristics are related to this variation? We
found that land use disturbance at the local scale did not interact with wetland-lake TP
relationships but regional agriculture was negatively related to wetland-TP slopes. In
regions with low proportions of agriculture wetland-TP slopes were positive, and in
regions with high proportions of agriculture wetland-TP slopes were negative. In
addition, lake TP and water color were related to similar local lake and catchment
characteristics but different regional characteristics. Regional agriculture explained 22%
of the variability in TP and regional groundwater explained 4% of the variability in water
color. Regional groundwater was negatively related to wetland-water color slopes such
that wetland effects were weak in regions with high groundwater. These results
demonstrate that lake TP and water color are hierarchically controlled by landscape
variables at both local and regional scales and that the regional setting influences wetland
relationships with both lake TP and color.

ii

ACKNOWLEDGEMENTS
My thesis research and graduate education were supported by the Michigan State
University AGEP graduate assistantship, Graduate School, and Center for Water
Sciences. I would like to thank my committee members Mary Bremigan and Catherine
Yansa for their constructive comments and suggestions. I also would like to thank Mary
Bremigan and Kendra Spence Cheruvelil for their help in preparing my thesis research
for publication. I also would like to acknowledge Kendra Spence Cheruvelil, Ty Wagner,
and Steven Pierce (CSTAT) for passing on their statistical knowledge and resources. My
thanks goes out to Katherine Webster for sharing her unique and valuable perspective to
my research questions, and also for graciously sharing her knowledge about the six states
database I worked with. I would like to acknowledge the MSU limnology lab members
Bret Alger, Tom Alwin, Stacie Auvenshine, Joshua Booker, Katie Droscha, Geoff Horst,
Andrea Miehls, Dianna Miller, Emily Norton, Kevin Pangle, Kim Peters, Jeff White, and
Heidi Ziegenmeyer for their support and making the graduate school experience a lot of
fun! I also would like to thank my family especially my parents Esther Onaga and Ted
Fergus and brother Scott Fergus for their unwavering support, encouragement, and
reminders of the important things in life. Last but not least I would like to acknowledge
my advisor, Patricia Soranno, who has been an incredible mentor preparing me both
academically and professionally for a career in limnology. I have learned so much
through your gracious mentoring style and unique perspective, and I look forward to
future work together as I begin the PhD.

iii

PREFACE

The document was written in the form of a manuscript to be submitted to a peer reviewed
journal with the following authors: C. Emi Fergus, Patricia A. Soranno, Kendra Spence
Cheruvelil, and Mary T. Bremigan. In the text “we” is used rather than “I” to reflect the
contribution of the coauthors in preparing the manuscript. However, I have done all
statistical analyses, literature reviews, and writing of this thesis.

iv

TABLE OF CONTENTS

LIST OF TABLES………………………………………………………………………..vi
LIST OF FIGURES……………………………………………………………………...ix
CHAPTER 1: THE EFFECTS OF LAND USE DISTURBANCE AT LOCAL
AND REGIONAL SCALES ON WETLAND RELATIONSHIPS WITH
LAKE TOTAL PHOSPHORUS AND WATER COLOR………………………………..1
Introduction……………………………………………………………………….1
Methods…………………………………………………………………………...6
Results……………………………………………………………………………16
Discussion………………………………………………………………………. 21
APPENDIX……………………………………………………………………………...45
REFERENCES…………………………………………………………………………..54

v

LIST OF TABLES

Table 1. Candidate models to test hypothesized relationships between local and regional
landscape-context variables and lake total phosphorus (TP) and water color. Landscapecontext variables were grouped into prediction type classes: hydrogeomorphic variables
(HGM), wetlands (WET), human land use disturbance (LU), local disturbance interactive
effects (WETL* LUL), random local wetland effects across regions (WETL-RAN), and
cross-level interactive effects (WETL * XR). Variables excluded from the models are
denoted by “–”. Variables with the same hypothesized relationship as the previous model
are symbolized by “*”. “All-lakes” is the complete lake dataset, including all levels of
land use disturbance, and “Min-dist.” is the minimally-disturbed lake dataset (<5%
agriculture and urban land use). Regional landscape-context variables are italicized. …29
Table 2. Summary statistics of local and regional landscape-context variables used to
predict lake TP for all lakes (n = 1,790 in 23 regions) and for minimally-disturbed lakes
(n = 777 in 16 regions). Minimally-disturbed lakes were defined to have < 5% agriculture
and urban land use within a 500 m buffer surrounding lakes from the all-lakes dataset.
Regional variables are italicized. ………………………………………………………..31
Table 3. Summary statistics of local and regional landscape-context variables used to
predict lake water color for all lakes (n = 1,527 in 21 regions) and for minimallydisturbed lakes (n = 686 in 16 regions). Minimally-disturbed lakes were defined to have <
5% agriculture and urban land use within a 500 m buffer surrounding lakes from the alllakes dataset. Regional variables are italicized. …………………………………………32
Table 4. Unconditional multilevel mixed effect models for lake TP and color with no
2
predictor variables and random intercepts. Variance estimates for within-region (σ ) and
between-region (τ00) are provided. The intraclass correlation (%ICC or ρ) is the
estimated proportion of regional variation from the total variation, calculated as τ00/(τ00
2

+ σ ) and is presented as a percent..…………………………………………………...33
Table 5. Mixed effects models predicting lake TP for the all-lakes dataset (n = 1,790
lakes, N = 23 regions) using local and regional landscape-context variables and allowing
for random intercepts and random wetland slopes across regions. Candidate models were
developed using local hydrogeomorphic (HGM), wetland (WET), and land use
disturbance (LU) variables. Regional variables were added to the top-ranked local models
based on AIC. Models with random local wetland effects among regions are noted as
WETL-RAN. Local wetland interactions at local and regional scales are presented in
Models 2b and 5. Estimates of variance components are presented for within-region lake
2
TP (σ ), between-region intercept (τ00), and between-region wetland slopes (τ10). Models
were compared using Akaike Information Criteria (AIC) by taking the difference in AIC
vi

values with the lowest AIC model (Δi) and calculating Akaike weights (wi). Local
predictor variables were transformed to correct for normality and group-mean centered.
Regional variables were grand mean centered and are italicized. ………………………34
Table 6. Mixed effects models predicting water color for the all-lakes dataset (n = 1,527
lakes, N = 21 regions) using local and regional landscape-context variables and allowing
for random intercepts and random wetland slopes across regions. Candidate models were
developed using local hydrogeomorphic (HGM), wetland (WET), and land use
disturbance (LU) variables. Regional variables were added to the top-ranked local models
based on AIC. Models with random local wetland effects among regions are noted as
WETL-RAN. Local wetland interactions at local and regional scales are presented in
Models 2b and 5. Estimates of variance components are presented for within region color
2

(σ ), between region intercept (τ00), and between region wetland slopes (τ10). Models
were compared using Akaike Information Criteria (AIC) by taking the difference in AIC
values with the lowest AIC model (Δi) and calculating Akaike weights (wi). Local
predictor variables were transformed to correct for normality and group-mean centered.
Regional variables were grand mean centered and are italicized. ………………………36
Table 7. Mixed effects models predicting lake TP for the minimally-disturbed lake
dataset (n = 777 lakes, N = 16 regions) using local and regional landscape-context
variables and allowing for random intercepts and random wetland slopes across regions.
Candidate models were developed using local hydrogeomorphic (HGM) and wetland
(WET) variables. Regional variables were added to the top-ranked local models based on
AIC. Models with random local wetland effects among regions are noted as WETL-RAN.
Wetland interactions at the local scale are presented in Model 2b. Estimates of variance
2

components are presented for within region lake TP (σ ), between region intercept (τ00),
and between region wetland slopes (τ10). Models were compared using Akaike
Information Criteria (AIC) by taking the difference in AIC values with the lowest AIC
model (Δi) and calculating Akaike weights (wi). Local predictor variables were
transformed to correct for normality and group-mean centered. Regional variables were
grand mean centered and are italicized. …………………………………………………38
Table 8. Mixed effects models predicting water color for the minimally-disturbed lake
dataset (n = 686 lakes, N = 16 regions) using local and regional landscape-context
variables and allowing for random intercepts and random wetland slopes across regions.
Candidate models were developed using local hydrogeomorphic (HGM) and wetland
(WET) variables. Regional variables were added to the top-ranked local models based on
AIC. Models with random local wetland effects among regions are noted as WETL-RAN.
vii

Wetland interactions at the local scale are presented in Model 2b. Estimates of variance
2

components are presented for within region lake TP (σ ), between region intercept (τ00),
and between region wetland slopes (τ10). Models were compared using Akaike
Information Criteria (AIC) by taking the difference in AIC values with the lowest AIC
model (Δi) and calculating Akaike weights (wi). Local predictor variables were
transformed to correct for normality and group-mean centered. Regional variables were
grand mean centered and are italicized. …………………………………………………39
Table A1. Description of Ecological Drainage Unit (EDU) Codes. …………………….50
Table A2. Correlation values for regional hydrology, climate, and landscape variables for
the all-lake TP dataset (N = 23 regions). Correlation coefficients that were greater than
0.50 were considered to be significantly correlated (α = 0.05). …………………………51

viii

LIST OF FIGURES

Figure 1. Distribution of the lakes and the regions (Ecological Drainage Units – Higgins
et al. 2005) used in this study. Minimally-disturbed lakes were defined as having less
than 5% agriculture/urban land use within 500m buffer surrounding the lake. …………40
Figure 2. The regional deviation of lake response variables from the grand mean – best
linear unbiased predictors (Blups) for lake TP (2a) and water color (2b). Error bars are
the standard error in the Blup estimate. Region codes are shown in Figures 4 – 5 are
described in Table A1. in the appendix. ………………………………………………...41
Figure 3a. The mapped regional deviation of lake response variables from the grand mean
for lake TP. ………………………………………………………………………………42
Figure 3b. The mapped regional deviation of lake response variables from the grand
mean for water color. ……………………………………………………………………43
Figure 4. Regional deviation of lake response variables (Blups – best linear unbiased
predictors) versus regional landscape-context variables. 4a) TP Blup vs. regional
agriculture (arcsine square root transformed); 4b) Color Blup vs. regional baseflow index
(% groundwater contribution). …………………………………………………………..44
Figure A1. Percent agriculture in the region (Ecological Drainage Unit). ……...………52
Figure A2. Percent regional groundwater contribution (baseflow index). ………..…….53

ix

CHAPTER 1: THE EFFECTS OF LAND USE DISTURBANCE AT LOCAL AND
REGIONAL SCALES ON WETLAND RELATIONSHIPS WITH LAKE TOTAL
PHOSPHORUS AND WATER COLOR
Introduction
It is becoming widely recognized that lake trophic status is more completely
defined by both phosphorus and dissolved organic carbon (DOC) (Williamson et al.
1999) than by phosphorus alone. Phosphorus promotes primary productivity especially at
moderate to high levels of total phosphorus (TP); and the humic proportion of DOC
(water color) limits algae growth by restricting light availability, especially at low levels
of TP (Karlsson et al. 2009). The majority of both substances originate from the
landscape surrounding lakes; and because TP and DOC have been shown to be positively
correlated to one another (Detenbeck et al. 1993; Dillon and Molot 1997), they may be
controlled by similar landscape features. However, few studies have examined the
landscape controls of lake TP and DOC together (but see Webster et al. 2008).
Empirical cross-sectional studies of many lakes have found both TP and DOC to
be related to natural hydrogeomorphic (HGM) variables, including lake and catchment
morphometry and natural land cover in the surrounding landscape. Some of these
variables have similar effects on lake TP and DOC. For example, lake depth is negatively
correlated with lake TP and DOC (Rasmussen et al. 1989, Webster et al. 2008), and
drainage ratio (the catchment area to lake area ratio) is positively correlated with lake TP
and DOC (Prepas et al. 2001, Mullholland 2003). However, landscape studies are often
conducted within single geographic regions or ignore how regional heterogeneity may
affect these relationships (but see Xenopoulos et al. 2003). Because HGM variables are

1

controlled by a hierarchy of factors, it is likely that there are both regional and local
drivers of variation in such relationships (Soranno et al. 2010).
For example, we know that local wetlands are considered to be influential
landscape modifiers of lake and stream water chemistry because 1) they are active sites of
biogeochemical processes and 2) their location at the interface between uplands and open
water allows them to regulate nutrient and material transport (Mitsch and Gosselink
2007). But wetlands are found in diverse regional HGM settings and these differences
may affect wetland relationships with surface water (Cole et al. 2002). Wetlands
positively affect lake carbon measures such as DOC and color (Gergel et al. 1999;
Canham et al. 2004); a result that has been found across a wide range of different
regional settings (Xenopoulos et al. 2003). However, wetland effects on lake TP are more
variable, with wetlands being either positively or negatively correlated to lake TP
depending on the study (Detenbeck et al. 1993, Devito et al. 2000, Diebel et al. 2009). In
fact, few studies have considered regional patterns in wetland TP relationships.
One way to reconcile the conflicting relationships between wetlands and lake TP
is to consider human land use disturbance such as agriculture. In north temperate regions
with minimal agriculture and urban disturbance, wetlands have been associated with
increased lake and stream TP concentrations suggesting that these wetlands are sources of
phosphorus to lakes (Dillon and Molot 1997, Devito et al. 2000). However, in moderate
to highly disturbed settings, wetlands have been associated with decreased lake and
stream TP concentrations (Johnston et al. 1990, Detenbeck et al. 1993, Weller et al. 1996,
Diebel et al. 2009). Wetlands in these latter settings may act as nutrient modifiers by
intercepting some of the excess phosphorus from agriculture and urban activities from

2

receiving surface waters. Together, these studies suggest that wetlands interact with land
use to affect lake or stream TP concentrations. However, given that nutrient and material
transport between terrestrial and aquatic systems may operate at different spatial scales
(Gergel et al. 1999, Sliva and Williams 2001), it is not known at what spatial scale these
interactions occur.
Land use disturbance may also affect lake and stream carbon concentrations.
Studies have shown land use activities that disrupt surface and soil hydrology, such as
logging, are associated with increased DOC concentrations. For example, clearcut
logging was related to higher lake DOC concentrations in comparison to lakes from
reference areas (Carignan et al. 2000). Agriculture land use in the catchment has not been
strongly related to stream DOC, but Wilson and Xenopoulos (2008) indicated that land
use activities that alter soil moisture or flow paths may affect DOC. In summary, few
studies have thoroughly explored how land use disturbance such as agriculture affects
either DOC or water color.
Most landscape studies have focused on catchment level drivers of water
chemistry and have been restricted to single geographic regions (Detenbeck et al. 1993,
Weller et al. 1996). However, freshwater systems are hierarchically structured such that
lakes are nested within local catchments and local catchments are nested within broader
regional boundaries and these differences can lead to naturally varying water chemistry.
Studies have shown that lake TP and DOC concentrations are different between regions
and that broad-scale features such as topography and climate may explain these
differences (D’Arcy and Carignan 1997, Xenopoulos et al. 2003, Sobek et al. 2007).
This regional variation may obscure the functional role of wetlands on lake nutrient and

3

carbon dynamics, especially if the among-region TP and DOC variation is greater than
variation accounted for by wetlands at the local scale and if these two scales are not
explicitly modeled. Additionally, wetland functional role in nutrient and carbon transport
may be different among regions due to characteristics in the regional HGM setting. It has
been shown that similar wetland types can have very different hydrologic characteristics
among regions. For example, the frequency of wetland inundation, which affects the
chemical environment for nutrient and carbon cycling, differs among regions and may be
attributed to regional variation in climate and soil characteristics (Cole et al. 2002). In
sum, the relative importance of local and regional landscape features on lake TP and
DOC is not well documented or understood.
One obstacle to assessing multi-scaled landscape effects on lake TP and DOC is
the limitation of common statistical approaches. Classical regression models are often
inappropriate for regional analyses because they either 1) pool observations across
regions (i.e. ignore regional variation) and violate assumptions of statistical independence
or 2) separate observations into regional groups which can lead to inflated variance
estimates within each region (Gelman and Hill 2007). Statistical approaches such as
multilevel mixed effects models can account for hierarchical drivers of water chemistry
by modeling variation in lake TP and color attributed to both local-catchment scale and
broad-region scale features. Mixed effects models have not been widely used in
freshwater- land use/land cover studies (but see Taranu and Gregory-Eaves 2008) but are
highly relevant to the hierarchical structure of freshwater ecosystems and ecological
processes (Wagner et al. 2006).

4

In this study, we used multilevel mixed-effects models to explore wetland
relationships with lake TP and water color and to consider both land use/land cover and
other HGM features at both local and regional scales. Our specific research questions
were, 1) Does local human land use disturbance affect wetland relationships with lake TP
and water color?, and 2) Are wetland relationships with lake TP and water color different
across regions, and if so, which regional characteristics are related to this variation? To
address these questions we compiled lake water chemistry and landscape data for lakes
distributed across twenty-three regions (see below for how we define region) within
several states (Wisconsin, Michigan, New Hampshire, and Maine).
We took an information theoretic approach to examine our hypotheses by
comparing six candidate models based on hypothesized mechanisms and ranking them by
best fit given the data. We developed candidate models that included both direct effects
and interactive effects between landscape-context variables and lake-response variables
(TP and water color) (Table 1). The citations for our candidate models are not an
exhaustive review but are examples of the hypothesized landscape-context relationships.
In addition, we hypothesized that lake TP and water color concentrations are different
among regions and thus, should be modeled as such. The five candidate models are as
follows. Local or within-region lake TP and color variation was modeled using localscale landscape variables including natural HGM features and wetlands (Model 1),
additive human land use disturbance (Model 2a), and land use interactions with wetlands
(Model 2b). Regional-scale models were added to the local models to account for
among-region TP and color variation (Model 3). We tested whether local wetland effects

5

(i.e. the slope parameters for wetlands) were different among regions (Model 4) and if
regional-scale features could explain these differences (Model 5).

Methods
Study lakes and lake data
Landscape, hydrology, and lake water chemistry data were compiled for 1,790
lakes in Maine, New Hampshire, Michigan, and Wisconsin, USA from existing databases
collected and obtained from state agencies (Figure 1). We included lakes with surface
2

area greater than 0.01 km and maximum depths greater than 2 m. See Webster et al.
(2008) for a complete description of the compiled database.
Water chemistry measurements were collected by state agencies with the majority
of sampling occurring from 1990-2003 following standardized field and laboratory
procedures. We restricted water chemistry measurements to single samples collected
from the epilimnion from June – September when summer stratification is likely to occur.
Lake TP was measured using colorimetric analyses with persulfate digestion. Lake water
color was measured using visual comparators in platinum cobalt units (PCU). Lakes had
either true color estimates, measured from filtered samples, or apparent color estimates
(unfiltered). To correct for the positive bias attributed to apparent color measures,
Webster et al. (2008) developed a regression equation to convert apparent color to true
color.
Analyses were conducted on two datasets (Figure 1). The first dataset (“alllakes”) included 1,790 lakes for TP and 1,527 lakes for color. The second dataset
(“minimally-disturbed”) included 777 lakes for TP and 686 lakes for color. Lakes were
6

considered minimally disturbed if they had less than 5% agriculture or urban land use
within a 500 m buffer surrounding the lake. Although there are other landscape
disturbances that affect lakes, we use agriculture and urban land use as an indicator of
many of the major human impacts in the catchment that influence water chemistry
(Morrice et al. 2007).
Landscape-context variables defined as hydrogeomorphic (HGM), wetland, and
land use variables were collected from existing state databases and geographic
information systems (GIS) data layers. Landscape-context variables that were considered
“natural” were lake depth, drainage ratio (catchment area to lake area ratio), baseflow
(groundwater contribution), and wetlands; and landscape-context variables that
represented “disturbance” were agriculture and urban land use. We obtained land use /
cover (LULC) and wetland data from the 1992 National Land Cover Dataset
(http://landcover.usgs.gov/natllandcover.php) and included categories for upland forest,
wetland (forested, emergent, and total), agriculture, and urban land. Baseflow, the
proportion of streamflow from ground-water discharge, was quantified using a baseflow
index based on United States Geological Survey (USGS) stream gage data (Wolock
2003).
Measures of lake depth (maximum and mean) and total catchment area were
compiled from state databases. Lake surface area was quantified using state GIS data at
1:24,000 resolution. We calculated lake drainage ratio by dividing the catchment area by
lake surface area, which is a morphometric measure related to water residence time and
allochthonous sources of material to recipient waters.

7

All of the landscape-context variables except morphometry metrics were
quantified at two spatial scales: local and regional. Local variables were quantified within
a 500 m buffer surrounding lakes. The 500 m lake buffer was intended to simulate the
local lake catchment area, and correlation analyses using a subset of these study lakes
showed that LULC proportions in the 500 m lake buffer were highly correlated (r > 0.50)
with proportions in the local catchment (Webster et al. 2008). In addition, land cover
quantified in the 500 m lake buffer was just as good a predictor of lake DOC as land
cover measured in the whole catchment for several Wisconsin lakes (Gergel et al. 1999).
Thus, landscape context features measured at the 500 m buffer scale can capture
ecologically relevant features at the local scale that are related to lake nutrient and carbon
dynamics. We used only one measure for wetland cover to simplify the model building
process. We tested univariate relationships between TP and color and the three wetland
measures and found that forested wetlands had the strongest relationship with both TP
and color in comparison to total wetland and emergent wetland measures. In addition
forested wetland land cover was highly correlated with all other wetland variables
(Fergus, unpublished data) and thus we used only forested wetlands in all models.
At the regional scale we quantified the proportion of agriculture and baseflow
within regions defined as Ecological Drainage Units (EDU) (Higgins et al. 2005).
Ecological Drainage Units are regionalization schemes developed to classify freshwater
ecosystems and aquatic biodiversity using topographic, hydrologic, and climatic features.
These units have been shown to capture more variation in lake water chemistry variables
on average than other regionalization schemes (Cheruvelil et al. 2008). Background on
the landscape-contextual similarities and differences among the Midwest and New

8

England regions of the U.S. is provided in the Appendix. The regional-scale features used
in the models were narrowed down to agriculture for TP and baseflow for color based on
preliminary correlation analyses.
Statistical analyses
We fit the five candidate mixed effects models (Table 1) for each of the datasets
and for each response variable, except where noted below. Models 2a and 2b, which
included local agriculture and urban predictor variables and their interactions with
wetlands, were excluded from the minimally-disturbed dataset analysis because
agriculture and urban development were not hypothesized to account for variation in TP
and color under such low proportions in the landscape. However, models 2a and 2b were
included in the all-lake dataset analysis.
Lake TP and water color were natural log transformed, lake and catchment
morphometry variables were natural log transformed, and land cover/use proportions
were arc sine square root transformed to meet statistical assumptions of normality. Local
predictor variables were group mean centered (Xij – X .j) to control for among region
differences and reduce correlations with regional predictors (Enders and Tofighi 2007).
Regional variables were grand mean centered by subtracting individual observations from
the overall mean value (Xij – X ..). Multilevel mixed effects analyses were performed
using SAS MIXED procedure (SAS, Cary, NC).
Mixed effects model building process
We describe the model-building process using TP as the response variable but the
exact steps were applied to color as the response variable.

9

Regional variation in lake TP: unconditional models with random intercepts –
Prior to testing our specific research questions, we determined whether there was
significant regional variation in lake TP by developing unconditional models that had no
predictor variables but allowed for TP to vary by region (random intercept). This
preliminary analysis partitioned the total variation in lake TP into two variance
2

components: σ within region (local variation) and τ00 among region (regional variation).
The unconditional model to predict lake TP is:
Yij = β0j + rij ,
β0j = γ00 + u0j ,

(1)
2

where rij ~ N (0,σ ) and u0j ~N(0, τ00)
In this model Yij is lake TP for lake i in region j; β0j is the mean region TP or
random intercept; and rij is the individual error for lake i in region j. The intercept (β0j) is
the sum of the overall mean TP (γ00) across all regions and u0j, a random regional error
for region j (regional residual from the grand mean intercept), where u0j ~N(0, τ00) and
τ00 represents the among region variability in lake TP. Within region lake TP variation is
2

represented by σ .
From this unconditional model we can calculate the intraclass correlation
coefficient (ICC, ρ) which measures the proportion of TP variance that is among regions:
2

ρ = τ00 / (τ00 + σ ).

(2)

10

The intraclass correlation coefficient (ICC) provides an estimate of the extent to
which lakes in a region are correlated to one another with regard to lake TP
concentrations. The ICC value was used to evaluate if a mixed effect model approach
was needed or if regional variation was sufficiently low that ordinary least squares
regression models were appropriate. The unconditional models indicated that there were
significant in TP and color among regions and all subsequent models were built as
mixed-effect models allowing for random intercepts.
Local TP variation models with local landscape-context predictors [Models 1-2b]
– To address research question 1 we modeled within-region (local) lake TP variation
using the a priori local landscape-context variables hypothesized to affect lake TP and
color (Models 1-2b in Table 1). Local predictor variables that were correlated to one
another (r > 0.5) were not included in the same model to avoid problems of collinearity.
Local landscape variables were treated as fixed effects across regions. Model parameters
were derived using full maximum likelihood estimation. The premise for this approach is
to maximize the likelihood of the parameters given the data as opposed to the ordinary
least squares method of minimizing the model residual error.
A local model with a single fixed-effect predictor variable (lake depth) for lake
TP is presented below:
Yij = β0j + β1j (Depthij) + rij
β0j = γ00 + u0j
β1j = γ10
2

where rij ~ N (0,σ ) and u0j ~N(0, τ00)
11

(3)

In this model Yij is lake TP for lake i in region j; β0j is the mean region TP or
random intercept; β1j is the fixed effect of lake depth on lake TP; and rij is the individual
2

2

error for lake i in region j, where rij ~ N (0, σ ) and σ represents the within-region lake
TP variation. β1j is the fixed average lake depth-TP regression slope across regions (γ10).
More complex local models can be developed from this simplified model by including
other local landscape-context predictors.
2

Local predictor variables that reduced the within-region variation (σ ) were
retained in the model building process. The proportion of variance explained by the local
landscape-context variables was calculated as:
2

2

2

Varlocal = (σ unconditional - σ current model )/ σ unconditional,
(4)
where Varlocal is the local TP or water color variation explained by the local
2

model; σ unconditional is the within-region variation in the unconditional random
2

intercept model, and σ current model is the within-region variation in the models
conditioned with predictor variables.
Local and regional TP variation modeled with local and regional landscapecontext predictors [Model 3] – Regional landscape-context variables were added to the
best ranked local model to explain among-region TP variation. An example of a model
including local and regional variables to predict TP is illustrated below using lake depth
at the local scale and the proportion of agriculture at the regional scale:

12

Yij = β0j + β1j (Depthij) + rij ,

(5)

β0j = γ00 + γ01(Regional Agriculture)j + u0j ,
β1j = γ10
2

where rij ~ N (0,σ ) and u0j ~N(0, τ00)
In this model the random intercept (β0j) is a function of the overall mean (γ00) and
the proportion of agriculture in the region (γ01). The variance τ00 is the regional level
variance after controlling for regional agriculture.
Many of the variables quantified at the regional scale were strongly correlated to
one another (Table A2) so we did not include more than one regional variable at a time
into any model. The regional variation explained by the addition of the regional
landscape-context predictor was calculated as:
Varregion = (τ00 unconditional - τ00 conditional model)/ τ00 unconditional, (6)
where Varregion is the regional TP variation explained by the addition of the regional
variable to the model; τ00 unconditional is the between-region variation in the
unconditional random intercept model, and τ00 conditional model is the between-region
variation in the model conditioned with predictor variables.
Regional wetland effects variation modeled with local & regional predictors and
random wetland effects [Model 4] – To address research question 2, local wetlands were
treated as random effects to determine whether wetland relationships with lake TP were

13

significantly different among regions (i.e. random wetland slopes). A model with local
and regional landscape-context variables and random wetland slopes is presented below:
Yij = β0j + β1j(Depthij) + β2j(Wetij) + rij ,

(7)

β0j = γ00 + γ01(Regional Agriculture)j + u0j ,
β1j = γ10,
β2j = γ20 + u2j ,
⎛ u 0j
where rij ~ N (0,σ ) and ⎜
⎜u
⎜ 2j
⎝
2

⎞
⎟~N
⎟
⎟
⎠

⎡⎛ 0 ⎞ ⎛τ τ ⎞⎤
⎢⎜ ⎟, ⎜ 00 01 ⎟⎥
⎢⎜ 0 ⎟ ⎜τ τ ⎟⎥
⎣⎝ ⎠ ⎝ 10 11 ⎠⎦

In this model, lake TP is a function of lake depth (fixed effect) and random local
wetland effects across regions. Both the intercept and wetland slope are allowed to vary
among regions by including error terms u0j and u2j where u0j is the regional intercept
error for region j and τ00 represents the among-region variability in lake TP after
controlling for regional agriculture; u2j is the regional error to the slope associated with
region j and τ11 represents the among-region wetland effect variability; and τ01 is the
covariance between u0j and u2j.
These models were evaluated using an alpha of 0.05 to determine whether
wetlands should be treated as random effects or fixed effects across all regions.
Regional wetland effects variation modeled with random wetland effects and
cross-scale interactions [Model 5] – After determining whether wetland effects were
variable across regions, we attempted to explain this variation by including interaction
14

terms between local wetlands and regional landscape-context variables referred to as
cross-scale interactions. Building on equation 7, we can include an interaction term
between local wetlands and a regional-scale variable such as agriculture which yields the
following example model:
Yij = β0j + β1j(Depthij) + β2j( Wetij) + β3j(Wetij x Regional Agric.) rij ,

(8)

β0j = γ00 + γ01(Regional Agric.)j + u0j ,
β1j = γ10,
β2j = γ20 + u2j ,
β3j = γ21 ,
⎛ u 0j
where rij ~ N (0,σ ) and ⎜
⎜u
⎜ 2j
⎝
2

⎞
⎟~ N
⎟
⎟
⎠

⎡⎛ 0 ⎞ ⎛τ τ ⎞⎤
⎢⎜ ⎟, ⎜ 00 01 ⎟⎥
⎢⎜ 0 ⎟ ⎜τ τ ⎟⎥
⎣⎝ ⎠ ⎝ 10 11 ⎠⎦

In this model lake TP is a function of the overall intercept (γ00), the main effect of
regional agriculture (γ01), fixed effect of lake depth (γ10), random effect of wetlands
(γ20), and the cross- scale interaction between local wetlands and regional agriculture
(γ21). A decrease in τ11, the variation in wetland slope among regions, indicates that the
regional interaction explained differences in wetland effects among regions.
Model evaluation: information criteria theory – All candidate models were
evaluated and ranked based on information criteria theory. For each model, Akaike
Information Criteria values (AIC) were calculated with lower values indicating better
model fits for the data. Models with the lowest AIC value were compared to other
15

candidate models in each dataset based on the differences in AIC values (Δi) and Akaike
weights. The Akaike weights (wi) quantify the evidence for empirical support of that
model (Anderson and Burnham 2002).

Results
Lake water chemistry and local and regional landscape context variables
There was a wide range of lake TP and color concentrations across the study
-1

lakes. Lake TP concentration ranged from 1 to 193 μgL , water color concentration
ranged from 1 to 140 PtCo (Tables 2 and 3). In the minimally-disturbed datasets, lake TP
-1

did not exceed 30 μgL , but water color ranged from 1 to 140 PtCo. Local landscapecontext variables including lake depth, drainage ratio, and mean baseflow also captured a
wide range of values across all datasets (Tables 2 and 3). Local LULC surrounding lakes
was highly variable with agriculture ranging from 0 to 93% in the all-lake TP dataset.
Most of the lakes were dominated by forest cover (median 78%) with low agriculture
(median 4%) and low urban (median 1%) land use coverage.
Lakes were contained within twenty-three regions composed of diverse climatic,
hydrologic, and landscape characteristics. Regional landscape characteristics were similar
for the all-lake and minimally-disturbed lake datasets (Table 2 and 3). Regional
agriculture and urban land use in the all-lake TP dataset ranged from 1 - 78% and 0 - 20%
respectfully and were greater than those in the minimally-disturbed TP dataset (1- 61%
for agriculture and 0 - 6% for urban). The HGM variables among regions varied as well:

16

baseflow ranged from 46 to 78% and wetlands ranged from 1 to 35% in the all-lake TP
dataset and had similar values for the minimally-disturbed TP dataset.
Regional variation in lake TP: unconditional models with random intercepts
The unconditional models showed that there was significant variation in lake TP
and water color among regions (Table 4). This result means that lakes within the same
region have more similar lake TP and color concentrations compared to lakes from other
regions. The all-lake TP dataset had the greatest TP regional variation with 31% of the
total variation in TP occurring at the regional scale. For the minimally-disturbed lake TP
dataset, 16% of the total TP variation occurs at the regional scale and 84% of the
variation occurs at the local scale. This indicates that within-region differences may
account for more TP variation for minimally-disturbed lakes.
For the all-lake color dataset, about 13% of the total color variation occurs at the
regional scale and the remaining 87% variation occurs at the local scale (Table 4). For the
minimally-disturbed water color dataset, 20% of the total water color variation occurs at
the regional scale and 80% of the variation occurs at the local scale (Table 4).
The above unconditional models for all four datasets showed that a significant
proportion of TP and color variation occurs at the regional scale and so to treat lakes as
independent observations would violate assumptions of independence. Thus, all
remaining analyses were conducted using multilevel mixed effects models in which we
account for this regional variation in lake TP and water color by allowing for random
intercepts among regions.
Local TP variation models with local landscape-context predictors [Models 1-2b]

17

Across the local candidate models, the a priori hypothesized HGM variables (lake
depth, drainage ratio, and baseflow) had similar relationships with TP and color. Total
phosphorus and color relationships were both negatively related to lake depth, positively
related to drainage ratio, and negatively related to baseflow (Tables 5 – 8) as were
hypothesized in Table 1. Mean lake depth was used in all models except the all-lakes TP
dataset in which maximum depth was more strongly related to TP.
The relationships between wetlands and lake chemistry differed between TP and
color and across datasets. As expected (Table 1), local wetlands were positively related to
lake color for both the all-lake and minimally-disturbed datasets (Tables 6 and 8).
However, in the all-lake TP dataset, wetland effects were not significant across candidate
models (Table 5). In the minimally-disturbed lake TP dataset, wetlands were related to
TP (α < 0.1) such that higher TP occurred in lakes with increased wetlands around the
lake (Table 7), which matched our expectation.
The local candidate models that had the best fit given the data were different for
TP and water color for the all-lake datasets. Of the local candidate models evaluated for
the all-lake TP dataset, Model 2a had the lowest AIC value relative to the other candidate
models (Table 5). Including agriculture and urban development as positive additive
effects explained 21% of the TP variation at the local scale (within-region variation). For
the all-lake color dataset, local agriculture and urban land use in Model 2a were not
significantly related to water color (Table 6). Local models that included interaction
terms between wetlands and land use did not significantly improve the model fit for TP
nor water color (Model 2b in Tables 5 and 6).

18

In the all-lake TP datasets the lack of a wetland–TP relationship suggests that
local land use disturbance may obscure wetland effects and modeling wetland-land use
interactions at different spatial scales may explain wetland relationships with TP.
Local and regional TP variation modeled with local and regional landscape-context
predictors [Model 3]
After accounting for local landscape context variables, there was significant
regional lake TP and water color variation remaining. We calculated the deviation of
regional mean lake TP and water color from the grand mean, referred to as best linear
unbiased predictors (Blups), and geographically portrayed Blup values in Figures 2 and 3.
On average, regions in Michigan and Wisconsin had higher TP compared to lakes in
Maine and New Hampshire (Figure 4). Regional water color had less obvious regional
spatial patterns (Figure 5).
The addition of regional landscape context variables to the best local models
reduced among-region variance (τ00) in the TP and color models and improved the model
fit for all datasets (Model 3 in Tables 5 – 8). In the all-lake TP dataset, the proportion of
agriculture within the region was positively associated with region-average TP
concentrations and explained about 22% of the regional TP variation (Table 5). Lakes
within regions with high agriculture tended to have higher TP concentrations compared to
lakes from regions with low agriculture (Figure 6).
In the all-lake color dataset, the proportion of baseflow in the region accounted
for 4% of the regional color variation (Table 6). Regions with high baseflow tended to
have lakes with lower color concentrations compared to regions with low baseflow
(Figure 7).
19

The same regional-scale variables accounted for lake TP and color regional
differences for both the minimally-disturbed and all-lake datasets. In the minimallydisturbed TP dataset, regional agriculture accounted for 8% of the regional TP variation
and was associated with increased regional average TP concentrations (Table 7). In the
minimally-disturbed color dataset, the percent baseflow at the regional scale explained
14% of the regional water color variation and was negatively related to regional average
color concentrations (Table 8). Overall, including regional-scale variables to the model
improved the model fit (lowered AIC values) compared to models with only local
variables.
Regional wetland effects variation modeled with local & regional predictors and random
wetland effects [Model 4]
Wetland effects were different among regions for both TP and color, and allowing
for random wetland slopes among regions (Model 4) lowered AIC values and improved
model fit compared to models treating wetlands as fixed effects for the all-lakes datasets.
In the all-lakes TP dataset (Table 5), variation in the random wetland effects was
significantly different from zero (α <0.10) and lowered the AIC values by 6 units from
Model 3 which kept wetland effects fixed across regions. In the all-lakes color dataset
(Table 6), variation in mean wetland regression effects was significantly different from
zero (α <0.10) and lowered the AIC values by 11 units from the model which kept
wetland effects fixed across regions. However, in the minimally-disturbed TP and color
datasets, including random slope terms for wetlands did not improve the model fit,
indicating that wetland effects on TP and color in catchments dominated by natural land
cover (>95%) were not highly different across regions (Model 4 in Tables 7 and 8).

20

Regional wetland effects variation modeled with random wetland effects and cross-scale
interactions [Model 5]
In the all-lakes datasets, among-region differences in wetland effects on TP and
color were explained by regional-scale landscape-context variables. Cross-scale
interaction terms between local wetlands and regional landscape context variables
explained variation in wetlands slopes for TP and color (Model 5 in Tables 5 and 6).
Regional agriculture explained differences in wetland-TP relationships across regions,
such that wetland effects on lake TP were greater in regions with low proportions of
agriculture than in regions with high proportions of agriculture. Regional baseflow
explained differences in wetland-color relationships, such that wetland effects on color
were greater in regions with low baseflow than in regions with high baseflow.

Discussion
We found that when examining the relationship between wetlands and lake
chemistry, we must consider multiple spatial scales as well as human disturbances, both
of which are closely intertwined. We have two main conclusions. First, we found that
lake TP and color were controlled by similar lake and catchment features but different
regional variables. Differences in regional lake TP were related to regional land use,
whereas differences in regional water color were related to regional ground water
contribution (baseflow). Second, we found that across regions, local wetland effects on
TP and color were different and related to different factors. Local wetlands interacted
with regional land use to affect lake TP concentrations such that increasing local
wetlands decreased agricultural effects on lake TP concentrations or vice versa. For water

21

color, wetlands were positively related to color and differences in the region-specific
slopes were attributed to regional groundwater contribution.
Our mixed-effect analyses indicated that both lake TP and color variation were
greater within regions than between regions highlighting the importance of local lake and
catchment features to predict lake TP and color. Other studies have supported this
finding. They have demonstrated that within regions, TP variation is high and can
partially be accounted for by considering lake depth (Omernik et al. 1991) and local
catchment land use (Wickham et al. 2005). Nevertheless, regional differences in lake
phosphorus and color were large enough that across-region variation should not be
ignored.
Local landscape variables that predict within-region variation in TP and color
were similar to past studies and matched our expectations (Table 1). Lake depth was
negatively associated with both lake TP and color (Rasmussen et al. 1989). Drainage
ratio was positively associated with lake TP and color suggesting that small lakes with
large catchments have increased TP and color concentrations, which supports other
studies (Rasmussen et al. 1989, Kortelainen 1993). Webster et al. (2008) found that
drainage ratio is strongly correlated to water residence time and thus provides the
mechanisms for this relationship (i.e. internal controls of lake TP and water color). Lakes
with small drainage ratios are likely to have long water residence time which has been
related to decreased internal lake TP loading (D’Arcy and Carignan 1997) and increased
photodegredation of color (Molot and Dillon 1997).
We found that lake TP and color exhibit regional patterns, which is supported by
other cross-region studies for lake TP (Omernik et al. 1991, Rohm et al. 1995) and lake

22

DOC (Xenopoulos et al. 2003, Sobek et al. 2007). Our measure of region (EDU) captured
TP variation (as found by Cheruvelil et al. 2008 for Michigan) and it also captured
significant water color variation. However, we did not test EDU against different
regionalization frameworks and there may be other region measures that capture even
more lake water color variation. This is an area that warrants further study.
Regulation of lake TP and color by regional factors and the role of disturbance
Lake phosphorus and water color are hierarchically regulated by landscapecontext features at local and regional scales. Numerous studies have examined how lake
nutrients and organic carbon are affected by local features but few have explicitly
examined how they are affected by regional variables. Using mixed effect model
analysis, we found that regional variables consistently improved the overall model fit and
explained regional differences in lake TP and water color.
An important outcome of our study is the potential role of regional disturbance on
lake TP. It is well documented that agriculture exports phosphorus to surface waters at
local scales within a catchment or in near-lake buffers (Hunsaker et al. 1992). But our
results imply that agricultural activities within a region affect lake phosphorus
concentrations as well. Agriculture and urban development have been related to higher
regional lake TP concentrations compared to regions with low amount of agriculture and
urban land use within the Northeastern USA (Rohm et al. 1995). These relationships
indicate that human land use development, particularly agriculture, have diffuse, far
reaching effects on surface water chemistry such that land use within a region may even
affect lakes buffered by minimally disturbed local catchments. One potential explanation
for this pattern is that agriculture-dominant regions may have more nutrient rich soils as

23

opposed to regions with low amounts of agriculture. However, we did not have soils data
available to tease apart the relative influence of soil composition versus agricultural
activity for this study, and this issue warrants further study. Regional agricultural effects
on lake water chemistry have important implications for identifying reference lakes used
to set nutrient criteria. For example, identifying reference lakes simply by land use
disturbance at the local scale may be misleading because it does not take into account
regional land use disturbance that may be having an effect on lake nutrient levels.
Although water color was not related to regional disturbance in our study, it was
related to regional natural features such as the percent of groundwater contribution
(baseflow). In the all-lake dataset, the proportion of groundwater contribution in the
region was negatively related to water color. Negative relationships between baseflow
and carbon measures have been observed in studies conducted at the local scale (Jordan
et al. 1997). In our study regions, areas with low baseflow are characterized as having
clay-rich glacial deposits which promote surface runoff contribution (Wolfson 2009) and
this runoff could carry greater humic carbon concentrations (Thierfelder 1999).
Conversely, lakes in regions with high baseflow likely receive large groundwater input
which is low in organic carbon (Rasmussen et al. 1989). Regional baseflow may also be
indicative of other regional characteristics such as geology and vegetation cover classes
that have clear mechanistic relationships with carbon transport to lakes. For example, in
Michigan, areas associated with high groundwater are characterized as having sand and
gravel substrate (Lusch 2009) and less organic rich soils (Schatzel and Isard 1991). In
addition, we found a negative correlation between regional baseflow and regional upland

24

forest cover (r = -0.51) across all states. Together these findings indicate that high
baseflow regions may have reduced allochthonous humic carbon sources to lakes.
Wetland effects on lake TP and water color and the role of disturbance
Across regions, the local effects of wetlands on lake TP were related to regional
agriculture. In regions with high amounts of agriculture, wetlands surrounding lakes were
associated with decreased lake TP concentrations, indicating that wetlands may reduce
nutrient loading from agriculture at a regional scale. Given the nature of the analysis we
are not able to determine the specific mechanisms underlying the interaction term so we
can interpret the finding as either a result of local wetlands decreasing the regional
agricultural effects on lake TP or conversely regional agriculture decreasing the local
wetland effects on lake TP.
We were surprised to not find support for local-scale land use interactions with
wetland-lake TP relationships. The lack of a local agriculture-wetland interaction may be
explained by considering the coarse wetland metric that we used. Quantifying wetlands
within the 500 m lake buffer ignores wetland spatial configuration and connectivity in the
catchment, and landscape metrics that explicitly capture wetland connectivity have been
shown to be better predictors of lake nutrients than non-spatial metrics (Devito et al.
2000). For wetlands to interact or affect land use disturbance effects on lake phosphorus,
it is expected that they must be located in-between nutrient sources (agriculture or urban
land use) and sinks (lakes) in the catchment to intercept nutrient runoff. The wetland
proportion metric used in this study quantifies all wetlands in the surrounding landscape
near the lake, but not necessarily the wetlands that are most hydrologically connected to
the lake. Thus, wetlands were weakly associated with lake TP concentrations in our

25

dataset that included catchments with the full range of disturbance. We were not able to
quantify wetland spatial metrics for our study lakes but we expect that these spatial
characteristics are important to capture wetland effects on lake TP at local scales.
In the minimally-disturbed dataset, local wetland-TP relationships were generally
fixed across regions and were positively associated with lake TP. These findings are in
support of Devito et al. (2000) which indicated that forested wetlands may be phosphorus
sources to lakes in settings with minimal land use disturbance.
Local wetlands were positively related to water color and the magnitudes of these
relationships were different among regions (Table 6). Past studies link wetlands to
increased DOC concentrations in lakes and streams and support the finding that wetlands
act as carbon sources (Gergel et al. 1999, Xenopoulos et al. 2003). Land use disturbance
was not related to water color, which is consistent with other studies (Trebitz et al. 2007,
Wilson and Xenopoulos 2008). Unlike wetland-TP relationships, land use disturbance did
not interact with wetland effects on lake water color. Rather wetland effects were
different across regions due to differences in regional baseflow. Local wetlands
negatively interacted with regional baseflow such that one reduced the effect of the other
on lake water color. Greater regional baseflow or groundwater contribution could be
associated with less surface water runoff, reducing wetland effects on humic carbon
transport (Jordan et al. 1997). However, it is challenging to link regional hydrology to
specific wetland function when it is not known how regional hydrology affects local
wetland hydrology and carbon transport.
Limitations of analysis

26

Similar to landscape studies conducted at the local scale, regional studies suffer
from problems of multi-collinearity among landscape variables (King et al. 2005).
Several regional scale landscape context variables were highly correlated to one another
(Table A2), and we included only one regional variable in a model at a given time. In our
study, regional agriculture was negatively correlated to regional runoff (r = -0.8) such
that we cannot distinguish agricultural effects from runoff effects on lake TP. However,
agricultural mechanisms related to phosphorus transport to lake bodies at the local scale
are well supported in the literature and likely could persist at the regional scale. We tried
to avoid collinearity problems by keeping our landscape models simple and grounded in
ecological theory, but it is important to consider these spatial limitations when
conducting multi-level landscape studies.
Multi-level landscape studies are empirical in nature and lack the fine spatial and
temporal information integral to process-based mechanistic models. This is both a
weakness and strength of the empirical landscape approach (Ahlgren et al. 1988). The
landscape-context variables used in this analysis were fixed in time and were of coarse
spatial resolution. It is known that seasonal variability in features such as precipitation
can affect freshwater hydrology and thus nutrient and carbon processes (Cole et al. 2002),
and phosphorus and dissolved organic carbon dynamics can exhibit spatial and temporal
variability (Johnston et al. 2001, Wilson and Xenopoulos 2008). However, it can be
challenging to acquire both high frequency and high resolution data for a single lake
system let alone multiple lakes across different regions. In addition, we can only
hypothesize about the ecological processes and internal lake mechanisms that propagate
the lake TP and color patterns observed. However, without broadening the study scope to

27

include local and regional characteristics, we would be missing important ecological
relationships that hierarchically link systems together. Regional-scale studies are
important steps to identify broad-scale patterns, and they are a starting point to develop
more in-depth studies to better understand how the regional landscape context may
constrain finer-scale ecological processes. In this way multi-level landscape studies can
complement process based models and improve upon their efficiency (Strayer et al.
2003).
Implications
In this study we have shown that lake TP and color are controlled by variables at
both local and regional scales. Wetlands are important land cover features that affect lake
TP and color differently and interact with land use at broad spatial scales. Considering
multi-scaled landscape controls on lake TP and color is important to be able to predict
water chemistry and set realistic nutrient guidelines. In addition actively incorporating
multi-scaled interactions among landscape variables will help elucidate mechanisms for
TP and color transport between ecosystems. A multilevel mixed-effect model approach is
a useful analytic technique to better understand how ecological processes operate across
multiple spatial scales which has been identified as an important research challenge in
landscape and ecosystem ecology.

28

29

TP

Lake
response
*

*

3. & 4.

–
–
–
–

–
–
–
–

Agr.R * Wet.L

WETL-RAN

AgricultureRegional

Wet.L * Urb.L

–

Wet.L * Agr.L

4.

–

–

–

Neg. interact.
Neg. interact.

*

Pos.
–

–

Urban

5.

*

*

*

*

*

–

n.s.

Agriculture

Wetland

Neg.
Pos.

*

3.

Baseflow

*

Drainage ratio
Pos.

2.

Local
2b. HGM +
WET + LU +
(WETL * LUL)

4.

2a. HGM +
WET + LU

Lake depth
Neg.

1. HGM +
WET

Local

1.

Landscape-context
variable

Scale of predictor variables:
Scale of interaction:
Candidate models:

–

–

Pos.

*

*

*

*

*

*

*

*

3. Top
local +
regional

–

Random

*

*

*

*

*

*

*

*

*

4. Top local
+ regional +
WETLRAN

Neg. interact.

*

*

*

*

*

*

*

*

*

(WETL * XR)
*

Regional
5. Top local +
regional +
WETL-RAN

Local and Regional

wetland effects across regions (WETL-RAN), and cross-level interactive effects (WETL * XR). Variables excluded from the models
are denoted by “–”. Variables with the same hypothesized relationship as the previous model are symbolized by “*”. “All-lakes” is the
complete lake dataset, including all levels of land use disturbance, and “Min-dist.” is the minimally-disturbed lake dataset (<5%
agriculture and urban land use). Regional landscape-context variables are italicized.

Table 1. Candidate models to test hypothesized relationships between local and regional landscape-context variables and lake total
phosphorus (TP) and water color. Landscape-context variables were grouped into prediction type classes: hydrogeomorphic variables
(HGM), wetlands (WET), human land use disturbance (LU), local disturbance interactive effects (WETL* LUL), random local

30

6.

5.

4.

3.

2.

1.

*

6.

Rasmussen et al. 1989

Soranno et al. 1996

Detenbeck et al. 1993

Devito et al. 2000

Prepas et al. 2001

Canfield and Bachmann 1981

BaseflowR * Wet.L

WET.L-RAN

BaseflowRegional

Urb.L * Wet.L

Urban
Agr.L * Wet.L

Agriculture

Wetland

–
–
–
–

–
–
–
–

–
–

Pos.
–

Neg.
n.s.
–
–

*

8.

Baseflow

*

Drainage ratio
Pos.

*

7.

Lake depth

Color
Neg.

6.

9.

2a. HGM +
WET + LU

Landscape-context
variable

1. HGM + WET

Local

Lake
response

Scale of predictor variables:
Scale of interaction:
Candidate models:

Table 1. (cont’d)

–

–

–

Neg.

–

–
–

–

*

*

*

*

3. Top
local +
regional

Xenopoulos et al. 2003

Wilson and Xenopoulos 2008

Jordan et al. 1997

D’Arcy and Carignan 1997

10.

9.

8.

7.

–

–

n.s.

n.s.
n.s.

*

*

*

*

*

Local
2b. HGM +
WET + LU +
(WETL * LUL)

10.

Random
–

*

–

–
–

–

*

*

*

*

Neg. interact.

*

*

–

–
–

–

*

*

*

(WETL * XR)
*

Local and Regional
Regional
4. Top local 5. Top local +
+ regional +
regional +
WETL-RAN WETL-RAN

31

WET
LU

HGM

WET
LU

HGM

response

Prediction
type
Lake

Regional

Local

Scale
-1

2

Regional area (km )
% Baseflow
% Forest
% Wetland
% Agriculture
% Urban

Catchment area (ha)
Lake area (ha)
Drainage ratio
Mean depth (m)
Max. depth (m)
% Baseflow
% Forest
% Wetland
% Agriculture
% Urban

TP (μg L )
Color (PtCo)

Landscape-context
variable

53
76
3
7
2

15,250

734
51
13
4
10
54
78
3
4
1

14

8
20
6
19
3

11,524

19,172
681
191
3
9
10
22
9
17
0.9

19

46 – 78
9 – 86
1 – 35
1 – 78
0 – 20

2,830 – 48,950

9 – 197,357
1 – 18,043
1 – 2,800
0.3 – 26
2 – 68
34 – 89
4 – 100
0 – 66
0 – 90
0 – 93

1 – 140

All lakes (n = 1,790; N = 23)
Median
SD
Range
11
14
1 – 193

53
76
3
7
1

19,670

507
47
12
4
10
55
89
3
1
0.04

14

8
12
7
9
2

12,451

18,560
581
72
3
10
10
13
10
1
0.6

18

46 – 78
23 – 86
2 – 35
1 – 61
0–6

3,787 – 48,950

9 – 185,443
1 – 6,724
1 – 856
1 – 26
2 – 68
43 – 88
9 – 100
0 – 66
0–5
0–4

1 – 140

Minimally disturbed lakes (n = 777; N = 16)
Median
SD
Range
9
6
1 – 30

Table 2. Summary statistics of local and regional landscape-context variables used to predict lake TP for all lakes (n = 1,790 in 23
regions) and for minimally-disturbed lakes (n = 777 in 16 regions). Minimally-disturbed lakes were defined to have < 5% agriculture
and urban land use within a 500 m buffer surrounding lakes from the all-lakes dataset. Regional variables are italicized.

32

WET
LU

HGM

WET
LU

HGM

response

Prediction
type
Lake

Regional

Local

Scale
-1

2

Regional area (km )
% Baseflow
% Forest
% Wetland
% Agriculture
% Urban

Catchment area (ha)
Lake area (ha)
Drainage ratio
Mean depth (m)
Max. depth (m)
% Baseflow
% Forest
% Wetland
% Agriculture
% Urban

TP (μg L )
Color (PtCo)

Variable

12,020

19,474
718
200
3
9
10
21
9
15
9

19

8
18
6
17
3

15,250

710
46
14
4
9
53
80
3
4
4

14

55
70
6
13
3

46 – 78
9 – 86
2 – 35
1 – 78
0 – 20

2,830 – 48,950

9 – 197,357
1 – 18,043
1 – 2,800
1 – 26
2 – 68
36 – 89
4 – 100
0 – 66
0 – 88
0 – 93

1 – 140

All lakes (n = 1,527; N = 21)
Median
SD
Range
10
11
1 – 155

54
75
6
7
2

19,670

488
41
12
4
10
51
90
3
1
0.04

14

7
10
7
7
2

12,451

19,208
604
73
4
10
9
13
10
1
0.6

18

46 – 78
23 – 86
2 – 35
1 – 61
0–6

3,787 – 48,950

9 – 185,443
1 – 6,724
1 – 856
1 – 26
2 – 68
44 – 88
9 – 100
0 – 66
0–5
0–4

1 – 140

Minimally disturbed lakes (n = 686; N = 16)
Median
SD
Range
9
6
1 – 30

Table 3. Summary statistics of local and regional landscape-context variables used to predict lake water color for all lakes (n = 1,527
in 21 regions) and for minimally-disturbed lakes (n = 686 in 16 regions). Minimally-disturbed lakes were defined to have < 5%
agriculture and urban land use within a 500 m buffer surrounding lakes from the all-lakes dataset. Regional variables are italicized.

Table 4. Unconditional multilevel mixed effect models for lake TP and color with no
2
predictor variables and random intercepts. Variance estimates for within-region (σ ) and
between-region (τ00) are provided. The intraclass correlation (%ICC or ρ) is the
estimated proportion of regional variation from the total variation, calculated as τ00/(τ00
2

+ σ ) and is presented as a percent.
Lake dataset

Estimate

σ

All-lake TP

Lake and
(region) N
1,790 (23)

2.66

All-lake color

1,527 (21)

Minimally-disturbed TP
Minimally-disturbed color

2

τ00

% ICC

0.42

0.19

31

2.55

0.67

0.10

13

777 (16)

2.28

0.35

0.07

16

686 (16)

2.54

0.62

0.16

20

33

34

Regional

Interaction

WET
LU

HGM

Agr.R * Wet.L
Intercept

AgricultureR

Wet.L * Urb.L

Max. depth
Drainage ratio
Baseflow
Wetland
Agriculture
Urban
Wet.L * Agr.L

2a. HGM +
WET + LU

Local

-0.43**
0.06**
-0.008*
n.s.
0.76**
0.42**
–

Scale of predictor variables:
Scale of interaction:
Candidate models: 1. HGM +
WET
Landscape-context
variable
Category

–

–

–
2.69**

–

–

-0.43**
0.08**
-0.008*
n.s.
–
–
–

–
2.68**

2.69**

–

–

n.s.

-0.43**
0.06**
-0.007*
n.s.
0.74**
0.46**
n.s.

Local
2b. HGM + WET
+ LU +
(WETL * LUL)

2.44**

–

1.20**

–

-0.43**
0.06**
-0.008*
n.s.
0.76**
0.42**
–

3. Top
local +
regional

2.44**

–

1.20**

–

-0.43**
0.06**
-0.008*
n.s.
0.72**
0.38**
–

L

2.44**

-1.59**

1.20**

–

-0.43**
0.06**
-0.008*
n.s.
0.70**
0.38**
–

R

Local and Regional
Regional
4. Top local 5. Top local +
+ regional +
regional +
WETLWETL-RAN
RAN
(WET * LU )

calculating Akaike weights (wi). Local predictor variables were transformed to correct for normality and group-mean centered.
Regional variables were grand mean centered and are italicized.

compared using Akaike Information Criteria (AIC) by taking the difference in AIC values with the lowest AIC model (Δi) and

Table 5. Mixed effects models predicting lake TP for the all-lakes dataset (n = 1,790 lakes, N = 23 regions) using local and regional
landscape-context variables and allowing for random intercepts and random wetland slopes across regions. Candidate models were
developed using local hydrogeomorphic (HGM), wetland (WET), and land use disturbance (LU) variables. Regional variables were
added to the top-ranked local models based on AIC. Models with random local wetland effects among regions are noted as WETLRAN. Local wetland interactions at local and regional scales are presented in Models 2b and 5. Estimates of variance components are
2
presented for within-region lake TP (σ ), between-region intercept (τ00), and between-region wetland slopes (τ10). Models were

35

AIC
Δi

Model
selection

wi
p-values: <0.1†, <0.05*, <0.001**

Within region
Between region

Variation
explained

τ11

00

Model results:
2
Variance
σ
components
τ
0.29**

2a. HGM +
WET + LU

Local

–

2954
34
<0.001

3085
165
<0.001

21%
–

–

0.21*

18%
–

0.22*

0.31**

Scale of predictor variables:
Scale of interaction:
Candidate models: 1. HGM +
WET
Landscape-context
variable
Category

Table 5. (cont’d)

<0.001

2956
36

21%
–

–

0.22*

0.29**

Local
2b. HGM + WET
+ LU +
(WETL * LUL)

0.004

2931
11

21%
22%

–

0.06*

0.29**

3. Top
local +
regional

0.05

2926
6

22%
22%

0.22†

0.06*

0.29**

n.s.

0.95

2920
0

22%
24%

0.07

0.06*

0.29**

LUR)

L

Local and Regional
Regional
4. Top local 5. Top local +
+ regional +
regional +
WETLWETL-RAN
RAN
(WET *

36

Regional

Interaction

WET
LU

HGM

Category

BaseflowR * Wet.L
Intercept

BaseflowR

Wet.L * Urb.L

Mean Depth
Drainage ratio
Baseflow
Wetland
Agriculture
Urban
Wet.L * Agr.L

Landscape-context
variables

Scale of predictor variables:
Scale of interaction:
Candidate models:

-0.47**
0.20**
-0.02**
1.46**
n.s.
n.s.
–

2a. HGM +
WET + LU

–

–

–
2.55 **

–

–

-0.47**
0.20**
-0.02**
1.45**
–
–
–

1. HGM +
WET

Local

–
2.55**

2.56**

–

–

n.s.

-0.47**
0.20**
-0.02*
1.46**
n.s.
n.s.
n.s.

Local
2b. HGM +
WET + LU +
(WETL * LUL)

2.59**

–

-0.02*

–

-0.47**
0.20**
-0.02**
1.45**
–
–
–

3. Top
local +
regional

2.59**

–

-0.02*

–

-0.46**
0.20**
-0.02**
1.38**
–
–
–

2.59**

-0.06*

-0.02*

–

(WETL * LUR)
-0.46**
0.20**
-0.02**
1.63**
–
–
–

Local and Regional
Regional
4. Top local
5. Top local +
+ region +
regional +
WETL-RAN WETL-RAN +

calculating Akaike weights (wi). Local predictor variables were transformed to correct for normality and group-mean centered.
Regional variables were grand mean centered and are italicized.

compared using Akaike Information Criteria (AIC) by taking the difference in AIC values with the lowest AIC model (Δi) and

Table 6. Mixed effects models predicting water color for the all-lakes dataset (n = 1,527 lakes, N = 21 regions) using local and
regional landscape-context variables and allowing for random intercepts and wetland slopes across regions. Candidate models were
developed using local hydrogeomorphic (HGM), wetland (WET), and land use disturbance (LU) variables. Regional variables were
added to the top-ranked local models based on AIC. Models with random local wetland effects among regions are noted as WETLRAN. Local wetland interactions at local and regional scales are presented in Models 2b and 5. Estimates of variance components are
2
presented for within region color (σ ), between region intercept (τ00), and between region wetland slopes (τ10). Models were

37

AIC
Δi

Model
selection

wi
p-values: <0.1†, <0.05*, <0.001**

Within region
Between region

Variation
explained

τ11

components τ
00

Model results:
2
Variance
σ

Category

Landscape-context
variables

Scale of predictor variables:
Scale of interaction:
Candidate models:

Table 6. (cont’d)

3082
21
<0.001

3081
20
<0.001

32%
–

–

–
32%
–

0.12*

0.42**

2a. HGM +
WET + LU

0.12*

0.42**

1. HGM +
WET

Local

<0.001

3085
24

32%
–

–

0.12*

0.42**

LUL)

Local
2b. HGM +
WET + LU +
(WETL *

<0.001

3077
16

32%
4%

–

0.07*

0.42**

3. Top
local +
regional

0.08

3066
5

33%
4%

0.50†

0.07*

0.42**

4. Top local
+ region +
WETL-RAN

n.s.

0.92

3061
0

33%
4%

0.06

0.07*

0.42**

(WETL * LUR)

Regional
5. Top local +
regional +
WETL-RAN +

Local and Regional

Table 7. Mixed effects models predicting lake TP for the minimally-disturbed lake
dataset (n = 777 lakes, N = 16 regions) using local and regional landscape-context
variables and allowing for random intercepts and random wetland slopes across regions.
Candidate models were developed using local hydrogeomorphic (HGM) and wetland
(WET) variables. Regional variables were added to the top-ranked local models based on
AIC. Models with random local wetland effects among regions are noted as WETL-RAN.
Wetland interactions at the local scale are presented in Model 2b. Estimates of variance
2
components are presented for within region lake TP (σ ), between region intercept (τ00),
and between region wetland slopes (τ10). Models were compared using Akaike
Information Criteria (AIC) by taking the difference in AIC values with the lowest AIC
model (Δi) and calculating Akaike weights (wi). Local predictor variables were
transformed to correct for normality and group-mean centered. Regional variables were
grand mean centered and are italicized.
Scale of predictor variables:
Candidate models:
Category
HGM
WET
Regional
Model
results:
Variance
components

Local
1. HGM +
WET

Landscapecontext variables
Mean depth
Drainage ratio
Wetland
AgricultureR
Intercept

-0.47**
0.05*
0.27†
–

Local and Regional
3. Top local + 4. Top local
regional
+ regional +
WETLRAN
-0.47**
-0.47**
0.05*
0.05*
0.27†
n.s.
0.95*
0.95**

2.28**

2.19**

2.19**

2

0.25**

0.25**

0.25**

τ00

0.07*

0.03*

0.03*

–

–

σ

τ11

n.s.

0.10

Variation
explained

Within region
Between region

23%
–

23%
8%

23%
8%

Model
selection

AIC
Δi

1189
7

1182
0

1183
1

0.02

0.60

0.38

wi
p-values: <0.1†, <0.05*, <0.001**

38

Table 8. Mixed effects models predicting water color for the minimally-disturbed lake
dataset (n = 686 lakes, N = 16 regions) using local and regional landscape-context
variables and allowing for random intercepts and random wetland slopes across regions.
Candidate models were developed using local hydrogeomorphic (HGM) and wetland
(WET) variables. Regional variables were added to the top-ranked local models based on
AIC. Models with random local wetland effects among regions are noted as WETL-RAN.
Wetland interactions at the local scale are presented in Model 2b. Estimates of variance
2
components are presented for within region lake TP (σ ), between region intercept (τ00),
and between region wetland slopes (τ10). Models were compared using Akaike
Information Criteria (AIC) by taking the difference in AIC values with the lowest AIC
model (Δi) and calculating Akaike weights (wi). Local predictor variables were
transformed to correct for normality and group-mean centered. Regional variables were
grand mean centered and are italicized.
Scale of predictor variables:
Local
Candidate models: 1. HGM +
WET
Landscape-context
variables
Category
HGM
Mean Depth
-0.50**
Drainage ratio
0.23**
Baseflow
-0.01†
WET
Wetland
1.26**
Regional
–
BaseflowR
Intercept
2.54**
Model
results:
2
Variance
0.39**
σ
components

τ00
τ11

Local and Regional
3. Top local 4. Top local +
+ regional
region +
WETL-RAN
-0.50**
0.23**
-0.01†
1.26**
-0.04**

-0.50**
0.23**
-0.01†
1.26**
-0.04**

2.63**

2.63**

0.40**

0.39**

0.18*

0.05*

0.05*

–

–

n.s.

0.17

Variation
explained

Within region
Between region

29%
–

29%
14%

28%
14%

Model
selection

AIC
Δi

1360
9

1351
0

1352
1

0.01

0.65

0.34

wi
p-values: <0.1†, <0.05*, <0.001**

39

40
State boundaries

Figure 1. Distribution of the lakes and the regions (Ecological Drainage Units – Higgins et al. 2005) used in this study. Minimallydisturbed lakes were defined as having less than 5% agriculture/urban land use within 500m buffer surrounding the lake.

Kilometers

Study lakes
Minimally disturbed
Ecological Drainage Unit

Kilometers

Regional deviation of TP
from grand mean

2a.

1.5

1

0.5

0

-0.5

-1

Region

Regional deviation of color
from grand mean

2b.

1

0.5

0

-0.5

-1

Region
Figure 2. The regional deviation of lake response variables from the grand mean – best
linear unbiased predictors (Blups) for lake TP (2a) and water color (2b). Error bars are the
standard error in the Blup estimate. Region codes are shown in Figures 4 – 5 are
described in Table A1. in the appendix.
41

42

Figure 3a. The mapped regional deviation of lake response variables from the grand mean for lake TP.

Kilometers

-0.57 – -0.36
-0.36 – 0.05
0.05 – 0.26
0.26 – 0.53

Regional deviation of
lake TP from grand mean

Kilometers

43
-0.08 – 0.07
0.07 – 0.29

Figure 3b. The mapped regional deviation of lake response variables from the grand mean for water color.

Kilometers

-0.27 – -0.19
-0.19 – -0.08

Regional deviation of water
color from grand mean

Kilometers

1.0
0.5
0.0
-0.5

Regional deviation of TP from grand mean

y = -0.65 + 1.15x
p < 0.0001

0.2

0.4

0.6

0.8

1.0

-0.4

-0.2

0.0

0.2

0.4

y = 0.98 – 0.02x
p < 0.05

-0.6

Regional deviation of color from grand mean

Regional agriculture

45

50

55

60

65

70

75

Regional baseflow

Figure 4. Regional deviation of lake response variables (Blups – best linear unbiased
predictors) versus regional landscape-context variables. 4a) TP Blup vs. regional
agriculture (arcsine square root transformed); 4b) Color Blup vs. regional baseflow index
(% groundwater contribution).
44

APPENDIX

45

Regional background of study lakes
The study lakes are distributed in Wisconsin and Michigan in the Midwest and New
Hampshire and Maine in the Northeast of the United States. These continental regions
have a shared glacial history from the last ice age and thus have similar freshwater
landscapes, rich in lakes, streams, wetlands and groundwater. But there are striking
differences in topography and human land use intensity. We highlight some of the
similarities and differences found in the study areas to aid in the interpretation of the
regional analyses.
Glacial history
The last glaciation ended around 10,000 years ago in the Northern U.S. and covered
most of Wisconsin, and all of Michigan, New Hampshire and Maine. As the glaciers
retreated and melted away, they shaped the terrestrial and freshwater landscapes that we
see today. Glacial modified landscape features include: moraines – where sediment was
piled up at the edge of glacial lobes; drumlins that form elongated hills; and eskers or
sand and gravel ridges (Blewett et al. 2009). Inland lakes were formed through glaciers
physically scouring the landscape, depositing moraine sediment and damming river flow,
or leaving behind blocks of ice buried in glacial drift that later formed kettle lakes
(Wolfson 2009). Water under high pressure beneath glaciers created long, deep “tunnel
channel lakes” found in northwestern MI (Wolfson 2009). Many inland lakes are found in
areas where two glacial lobes intersected (interlobate zones) (Blewett et al. 2009). Other
glacial activity has lead to lake-free landscapes such as where glacial lakes have since
dried up and left behind flat plains, fine-textured soils and wetlands (Blewett et al. 2009).
An example of this type of landscape is found in around the Saginaw Bay region of
46

Michigan. In the Driftless region of Wisconsin, which was never glaciated by any of the
ice advances, there are few natural lakes because “millions of years of stream incision
have removed the necessary dams and integrated the drainage” of the area (Mickelson
1997).
Regional topography
The topographic elevation and relief greatly differs between the Midwest and
Northeastern U.S. In Michigan and Wisconsin, the land can be characterized as fairly flat
to hilly. In Maine and New Hampshire there are diverse topographic features and changes
in elevation. Most of Maine is upland, and New Hampshire is covered by the White
Mountain section, which includes Mount Washington, the highest peak in the
Northeastern U.S. with an elevation of 6,288 ft above sea level (Olcott 1995). Lowlands
reside along the Atlantic seaboard of Maine and New Hampshire.
Groundwater hydrology
Groundwater characteristics among regions are varied and can be related to the
glacial activity. Sand and gravel sediment deposited by glaciers (outwash) “readily
receive, store, transmit, and discharge water” (Olcott 1995) and supply much of the
baseflow of streams. In our study extent, we find regions with high baseflow in northern
Michigan (Figure A2). Regions with low baseflow include southeastern Wisconsin, the
Thumb region of Michigan, western New Hampshire, and northern Maine. The Thumb
region is characterized as having glacial-lake deposits consisting of clay, silt, and fine
sand which make for poor surficial aquifer systems. Glacial till in the Northeast transmits
water slowly and results in low groundwater yields. But outwash deposits in the

47

Northeast valleys yield high groundwater, which were not captured in the regional scale
of the baseflow map presented in Figure A2.
Climate and surface hydrology
The contrasting topography of the Midwest from the Northeastern U.S. influences
differences in climate and hydrology among regions. In regions with high elevation in
New Hampshire and Maine there is greater precipitation than in regions with low
elevation (Olcott 1995). The higher precipitation coupled with the steep slope results in
greater regional runoff in areas with high elevation in the Northeast. About half of the
annual precipitation becomes runoff to streams and rivers in New Hampshire and Maine
(Olcott 1995).
In the Midwest, topographic differences are moderate and climatic patterns are more
strongly influenced by other features in the landscape such as the Great Lakes. For
Michigan, precipitation decreases as we move from west to east in the state. For
Wisconsin, precipitation is greater in the west central part of the state (PRISM
http://www.prism.oregonstate.edu). Mean annual runoff in streams and rivers reflects the
annual precipitation and in no areas does runoff exceed precipitation values (Olcott
1992).
Vegetation cover
Climate and other physical factors have influenced land cover patterns in the study
regions. Temperature and precipitation gradients in the Midwest, referred to as the
Wisconsin and Michigan tension zones, have resulted in two distinct vegetation
compositions. In the south, summers are warmer and longer than in the north and
promote drought tolerant oak and sugar maple forests. North of the tension zone, the

48

climate is cooler and there is greater moisture in the soils, which allows for mixed
conifer-hardwood forests to dominate.
In the Northeastern USA, forest composition also is influenced by north-south
temperature gradients but also microclimate differences due to topographic elevation and
relief differences. In Maine, deciduous-dominated forests are found in the south and
spruce-fir or northern hardwood forests are found in the north (Dept. Conservation Maine
http://www.maine.gov/doc/nrimc/mnap/features/eco_whitepinemixforest.htm).
Human land use
In the Midwest, the dominant human land use activity in the surrounding landscape is
agriculture. Productive soils found in the Midwest are attributed to glacial activity that
has significantly altered the surficial landscape. During the ice age, “continental glaciers
planed off soil and weathered bedrock and redeposited these materials as thin mantel of
glacial debris over the bedrock surface” (Olcott 1995). The smaller sediments and silt
produced by the ice-activity, readily breaks down chemically and releases nutrients when
exposed to water and temperature changes. Agricultural activity also benefits from clays
produced by weathering processes which retain moisture and promote plant growth. The
percent of regional agriculture is greatest in central and southern WI and MI (Figure A1).
In the Northeast, there is a small percentage of agriculture around rivers and valleys and
small urban centers, but these regions are mostly forested with minimal human landscape
modifications in comparison to the Midwest. Surface glacial debris has resulted in rolling
flat topography in the lowland area and partially filled valleys between highlands,
moderating the sharp topographic relief of the preglacial landscape of the Northeast
(Olcott 1995).

49

Table A1. Description of Ecological Drainage Unit (EDU) Codes.
Description

EDU Code
BPU

Bayfield Peninsula and Uplands

CBR

Chippewa-Black River

CUP

Central Upper Peninsula

D

Driftless

ECW

East Central Wisconsin

EUP

Eastern Upper Peninsula

LSC

Lower St. Croix – Downeast Maine Coastal

MC

Middle Connecticut

MHS

Northern Lake Michigan, Lake Huron, and Straits of Mackinac

MIP

Southeast Michigan Interlobate and Lake Plain

PKA

Penobscot – Kennebec – Androscroggin

RR

Rock River

SB

Saginaw Bay

SCR
SD

St. Croix River
Southern Driftless

SLM

Southeast Lake Michigan

SMC

Saco – Merrimack – Charles

UC

Upper Connecticut

ULR

Upper Illinois River

USJ

Upper St. John – Aroostook

WMD

Western Lake Michigan and Door Peninsula

WPK

Western Upper Peninsula and Keweenaw Peninsula

WR

Wisconsin River

50

51

Baseflow
Runoff
Precipitation
Wetland
Forest
Agriculture
Urban

Baseflow
1
-0.57
-0.68
0.67
-0.51
0.31
-0.18

Runoff
–
1
0.85
-0.48
0.88
-0.82
-0.04

Precipitation
–
–
1
-0.69
0.80
-0.63
-0.01

Wetland
–
–
–
1
-0.36
0.08
-0.32

Forest
–
–
–
–
1
-0.94
-0.22

Agriculture
–
–
–
–
–
1
0.22

Urban
–
–
–
–
–
–
1

Table A2. Correlation values for regional hydrology, climate, and landscape variables for the all-lake TP dataset (N = 23 regions).
Correlation coefficients that were greater than 0.50 were considered to be significantly correlated (α = 0.05).

52
Figure A1. Percent agriculture in the region (Ecological Drainage Unit).

Kilometers

1–9%
9 – 16 %
16 – 45 %
45 – 78 %

Percent agriculture
in region

Kilometers

53
Figure A2. Percent regional groundwater contribution (baseflow index).

Kilometers

66 – 78 %

57 – 66 %

50 – 57 %

46 – 50 %

Percent baseflow
in region

Kilometers

REFERENCES

54

REFERENCES
Ahlgren, I., T. Frisk, and L. Kamp-Nielsen. 1988. Empirical and theoretical models of
phosphorus loading, retention and concentration vs. lake trophic state. Hydrobiol.
170: 285–303.
Anderson, D.R. and K.P. Burnham. 2002. Avoiding pitfalls when using informationtheoretic methods. J. Wildl. Manage. 66: 912–918.
Blewett, W.L, D.P. Lusch, and R.J. Schaetzl. 2009. The physical landscape: a glacial
legacy, p. 249–273. In R. Schaetzl, J. Darden, and D. Brandt [eds.], Michigan
Geography and Geology. Pearson, New York, New York.
Canfield, D.E. and R.W. Bachmann. 1981. Prediction of total phosphorus concentrations,
chlorophyll a, and Secchi depths in natural and artificial lakes. Can. J. Fish. Aquat.
Sci. 38: 414–423.
Canham, C.D., M.L. Pace, M.J. Papaik, A.G.B. Primack, K.M. Roy, R.J. Maranger, R.P.
Curran, and D.M. Spada. 2004. A spatially explicit watershed-scale analysis of
dissolved organic carbon in Adirondack lakes. Ecological Applications. 14: 839–854.
Carignan, R., P. D’Arcy, and S. Lamontagne. 2000. Comparative impacts of fire and
forest harvesting on water quality in Boreal Shield lakes. Can. J. Fish. Aquat. Sci. 57:
105–117.
Cheruvelil, K.S., P.A. Soranno, M.T. Bremigan, T. Wagner, and S.L. Martin. 2008.
Grouping lakes for water quality assessment and monitoring: the roles of
regionalization and spatial scale. Environ. Manage. 41: 425–440, doi:
10.1007/s00267-007-9045-7
Cole, C.A., R.P. Brooks, P.W. Shaffer, and M.E. Kentula. 2002. Comparison of
hydrology of wetlands in Pennsylvania and Oregon(USA) as an indicator of
transferability of hydrogeomorphic (HGM) functional models between regions.
Environmental Management. 30: 265–278, doi: 10.1007/s00267-001-0055-6
D’Arcy, P. and R. Carignan. 1997. Influence of catchment topography on water
chemistry in southeastern Quebec Shield lakes. Can. J. Fish. Aquat. Sci. 54: 2215–
2227.
Detenbeck, N.E., C.A. Johnston, & G.J. Niemi. 1993. Wetland effects on lake water
quality in the Minneapolis/St. Paul metropolitan area. Landscape Ecology. 8: 39–61.
Devito, K.J., I.F. Creed, R.L. Rothwell, & E.E. Prepas. 2000. Landscape controls on
phosphorus loading to boreal lakes: implications for the potential impacts of forest
harvesting. Can. J. Fish. Aquat. Sci. 57: 1977–1984.
55

Diebel, M. W., J.T. Maxted, D.M. Robsertson, S. Han, and M.J. Vander Zanden. 2009.
Landscape planning for agricultural nonpoint source pollution reduction III: assessing
phosphorus and sediment reduction potential. Environ. Manage. 43: 69–83.
Dillon, P.J. and L.A. Molot. 1997. Effect of landscape form on export of dissolved
organic carbon, iron, and phosphorus from forested stream catchments. Water
Resources Research. 33(11): 2591–2600.
Enders, C.K. and D. Tofighi. 2007. Centering predictor variables in cross-sectional
multilevel models: a new look at an old issue. Psychological Methods. 12: 121–138.
Gelman, A. and J. Hill. 2007. Data Analysis Using Regression and
Multilevel/Hierarchical Models. Cambridge University Press, New York, New York.
Gergel, S.E., M.G. Turner, & T.K. Kratz. 1999. Dissolved organic carbon as an indicator
of the scale of watershed influence on lakes and rivers. Ecological Applications. 9(4):
1377-1390.
Higgins, J.V., M.T. Bryer, M.L Khoury, and T.W. Fitzhugh. 2005. A freshwater
classification approach for biodiversity conservation planning. Conservation Biology.
19(2): 432-445.
Hunsaker, C.T., D.A. Levine, S.P. Timmins, B.L Jackson, and R.V. O’Neill. 1992.
Landscape characterization for assessing regional water quality. In Ecological
Indicators. ed. D.H. McKenzie, D.E. Hyatt, and V.J. McDonald, pp. 997-1006. New
York: Elsevier.
Johnston, C.A., N.E. Detenbeck, and G.J. Niemi. 1990. The cumulative effect of
wetlands on stream water quality and quantity: a landscape approach.
Biogeochemistry. 10: 105–141.
Johnston, C.A., S.D. Bridgham, and J.P. Schubauer-Berigan. 2001. Nutrient dynamics in
relation to geomorphology of riverine wetlands. Soil Sci. Soc. Am. J. 65: 557–577.
Jordan, T.E., D.L. Correll, and D.E. Weller. 1997. Relating nutrient discharges from
watersheds to land use and streamflow variability. Water Resources Research. 33:
2579–2590.
Karlsson, J., P. Bystrom, J. Ask, P. Ask. L. Persson, and M. Jansson. 2009. Light
limitation of nutrient-poor lake ecosystems. Nature. 460: 506–509, doi:
10.1038/nature08179
King, R.S., M.E. Baker, D.F. Whigham, D.E. Weller, T.E. Jordan, P.F. Kazyak, and M.K.
Hurd. 2005. Spatial considerations for linking watershed land cover to ecological
indicators in streams. Ecological Applications. 15: 137–152.
56

Kortelainen, P. 1993. Content of total organic carbon in Finnish lakes and its relationship
to catchment characteristics. Can. J. Aquati. Sci. 50: 1477–1483.
Lusch, D.P. 2009. Groundwater and karst, p. 223 –248. In R. Schaetzl, J. Darden, and D.
Brandt [eds.] , Michigan Geography and Geology. Pearson, New York, New York.
Mickelson, D.M. 1997. Wisconsin’s glacial landscapes, p. 35 – 48. In R.C. Ostergren and
T.R.Vale [eds]. Wisconsin Land and Life. University of Wisconsin Press.
th

Mitsch, W.J. and J.G. Gosselink. 2007. Wetlands. 4 edition. Wiley. Hoboken, NJ.
Molot, L.A. and P.J. Dillon. 1997. Colour – mass balances and colour – dissolved organic
carbon relationships in lakes and streams in central Ontario. Can. J. Fish. Aquat. Sci.
54: 2789–2795.
Morrice, J.A., N.P. Danz, R.R. Regal, J.R. Kelly, G.J. Niemi, E.D. Reavie, T.
Hollenhorst, R.P. Axler, A.S. Trebitz, A.M. Cotter, and G.S. Peterson. 2007. Human
influences on water quality in Great Lakes coastal wetlands. Environ. Management.
doi: 10.1007/s00267-007-9055-5
Mullholland, P.J. 2003. Large-scale patterns in dissolved organic carbon concentration,
flux, and sources, p. 139–159. In S.E.G. Findlay and R.L. Sinsabaugh [eds.], Aquatic
Ecosystems: Interactivity of Dissolved Organic Matter. Elsevier Science.
Olcott, P.G. 1992. Ground water atlas of the United States: Iowa, Michigan, Minnesota,
Wisconsin. In Ground Water Atlas of the United States by J.A. Miller. U.S.
Geological Survey. HA 730-J, Available from
http://pubs.usgs.gov/ha/ha730/ch_a/index.html. Accessed June 2010.
Olcott, P.G. 1995. Ground water atlas of the United States: Connecticut, Maine,
Massachusetts, New Hampshire, New York, Rhode Island, Vermont. In Ground
Water Atlas of the United States by J.A. Miller. U.S. Geological Survey. HA 730-M,
Available from http://pubs.usgs.gov/ha/ha730/ch_a/index.html. Accessed June 2010.
Omernik, J.M, C.M. Rohm, R.A. Lillie, and N. Mesner. 1991. Usefulness of natural
regions for lake management: analysis of variation among lakes in northwestern
Wisconsin, USA. Environ. Management. 15: 281–293.
Prepas, E.E., D. Planas, J.J. Gibson, D.H. Vitt, T.D. Prowse, W.P. Dinsmore, L.A.
Halsey, P.M. McEachern, S. Paquet, G.J. Scrimgeour, W.M. Tonn, C.A. Paszkowski,
and K. Wolfstein. 2001. Landscape variables influencing nutrients and phytoplankton
communities in Boreal Plain lakes of northern Alberta: a comparison of wetland- and
upland-dominated catchments. Can. J. Fish. Aquat. Sci. 58: 1286–1299.

57

Rasmussen, J.B., L. Godbout, M. Schallenberg. 1989. The humic content of lake water
and its relationship to watershed and lake morphometry. Limnol. Oceanogr. 34: 13361343.
Rohm, C. M., J. M Omernik, and C.W. Kiilsgaard. 1995. Regional patterns of total
phosphorus in lakes of the northeastern United States. Lake and Reservoir
Management. 11: 1–14.
SAS, Cary, NC
Schatzel, R.J. and S.A. Isard. 1991. The distribution of spodosol soils in southern
Michigan: a climatic interpretation. Annals of the Association of American
Geographers. 81: 425–442.
Sliva, L. and D.D. Williams. 2001. Buffer zone versus whole catchment approaches to
studying land use impact on river water quality. Water Resources. 35: 3462–3472.
Sobek, S., L.J. Tranvik, Y.T. Prairie, P.Korelainen, and J.J. Cole. 2007. Patterns and
regulation of dissolved organic carbon: an analysis of 7,500 widely distributed lakes.
Limnol. Oceanogr. 52(3): 1208–1219.
Soranno, P.A., S.L. Hubler, S.R. Carpenter, and R.C. Lathrop. 1996. Phosphorus loads to
surface waters: a simple model to account for spatial pattern of land use. Ecol. App.
6: 865–878.
Soranno, P.A., K.S. Cheruvelil, K.E. Webster, M.T. Bremigan, T. Wagner, and C.A.
Stow. 2010. Using landscape limnology to classify freshwater ecosystems for multiecosystem management and conservation. Bioscience. 60: 440 –454.
Strayer, D.L, R.E. Beighley, L.C. Thompson, S. Brooks, C. Nilsson, G. Pinay, and R.J.
Naiman. 2003. Effects of land cover on stream ecosystems: roles of empirical models
and scaling issues. Ecosystems. 6: 407–423.
Taranu, Z.E. and I. Gregory-Eaves. 2008. Quantifying relationships among phosphorus,
agriculture, and lake depth at an inter-regional scale. Ecosystems. doi:
10.1007/s10021-008-9153-0
Thierfelder, T. 1999. The role of catchment hydrology in the characterization of water
quality in glacial/boreal lakes. J. Hydrology. 216: 1–16.
Trebitz, A.S., J.C. Brazner, A.M. Cotter, M.L. Knuth, J.A. Morrice, G.S. Peterson, M.E.
Sierszen, J.A. Thompson, and J.R. Kelly. 2007. Water quality in Great Lakes coastal
wetlands: basin-wide patterns and responses to an anthropogenic disturbance
gradient. J. Great Lakes Res. 33: 67–85.

58

Wagner, T., D.B. Hayes, and M.T. Bremigan. 2006. Accounting for multilevel data
structures in fisheries data using mixed models. Fisheries. 31: 180–187.
Webster, K.E., P.A. Soranno, K.S. Cheruvelil, M.T. Bremigan, J.A. Downing, P.D. Vaux,
T.R. Asplund, L.C. Bacon, and J. Connor. 2008. An empirical evaluation of the
nutrient-color paradigm for lakes. Limnology and Oceanography. 53(3): 1137–1148.
Weller, C.M., M.C. Watzin, and D. Wang. 1996. Role of wetlands in reducing
phosphorus loading to surface water in eight watersheds in the Lake Champlain
Basin. Environmental Management. 20(5): 731–739.
Wickham, J.D., K.H Riitters, T.G. Wade, and K.B. Jones. 2005. Evaluating the relative
roles of ecological regions and land-cover composition for guiding establishment of
nutrient criteria. Landscape Ecology. 20:791–798.
Williamson, C.E., D.P. Morris, M.L. Pace, and O.G. Olson. 1999. Dissolved organic
carbon and nutrients as regulators of lake ecosystems: resurrection of a more
integrated paradigm. Limnology and Oceanography. 44(3): 795–803.
Wilson, H.F. and M.A. Xenopoulos. 2008. Ecosystem and seasonal control of stream
dissolved organic carbon along a gradient of land use. Ecosystems. 11: 555–568
Wolfson, L. 2009. Surface water resources, p. 206 –222. In R. Schaetzl, J. Darden, and D.
Brandt [eds.] , Michigan Geography and Geology. Pearson, New York, New York.
Wolock, D.M. 2003. Estimated mean annual natural ground-water recharge estimates in
the conterminous United States: U.S. Geological Survey Open-File Report 03-311,
digital dataset. Available from
http://water.usgs.gov/GIS/metadata/usgswrd/XML/rech48grd.xml . Accessed August
2003.
Xenopoulos, M.A., D.M. Lodge, J. Frentress, T.A. Kreps, S.D. Bridgham, E.Grossman,
& C. J. Jackson. 2003. Regional comparisons of watershed determinants of dissolved
organic carbon in temperate lakes from the Upper Great Lakes region and selected
regions globally. Limnology and Oceanography. 48: 2321–2334.

59