A GIS-Based Approach for Determining Potential Runoff Coefficient and Runoff Depth for the Indian River Lagoon, Florida, USA

Abstract

The Indian River Lagoon system (IRL), spanning ~40% of Florida’s east coast, is one of the nation’s biggest and most biodiverse estuaries. In 2011, a super algal bloom event occurred in the IRL with total nitrogen and phosphorus levels that exceeded historical levels. Scientists suspect that nonpoint source pollution through surface runoff may have had a significant impact on the recent recurring algal blooms. Digital Elevation Model, land cover/land use, and soil data were used to calculate a runoff coefficient for the IRL drainage basin. Rainfall data were used to calculate runoff depth for the study area between the years of 2006–2016. When the monthly runoff depth data for 2011 were compared to a previous study on the 2011 super algal bloom in the lagoon, areas with high runoff visually matched the areas with higher chlorophyll a concentrations. Land development was a significant variable for determining runoff depth (p < 0.0001), and although used to derive runoff depths, the influence of precipitation was marginally significant (p = 0.06). Significant spatial autocorrelation indicated local trends between land development and runoff depth (p < 0.0001). Outputs will aid with decisions on stormwater management to more sustainable land development planning.

1. Introduction

Algae blooms within coastal estuarine systems have been a threat to vital key ecosystem components causing the degradation of ecological integrity. With non-point source pollution being a primary concern, using geographic information system (GIS) approaches to assess the impacts is effective for stormwater management. Therefore, with the use of land use/land cover (LC/LU), soil, and elevation data, the Potential Runoff Coefficient (PRC) and runoff depth were calculated for the IRL and Halifax River watershed. The analysis consisted of manipulating the geospatial data to derive the potential runoff coefficients and runoff depths.

Considering the contributing factors of surface runoff, the overall goal of the study is to estimate the quantities of runoff within the Indian River Lagoon (IRL) watershed based on a method that encompasses those parameters. The findings can also address whether such method and similar approaches can indicate locations of algae blooms, and aid in stormwater/watershed management. The objectives of this study are listed respectively; Objective 1: to calculate the potential runoff coefficients within the IRL watershed. The values will be based on the satellite image classification and validation for land cover/land use, elevation data, and soil data of the study area. Objective 2: to calculate the runoff depth of the IRL watershed over an eleven-year duration (2006–2011) using the derived value of the runoff coefficients and rainfall data provided by National Oceanic and Atmospheric Administration National Weather Service (NOAA NWS) River Forecast Centers (RFCs) collected from the Hydrologic Rainfall Analysis Project (HRAP). The outcome will represent the actual quantity of rainfall that was converted to runoff for the year. Objective 3: to visually assess if there is a geographic correlation of surface runoff and algae concentrations during months of the 2011 super algal bloom. The finished products can aid in gaining coastal resilience to help adapt to storms, flooding events, and parameters can be used to determine suitability for stormwater parks and infrastructure. The data acquired from the public GIS databases include ground-truthed information and remotely sensed data which were carefully interpreted and validated by professionals.

1.1 The Indian River Lagoon system

The Indian River Lagoon (IRL), spanning ~40% of Florida’s east coast, is one of the nation’s biggest and most biodiverse estuaries. The IRL consists of barrier islands separating its water from Atlantic Ocean [1]. The exchange of the IRL water with the ocean occurs naturally at Ponce De Leon Inlet in New Smyrna Beach, and Jupiter Inlet near West Palm Beach. The other man-made inlets include Sebastian Inlet, Fort Pierce Inlet, Port Canaveral, and St. Lucie inlet. The estuary stretches 251 km along the east coast of Florida with numerous tributaries [2]. The IRL system is made up of three sub lagoons that include the Mosquito Lagoon, which is in the northern section, the Banana River, and the IRL (Figure 1). The natural sources of freshwater for the IRL include Crane Creek (Melbourne, FL), Eau Gallie River, St. Lucie River, St. Sebastian River, and Turkey Creek. A secondary natural source of freshwater in the IRL is the Tomoka River which is located west of the lagoon running north connecting to the Halifax River then eventually the Mosquito Lagoon. Although the Tomoka River is not directly connected to the IRL or in its watershed, the Halifax River (Figure 1) is partially connected to the northern lagoon at Ponce Inlet and therefore its watershed is included in this study. The IRL and Halifax River watershed contains ~40 cities. The developed urban land comprises impervious surfaces and residential communities that primarily contain turf grass.

Figure 1.

A map of Indian River Lagoon and Halifax River, Florida. The Indian River Lagoon is composed of three waterbodies: the Mosquito Lagoon, Indian River, and the Banana Lagoon. The inset map provides a reference for the location of the lagoon in Florida.

In the summer of 2011, a super algal bloom event occurred in the IRL which reached a high biovolume of dinoflagellate Pyrodinium bahamense var. bahamense (33.9 × 106 μm³ mL⁻¹) with mean chlorophyll a concentrations (6.2–16.4 μg/L) that positively correlated with total nitrogen and total phosphorus levels that exceeded historical levels in various locations [3]. Following the massive algae bloom, there have been recurrent blooms consisting of green macroalgae such as Chaetomorpha sp. since 2013 [4, 5]. As a result of the 2011 super algal bloom, the coverage of seagrass within the IRL drastically declined from the loss of photosynthetic light by the surface algae [6]. Although fluctuations in seagrass bed percent cover in the lagoon have been understood as a part of a natural cycle of decline and recovery as seagrass abundance, scientists suspect that nonpoint source pollution via surface runoff may have had a significant impact on the recent recurring algal blooms in the lagoon [7, 8].

1.2 Surface runoff and runoff models

Surface runoff is water from rain or snowmelt that travels over the land before entering nearby waterbodies. Stormwater flowing across surrounding land transports various pollutants, and ultimately contributes to non-point source pollution. Surface runoff negatively affects many aquatic ecosystems as the runoff transports pollutants and other substances into waterbodies, which can alter turbidity, phosphorus and nitrogen concentrations, and organic matter content in receiving waterbodies [9]. The effects of surface runoff can also be intensified by climate change in specific regions that may have highly developed land and altered hydrology from the addition of artificial stormwater structures that modify the flow of water [10]. Human activities have been shown to have a stronger impact on runoff than climate change, but both stressors significantly impacts runoff quantities [10].

Hypothetical land cover change scenarios in a simulated hydrological study within the Lavrinha watershed in Minas Gerais State, Brazil showed that deforestation in the Atlantic Forest biome would lead to increases in soil moisture (5%), runoff (22%), and decreases in runoff interception (71%) from the loss of roots and extensive rhizomes [11]. Impervious surfaces in urban watersheds can influence the biogeochemical processes, organism’s abundance, stress, and vulnerability from heated surface runoff during hot summers [12]. Incorporation of the contributing factors such as vegetative cover that may enhance or influence the effects is effective in hydrological modeling for determining the amount of runoff.

The characteristics of the land surrounding waterbodies affect the amount of surface runoff. During the process of rainfall becoming runoff, various characteristics of the land’s surface, such as land use, soil type, and topography, will heavily impact the quantity of runoff [13]. Vegetative cover of the surrounding land can potentially act as a buffer for aquatic systems receiving runoff [14]. During rainfall events, impervious surfaces such as roads, parking lots, and other pavements increase runoff due to extremely limited infiltration into the ground. Areas with 75–100% imperviousness can yield runoff that represents up to 55% of any rainfall [15].

The potential runoff coefficient (PRC) represents the portion of rain that becomes surface runoff during a rain event, and it is determined by the land use, soil texture, and slope [16]. The potential runoff coefficient was derived from methods of developing a unit hydrograph (UH) for specific depths of rainfall. The hydrograph provided the assumption that discharge at any time is proportional to the volume runoff, and the temporal factors for a given duration are constant [17]. Runoff coefficients have been widely utilized in the hydrological modeling along with other computational factors for research in flood frequency, flood prediction, and storm management [18, 19, 20]. Hydrologic simulation model software’s have been developed using spatial data and GIS [5]. Another method that includes runoff coefficients to hydrological modeling is the runoff curve number (CN) method created by the United States Department of Agriculture Natural Resources Conservation Service. Unlike the potential runoff coefficient, the CN can be calculated for each watershed and encompass the potential maximum retention of the soil over a given period of time.

PRC can be determined for land surfaces with different characteristics along with the quantity of runoff known as the “runoff depth.” Given the quantity of runoff being influenced by the determining conditions, spatial variation in PRC can be estimated for a specific time duration within a given estuarine drainage area. In order to demonstrate this, runoff coefficients and runoff depths were calculated using geographic information systems (GIS) for the drainage basin of the Indian River Lagoon (IRL), Florida, which recently had recurring severe algal blooms. Nonpoint source pollution from surface runoff may have had been a cause for the recurring algal blooms in the lagoon [7, 8]. Use of the spatially contiguous PRC across an area of interest provides additional resources and information for stormwater research within a coastal watershed. Runoff coefficients of a watershed along with other information can be utilized for analytical processes to gain further insight of stormwater dynamics on local and regional scales.

Since 2011, IRL experienced severe algal bloom events; and non-point source pollution through surface runoff is suspected to be one of the causes for the algal blooms. The goal of this study is to calculate the spatially contiguous PRC and runoff depth for the drainage basin of IRL and the connected estuary, the Halifax River, Florida for an eleven-year period (2006–2016) in order to determine which areas and factors contribute to the runoff. The 2011 monthly runoff depth of the draining areas was compared with the 2011 monthly algal bloom maps of a previous study in order to see any visible correspondence between the locations of algal bloom initiation and the locations with high runoff depth values.

2. Data for model components

The procedure to derive the PRCs and runoff depths for the IRL consisted of processing satellite imagery to derive the land cover and land use (LC/LU), collecting the soil textures throughout study area, and calculating slope using terrain elevation data within the watershed.

2.1 Land cover/land use

Land cover and land use (LC/LU) is one of the factors for calculating PRC. The LC/LU was derived by classifying the European Space Agency (ESA) Sentinel 2 Level 1C 10 m satellite imagery from November of 2016. Four images were downloaded from the ESA Sentinel Scientific Data Hub website to encompass the elongated watershed of the IRL (https://scihub.copernicus.eu/dhus/#/home). The images were preprocessed with the ESA Sentinel Application Platform (SNAP) remote sensing software along with the Sentinel 2 toolbox. Before classifying LC/LU of the study area, an atmospheric correction was applied to the images using the Sen2cor 2.3.2 plugin within ESA SNAP to eliminate the effects of water vapor, aerosols, and cirrus clouds when utilizing spectral reflectance data. The preprocessing output of the Sentinel 2 data produces Sentinel Level-2A data which includes values that represent the radiation at the bottom of the atmosphere (BOA). Once the BOA output was produced, the four images were mosaicked to produce a continuous raster image of the IRL. By applying a supervised maximum likelihood classification in ENVI 5.4, the images were classified into five categories; forest, grass, bare soil, crop, and impervious.

2.2 Slope

A digital elevation model (DEM) from the United States Geological Survey National Elevation Dataset (USGS NED) was used to generate terrain slope at a spatial resolution of 10 m (http://viewer.nationalmap.gov/basic/?howTo=true) within ESRI ArcMap 10.5 (380 New York Street, Redlands, CA 92373-8100). The elevation values were collected with Interferometric Synthetic Aperture Radar, and referenced to the North American Vertical Datum of 1988 (NAVD 88). A preliminary analysis of the DEM was performed to fill in the low areas or “sinks” that are considered to be errors so that modeled runoff would flow smoothly across the land’s surface. The filled output map was used to create the slope for the areas surrounding IRL. The percent slope was classified into three classes due to the low elevation throughout Florida.

2.3 Soil

The soil data used in the analysis were obtained from the Web Soil Survey (WSS) (http://websoilsurvey.sc.egov.usda.gov/App/WebSoilSurvey.aspx). The WSS is operated by the United States Department of Agriculture Natural Resources Conservation Service (NRCS) and contains geospatial data and information produced by the National Cooperative Soil Survey. The NRCS soil data are produced from soil samples collected from NRCS State Soil Scientist for counties throughout the United States and are available in tabular and geospatial data. The spatial data are provided in the Geographic Coordinate System and World Geodetic System of 1984 datum (GCS_WGS_84). The data for soil classification were acquired for six Florida coastal counties: Volusia, Brevard, Indian River, St. Lucie, Martin, and Palm Beach. Tabular information for the soil texture was extracted from the Web Soil Survey Microsoft Access Database file and imported into the ArcMap 10.5 software. The data contained a variety of different soil names for classification: muck, Myakka fine sand, and Turnbull muck, that are all used for determining the slope constant along with LC/LU.

2.4 Precipitation for runoff depth

The runoff depth represents the amount of rainfall that is converted into runoff [16]. Therefore, rainfall data for the Halifax River and IRL watershed were collected to calculate the runoff depth using the runoff coefficients. The data were acquired from the National Oceanic and Atmospheric Administration (NOAA) National Weather Service (NWS) River Forecast Center (RFC) website. The data were downloaded in the ArcGIS shapefile format as point data with a projection of the Hydrologic Rainfall Analysis Projection (HRAP) grid coordinate system that has a North Pole Stereographic projection, and a grid resolution of 4762.5 m (https://water.weather.gov/precip/download.php). The rainfall data are acquired in a multi-sensor process that uses radar and rain gauge to estimate the precipitation. After extracting the point data, the shapefiles were converted into raster data. Due to estimation of multi-sensor collected data, the data are first stored in a binary file format called XMRG. This file is then read into the HRAP grid coordinate system through the NWSRFS Operational Forecast System using the NEXRAD Mean Areal Precipitation Preprocessor (MAPX) to associate grid points from XMRG data to represent the hourly average precipitation for each area [21].

3. Data analysis

3.1 LC/LU accuracy assessment

Before assessing the PRCs for the study area, the LC/LU image was tested for its accuracy. The accuracy assessment test consisted of collecting 600 referenced points using a stratified random method that randomly assigns points in each class. A 2016 Digital Globe basemap in ArcMap 10.5 of a higher spatial resolution (0.62 m) was utilized to visually interpret the land cover for each reference point. The output table consisted of a confusion matrix that displays the error of omission, the error of commission per class, and overall accuracy ranging from 0 to 1. Another test for accuracy of the LC/LU classification image included calculating the Cohen’s kappa (K) coefficient [22].

3.2 PRC

Potential runoff coefficient (PRC) values were derived to represent ratio of the rainfall that would convert to surface runoff per pixel. The PRC for the IRL area was determined by combining the soil texture, LC/LU, and the slope data. The PRC is calculated from a linear relationship between the runoff coefficients and slope, which is shown in Eq. (1) [23].

$$\text{𝐶}=C_0+\left(1-C_0\right)\frac{S}{s+s_o}$$E1

𝐶 is the PRC, S (%) is the slope of the land surface, $C_0$ is the PRC for the near zero slope in reference to the first row of every land use class in PRC values for different land use, slope, and soil texture published in [23] which is sourced from published material [24, 25, 26, 27, 28]. $s_o$represents the slope constant for different land use and soil textures that were empirically derived over a collection of studies. Following reclassification, classes for both parameters were assigned arbitrary weighted values and the soil and LC/LU values for the data were multiplied in the ArcMap 10.5 “Raster Calculator” tool. The arbitrary values were assigned to the classes to conveniently identify each combination of LC/LU and soil texture per pixel from the products. The products of the combinations helped derive the $C_0$and $s_0$ per pixel in the image. The products for the variables were also used to calculate the PRC (Eq. (1)) using the Raster Calculator Tool.

3.3 Runoff depth

The total precipitation values were collected for eleven years (2006–2016), and imported into ArcMap 10.5 to be interpolated. The precipitation values (in.) for each of the years were interpolated using the Kriging method with a spherical semivariogram model. The method assumes that the values are more related when in close proximity, and the spatial autocorrelation decreases with distance. After the precipitation was interpolated for each year, the data were multiplied (cell-by-cell) by the PRC raster of the corresponding year using the raster calculator tool provided in ArcMap 10.5 toolboxes. The output of the images provided the runoff depth (centimeters) for each year, and the average runoff depth for the eleven-year period (2006–2016) was calculated per pixel (10 m). The outputs of this image can delineate potential sources of runoff for inland waterbodies that may be connected to the lagoon through a network of drainage systems.

Concentrations of chlorophyll a in the IRL during the 2011 super algal bloom were compared to runoff depth of surrounding areas. Kamerosky et al. [29] estimated and mapped the Chi a concentrations using the Medium Resolution Imaging Spectrometer (MERIS) platform aboard the European Space Agency (ESA) Environmental Satellite (ENVISAT) and calculated Normalized Difference Chlorophyll Index (NDCI) [29, 30].

3.4 Linear regression between LDI and runoff depth

In order to meet proper data conditions for linear regression analysis in ArcMap 10.5, the raster images were sampled into vector data as a point feature class. Land development intensity (LDI) data was collected from the Florida Department of Environmental Protection (FDEP) Geospatial Open Data Site (http://geodata.dep.state.fl.us/). The LDI serves as a human disturbance gradient that incorporates land use and energy used per unit area [31]. It is used in watershed modeling to delineate human-dominated areas, and to scale the human induced impacts on physiological, biological, and chemical processes. A total of 600 points were randomly placed within the Halifax River and IRL watershed via “Create Random Points” tool. Points that were placed over large waterbodies were deleted, leaving 528 sample points left for the analysis. Values from the LDI and runoff depth were extracted to the points.

To adequately assess the relationship between urbanized areas of intense impervious coverage and surface runoff, an ordinary least squares (OLS) regression analysis and geographically weighted regression (GWR) was performed in ArcGIS 10.5. These regression analyses use bandwidth methods to find the optimal sampling distances between data points, adding a geospatial component to regression analysis. The OLS regression is designed as a “Global Model” with an assumption that the explanatory and dependent variables have global trends over a particular study area. In simplified context, it is assumed that the data are continuous throughout the area therefore being “stationary” data.

For the OLS analysis, the Jarque-Bera statistic tests for model bias that can arise from nonstationary data, misspecification of independent variables, and skewed residuals [32]. Due to the positively skewed data for LDI and runoff depth values, a logarithmic transformation was applied to data to ensure a normal distribution of the datasets while making the variance independent of the mean. The Koenker’s Studentized Breusch-Pagan (Koenker BP) statistic tests for nonstationary with a null hypothesis that the dependent and independent variables have a consistent relationship in geographic space, thus being stationary [33]. A rejected null hypothesis of this test indicates that there are local trends between the variables within the study area.

Presence of significant spatial autocorrelation using the Global Moran’s Index (Moran’s I) is based on the assumption of stationary data. In this case there will be clustering of standard residuals from heteroscedasticity, thus indicating a local model such as GWR is more appropriate. Therefore, the standard residuals produced from both regression analyses were tested for significant clustering using the Moran’s I test. On the other hand, a GWR is a “nonstationary” model that accounts for the local trends in relationships between the variables. In OLS analysis, LDI data from the FDEP and 11-year mean precipitation were used as the independent variable, and the 11-year mean runoff depth as dependent variable. The GWR only included the LDI as independent, and runoff depth as dependent variable due to collinear relationships with rainfall within clustered locations within the study area.

4. Results

4.1 LC/LU classification

There are six LC/LU classes delineated from the supervised classification. The land that mostly consists of agriculture occurs in the southern section of the watershed. The overall accuracy of the LC/LU classification image was 0.82, and the lowest accuracies were in the impervious (User accuracy of 0.65) and bare soil (0.53) classes. This may have been due to the spectral similarity between bare sand along the coast and impervious surfaces such as roof tops. The reference points that appeared to exist in unhealthy brown vegetation were misclassified as bare soil. The kappa coefficient for the LC/LU image was 0.77, with an overall average of 0.82.

4.2 Slope

The slope for the study area was assigned a quantile classification to exclude the effects from outliers in the digital elevation map. The elevation in the state of Florida is relatively flat with an average slope of ~0.47 m per pixel. The areas of high percent slope are manmade structure such as buildings, walls, or homes in developed areas. Some manmade structures with unusually high slopes were identified as the outliers. The other cities have slopes ranging from 0.50 to 1.79 average percent.

4.3 Soil

The soil texture classification for central east Florida consists of mostly fine with Myakka Fine Sand as a native soil, covering more than 1.5 million acres of land, and labeled as the Florida Official State Soil [34]. The soil data were reclassified into the 12 different textures within the USDA Soil Texture Triangle to accurately implement the values: sand, loamy sand, sandy loam, silt loam, sandy clay loam, silty clay loam, sandy clay, silty clay, and clay. The total study area was composed of 67.7% sand, 4.5% loamy sand, and 8.7% silty clay.

4.4 PRC

The PRC values range from 3 to 100% (Figure 2). The PRCs are higher in runoff values in developed areas that are in close proximity to the coastal waterbodies of the IRL and Halifax River. The spatial resolution (10 m) of the image shows a detailed delineation of the manmade infrastructure within urban coastal communities such as roads, buildings, homes, and airports.

Figure 2.

Map displaying the potential runoff coefficients (PRC; %) for the Indian River Lagoon and Halifax River watersheds, FL.

4.5 Precipitation

The precipitation data in the IRL watershed from 2006 to 2016 were divided into four quarters with each quarter representing the average of three-month intervals: January–March, April–June, July–September, and October–December. Although the quarterly intervals do not accurately align with seasons, the data are segmented to show the temporal shifts of the rainfall pattern in this area. The IRL watershed precipitation is usually the lowest within the first quarter averaging ~5.48 cm. As the seasonal rainfall increases in spring and summer moving from 10.03 cm in second quarter to 15.55 cm in the third quarter. The rainfall decreases towards the end of the year with an average of 5.58 cm.

4.6 Runoff depth

The runoff depth values for the IRL ranges from 2.51 to 141.48 cm. The monthly runoff depth was calculated for 2011, the year of the super algal bloom in the IRL, to serve as potential explanation for the contribution of high surface runoff to locations of the algal blooms (Figures 3 and 4). The average runoff per sub-basin (Figures 3 and 4) was compared to the chlorophyll a concentrations quantified from European Space Agency’s Medium Resolution Imaging Spectrometer (MERIS) for 2011 [29] (Figures 5 and 6). The maps were all assigned the same symbology to aid easier depictions of changes in quantities, and for comparison between months.

Figure 3.

The monthly mean runoff depth per sub-basin in the Indian River Lagoon from January 2011 to June 2011. The values increase form dark green, to warmer colors.

Figure 4.

The monthly mean runoff depth per sub-basins in the Indian River Lagoon from July 2011 to December 2011. The values increase form dark green, to warmer colors.

Figure 5.

Indian River Lagoon chlorophyll a concentrations for spring 2011. The concentrations are estimated using medium resolution imaging spectrometer normalized difference chlorophyll index (image source: [29]).

Figure 6.

Indian River Lagoon chlorophyll a concentrations for 2011 (September–December 2011). The concentrations are estimated using medium resolution imaging spectrometer normalized difference chlorophyll index (image source: [29]).

4.7 Ordinary least squares regression

For the OLS model, LDI is statistically significant (p < 0.0001) for and robust probability. Precipitation determines the runoff depth, however, appears to be also significant at the 5%, but not as significant according to the robust p value (p = 0.05, robust p = 0.06). The variance inflation factor (VIF) tests for the redundancy amongst the explanatory variables that are added to the model. If two or more explanatory variables tell the same story because they are linearly related, the error variances are inflated, and the resulting multicollinearity produces a higher VIF. Studies suggest that accepting VIFs fewer than 7.5 or 10 is the rule of thumb for determining if there is multicollinearity within a dataset [35].

Within the OLS diagnostic results, statistical values provide information describing the performance of the model along with indicators for choosing an alternative model to adequately address the overall question (Table 1). The Akaike’s information criterion (AIC) measures the overall model performance, which can be used in comparison to other regression analyses [36]. The multiple R² explains how much the independent variables explain the variation in the dependent variable. In relation to the multiple R² the Adjusted R² accounts for the model complexity. The multiple R² the Adjusted R² for this tests shows a small R² between the variables (R² = 0.15). The OLS regression also tests for the model significance with the Joint F-statistic and Joint Wald statistic to support the significance of R² values (Table 2).

Variable	Coefficient	Standard error	t-statistic	Probability	Robust SE	Robust probability	VIF
Intercept	1.032	0.531	1.943	0.053	0.571	0.071
Precipitation	0.009	0.005	1.967	0.050	0.005	0.060	1.001
LDI	0.692	0.071	9.780	<0.0001	0.085	<0.0001	1.001

Table 1.

A table of the ordinary least squares model variables.

Number of observations:	564	AICc	1573.001
Multiple R-squared	0.151	Adjusted R-squared	0.147
Joint Wald statistic	67.250	Prob(>chi-squared), (2) df:	<0.0001
Koenker (BP) statistic	10.963	Prob(>chi-squared), (2) df:	0.004
Jarque-Bera statistic	5.500	Prob(>chi-squared), (2) df:	0.064

Table 2.

Statistical diagnostic results from the ordinary least squares regression.

The Koenker’s BP statistic tests for nonstationary and heteroscedasticity. The null hypothesis is that the dependent and independent variables have a consistent relationship in geographic space, thus being stationary [33]. The Koenker’s BP statistic shows significant existence of nonstationary trends between runoff depth and LDI (p = 0.004). Therefore, the model significance was interpreted based on the Joint Wald statistic (p < 0.0001) which also indicates that the relationship was statistically significant. However, the overall measure of how well the explanatory variables explained the variation in the runoff depth from the OLS analysis was relatively small (R² = 0.15). The Jarque-Bera statistic tests for model bias that can arise form nonstationary data, misspecification of independent variables, and skewed residuals [30]. In this case, the Jarque-Bera statistic shows no significant model bias (p = 0.064). A Global Moran’s Index was performed on the residuals of the output file to test for the assumption of no spatial autocorrelation or clustering in the data. The Global Moran’s Index showed statistically significant clustering rejecting the null hypothesis that the data are randomly distributed spatially within a global assumption (Moran’s I = 0.07, p < 0.0001). Therefore, the OLS results should not be used to adequately interpret relationship between the explanatory variables and runoff depth.

4.8 Geographically weighted regression

Due to the detection of nonstationary and/or heteroscedasticity in the datasets, a GWR was used to adequately assess the relationship. The Global Moran’s Index and Anselin Local Moran’s Index was performed to test for clustering and local patterns of spatial association [37]. The Global Moran’s I indicated that the standard residuals produced from the GWR were significantly dispersed indicating the absence of spatial autocorrelation (Moran’s I = −0.025, p = 0.111). The algorithm calculated an index for every feature, and 96% of the local p-values were not significant (p > 0.05). The validation for choosing the GWR also can be justified by the smaller AICc produced (AICc = 1522.83) compared to the OLS (1573). The R² for the GWR increased (R² = 0.35) with a lower adjusted R² of 0.26. As previously stated, the GWR accounts for nonstationary data that contain local trends for the relationship between the variables. Local trends within the dataset relationships were inevitable due to the complexity of different LC/LU within urban communities. The locally weighted regression coefficients can be seen on the coefficient raster produced by the GWR analysis (Figure 8). The coefficients show that LDI influences runoff in locations with more impervious surfaces and higher runoff depths, and forested land cover that consists of low LDI values and low runoff. The coefficients increase from blue to yellow to red, indicating higher relationships between the two variables.

Table 3 shows statistical values generated from the model according to the optimal sampling distance for nearest neighbors (bandwidth). The sampling kernel type for the GWR was fixed, and therefore provides the bandwidth in meters. The Residual Squares is the sum of squared residuals that represent the distance between the observed and estimated values. Therefore, the data are more related when this value is smaller. With a strong influence from bandwidth, the Effective Number is a measure of the complexity of the model that is used to calculate other variables within the GWR model, and it is useful when compared to other models. The sigma is the estimated standard deviation for the residual sum of squares, which shows that the standard deviation of the observed values for runoff depth were relatively close to the predicted values calculated for the regression model (σ = 0.90).

Variable name	Values
Bandwidth	11228.83 m
Residual squares	406.15
Effective number	66.62
Sigma	0.90
AICc	1522.83
R-squared	0.35
Adjusted R-squared	0.26
Dependent field	0.00
Explanatory field	1.00

Table 3.

Results from the geographic weighted regression analysis.

5. Discussion

5.1 Spatially continuous PRC of the IRL and Halifax River watersheds

The goal of this research is to calculate spatially continuous potential runoff coefficient (PRC) and runoff depth. In order to demonstrate how spatial and temporal variation in PRC can be estimated within a given estuarine drainage area, this study calculated PRC as a proportion of rainfall becoming surface runoff and calculated the runoff depth that is the amount of rainfall converted to runoff.

The average PRCs increase from forest, grass, agriculture, bare soil, and to impervious. Ideally, forested areas would have the highest interception of precipitation because of high percent cover of vegetation. Forested areas also have an increase in absorption from the abundance of extensive rhizome systems in the substrate. Areas of grass may have high percent cover of vegetation, but the interception of storm water may not as efficient due to the small biomass of plants. The classification image shows that forest (27.4%) and grass (24.7%) are the most dominant land covers within the IRL watershed. Forests mostly cover the northern section of the IRL watershed and the Halifax River watershed westward of coastal cities. The higher PRC values are located along the Halifax River and IRL in more developed urban communities such as Daytona Beach, Melbourne, and Palm Bay, Florida. Throughout the state of Florida, St. Augustine grass (Stenotaphrum secundatum [Walt.] Kuntze) is a popular turf grass used for urban lawns. This rhizome structure of this grass is dense, but relatively short in length which increases yields in runoff and shoreline erosion. Regardless of specific lawn grass, runoff coefficients are higher than forest cover. The recommended runoff coefficient value table for Georgia Stormwater Management also shows different values for grass covered lands based on the soil texture and slope [38]. However, there is only one value for forested areas despite the slope and texture. Although different from forests PRCs used in this study, a change in land cover can impact runoff yields particularly in areas of dense vegetation.

Impervious surfaces make up 15.3% of the study area much of which is located along the coastlines. Based on the National Atlas of the United States Spatial Data collected from the Florida Geographic Data Library (FGDL), there is a total of 40 cities within the IRL watershed. For this study, the ten coastal communities with the highest cover of impervious surface were included: Palm Bay, Port St. Lucie, Melbourne, West Melbourne, Daytona Beach, Port Orange, Ormond Beach, New Smyrna Beach, Titusville, and Fort Pierce. The cities of the most impervious surfaces are Palm Bay with ~50,554 acres of impervious surfaces, and Port St. Lucie with ~73,959 acres of impervious cover.

5.2 Temporal variation in runoff depth

The runoff depth varies with changes in LC/LU and intensity of precipitation. The estimated average runoff depth for the IRL ranges from 2.5–141.5 cm for the 11-year interval. The runoff depth throughout the study area fluctuates among the years (Figure 7), due to the changes in precipitation. With PRC values, the areas with potential nonpoint source pollution can be used as target locations for management or mitigation. Runoff depth values above the 11-year mean varied across the area and amongst the years. Runoff deviation from the mean indicated that heavy runoff depth values above the 11-year mean are years of 2008, 2014, 2015, and 2016. Causes for runoff differences can be contributed to fluctuations in climatic and annual weather patterns for rainfalls. Climatic and temporal trends have been related to changes in the IRL water quality such as El Niño years linked to declines in salinity levels in 1997, and the extended period of La Niña drought events that persisted in autumn of 2006 and summer 2007 [3]. The precipitation data for the IRL were low quantities in 2006 and 2007, and a gradual increase to 2016. The runoff depth appeared to be above the 11-year mean (35.89 cm) throughout the watershed for the years of 2014 and 2016. These are the years of strong El Niño events during which recurrent algal blooms occurred in IRL. La Niña events in Florida have shown nitrate levels to higher in ground water than in streams, which results that nutrients in aquifers accumulates from fertilizer, septic tank effluent, and animal wastes [39, 40]. The average runoff depth for 2006 was the lowest of the years, with 2016 being the highest.

Figure 7.

The mean of the total runoff depth for each year with standard error bars for standard deviation. The dotted line represents the average 11-year runoff depth (μ = 35.89 cm) (created in Microsoft Excel 2010).

Nutrient loads vary with different land use. For example, golf courses are suspected to be a major contribution to nutrient loading in waterways aside from agricultural lands. Recorded nitrate and phosphorus concentrations significantly increased at the outflow locations from the inflow concentrations for the Morris Williams Municipal Golf Course in Austin, TX [41]. Above average runoff depths in such locations can be monitored as an indicator for early warnings of algal blooms.

5.3 Linear regression between runoff depth, precipitation, and land development

Developed land within the IRL watershed contains impervious surfaces consisting of roads, parking lots, and also vegetated lots that are highly altered by human development. Precipitation undoubtedly contributes to runoff quantities, but LC/LU, and development can influence runoff yields. The regression analyses were used to test the relationship between runoff and development, as well as between runoff and precipitation. The OLS regression (Tables 1 and 2), the test to analyze if precipitation is an important factor for determining areas and timing with high runoff contribution, could not be adequately assessed due to spatial autocorrelation (Moran’s I = 0.07, p < 0.0001). However, the results appear to be marginally significant at the 95% confidence interval (p = 0.05; robust p = 0.06). The LDI showed to be a significant (p < 0.0001) variable at explaining a significant amount of variation in runoff depth. Presence of significant spatial autocorrelation using the Global Moran’s I is based on the assumption of stationary data. In this case there will be clustering of standard residuals from heteroscedasticity, thus indicating a local model such as GWR is more appropriate. The Global Moran’s Index indicated no spatial autocorrelation with a negative index and rejecting the null hypothesis at with 95% confidence (Global Moran’s I = −0.025, p = 0.111).

Although the runoff depth was determined by precipitation, LC/LU can have a higher impact on the quantities of runoff. Empty grass lots within urban communities can have compacted soil from earlier construction activity which may decrease infiltration rates up to 70% in the central Florida region [42]. The runoff coefficients for the agricultural land surfaces include the effects of compacted soil from heavy machinery. Based on the coefficient raster generated from the GWR analysis, LDI values for forested and impervious areas may account for most of the linear relationship between development and runoff (Figure 8). The local trends between rainfall and runoff on smaller time intervals may have a strong linear relationship. However, the mean rainfall values may have reduced the weights in local rainfall-runoff relationships. This outcome also can be noticed within the 11-year mean runoff OLS regression from the existence of local relationships between runoff and the independent variables.

Figure 8.

The map above shows the geographically weighted regression locally weighted coefficients throughout the study area.

Precipitation estimates along the 251-kilometer IRL estuary and ~ 35-kilometer Halifax River varies within locations, with changing LC/LU as a factor. Areas with lower rainfall can have a higher runoff yield than areas with higher precipitation over forested areas that were assigned lower runoff coefficients. As a result, areas consisting of more human disturbance have a linear relationship with more runoff. Urban communities are often primary targets for some studies to analyze rainfall-runoff by enhancing methods to estimate the DCIA in developed catchments [43]. Based on the global trends of urbanization within coastal areas, stronger rainfall-runoff relationships have positive correlations with the increase of impervious cover percentages for urban zones within separate countries [44]. The purpose of choosing LDI was to indicate the contributions of surface runoff from vegetated lands affected by urban development. To further explain this relationship, future research should assess stormwater runoff using the impervious percentage images created by the USGS. The percentage of imperviousness can also be compared to increases in runoff depth.

5.4 Runoff depth during the 2011 super algal bloom

The IRL ecosystem recently suffered from a recurrence of algae blooms since 2011 which are heavily influenced by anthropogenic stressors within its watershed, such as surrounding developed land with possible higher surface runoff [8]. Based on a visual comparison; the runoff depth was higher prior to the algal bloom events between 2011 and 2016 particularly near the areas of recorded high Chlorophyll a concentrations (Figures 3–6). It is important to note that the monthly runoff reflects precipitation estimates collected at the end of the month. Therefore, the runoff depth map for March should be visually compared to the chlorophyll a concentrations in April 2011.

Although there is no available MERIS NDCI calculations collected throughout the summer, there was an increase of runoff to 10–18 cm in May and June for the Banana River (Figures 5 and 6). The increase may explain the 48.62 μg/L spike in chlorophyll a May 2011 from April 2011. Based on the estimated concentrations by MERIS NDCI and water quality samples from SJRWMD, the algal bloom became higher on the 14 September 2011 with the highest concentrations in October [29]. The runoff depth for October 2011 showed the high runoff with values above 15 cm in the southern IRL and Northern IRL. Subsequently, concentrations of chlorophyll a gradually decreased throughout October, and further dropped in November and December. Results also indicate that there were also smaller contributions of runoff during those months with decreasing trends. Visual comparisons between the chlorophyll a and runoff depth indicate that there may be a correlation positive association between the two variables for 2011. However, further analysis including statistical measures should be performed to assess the relationship.

5.5 Implication for coastal water management

Delineation of PRCs and runoff depths can provide a geographic depiction for assessing lands of interest to implement sustainable developmental designs and structures. In this study, runoff coefficients were calculated for each pixel regardless of surrounding pixel values. Therefore, computational methods used in this study to determine runoff depth were not assessed using methods to incorporate Directly Connected Impervious Areas (DCIA). DCIAs are areas that are considered to be hydraulically connected to the conveyance system according to the Southwest Florida Water Management District Resource Regulation Technical Guide [45]. Other study estimated runoff volumes in various sub-basins of rivers and tributaries within the IRL watershed using DCIA and non-DCIA methods [46]. Runoff from such studies use a measure called the Soil Conservation Survey Runoff Curve Number. The overall runoff was calculated from the sum of the DCIA and non-DCIA runoff, while runoff coefficients were derived by dividing the generated runoff by the total rainfall for the stations which was listed as “C values”. While this approach can be used to determine Total Maximum Daily Loads (TMDLs) for nutrients, PRCs from this study can be used to emphasize exact locations within the watershed that are suitable for LC/LU management practices.

Direct surface runoff into waterbodies can be significantly affected by impervious surfaces in close proximity. In other scenarios, runoff from developed lands may travel through vegetation before entering into a waterbody. The harmful effects of surface runoff from and on urban communities call for a need of more stringent regulations, and more efficient coastal urban planning and management. This approach of stormwater management provides a long term adaptation plan to be proactive to the future impacts from climate change. Delineation of potential runoff coefficients and runoff depths can provide a geographic depiction for assessing lands of interest to implement sustainable developmental designs and structures. Mean runoff depth and runoff coefficient values can be used to determine areas of high runoff to apply green infrastructure within a watershed. “Green Infrastructure” is the practice of utilizing natural vegetated areas for runoff treatment by mimicking natural stormwater flow paths, and is composed of many low impact development (LID) designs [47]. Developed land within the IRL watershed contains impervious surfaces consisting of roads, parking lots, and also vegetated lots that are highly altered by human development.

Precipitation undoubtedly contributes to runoff quantities, but LC/LU, and development can influence runoff yields. The regression analyses were used to test the relationship between runoff and development, as well as between runoff and precipitation. The OLS regression, the test to analyze if precipitation is an important factor for determining areas and timing with high runoff contribution, could not be adequately assessed due to spatial autocorrelation (Moran’s I = 0.07, p < 0.0001). However, the results appear to be marginally significant at the 95% confidence interval (p = 0.05; robust p = 0.06). The LDI showed to be a significant (p < 0.0001) variable at explaining a significant amount of variation in runoff depth. Presence of significant spatial autocorrelation using the Global Moran’s I is based on the assumption of stationary data. In this case there will be clustering of standard residuals from heteroscedasticity, thus indicating a local model such as GWR is more appropriate. The Global Moran’s Index indicated no spatial autocorrelation with a negative index and rejecting the null hypothesis at with 95% confidence (Global Moran’s I = −0.025, p = 0.111).

Although the runoff depth was determined by precipitation, LC/LU can have a higher impact on the quantities of runoff. Empty grass lots within urban communities can have compacted soil from earlier construction activity which may decrease infiltration rates up to 70% in the central Florida region [42]. The runoff coefficients for the agricultural land surfaces include the effects of compacted soil from heavy machinery. Based on the coefficient raster generated from the GWR analysis, LDI values for forested and impervious areas may account for most of the linear relationship between development and runoff. The local trends between rainfall and runoff on smaller time intervals may have a strong linear relationship. However, the mean rainfall values may have reduced the weights in local rainfall-runoff relationships. This outcome also can be noticed within the 11-year mean runoff OLS regression from the existence of local relationships between runoff and the independent variables.

Mitigation of stormwater runoff often includes developing more sustainable development strategies within urban communities. An aggregation of developmental designs for combating runoff in an urban community showed reductions in runoff, and also projected expansions in bare soil, impervious cover, and soil alteration will lead to higher runoff volumes [48]. As with vegetated and undeveloped surfaces within riverine systems, the hydrological changes in volume and base flows are reduced with this hybrid design. In reference to GIS approaches to correlating surface runoff to LCLU, urban areas and bare land also corresponded to the degradation of stormwater quality [49].

6. Conclusion

The PRCs for the IRL were applied to land surfaces based on soil, land cover, and slope. These coefficients were used as ratios to determine the runoff depth per pixel within the IRL using precipitation data. After calculating the runoff depth for an 11-year period (2006–2016), it was found that the recent years (2014, 2016) were above the average 11-year runoff matched years of strong El Niño. The runoff deviation from the 11-year mean was also calculated per pixel for each year and highlighted higher runoff quantities closer to the shore of the IRL within the watershed. It is well known that impervious surfaces decrease infiltration, thus increasing runoff yields. Even with vegetated landscape, highly developed land can have poor infiltration from compact soil. The linear regression analysis showed that land development has a significant relationship with runoff depth, and there are local trends between the variables. During the 2011 super algal bloom, the months of March and April 2011 showed increases in runoff, which matched the areas with higher chlorophyll a mapped with MERIS in the Mosquito Lagoon in the Northern IRL [29]. In October 2011, extremely high concentrations were detected from MERIS and sampled from St. Johns River Water Management District; this research also calculated high runoff depth concentrations, delineated in the IRL watershed for October 2011. Based on these analyses, the output of this research can possibly delineate areas within the coastal communities that experience higher runoff, and help locate more suitable areas for stormwater parks, green infrastructure, and sustainable stormwater structures. Future research can include using the indices such as LDI to further correct the runoff coefficients for a particular watershed. PRCs can be applied to other watersheds of coastal ecosystems for as a visual reference, or used as a parameter for more advanced hydrologic modeling.

Acknowledgments

This publication was made possible by the National Oceanic and Atmospheric Administration, Office of Education Educational Partnership Program award (NA16SEC4810009). Its contents are solely the responsibility of the award recipient and do not necessarily represent the official views of the U.S. Department of Commerce, National Oceanic and Atmospheric Administration.

Any opinions, findings, conclusions, or recommendations expressed in this publication are those of the author(s) and do not necessarily reflect the view of the U.S. Department of Commerce, National Oceanic and Atmospheric Administration.

Conflict of interest

There is no conflict of interest.