Assessment of suitable habitat for Phragmites australis ( common reed ) in the Great Lakes coastal zone

In the Laurentian Great Lakes, the invasive form of Phragmites australis (common reed) poses a threat to highly productive coastal wetlands and shorelines by forming impenetrable stands that outcompete native plants. Large, dominant stands can derail efforts to restore wetland ecosystems degraded by other stressors. To be proactive, landscape-level management of Phragmites requires information on the current spatial distribution of the species and a characterization of areas suitable for future colonization. Using a recent basin-scale map of this invasive plant’s distribution in the U.S. coastal zone of the Great Lakes, environmental data (e.g., soils, nutrients, disturbance, climate, topography), and climate predictions, we performed analyses of current and predicted suitable coastal habitat using boosted regression trees, a type of species distribution modeling. We also investigated differential influences of environmental variables in the upper lakes (Lakes Superior, Michigan, and Huron) and lower lakes (Lakes St. Clair, Erie, and Ontario). Basin-wide results showed that the coastal areas most vulnerable to Phragmites expansion were in close proximity to developed lands and had minimal topographic relief, poorly drained soils, and dense road networks. Elevated nutrients and proximity to agriculture also influenced the distribution of Phragmites. Climate predictions indicated an increase in suitable habitat in coastal Lakes Huron and Michigan in particular. The results of this study, combined with a publicly available online decision support tool, will enable resource managers and restoration practitioners to target and prioritize Phragmites control efforts in the Great Lakes coastal zone.


Introduction
Phragmites australis (Cav.)Trin.ex Steud.(henceforth Phragmites) is a highly invasive wetland plant growing to five meters in height and forming impenetrable stands that outcompete native plants and quickly degrade large areas of highly productive coastal wetlands and shorelines (Minchinton et al. 2006).Although there are native strains endemic to North America (Saltonstall 2002;Saltonstall et al. 2004), an invasive Eurasian strain began invading 150 years ago, threatens to degrade beaches and recreation areas, and alters natural habitat for fish and wildlife species (Weinstein and Balletto 1999;Meyerson et al. 2000).Current control methods include combinations of cutting, flooding, burning, and herbicide application, but treatments must be repeated annually or biannually to maintain benefits (Chambers 1997;Moreira et al. 1999;Ailstock et al. 2001;Warren et al. 2001;Carlson et al. 2009).Therefore, the cost to control established Phragmites stands is enormous, and sustained impacts are difficult to achieve.The need for improved strategies for limiting its spread (e.g., early identification of areas needing treatment) and managing existing populations is immense.
The attention that invasive species like Phragmites receive from the research and management community is growing.For example, the U.S. Environmental Protection Agency (EPA) recognizes invasive species as a significant stressor to aquatic ecosystems (EPA 2010).Invasion by Phragmites is a large-scale issue that requires a broad, coordinated approach.Landscapescale analyses are needed to help land managers target existing populations and identify areas most vulnerable to future invasions in order to limit spread, especially those areas being ecologically restored or those likely to experience the greatest change in regional climate.
Phragmites expansion has been particularly problematic in coastal areas surrounding the Laurentian Great Lakes in North America.Evidence suggests that invasion was assisted by seedbank inoculation and subsequent emergence on bottomlands exposed during the low lake levels that followed the 1997 high level (Wilcox 2012).Other researchers have found that an increase in the genetic diversity in a given area promotes seed viability and has led to expansion (McCormick et al. 2010;Kettenring et al. 2011).New research has shown that a second strain of invasive Phragmites is present in North America and may further enhance the genetic diversity needed for rapid expansion of Phragmites (Meyerson and Cronin 2013), underscoring the need for early-identification and rapid-response tools.
In addition to the aforementioned biotic factors, abiotic variables play a major role in Phragmites expansion.At site-specific scales, invasion success of Phragmites has been related to soil conditions (King et al. 2007;Tulbure and Johnston 2010), shoreline development and nitrogen pollution (King et al. 2007), disturbance (Kettenring et al. 2011), and site hydrology (Wilcox et al. 2003; Bart et al. 2006;Welch et al. 2006;Tulbure et al. 2007;Tulbure and Johnston 2010;Wilcox 2012).To our knowledge, however, associations between Phragmites distribution and environmental conditions have not been examined previously with continuous data across a large watershed such as the Great Lakes.Establishing such relationships will provide a better understanding of the invasion ecology of this species and will establish a foundation for mechanistic modeling that highlights the relative vulnerability of areas not yet impacted by this invasive species.A modeling environment also allows for exploration of the potential changes in Phragmites distribution under various temperature and precipitation scenarios related to global climate change.
For the purpose of developing a prioritized control strategy for Phragmites, estimating habitat suitability provides a means of assessing vulnerability to future invasions and expansion of existing Phragmites stands.Using species distribution modeling (SDM), we explored the following questions: 1) what environmental conditions are tied most strongly with Phragmites distribution in the Great Lakes coastal zone, 2) how do relationships between environmental variables and Phragmites distribution vary across the study area, 3) what areas are most vulnerable to Phragmites invasion and expansion, and 4) how could habitat suitability on the landscape change based on predicted climate changes?Utilizing remotely sensed non-native Phragmites distribution in the Great Lakes coastal zone that was mapped by Bourgeau-Chavez et al. 2013, our specific objectives were to examine associations with environmental conditions using SDM, assess areas vulnerable to invasion, and provide an online tool that will enable resource managers and restoration practitioners to target and prioritize control efforts.

Methods
To achieve these objectives, we used SDM, which estimates the ecological niche of a species based on overlap between species observations and a suite of environmental conditions (e.g., soil hydrologic type, road density, precipitation) in the same physical space.Specifically, we used a boosted regression tree (BRT) approach that is consistently ranked as one of the best performing SDMs (Elith et al. 2006;Elith and Graham 2009;Stohlgren et al. 2010).Our Phragmites occurrence data were obtained from a recent, basin-scale mapping effort performed from radar images collected between 2008 and 2010 within a 10-km inland buffer of the Great Lakes shoreline (Borgeau-Chavez et al. 2013).Our SDM followed the same spatial extent.To ascertain how future climate changes might affect Phragmites invasion success in the existing landscape (i.e., bottomlands potentially exposed by low water levels were not included), we parameterized the BRT model with modeled projections for precipitation and air temperature to obtain a projected vulnerability map for the year 2050.

Study area
The spatial extent of our species distribution model represented a 10-km inland buffer of the U.S. Great Lakes shoreline, included islands where available (e.g., Les Cheneaux Islands in Lake Huron, West Sister Island in Lake Erie) but excluded exposed bottomlands due to data limitations, and matched that of the Phragmites mapping effort by Bourgeau-Chavez et al. (2013) (Figure 1).We chose two spatial extents for investigating relationships between environmental conditions and Phragmites distribution: 1) basinwide and 2) upper and lower basins.For the latter, we developed separate models for the upper lakes and lower lakes.These two extents

Conceptual model and variable selection
To guide our SDM approach, we developed a conceptual causal model (Figure 2) that described our hypotheses relating environmental variables to Phragmites presence based on theoretical and empirical relationships.Our conceptual model was not necessarily comprehensive but rather included variables that could be obtained from the processing of pre-existing geospatial datasets.
In our conceptual model, the left half of the diagram represents variables considered to be natural drivers (e.g., climate, topography, soils, nutrients) (Figure 2).Toward the right of the diagram are anthropogenic drivers (e.g., nutrients again, urban development, agriculture, road density, population density, land-cover change), all of which influence the disturbance regime, also included in the diagram.Our implementation of the conceptual model into an SDM modeling environment allowed us to test differences in effects of environmental variables at two scales: basinwide and at a sub-basin scale representing the upper lakes (Lake Superior, Lake Michigan, and Lake Huron) and lower lakes (Lake St. Clair, Lake Erie, Lake Ontario).For example, we suspected that greater development, greater agriculture, and lesser soil permeability in the lower lakes may lead to differences in the way that environmental variables, either individually or interacting, affect Phragmites distribution.
Our response variable (Phragmites presence/ absence) resulted from the mapping methodology of Bourgeau-Chavez et al. (2013), which detected mature, monotypic stands of invasive Phragmites (0.2 ha in size or larger) and did not detect the native strain.To produce the Phragmites map, Bourgeau-Chavez et al. ( 2013) relied primarily on  2013).The C-CAP/CDL filter retained points in the data set that fell within classifications that were wetland-related and, therefore, likely to support Phragmites (i.e., C-CAP: grassland herbaceous, palustrine forested wetland, palustrine scrubshrub wetland, open water, palustrine aquatic bed, and unconsolidated shore; CDL: open water, shrubland, grassland herbaceous, and herbaceous wetlands).Therefore, we avoided using land cover directly as an independent environmental variable because it was assumed to act as an artificially strong predictor of Phragmites.
Existing geospatial data that spanned the entire study area were analyzed in a GIS (ESRI ArcGIS Desktop version 9.3 and 10.0) to determine the independent variables to include in our models.Variables were considered that either directly described environmental conditions pertinent to Phragmites distribution or that could act as surrogates where equivalent data were either unavailable or difficult to define in a spatially explicit manner.Consideration was given to 1) the watershed delivery processes that impact nutrient, sediment, and propagule dynamics; 2) phenomena that may influence Phragmites establishment based on landscape context; and 3) variables that affect Phragmites growth locally.When multiple environmental variables provided similar information, we chose the dataset that reported the greatest accuracy in the metadata and had the better spatial precision.An assessment of colinearity identified cases of quantifiable redundancy within the data set, and one variable from pairs with a Pearson correlation coefficient greater than 0.8 was eliminated prior to preliminary species distribution modeling.Variables with negligible influence on habitat suitability estimates during preliminary model parameterization were identified and eliminated prior to final model fits.The resulting set of 15 variables was used for the remainder of our analyses (Table 1).
As implemented in the SDM, our environmental variables fell into one of six categories: topography, disturbance, ecoregion, soils, nutrients, and climate (Table 1).All variables included native data taken directly from the source with the exception of the following derived variables: road density  (NOAA 1995-present).Nutrient data were comprised of watershed-scale modeled nitrogen (N-Conc) and phosphorus (P-Conc) concentrations in stream water (Robertson and Saad 2011) and were included as a proxy for soil nutrients, which were not available as a geospatial dataset on a large scale.Topographic  roughness was calculated as a 50-cell radius standard deviation of elevation from the one-arcsecond U.S. Geological Survey National Elevation Dataset (Gesch et al. 2002).We used three climate variables from the WORLDCLIM dataset (Hijmans et al. 2005): temperature of the coldest quarteryear (WinTemp), temperature of the driest quarteryear (DryTemp), and precipitation of the coldest quarter-year (WinPrecip).

Species distribution modeling
The BRT method of SDM is a machine learning technique that accommodates nonlinear relationships among categorical and numeric variables and automatically detects and accounts for interactions among predictors (Elith et al. 2008).This method requires both presence and absence observations but can handle datasets containing missing data.The BRT modeling combines two approaches: classification and regression trees, and boosting (De'ath 2007;Elith et al. 2008).As commonly implemented in BRT for binary response variables such as a species' presence or absence, classification and regression trees are based on decision trees that estimate the mean value for a subset of the data through a logistic regression approach (Breiman et al. 1984).Boosting refers to the combination of numerous models, which may be weak predictors individually, to produce a robust final model (Friedman 2002).The BRT analysis produces outputs on a continuous scale that represent the probability of occurrence of a response variable, in this case, Phragmites presence.When paired with a threshold analysis, this output represents the magnitude of vulnerability to future Phragmites invasion, a basin-wide issue of concern to private and public land management entities at all jurisdictional levels.
We used BRT methods to 1) produce a spatially continuous estimate of relative habitat suitability that assumes full dispersal and 2) generate a vulnerability assessment representing the potential to support Phragmites growth under current environmental conditions.BRT was implemented, interpreted, and post-processed using R version 2.14.0 (Pinheiro et al. 2011) and three R packages: DISMO (Hijmans and Elith 2013), GBM (Ridgeway 2010), and SDMTools (VanDerWal et al. 2012).In addition to suitability values provided as a geospatial dataset, BRT output included the relative percent influence of each environmental variable on Phragmites presence.Marginal response curves were used to interpret apparent trends between environmental variables and Phragmites distribution, recognizing the limitations of such analysis given that the mean values of all other variables are held constant in producing these curves.
BRT automatically models interactions between predictor variables due to the dependence of the model response to a given variable on the values of other predictors at higher tree levels (Elith et al. 2008).By default, the fifteen strongest twoway interactions were plotted as three-dimensional partial dependence plots (x-axis: first variable, y-axis: second variable, z-axis: marginal effect) to elucidate effects that interactions between variables had on Phragmites distribution.
Model fit assessment followed Elith et al. (2008) and included calculation of percent deviance explained and reduction in the area under the receiver operating characteristic curve (AUC).Models were trained on half of the available data for the basin-wide, upper basin, and lower basin spatial extents.Model performance was assessed on the remaining half of each dataset.Model accuracy (i.e., how often the model correctly predicted Phragmites presence and absence) was assessed using a threshold approach described in further detail below.
Since Phragmites prevalence was higher in the lower basin than the upper (Bourgeau-Chavez et al. 2013), we parameterized each of the three models (basin-wide, upper, lower) independently to obtain the best model fit.For each model, we optimized four parameters in R: number of trees, learning rate, tree complexity, and bag fraction.The optimal number of trees, the learning rate, and tree complexity were assessed using the gbm.step function on randomly selected 225,000point subsamples of our basin-wide, upper, and lower basin datasets.Subsamples were used during parameterization for computational efficiency because our full dataset exceeded the computational limits of our hardware and software.Following Elith et al. (2008), the gbm.step function used a ten-fold cross validation procedure to determine an optimal number of trees, and we investigated the effects of varying bag fraction and variable order using the gbm.fixed function.Final model training was completed using a bag fraction set to 1.0.This eliminated the stochastic component of the gradient boosting algorithm and was justified by preliminary model runs using an independently withheld subset for testing where a small but incremental increase in model performance was documented as the bag fraction was increased from 0.5 to 1.0.We also investigated the effects of reordering environmental variables and found that changing the variable order had a minor impact on model performance when compared to varying bag fraction values.Therefore, we maintained an arbitrary but consistent variable order in all models.
To infer predicted Phragmites presence or absence in the landscape, we chose an approach that optimized a presence-absence threshold (Liu et al. 2005).A consequence of working with large presence-absence datasets is a bias in model output toward the more commonly observed outcome, in this case, Phragmites absence (Cramer 1999;Jiménez-Valverde and Lobo 2006).In the case of very low prevalence samples, this effect produces outputs where a high proportion of presence observations will have low probability of occurrence scores.Although models from unbalanced samples can accurately depict relative habitat suitability within their domain, the numeric values of the output are counterintuitive to many users, and output from models with different prevalence rates are not directly comparable.To deal with this issue, we used a threshold approach determined in R using the SDMTools package (VanDerWal et al. 2012).By adopting a threshold that balanced the rates of successfully predicting presence and absence, we objectively defined a probability value above and below which habitat suitability conditions were generally more suitable and less suitable, respectively.Establishing this threshold provided a means of developing an intuitive and interpretable classification scheme for land managers to assess relative vulnerability to Phragmites invasion.In addition, adoption of a threshold permitted the calculation of model percent accuracy and implementation of a habitat suitability index that allowed comparisons among models and provided a convenient means for estimating area of suitable habitat.
More and less suitable areas for Phragmites habitat were summarized basin-wide and by lake, upper and lower basins, and county for the entire study area and for the wetland-related areas determined by the C-CAP/CDL filter described above.Whereas Phragmites mapping included offshore areas, the SDM did not due to lack of environmental data for areas off shore.Since many practitioners may be interested in the distribution of Phragmites off shore, especially given the amount of fertile bottomlands that could be exposed during predicted reductions in lake levels, we included county-level summaries of Phragmites area in offshore areas that are available at http://cida.usgs.gov/glri/phragmites/.

Climate change analysis
After developing a habitat suitability model based on current environmental conditions, we explored how the distribution of suitable Phragmites habitat might change with future climate scenarios, as applied to the basin-wide model.The three available environmental variables that described climate parameters (DryTemp, WinTemp, WinPre-cip; see Table 1) were modified in the prediction of suitable habitat given the Intergovernmental Panel on Climate Change (IPCC) A2a future emissions scenario for 2050 (Nakićenović and Swart 2000), reflecting a substantial increase in greenhouse gas emissions.The source of the input data for these future climate variables was the Institut Pierre Simon Laplace des Sciences de l'Environnement Global Climate Model, version 4 (IPSL CM4) (Marti et al. 2010).The IPSL CM4 model was appropriate because it directly accounted for the presence of the Great Lakes (Marti et al. 2010) and was more sensitive to greenhouse gas increases than other models (Solomon et al. 2007).
Data output from the climate-change models were downscaled based on the Delta method (Ramirez-Villegas and Jarvis 2010), which uses a splining interpolation to increase model output resolution to that of 30 arc-seconds (~1 km).This resolution matched that of the three climate variables used to develop our BRT (Hijmans et al. 2005), thus allowing the direct substitution of a future climate scenario into our modeling framework.

Results
Our results showed that the areas most closely associated with Phragmites expansion in the Great Lakes coastal zone were those places where the influences of minimal topographic relief, close proximity to developed areas, poorly drained soils, and dense road networks promoted the colonization and persistence of Phragmites stands (Table 2).The relative influence of these variables on Phragmites presence in the upper lakes (Superior, Huron, and Michigan) and lower lakes (Ontario, Erie) were similar to basin-wide results (Figure 3).However, whereas watershed nitrogen concentration exerted a greater relative influence in the upper basin, proximity to agricultural land cover was more influential in the lower basin (Table 2).In addition, Phragmites tended to be more prevalent in regions that experienced milder winters with less precipitation, particularly in the upper lakes (Figure 4).

Landscape characteristics
The most influential environmental variable of Phragmites presence in the Great Lakes in our study was topographic roughness (TopoRough) (Table 2).The soil hydrologic group (SoilHydro) also influenced habitat suitability (Table 2).This categorical variable places soils in one of four groups based on runoff potential, with an additional three dual categories for soils that can be drained in order to lower the water table.Soils with a high runoff potential, slow infiltration rate, or high water table (Group D) corresponded with greater Phragmites prevalence, as did the three dual categories (Groups A/D, B/D, and C/D) (Figure 3).

Anthropogenic influences
All three models (basin-wide, upper, lower) suggested that closer proximity to developed areas (ProxDev) positively influenced Phragmites presence (Table 2, Figure 3).Based on the marginal response curves, the habitat suitability was very poor at zero distance to development but relatively high a short distance away from developed areas (Figure 3).Greater road density in the catchment also was associated with greater Phragmites presence, to a point.The structure of this relationship also differed between the models; the basin-wide and lower lakes models had relatively similar responses that did not show a distinct trend (Figure 3a and 3c), whereas the upper lakes model showed a peak in Phragmites habitat suitability at an intermediate level of road density (Figure 3b).
We observed differences in the relative influence of some anthropogenic environmental variables on Phragmites distribution between the upper and lower lakes as well.For example, nitrogen concentration (N-Conc) in the catchment had a much larger relative influence in the upper lakes than in the lower lakes (Table 2) even though the upper basin had lower catchment nitrogen values overall in the coastal zone (mean N-Conc upper = 3.007.43mg/L, mean N-Conc lower = 7.0920.7 mg/L).Proximity to agricultural land cover (ProxAg) had a greater relative influence in the lower lakes than in the upper lakes (Table 2), reflecting a greater prevalence of agriculture in the lower lakes (49.4%) than in the upper lakes (26.8%).
Although an interaction between N-Conc and ProxAg was expected based on our conceptual model (Figure 2), this interaction was not present in the fifteen strongest two-way interactions.Any interaction effects were accounted for automatically in the BRT method.Additionally, the Pearson correlation coefficients between these two predictors were -0.02, -0.10, and 0.14 for the   basin-wide, upper, and lower basin models, respectively, further indicating minimal influence of ProxAg on N-Conc.We note, however, that the N-Conc did not represent soil nutrient levels but rather modeled concentrations in stream water.

Effects of temperature and precipitation
The temperature and precipitation variables (WinTemp, DryTemp, WinPrecip) did not have substantial effects on Phragmites presence individually, but a few notable patterns emerged collectively.The relative influence of these three variables was higher in the upper lakes than the lower lakes (Table 2).Taken in sum, their influence contributed 16.4% of the gains seen in the upper lakes model as opposed to 9.1% seen in the lower lakes model.Of the three variables, the most influential in determining Phragmites abundance basin-wide was precipitation during the coldest quarter (WinPrecip) (Table 2).The mean temperature of the coldest quarter (WinTemp), however, showed a greater relative influence in the upper basin than did WinPrecip.The mean temperature of the driest quarter (DryTemp) was minimally influential on Phragmites abundance for all three models (Table 2).
An analysis of pair-wise variable interactions suggested additional complexity in the role of these variables.In the basin-wide model, the strongest interactions were detected between WinTemp and WinPrecip.Perspective plots of marginal response curves for two-way interactions showed that places that were warmer and had less precipitation during winter months were more suitable for Phragmites in the basin-wide model (Figure 4a).In contrast, the strongest interactions in the upper lakes model were detected between WinPrecip and DryTemp (Figure 4c).The strongest interaction among these variables in the lower lakes model was observed between WinPrecip and WinTemp (Figure 4b), although these effects were far less dramatic than in the upper-lakes model.

Influence of ecoregion
The level II ecoregion (EcoregII), an ecologically relevant characterization of the large-scale landscape, was moderately influential in the basinwide model but was not influential in the upperand lower-lakes models (Table 2), as the latter two were relatively homogenous in ecoregion classification.The upper basins contained three ecoregions (5.2: mixed wood shield, 8.1: mixed wood plains, and 8.2: central USA plains), whereas the lower basin was composed primarily of ecoregion 8.1 along with a small area of ecoregion 8.2 on western Lake Erie basin.Phragmites presence was particularly high in ecoregion 8.2, which encompasses Saginaw Bay and the western basin of Lake Erie.

Model performance and accuracy
Our sampling approach resulted in a low observed Phragmites frequency of occurrence (0.37 %), although a large total number of presence observations across the study area (18,549) provided a rich data set that was used to train and test models.BRT learning rates of 0.3, 0.2, and 0.3 for the basin-wide, upper lakes, and lower lakes data sets, respectively, resulted.The machinecalculated optimal number of boosted regression trees for these data sets was 2,900, 3,200, and 3,750, respectively.
Our model fit statistics (percent deviance explained >50 and area under the receiver operating characteristic curve >95) for all three models (basin-wide, upper, lower) suggested that our models performed well and accurately represented the habitat suitability for Phragmites in the Great Lakes coastal zone (Table 3).The threshold methodology we used to differentiate more and less suitable habitat resulted in the following overall accuracies: basin-wide 93.0%, upper basin 94.7%, and lower basin 90.1%.

Vulnerability to Phragmites expansion
Our models showed that some areas of the Great Lakes contained a preponderance of suitable habitat and also supported extensive Phragmites stands (Figure 5a).Other regions contained vast areas of suitable habitat that did not yet support widespread Phragmites colonization, representing areas more vulnerable to future Phragmites expansion (Figure 5b).Of particular note is the paucity of currently occupied areas as suitable habitat in Lakes Ontario and Superior (Table 4).
Based on our analysis of 2008-2010 data, 0.36% of the area in the coastal zone was occupied  4) in the lower basin (105.7 km 2 , 0.92%) was also proportionately greater than in the upper basin (76.9 km 2 , 0.20%).

Climate change predictions
The IPSL future climate predictions for 2050 suggested that winters will become warmer and wetter, and that the driest quarter-year will be warmer, on average (Table 5).The range of values for the average temperature of the driest quarter-year (DryTemp) will expand from that of baseline conditions.Additionally, interactions between winter precipitation (WinPrecip) and winter temperature (WinTemp) were distributed across a greater temperature range in the future climate scenario results and were shifted towards warmer temperatures and greater precipitation,  indicating a more extensive effect of this interaction on Phragmites suitability in the future (Figure 6).The basin-wide habitat suitability model parameterized with the IPSL-CM4-predicted climate variables (DryTemp, WinPrecip, and WinTemp) suggested that suitable Phragmites habitat is likely to become more prevalent under the IPCC A2a future climate scenario (Figure 7).The area of current (2008)(2009)(2010) suitable habitat (colonized and not yet colonized) amounted to 3557 km 2 , or 7.0% of the total coastal zone (50,550 km 2 ) (Table 6).Our future model parameterization suggested that 25.7% (12,975 km 2 ) of the coastal zone will provide suitable habitat for Phragmites.This increase in suitable habitat will be disproportionately distributed in the upper basin, which predictions suggest will increase from 3.8% (1472 km 2 ) to 27.1% (10,560 km 2 ) of the total area.A smaller increase from 18.1% (2085 km 2 ) to 21.0% (2415 km 2 ) in suitable habitat is expected in the lower basin (Table 6).Our future climate predictions suggest a slight decrease in suitable habitat in Lakes Erie and St. Clair; however, the uncertainty in the climate predictions is likely to be well within the error bounds of these estimates of future suitable habitat.

Discussion
Habitat suitability for Phragmites Using our improved understanding of the major large-scale environmental variables influencing Phragmites presence in the Great Lakes, we found that the places most vulnerable to invasion were disproportionately located in the lower two lakes and included major areas along the western and central Lake Erie coastal zone and Lake St. Clair delta.However, climate predictions suggested a substantial increase in susceptibility to Phragmites invasion by 2050 in the coastal landscape of the upper lakes.Presently, relatively small portions of the Lake Ontario and Lake Superior coastal zone were affected.In the upper Great Lakes, the primary expansion is likely to occur in the Green Bay area and eastern coastal zone of Lake Michigan, as well as within the remaining portions of Lake Huron's Saginaw Bay area that have not yet been invaded by Phragmites.In areas where large tracts of contiguous suitable habitat exist, Phragmites is often already well-established in the landscape, so unoccupied wetlands and riparian zones are at high risk for future Phragmites invasion due to proximity to current seed sources as well as suitable environmental conditions.

Primary environmental drivers of Phragmites prevalence
Since Phragmites is a wetland plant and therefore grows in areas of low topographic relief and on poorly drained soils, it is no surprise that mature stands persist in such places.Nonetheless, it is telling that the local topography has the greatest influence of the environmental variables tested.Despite the ability of Phragmites to translocate water, gasses, and nutrients to advantageous areas, thereby enhancing its invasion success (Granéli et al. 1992;Asaeda and Rajapakse 2006;Tulbure 2008), our findings suggest that complex terrain (high degree of topographic roughness) limits its prevalence.These topographic highs may represent geographic barriers that restrict its vegetative expansion and seed dispersal.Although our full-dispersal modeling approach does not take into account these potential limitations to dispersal, our approach is more conservative and, therefore, likely to be more useful from an invasive species management perspective (Bateman et al. 2013).
Our results also showed that closer proximity to developed areas in the landscape increased the prevalence of mature Phragmites stands.If developed areas are taken as a proxy for disturbance, this result supports, on a large scale, the findings of Silliman and Bertness (2004) and King et al. (2007), who found similar effects of shoreline development, and McCormick et al. (2010) and Kettenring et al. (2011), who found that disturbance improves seed viability by increasing genetic diversity and, therefore, the probability of conspecific cross-pollination.Leblanc et al. (2010) showed, however, that development needs to coincide spatially and temporally with Phragmites dispersal phenology in order for disturbance to impact invasion success.
The influence of road networks as dispersal pathways for Phragmites has been documented by many (Jodoin et al. 2008;Lelong et al. 2007;Maheu-Giroux and de Blois 2007;Brisson et al. 2010;Leblanc et al. 2010).The creation of contiguous wetland habitat in the form of drainage ditches along roads, along with the use of de-icing salts and agricultural fertilizers, may favor Phragmites establishment (Jodoin et al. 2008).The high relative importance of road density observed in our study (Table 2) confirms the large-scale importance of roads as vectors of invasion.This is especially illuminating given that the map of Phragmites distribution used to develop the habitat suitability model consisted of mature stands greater than 0.2 ha and could not detect linear stands along roadways.Nonetheless, where road densities were greater, Phragmites was more prevalent.
Our finding that closer proximity to agricultural land use had a substantial influence on Phragmites presence suggests that soil nutrients may play a role in Phragmites invasion.With the exception of catchment nitrogen concentration in the upper lakes, however, we found that nutrient concentration in stream water was not highly influential in determining Phragmites abundance overall.Although other researchers (e.g., Minchinton and Bertness 2003;Saltonstall and Stevenson 2007) have observed that elevated soil nutrients improved the competitive ability of Phragmites by allowing for greater aboveground allocation of resources, we did not observe a substantial effect of nutrients levels in stream water on a basinwide scale.One potential explanation for this is that the nitrogen and phosphorus concentrations represented concen-trations in stream water at the catchment scale rather than soil values at a local scale and were derived from nutrient models (Robertson and Saad 2011).Our intent in using these modeled nutrient data was to investigate the possibility that nutrient levels in stream water could be used as a proxy for available soil nutrients, being that soils can be a major source of in-stream nutrients via subsurface flow pathways (Peterjohn and Correll 1984).However, the linkage between soil and stream nutrient content may be weak on a basin-wide scale.If soil nutrient data were available basin-wide, this clearly would be the preferred dataset.A second explanation may be that the nutrient concentrations observed in the majority of the basin are elevated beyond the threshold at which Phragmites holds a competitive advantage.This might explain why nitrogen concentration was more influential in the upper lakes than in the lower, as the upper lakes nutrient levels were much lower than in the lower lakes.

SDM limitations and uncertainties
Several limitations to SDM are worth mentioning.Ecological niche models assume conservation of the prior niche when a species invades a new area (Fitzpatrick and Weltzin 2005).Phragmites continues to expand and, therefore, has not yet realized its full niche (e.g., Tulbure and Johnston 2010).Although using a correlative model based on a fully realized niche to inform a vulnerability assessment may be preferred (Beaumont et al. 2009), such data are seldom available.Therefore, we used the current distribution of Phragmites in its unrealized niche to predict suitable habitat, an approach that has been used by many and has wide empirical and theoretical support [see review by Peterson (2003)].
The results are subject to inaccuracies of the input datasets as well, including the Phragmites distribution data that had an overall map accuracy of 87% (Bourgeau-Chavez et al. 2013).In addition, causal linkages could not be detected using the methods in this study.However, the correlations detected by this work are compelling and can be used to inform experimental studies incorporating mechanistic growth models.Finally, spatial and temporal limitations were inherent to the data used to generate the Phragmites map and environmental predictor datasets.For example, error can be introduced into large geospatial datasets that, by necessity, are compiled over multiple years.An acknowledgment of such limitations is needed in any such large-scale modeling exercise.
Lastly, while BRT directly accounts for nonlinearity and interactions between predictors, we assumed stationarity of the environmental variables and minimal spatial autocorrelation, which may have introduced error to the modeling process by increasing Type I errors (Hothorn el al. 2011;Crase et al. 2012;Dormann et al. 2007).We chose a 100-m by 100-m sampling grid to balance a reduction in spatial autocorrelation and, at the same time, allow for a substantial dataset that captured Phragmites presence and absence in a large landscape.Our qualitative analysis of residuals showed very little spatial patterning except along some narrow invading fronts of Phragmites, where high residuals suggested poor model performance.With one recent notable exception (Crase et al. 2012), spatial autocorrelation in BRT has received little attention.The methodology of Crase et al. (2012) that accounted for spatial autocorrelation in a mangrove dataset was published too recently to be incorporated into this study.However, they found that spatial autocorrelation was minimized at a distance of approximately five 25-m by 25m grid cells, which supports our assertion that our 100-m by 100-m grid cell sampling protocol likely produced a dataset with minimal spatial autocorrelation.Given that our training and testing datasets were equal and robust, our good model performance suggested minimal spatial autocorrelation, as well.

Differential responses between upper and lower lakes
Many of the environmental drivers of Phragmites distribution in the Great Lakes showed little difference between the upper and lower basins, but several notable differences were apparent.In addition to the differential response to nitrogen concentration observed in the upper lakes discussed above, the mean winter temperature (WinTemp) and precipitation (WinPrecip) also were more influential in the upper basin.Furthermore, these two variables showed the strongest two-way interaction (Figure 4a), indicating a preference of Phragmites for mild winters with minimal precipitation.The observation that these factors were more influential in the upper basin likely is due to the greater latitudinal range, and therefore climatic range, observed within the upper lakes than the lower (Figure 1).
In the lower lakes, close proximity to agricultural land cover was related to increased Phragmites presence, which supports the findings of several others (i.e., Jodoin et al. 2008;King et al. 2007;Chambers et al. 2008).Agricultural lands contribute nitrogen and other nutrients to the surrounding landscape.Because Phragmites has a competitive advantage when nitrogen levels are high (Minchinton and Bertness 2003;Kettenring et al. 2011), areas closer to agricultural lands could be more susceptible to Phragmites colonization.This phenomenon is likely to be more localized than watershed-scale nitrogen concentrations would suggest, which might explain why the nitrogen concentration in our study was not highly influential in determining Phragmites distribution in the lower lakes.The detected linkage between nitrogen concentration in stream water and Phragmites in the upper basin also suggests that nutrients play a role in Phragmites distribution.Why the lower and upper basins would differ in their specific response to nutrient-related predictors is an area requiring further research.
Of particular note is the paucity of Phragmites in Lakes Ontario and Superior.The mapping of Phragmites from radar did not detect mature stands of Phragmites in Lake Superior, although some smaller stands do exist and were detected in the field verification process, and very little Phragmites was observed in Lake Ontario (Bourgeau-Chavez et al. 2013).Hypotheses for the lack of Phragmites in these lake basins include 1) that the coastal zone of these lakes represents an unrealized niche and 2) that Typha has a competitive advantage over Phragmites in Lake Ontario, in particular, limiting the invasion success of the latter species (Tulbure and Johnston 2010).Our results, however, suggest that the suite of environmental conditions present in Lake Superior and Ontario is a major limiting factor to invasion.Almost certainly, Phragmites has not yet realized its full niche in the Great Lakes coastal zone, particularly when taken within the context of impending climate change, but its spread likely will be limited in areas such as Lake Superior and Ontario where environmental conditions are less favorable (Table 6).

Precipitation, temperature, and implications for climate change
The observation that the precipitation and temperature variables had less influence on Phragmites habitat suitability than most other variables included in the study suggests that landscape alteration could have a greater impact on Phragmites expansion than changes in climate.Nonetheless, our results also indicate that a major expansion of suitable habitat for Phragmites is expected in the landward coastal zone of Lakes Michigan and Huron by 2050 due to changes in climate, primarily due to drier and milder winters.Some areas along Lakes Ontario and Superior also are predicted to become more suitable for Phragmites despite the observation that the environmental conditions along these two lakes are not generally suitable at present.This highlights the need for early attention to controlling Phragmites suggested by many others (e.g., Kowalski and Wilcox 1999;Wilcox and Whillans 1999;Saltonstall and Stevenson 2007) before it can take hold and form extensive monocultures that reduce biodiversity (Weinstein and Balletto 1999;Keller 2000;Ailstock et al. 2001;Carlson et al. 2009), among other threats.The findings of Kettenring et al. (2011) and Meyerson and Cronin (2013) are particularly relevant in the context of changing climate and the resulting expansion of suitable habitat.Since greater diversity results in the production of more viable seeds, positive feedbacks between genetic diversity and population size may become more ubiquitous in the landscape, further underscoring the need for early detection and control.
Further complicating matters, the substitution of the 2050 climate predictions of the IPCC A2a climate scenario for the three climate variables in the SDM was implemented assuming that other variables would remain unchanged, which is unlikely to be the case; for example, climate changes could lead to changes in land cover, or socioeconomic systems could respond to climate change by altering land use.Additionally, altered lake levels almost certainly will result from climate change, and if large expanses of fertile unconsolidated sediment become exposed where Phragmites propagules are readily available, then expansion would accelerate (Tulbure and Johnston 2010;Wilcox 2012).This acceleration likely would be most prominent in areas with low bathymetric slope (e.g., western Lake Erie, Saginaw Bay in Lake Huron).Thus, the vulnerability of lakeward expansion of Phragmites will vary geographically.
Reduction or maintenance of suitable habitat for Phragmites is expected in Lakes Erie and St. Clair, which is somewhat encouraging since the coastal zones of these two lakes have been most affected per unit area by Phragmites expansion to date (Table 4).Predictions of the effects of climate change on suitable habitat for Phragmites, however, must be interpreted with caution for four main reasons.First, predictions for the climate variables are computed from the best models we have, but uncertainties are inherent in the models.Second, this niche-based climate modeling approach does not take into account adaptations or responses to changes in environmental conditions by Phragmites or its native and non-native competitors that may result from changes in climate (Morin and Thuiller 2009).Third, interactions between variables may occur in the future that were not accounted for in the baseline model.For example, warmer temperatures may lead to an expansion of agriculture in northern latitudes, thereby obscuring the true impacts of air temperature.Fourth, extrapolation was not attainable within this modeling framework for cases where values in the future climate dataset were outside the range of the baseline data used in model development.Rather, model response to the future climate data was limited to the maximum response set by baseline conditions, and thus may not accurately predict ecological response to a novel future climate (Williams and Jackson 2007).Nonetheless, our projections of future suitable habitat for Phragmites, when taken with these cautions in mind, can assist managers focus early detection efforts in areas identified as most vulnerable to invasion.Early identification allows for implementation of management measures before Phragmites plants get established completely and control becomes extremely difficult.

Conclusion
This study has demonstrated that the main drivers of Phragmites distribution in the Great Lakes coastal zone are topographic roughness, proximity to development, extent of soil drainage, and road density.Only nitrogen concentration in the upper basin and proximity to agriculture in the lower basin were differentially influential, but both point toward an influence of nutrients on Phragmites distribution.Our results suggest that climate changes will substantially increase habitat that is suitable for Phragmites, particularly in Lakes Michigan and Huron; however, climate factors had a lesser influence on Phragmites distribution than the aforementioned environmental variables, suggesting that climate change may have a lesser, albeit considerable, impact on landward invasion of Phragmites and expansion than would landscape modifications.
To make these data and model results available to the public and useful to managers, we developed an online decision-support tool (DST) (http://cida.usgs.gov/glri/phragmites/)that includes clickable, zoomable maps of 1) current Phragmites distribution (Bourgeau-Chavez et al. 2013); 2) corridor networks weighted by proximity to Phragmites and comprised of streams, wetlands, inland water bodies, and coastal corridors exposed during reduced lake levels; and 3) Phragmites habitat suitability.Users are able to view the data covering their area of interest, read metadata for the DST, and use a custom interface to download data.The DST serves as both a tool to explore maps online and disseminate data for individualized analysis.Natural resource managers dealing with Phragmites can use the map to identify potential sources of Phragmites invasion to sensitive areas (e.g., restoration or construction sites), where bare soils may need to be revegetated quickly to reduce the probability of Phragmites invasion.

Figure 1 .
Figure 1.Study area in white representing a 10-km buffer of the U.S. Great Lakes shoreline that was mapped by Bourgeau-Chavez et al. (2013) and subsequently analyzed using species distribution modeling.Sites A and B correspond to Figure 5 and Figure 7.

Figure 2 .
Figure 2. Conceptual causal model of environmental effects on Phragmites distribution in the coastal zone of the Great Lakes.

Figure 3 .
Figure 3. Response curves of the top five variables for the basin-wide (a), upper basin (b), and lower basin (c) models.Higher fitted function values represent relatively favorable Phragmites habitat conditions based on the values of a single variable.Decile rugs at the top of the plots show the sample distribution across the range of each variable.Curves are ordered from left to right in decreasing relative influence (see Table2).

Figure 4 .
Figure 4. Two-way interaction plots for variables representing precipitation of the coldest quarter (WinPrecip; winter precipitation), temperature of the coldest quarter (WinTemp; winter temperature), and temperature of the driest quarter (DryTemp; driest quarter-year temperature).Axes were scaled identically in order to compare differences in climate interactions for the three models.

Figure 5 .
Figure 5. Habitat suitability for two areas in the Great Lakes: near Bay City, MI on Saginaw Bay of Lake Huron (A) and near Muskegon, MI on Lake Michigan (B).

Figure 6 .
Figure 6.Two-way interactions between WinPrecip and WinTemp for the basin-wide model for baseline conditions (A) and the IPCC A2a future climate scenario (B).

Figure 7 .
Figure 7. Projected future (IPSL A2a climate predictions for 2050) habitat suitability for areas near Bay City, MI on Saginaw Bay of Lake Huron (A) and near Muskegon, MI on Lake Michigan (B).Climate models predict warmer and wetter winters in 2050, as well as warmer temperatures during the driest quarter-year.The model suggests that these conditions will generally support increased habitat suitability for Phragmites, with the trend more pronounced in the northern portion of the study area.

Table 1 .
Environmental variables used to model Phragmites occurrence in the coastal zone of the Great Lakes.

Table 2 .
Relative contributions of environmental variables to Phragmites distribution from the three model fits (basin-wide, upper basins, lower basins).Variable descriptions are provided in Table1.

Table 3
Predictive performance of the BRT models for the basin-wide, upper, and lower model parameterizations, as calculated from independent test data.Greater percent deviance explained and area under the receiver operating characteristic (ROC) curve values close to 1.0 indicate better model accuracy.

Table 4 .
Summary of total area (km 2 ) classified as more suitable or less suitable Phragmites habitat by the habitat suitability model.Total area supporting existing Phragmites stands is summarized by Great Lake, upper basins, lower basins, and entire basin.Wetland-Related Areas within NOAA C-CAP and USDA CDL filters

Table 5 .
Summary statistics of climate data used to develop habitat suitability BRT and estimate future (2050) habitat suitability under the IPCC A2a future climate scenario.Baseline values represent current climate normal, whereas IPSL values are predicted by the future climate model.Variable descriptions are presented in Table 2.

Table 6 .
Comparison of area (km 2 ) classified as suitable Phragmites habitat by the habitat suitability model for the current (2008-2010) baseline conditions and the future predicted IPCC A2a climate scenario for 2050.Total area is provided for reference.