Estimating River Sediment Discharge in the Upper Mississippi River Using Landsat Imagery

: With the decline of operational river gauges monitoring sediments, a viable means of quantifying sediment transport is needed. In this study, we address this issue by applying relationships between hydraulic geometry of river channels, water discharge, water-leaving surface reﬂectance (SR), and suspended sediment concentration (SSC) to quantify sediment discharge with the aid of space-based observations. We examined 5490 Landsat scenes to estimate water discharge, SSC, and sediment discharge for the period from 1984 to 2017 at nine gauging sites along the Upper Mississippi River. We used recent advances in remote sensing of ﬂuvial systems, such as automated river width extraction, Bayesian discharge inference with at-many-stations hydraulic geometry (AMHG), and SSC-SR regression models. With 621 Landsat scenes available from all the gauging sites, the results showed that the water discharge and SSC retrieval from Landsat imagery can yield reasonable sediment discharge estimates along the Upper Mississippi River. An overall relative bias of − 25.4, mean absolute error (MAE) of 6.24 × 10 4 tonne / day, relative root mean square error (RRMSE) of 1.21, and Nash–Sutcli ﬀ e E ﬃ ciency (NSE) of 0.49 were obtained for the sediment discharge estimation. Based on these statistical metrics, we identiﬁed three of the nine gauging sites (St. Louis, MO; Chester, IL; and Thebes, IL), which were in the downstream portion of the river, to be the best locations for estimating water and sediment discharge using Landsat imagery.


Introduction
Rivers are vital components of a drainage basin. They act as habitats for aquatic organisms, serve as drainage channels for surface water, and regulate the hydrological cycle. They provide water for domestic, industrial and agricultural use, power generation, and navigation [1,2]. Anthropogenic activities have increasingly altered rivers' physical and chemical properties [3][4][5][6]. These activities include agricultural practices, industrial operations, and land use and land cover changes, which all tend to negatively impact surface water quality through excessive point-or nonpoint-source loadings of pollutants [3,[7][8][9][10][11]. A major driver of this problem is human-induced sedimentation, which affects the hydrologic regime and water quality of a river and its drainage basin [12][13][14]. Human activities, such as agriculture, forest operations, mining, and urbanization, increase erosion and sediment transport within a drainage basin [15]. The nine study sites are located in four different states, namely Minnesota, Iowa, Illinois, and Missouri, with varying periods of available flow and sediment records (Table 1) Table 1. General information and periods of available flow and sediment data (both suspended sediment concentration (SSC) and sediment discharge) for the nine Upper Mississippi River study sites shown in Figure 1. Source: https://waterdata.usgs.gov/nwis/.   To estimate sediment discharge using remotely sensed data, we inspected the availability of remotely sensed imageries for three Landsat missions (5, 7, and 8) with different periods of operation in orbit ( Figure 2). Landsat 5 captured Earth images with a thematic mapper (TM) sensor from March 1984 to January 2013, when it was decommissioned. Landsat 7 started its mission in 1999 and is presently operational with its enhanced thematic mapper plus (ETM+) sensor. However, since May 2003, it has been adversely impacted by a failure of the scan line corrector. Landsat 8 was launched as the Landsat Data Continuity Mission with the Operational Land Imager (OLI). It began capturing Earth images in 2013 and remains operational.

Site
These Landsat missions share a common temporal resolution of 16-day repeat coverage, and 30-m spatial resolution for visible, near-infrared, and shortwave infrared bands. Multi-temporal Tier 1 Landsat scenes (pre-processed Landsat data at the highest available quality) were compiled using the Google Earth Engine platform [55], which allows a user to directly visualize, sort, and process imagery within the cloud platform. We identified 5490 Landsat scenes from March 1984 to September 2017 for which corresponding sediment data were available at the nine study sites in the Upper Mississippi River (Table 2). Some study sites share similar Landsat acquisition path/row numbers (e.g., Grafton, IL; Alton, IL; and St. Louis, MO, at 23/33), meaning that they had similar nominal scene centers.
To estimate sediment discharge using remotely sensed data, we inspected the availability of remotely sensed imageries for three Landsat missions (5, 7, and 8) with different periods of operation in orbit ( Figure 2). Landsat 5 captured Earth images with a thematic mapper (TM) sensor from March 1984 to January 2013, when it was decommissioned. Landsat 7 started its mission in 1999 and is presently operational with its enhanced thematic mapper plus (ETM+) sensor. However, since May 2003, it has been adversely impacted by a failure of the scan line corrector. Landsat 8 was launched as the Landsat Data Continuity Mission with the Operational Land Imager (OLI). It began capturing Earth images in 2013 and remains operational. These Landsat missions share a common temporal resolution of 16-day repeat coverage, and 30-m spatial resolution for visible, near-infrared, and shortwave infrared bands. Multi-temporal Tier 1 Landsat scenes (pre-processed Landsat data at the highest available quality) were compiled using the Google Earth Engine platform [55], which allows a user to directly visualize, sort, and process imagery within the cloud platform. We identified 5490 Landsat scenes from March 1984 to September 2017 for which corresponding sediment data were available at the nine study sites in the Upper Mississippi River (Table 2). Some study sites share similar Landsat acquisition path/row numbers (e.g., Grafton, IL; Alton, IL; and St. Louis, MO, at 23/33), meaning that they had similar nominal scene centers.

River Discharge Estimation using Landsat Imagery
River widths obtained from Landsat river masks were used as inputs for the river discharge retrievals. Prior to width extraction, the imagery collection for each site was filtered to obtain those with less than 10% cloud obstruction at the gauging site. The cross-sectional widths along the river channels were extracted using the RivWidthCloud algorithm in Google Earth Engine following Yang et al. [56]. RivWidthCloud automatically extracts the centerline and width transects from the river mask of the area of interest in the Landsat imagery. The algorithm also detects widths obstructed due to clouds, cloud shadows, and snow. The resulting river widths were visually inspected and filtered to retain only those having midpoints within a~1-km distance from their respective USGS gauging site (Figure 3), so as to best represent the reach of their site. In addition, widths within the standing confluence at an individual study site were excluded. The filtered river widths were subsequently used as inputs for the Bayesian-AMHG inference of river discharge. obstructed due to clouds, cloud shadows, and snow. The resulting river widths were visually inspected and filtered to retain only those having midpoints within a ~1-km distance from their respective USGS gauging site (Figure 3), so as to best represent the reach of their site. In addition, widths within the standing confluence at an individual study site were excluded. The filtered river widths were subsequently used as inputs for the Bayesian-AMHG inference of river discharge. The Bayesian-AMHG [37,39] was developed based upon the classic relationship between the hydraulic geometry of river channels and river discharge. Leopold and Maddock [57] posited a power-law statistical relationship between a river's width and discharge: where W is the river's width, Q is its discharge, and a and b are the classic AHG parameters. For an AMHG relationship, the AHG parameter a is substituted by the AMHG parameters, and the discharge equation becomes: where Wc and Qc are global AMHG parameters, and Wi and bi are, respectively, the river width and the classic AHG parameter for each river cross section [37,39]. Equation (2) is the basis for the inference of river discharge in the Bayesian-AMHG method. Bayesian inference requires both the likelihood of data and prior probability distributions of unknown parameters in order to obtain the full posterior distribution via Monte Carlo sampling [39]. Log-transforming Equation (2) yields a likelihood function and the prior probability distributions of the unknown parameters (bi, Wc, Qc, and Q). Based on a large number of acoustic doppler current profiler (ADCP) measurements at 10,081 USGS gauging stations, Hagemann et al. [39] obtained the likelihood function: The Bayesian-AMHG [37,39] was developed based upon the classic relationship between the hydraulic geometry of river channels and river discharge. Leopold and Maddock [57] posited a power-law statistical relationship between a river's width and discharge: where W is the river's width, Q is its discharge, and a and b are the classic AHG parameters. For an AMHG relationship, the AHG parameter a is substituted by the AMHG parameters, and the discharge equation becomes: where W c and Q c are global AMHG parameters, and W i and b i are, respectively, the river width and the classic AHG parameter for each river cross section [37,39]. Equation (2) is the basis for the inference of river discharge in the Bayesian-AMHG method. Bayesian inference requires both the likelihood of data and prior probability distributions of unknown parameters in order to obtain the full posterior distribution via Monte Carlo sampling [39]. Log-transforming Equation (2) yields a likelihood function and the prior probability distributions of the unknown parameters (b i , W c , Q c , and Q). Based on a large number of acoustic doppler current profiler (ADCP) measurements at 10,081 USGS gauging stations, Hagemann et al. [39] obtained the likelihood function: where ε is the AMHG error term with a standard deviation (SD) of 0.22. The prior distributions for Q and Q c both follow a truncated normal distribution with a coefficient of variation of 1. The empirical prior distribution of the AHG parameter b, which is normally distributed, is given by: where ε b is the random error with a standard deviation of 0.098. Further, the prior distribution of W c is log-normally distributed with its center at the mean width observed for the river reach, with a coefficient of variation of 0.01. From these likelihood and prior distributions, 1000 realizations of the posterior distribution were generated via Monte Carlo sampling, and a mean expected value of river discharge was obtained for each image [39]. We explored the sensitivity of the posterior mean value of Q to the prior distribution through the use of three different central tendency metrics. Following Feng et al. [38], we examined the potential of discharge estimation for gauged and ungauged settings. Three centers of distribution were analyzed, specifically, the median of daily discharge records Q median from the USGS database (gauged Remote Sens. 2020, 12, 2370 7 of 24 rivers) and two empirical water balance model mean discharges (Q wbm1 and Q wbm2 ) obtained from widths extracted for each study site (ungauged rivers). For the ungauged rivers, we developed an empirical power-law function to estimate mean discharge based solely on the widths obtained from Landsat imagery. Through reanalysis of widths and water balance model-based mean discharge for 34 independent rivers across the globe compiled by Hagemann et al. [39], we obtained this function: where Q wbm1 is the mean discharge from the first empirical water balance model and W is the mean width extracted for the study site (see Supplemental files, Figure S1). Q wbm2 is the mean discharge from the second empirical water balance model obtained by doubling Q wbm1 . This additional empirical metric was included to explore how doubling a center of distribution for Q prior affects the final discharge estimate.

Estimating Suspended Sediment Concentration
To estimate the SSC at the study sites, we fit the relationship between SSC and surface reflectance (SR) of Landsat images using a different set of images from those used in the discharge retrieval (see Supplemental files, Table S1). Using different image sets allowed an independent performance assessment of SSC and sediment discharge estimation. We identified and obtained the corresponding standard Landsat SR data, which were pre-processed using specialized atmospheric correction methods. The standard Landsat SR was generated from the Landsat Ecosystem Disturbance Adaptive Processing System (LEDAPS) for Landsats 5 and 7, and from the Surface Reflectance Code (LaSRC) for Landsat 8. Using Google Earth Engine, we applied a water mask for all Landsat imagery to identify for each study site a sampling area (diameter 240 m, except for study site 1 for which a radius of 120 m was used) that was at the center of the river reach, near the USGS gauge location, and clear of obstructions (e.g., bridges, ports, etc.) ( Figure 4). Mean surface reflectance within the sampling area was calculated for each available Landsat image.  Related studies reported a correlation between SSC and the ratio of red band (0.63-0.69 µm: Band 3 for Landsat 5 and 7, Band 4 for Landsat 8) to green band (0.52-0.60 µm: Band 2 for Landsat 5 and 7, Band 3 for Landsat 8) reflectance in Landsat imagery. This correlation was adopted for this study given the reasonable goodness of fit (R 2 > 0.5) [48][49][50][51]. Regression functions relating SSC and SR were developed using the in situ SSC measurements and the water-leaving red to green reflectance ratio for each of the specified sampling areas of the nine gauging stations, and for the Related studies reported a correlation between SSC and the ratio of red band (0.63-0.69 µm: Band 3 for Landsat 5 and 7, Band 4 for Landsat 8) to green band (0.52-0.60 µm: Band 2 for Landsat 5 and 7, Band 3 for Landsat 8) reflectance in Landsat imagery. This correlation was adopted for this study given the reasonable goodness of fit (R 2 > 0.5) [48][49][50][51]. Regression functions relating SSC and SR were developed using the in situ SSC measurements and the water-leaving red to green reflectance ratio for each of the specified sampling areas of the nine gauging stations, and for the Upper Mississippi River as a whole using the pooled data. A total of 1854 scenes with corresponding SSC data were used to develop the gauge-specific and the regional-scale regressive SSC-SR functions (for a detailed image count at each study site, see Supplement file, Table S1). These functions were then applied to a set of 621 Landsat images (1) that were successfully used for river discharge estimation, and (2) for which USGS sediment data were available to estimate the mean SSC from surface reflectance. We assume that the estimated mean SSCs are also representative of the channel depth as the in situ SSC measurement of the USGS follows the sampling techniques described by Edwards and Glysson [58]. Accordingly, this set of images with estimated Q and SSC was used for the final sediment discharge estimation.

Sediment Discharge Estimation and Performance Evaluation
The relationship between river sediment discharge, river discharge, and SSC is described by Porterfield [23] as: where Q s is the sediment discharge (t/d), Q is the river discharge (m 3 /s), SSC is the mean suspended sediment concentration (mg/L), and k is a coefficient based on the units of the river discharge that assumes a specific weight of 2.65 for sediment (0.0864 in SI units). We applied this relationship to obtain sediment discharges for each gauging station. Figure 5 summarizes the main steps to estimate river sediment discharge based on Landsat imagery. where Qs is the sediment discharge (t/d), Q is the river discharge (m 3 /s), SSC is the mean suspended sediment concentration (mg/L), and k is a coefficient based on the units of the river discharge that assumes a specific weight of 2.65 for sediment (0.0864 in SI units). We applied this relationship to obtain sediment discharges for each gauging station. Figure 5 summarizes the main steps to estimate river sediment discharge based on Landsat imagery. We evaluated the predictive performance of our approach to estimate river discharge, SSC, and sediment discharge for the study sites in the Upper Mississippi River with the set of 621 Landsat images from 1984 to 2017. Four statistical indicators, the Nash-Sutcliffe Efficiency (NSE), relative bias, mean absolute error (MAE), and the relative root mean square error (RRMSE), were calculated We evaluated the predictive performance of our approach to estimate river discharge, SSC, and sediment discharge for the study sites in the Upper Mississippi River with the set of 621 Landsat images from 1984 to 2017. Four statistical indicators, the Nash-Sutcliffe Efficiency (NSE), relative bias, mean absolute error (MAE), and the relative root mean square error (RRMSE), were calculated following Equations (7)-(10): where y i is the observed value (discharge, SSC, and sediment discharge), y i ' is the estimated value,ȳ is the average of the observed values, and N is the number of observations.

River Widths and Landsat Surface Reflectance Retrieval
Out of 5490 Landsat images available for the nine study sites (Table 2), 779 yielded river widths through the use of RivWidthCloud and post-processing techniques ( Table 3). The narrowest river channel, 121.8 m at Brooklyn Park, MN, was at the study site farthest upstream. The widest river channel, 947.6 m, was at Grafton, IL, toward the middle of the Upper Mississippi River. An advantage of automating the image-based width extraction is the consistency of width delineation from the river water mask, as seen in the similar coefficient of variation (CV) values. Such consistency reduces the uncertainty that would have been incurred by manually delineating river widths from remotely sensed images. From the same pool of 5490 images, those lacking either reflectance values in the sampling areas or corresponding sediment records from the USGS database were removed, leaving 2475 Landsat images for mean surface reflectance across the nine study sites. Lower mean surface reflectance was observed from the coastal aerosol, blue, and infrared bands rather than from green and red bands ( Figure 6). For Landsat 5, a higher mean reflectance in the red band (Band 3, 0.63-0.69 µm) was observed at St. Louis, MO, and Chester, IL. The remaining seven study sites had higher mean reflectance in the green band (Band 2, 0.52-0.60 µm), which can be associated with the presence of chlorophyll a and algae on the river surface due to lower turbulence in the water column [59][60][61]. On the other hand, the highest water-leaving reflectance was observed in the red band for both Landsat 7 (Band 3, 0.63-0.69 µm) and Landsat 8 (Band 4, 0.64-0.67 µm). A higher blue band reflectance indicates clear water while a lower infrared band reflectance confirms high absorption of water in the infrared spectrum. Further, a high water-leaving reflectance in the red band indicates a turbid river water column that brings up the red or brown suspended sediment areas visible on the water surface at the time of image capture. In remote sensing, the terms turbidity and SSC are used interchangeably because they tend to be highly correlated [46,62]. Note that the study sites in St. Louis, MO; Chester, IL; and Thebes, IL, have the highest water-leaving reflectance among the nine study sites.
Remote Sens. 2020, 12, x FOR PEER REVIEW 11 of 25 interchangeably because they tend to be highly correlated [46,62]. Note that the study sites in St. Louis, MO; Chester, IL; and Thebes, IL, have the highest water-leaving reflectance among the nine study sites.

Sensitivity analysis for Q prior
The dependence and sensitivity of the discharge estimates were demonstrated relative to the center of Q prior distribution. Overlaying the three centers of distribution on in situ discharge record density plots revealed which best estimated the peak of the density plot for three representative study sites (Figure 7). Qwbm1 is most adjacent to the peak of the discharge density plot for the study site at McGregor, IA ( Figure 7a). As expected, using this center for Q prior yielded a better discharge estimate than using either of the other two centers of distribution (Figure 7b). Similarly, Qmedian appears most adjacent to the peak of the in situ discharge density plot for the study site at Grafton, IL (Figure 7c), so using Qmedian as the center of distribution for Q prior led to better performance in estimating the discharge than using either Qwbm1 or Qwbm2 (Figure 7d). For the study site at Thebes, IL (Figure 7e,f), the use of Qwbm2 as the center for the Q prior distribution produced the best discharge estimate because this center is the most adjacent to the peak of the discharge density plot for this river segment. Using an underestimated center for the Q prior distribution tends to underestimate the discharge, and using an overestimated center for the Q prior will likely overestimate the discharge.

Sensitivity Analysis for Q Prior
The dependence and sensitivity of the discharge estimates were demonstrated relative to the center of Q prior distribution. Overlaying the three centers of distribution on in situ discharge record density plots revealed which best estimated the peak of the density plot for three representative study sites (Figure 7). Q wbm1 is most adjacent to the peak of the discharge density plot for the study site at McGregor, IA ( Figure 7a). As expected, using this center for Q prior yielded a better discharge estimate than using either of the other two centers of distribution (Figure 7b). Similarly, Q median appears most adjacent to the peak of the in situ discharge density plot for the study site at Grafton, IL (Figure 7c), so using Q median as the center of distribution for Q prior led to better performance in estimating the discharge than using either Q wbm1 or Q wbm2 (Figure 7d). For the study site at Thebes, IL (Figure 7e,f), the use of Q wbm2 as the center for the Q prior distribution produced the best discharge estimate because this center is the most adjacent to the peak of the discharge density plot for this river segment. Using an underestimated center for the Q prior distribution tends to underestimate the discharge, and using an overestimated center for the Q prior will likely overestimate the discharge. The performance of Bayesian-AMHG inference using the river widths obtained from the 779 Landsat images (filtered and post-processed) and three different Q prior centers is summarized in Table 4. Results varied across the nine study sites. As expected, the relative bias was found to be lower with the aid of USGS gauge records (for Qmedian) than the approximated mean discharge obtained from either empirical water balance model. An average relative bias of ±0.12 was obtained for all the study sites and found within the "satisfactory" (±0.25) category [63]. Large errors were observed with the use of Qwbm1 and Qwbm2 with frequent underestimation (negative bias) and overestimation (positive bias). These results are similar in terms of RRMSE. A lower RRMSE with an average of 0.49 was obtained using Qmedian, whereas the estimations with empirically approximated Q priors resulted in a higher average RRMSE (> 0.60). Further, the highest NSE was also obtained with Qmedian, averaging 0.36 for all the study sites. Among all sites, study site 5 at Grafton, IL, produced the highest MAE (>1850 m 3 /s) and RRMSE (>0.70) with negative NSEs regardless of which Q prior was used. Even with a rather large number of images (N = 107), this site posed a challenge for river width inference, possibly due to the large width (947.6 m) and standard deviation (98.9 m) ( Table 3). In relative terms, study sites 7, 8, and 9, the three farthest downstream sites, were the best sites for predicting discharge using either Qmedian or Qwbm2 as prior information. Despite the unexpected result that Qwbm2 improved the discharge estimates for these three sites, the estimates themselves were poor for the remaining sites, as indicated by a higher positive bias and RRMSE averaging 0.62 and 0.92, respectively.

Width-Based Discharge Estimates
The performance of Bayesian-AMHG inference using the river widths obtained from the 779 Landsat images (filtered and post-processed) and three different Q prior centers is summarized in Table 4. Results varied across the nine study sites. As expected, the relative bias was found to be lower with the aid of USGS gauge records (for Q median ) than the approximated mean discharge obtained from either empirical water balance model. An average relative bias of ±0.12 was obtained for all the study sites and found within the "satisfactory" (±0.25) category [63]. Large errors were observed with the use of Q wbm1 and Q wbm2 with frequent underestimation (negative bias) and overestimation (positive bias). These results are similar in terms of RRMSE. A lower RRMSE with an average of 0.49 was obtained using Q median, whereas the estimations with empirically approximated Q priors resulted in a higher average RRMSE (>0.60). Further, the highest NSE was also obtained with Q median , averaging 0.36 for all the study sites. Among all sites, study site 5 at Grafton, IL, produced the highest MAE (>1850 m 3 /s) and RRMSE (>0.70) with negative NSEs regardless of which Q prior was used. Even with a rather large number of images (N = 107), this site posed a challenge for river width inference, possibly due to the large width (947.6 m) and standard deviation (98.9 m) ( Table 3). In relative terms, study sites 7, 8, and 9, the three farthest downstream sites, were the best sites for predicting discharge using either Q median or Q wbm2 as prior information. Despite the unexpected result that Q wbm2 improved the discharge estimates for these three sites, the estimates themselves were poor for the remaining sites, as indicated by a higher positive bias and RRMSE averaging 0.62 and 0.92, respectively. Discharge estimates were generally improved by using Q median in the priors (Figure 8). The Bayesian-AMHG estimates and the in situ discharge records are agreeable (NSE ≥ 0.25) for all the study sites except McGregor, IA, and Grafton, IL. Further, the performance of the Bayesian-AMHG was the best for the three farthest downstream sites, which generally matched the measured discharge with an average relative bias of 0.08, MAE of 1591 m 3 /s, RRMSE of 0.37, and NSE of 0.62.
Remote Sens. 2020, 12, x FOR PEER REVIEW 14 of 25 Discharge estimates were generally improved by using Qmedian in the priors (Figure 8). The Bayesian-AMHG estimates and the in situ discharge records are agreeable (NSE ≥ 0.25) for all the study sites except McGregor, IA, and Grafton, IL. Further, the performance of the Bayesian-AMHG was the best for the three farthest downstream sites, which generally matched the measured discharge with an average relative bias of 0.08, MAE of 1591 m 3 /s, RRMSE of 0.37, and NSE of 0.62.

Gauge-specific SSC-SR Model
The regression models for SSC vs. SR are significant at a significance level of 0.05 for all study sites (Table 5). SSC tends to increase as the water-leaving reflectance from the red band exceeds that from the green band. The functions are similar for study sites 1, 2, 3, 4, and 6 with an average R² of 0.28. Meanwhile, study sites 5, 7, 8, and 9 have similar regression coefficients, with higher slopes and an average R 2 of 0.61. Of the nine sites, the three farthest downstream have the greatest R² for the SSC-SR regressions: 0.66 for St. Louis, MO; 0.64 for Chester, IL; and 0.54 for Thebes, IL.

Gauge-Specific SSC-SR Model
The regression models for SSC vs. SR are significant at a significance level of 0.05 for all study sites (Table 5). SSC tends to increase as the water-leaving reflectance from the red band exceeds that from the green band. The functions are similar for study sites 1, 2, 3, 4, and 6 with an average R 2 of 0.28. Meanwhile, study sites 5, 7, 8, and 9 have similar regression coefficients, with higher slopes and an average R 2 of 0.61. Of the nine sites, the three farthest downstream have the greatest R 2 for the SSC-SR regressions: 0.66 for St. Louis, MO; 0.64 for Chester, IL; and 0.54 for Thebes, IL.

Regional-Scale SSC-SR Model
The regional SSC-SR model is significant at a significance level of 0.05 for all missions: Landsat 5 (p-value < 0.001, R 2 = 0.44), Landsat 7 (p-value < 0.001, R 2 = 0.64), and Landsat 8 (p-value < 0.001, R 2 = 0.62). Landsats 7 and 8 exhibit similar patterns of clustering, while as previously discussed, the dispersion in the low SSC range in the Landsat 5 regression is due to the upstream Landsat 5-dominated sites. The regression model for the region with data from all three missions used was also significant at the 0.05 level (p-value < 0.001, R 2 = 0.50) and "acceptable" (Figure 9) [63]. 3.3.2. Regional-scale SSC-SR Model The regional SSC-SR model is significant at a significance level of 0.05 for all missions: Landsat 5 (p-value < 0.001, R 2 = 0.44), Landsat 7 (p-value < 0.001, R 2 = 0.64), and Landsat 8 (p-value < 0.001, R² = 0.62). Landsats 7 and 8 exhibit similar patterns of clustering, while as previously discussed, the dispersion in the low SSC range in the Landsat 5 regression is due to the upstream Landsat 5-dominated sites. The regression model for the region with data from all three missions used was also significant at the 0.05 level (p-value < 0.001, R² = 0.50) and "acceptable" (Figure 9) [63]. The estimation of SSC with the regional model results in a high positive bias, RRMSE, and mostly negative NSEs for most of the gauging sites ( Table 6). The only positive NSE was obtained for the study site at Alton, IL, with a relative bias of −0.13 and MAE of 77.9 mg/L. In comparison to the regional model, most gauge-specific models yielded a lower RRMSE and higher NSE. Interestingly, the SSC estimates for the three farthest downstream sites at St. Louis, MO; Chester, IL; and Thebes, IL, using their respective gauge-specific models, were all reasonable, with an SSC ranging from 2-1250 mg/L and NSE averaging 0.20. These results were supported by other statistical metric The estimation of SSC with the regional model results in a high positive bias, RRMSE, and mostly negative NSEs for most of the gauging sites ( Table 6). The only positive NSE was obtained for the study site at Alton, IL, with a relative bias of −0.13 and MAE of 77.9 mg/L. In comparison to the regional model, most gauge-specific models yielded a lower RRMSE and higher NSE. Interestingly, the SSC estimates for the three farthest downstream sites at St. Louis, MO; Chester, IL; and Thebes, IL, using their respective gauge-specific models, were all reasonable, with an SSC ranging from 2-1250 mg/L and NSE averaging 0.20. These results were supported by other statistical metric averages: a relative bias of 0.19, MAE of 116.3 mg/L, and RRMSE of 0.67.
Clearly, using gauge-specific models is better than using a regional model in estimating the SSC for the study sites. The regional model led to an overall MAE of 182.2 mg/L, RRMSE of 2.16, and a negative NSE of −3.06 (Figure 10). In contrast, the gauge-specific SSC-SR models yielded better estimations, with a lower MAE of 94.7 mg/L, RRMSE of 0.90, and a positive NSE of 0.25.
Remote Sens. 2020, 12, x FOR PEER REVIEW 17 of 25 Clearly, using gauge-specific models is better than using a regional model in estimating the SSC for the study sites. The regional model led to an overall MAE of 182.2 mg/L, RRMSE of 2.16, and a negative NSE of −3.06 (Figure 10). In contrast, the gauge-specific SSC-SR models yielded better estimations, with a lower MAE of 94.7 mg/L, RRMSE of 0.90, and a positive NSE of 0.25. Figure 10. Comparison of gauge-specific (left) and regional (right) SSC-SR models for the Upper Mississippi River. Note the different y-axis scales. The dashed line represents the 1:1 relationship.

Sediment Discharge Estimation
By combining the discharge estimates with the SSC estimates from the gauge-specific model, we calculated sediment discharge for the nine study sites for 621 Landsat images ( Figure 11). Interestingly, the site with the best sediment discharge estimates, in terms of RRMSE (0.88) and NSE (0.51), is 5, Grafton, IL. Estimated and observed sediment discharges are generally agreeable among the nine sites, following a similar trend over time. However, extreme events were noticeably overestimated for some images, especially for McGregor, IA; Clinton, IA; and Alton, IL. The second-best study site is 7, St. Louis, MO, with an RRMSE of 1.05 and NSE of 0.47. The estimation appeared problematic at study site 3, McGregor, IA, with a high RRMSE (2.44) and a negative NSE (−5.33). Figure 10. Comparison of gauge-specific (left) and regional (right) SSC-SR models for the Upper Mississippi River. Note the different y-axis scales. The dashed line represents the 1:1 relationship. Table 6. Performance of gauge-specific and regional SSC-SR models.

Site
Location Observed Range (mg/L) N

Relative Bias MAE (mg/L) RRMSE NSE
Gauge-Specific Regional Gauge-Specific Regional Gauge-Specific Regional Gauge-Specific Regional

Sediment Discharge Estimation
By combining the discharge estimates with the SSC estimates from the gauge-specific model, we calculated sediment discharge for the nine study sites for 621 Landsat images ( Figure 11). Interestingly, the site with the best sediment discharge estimates, in terms of RRMSE (0.88) and NSE (0.51), is 5, Grafton, IL. Estimated and observed sediment discharges are generally agreeable among the nine sites, following a similar trend over time. However, extreme events were noticeably overestimated for some images, especially for McGregor, IA; Clinton, IA; and Alton, IL. The second-best study site is 7, St. Louis, MO, with an RRMSE of 1.05 and NSE of 0.47. The estimation appeared problematic at study site 3, McGregor, IA, with a high RRMSE (2.44) and a negative NSE (−5.33).
Remote Sens. 2020, 12, x FOR PEER REVIEW 18 of 25 Figure 11. Estimated and observed sediment discharges for the individual study sites.
Combining the discharge and SSC retrieval from Landsat imagery for certain river segments can produce reasonable sediment discharge estimates, as confirmed by comparison with in situ field measurements ranging from 1.78 × 10 2 to 1.64 × 10 6 t/d. The relative biases and errors were reduced by using the discharge estimates for gauged settings and a gauge-specific SSC retrieval model. We  (Table 7). On the other hand, the estimations for the remaining sites were problematic with a higher relative bias (> 0.60) and lower NSE (< 0.20).  Combining the discharge and SSC retrieval from Landsat imagery for certain river segments can produce reasonable sediment discharge estimates, as confirmed by comparison with in situ field measurements ranging from 1.78 × 10 2 to 1.64 × 10 6 t/d. The relative biases and errors were reduced by using the discharge estimates for gauged settings and a gauge-specific SSC retrieval model. We obtained an average relative bias of 0.23, RRMSE of 0.95, and NSE of 0.40 in simulating the sediment discharge at the USGS gauging sites at Winona, MN; Grafton, IL; St. Louis, MO; Chester, IL; and Thebes, IL (Table 7). On the other hand, the estimations for the remaining sites were problematic with a higher relative bias (>0.60) and lower NSE (<0.20). The overall performance of the sediment discharge estimation is illustrated in Figure 12. With 621 Landsat scenes available across all study sites, we obtained a relative bias of −25.4, MAE of 6.24 × 10 4 t/d, RRMSE of 1.21, and NSE of 0.49. Moriasi et al. [63] recommended that model performance for sediment simulation have a relative bias ± 0.55 and NSE > 0.50 in order to be classified as satisfactory. Their model performance rating is aimed at monthly time step hydrologic simulations for a drainage basin, whereas the method described here is for independent events, at time intervals given by available Landsat imagery for the river reach.
Remote Sens. 2020, 12, x FOR PEER REVIEW 19 of 25 The overall performance of the sediment discharge estimation is illustrated in Figure 12. With 621 Landsat scenes available across all study sites, we obtained a relative bias of −25.4, MAE of 6.24 × 10 4 t/d, RRMSE of 1.21, and NSE of 0.49. Moriasi et al. [63] recommended that model performance for sediment simulation have a relative bias ± 0.55 and NSE > 0.50 in order to be classified as satisfactory. Their model performance rating is aimed at monthly time step hydrologic simulations for a drainage basin, whereas the method described here is for independent events, at time intervals given by available Landsat imagery for the river reach.

Discussion
Our results of estimating river discharge using Landsat imagery support the finding of Feng et al. [38] that the use of Bayesian-AMHG provides better discharge estimates for gauged rivers than for an ungauged setting for which Q priors are approximated from data reanalysis. We also observed that the closer the center of distribution for Q prior was to the peak of the discharge density plot, the more likely the discharge estimate would agree with the observed value. Hence, future studies should further explore how to best approximate the center of distribution of river discharge solely from remotely sensed data.
Our results also suggest that outputs from RivWidthCloud can be effectively used for Bayesian-AMHG discharge inference. The consistency of this automated width retrieval provides a substantial advantage over manual delineation of river widths from remotely sensed imagery. The sediment discharge estimates for the three study sites farthest downstream were "satisfactory" (±0.25 relative bias and NSE > 0.50) following the model performance evaluation criteria recommended by Moriasi et al. [63]. We posit that these results are associated with the morphologic

Discussion
Our results of estimating river discharge using Landsat imagery support the finding of Feng et al. [38] that the use of Bayesian-AMHG provides better discharge estimates for gauged rivers than for an ungauged setting for which Q priors are approximated from data reanalysis. We also observed that the closer the center of distribution for Q prior was to the peak of the discharge density plot, the more likely the discharge estimate would agree with the observed value. Hence, future studies should further explore how to best approximate the center of distribution of river discharge solely from remotely sensed data.
Our results also suggest that outputs from RivWidthCloud can be effectively used for Bayesian-AMHG discharge inference. The consistency of this automated width retrieval provides a substantial advantage over manual delineation of river widths from remotely sensed imagery. The sediment discharge estimates for the three study sites farthest downstream were "satisfactory" (±0.25 relative bias and NSE > 0.50) following the model performance evaluation criteria recommended by Moriasi et al. [63]. We posit that these results are associated with the morphologic difference between reaches of the study sites. The Rosgen Type I plan view classification of natural rivers broadly characterizes channels as relatively straight (class A), low sinuosity (class B), meandering (class C), braided (class D), anastomosed (class DA), and tortuously meandering (class E) [64]. Study sites 1-6 possess class D channel type features (see Figure 4): complex multiple or braided channels with wide eroding banks. In contrast, study sites 7, 8, and 9 have entrenched and stable channels, more characteristic of class A or B. Consequently, we consider that the selection of a representative reach for a river system is crucial for multi-river studies because the performance in estimating discharge from optical sensors varies widely from one reach to another.
We also observed different patterns in our developed SSC-SR retrieval models. Study sites 1, 2, 3, 4, and 6, which have nearly identical regression functions (Table 5), were likely influenced by the low SSC levels (2-324 mg/L) and by Landsat 5 being the predominant sensor platform used to acquire their images. The lower SSC levels at these study sites may have led to a weaker correlation with the water-leaving reflectance captured in the Landsat 5 images. Further, most of these sites have a higher mean green band than red band reflectance ( Figure 6), suggesting that the river segments with low SSC levels or turbidity do not fully illuminate the brownish sediment color in Landsat imagery. This is further confirmed with the regional SSC-SR retrieval model, which shows clearly that Landsat 5 data points differ from Landsats 7 and 8 data points, while the latter two closely match each other (Figure 9). The matching functions for study sites 5, 7, 8, and 9 may be a result of the dominance of the Landsat sensor used and to the proximity of the study sites. In most cases, these study sites were in the same Landsat image, as they share adjacent or similar Landsat acquisition paths/rows (path 24, row 33 and path 23, row 34; Table 2). It is worth noting that SSCs at these sites are higher (12.2-2340 mg/L) than the other five sites, suggesting that sites with lower SSC levels and captured by Landsat 5 likely have a lower coefficient of determination for the SSC-SR model.
Most published SSC-SR models were based on data from a single Landsat mission (e.g., Landsat 8). Pham et al. [51] presented SSC-SR models (R 2 = 0.75) using the red to green band ratio from Landsat 8, with 40 images and SSC levels ranging from 22.4-178 mg/L. A related study by Pereira et al. [49] using the red to green band reflectance ratio of 26 Landsat 8 images yielded an SSC-SR model with a similarly high coefficient of determination (R 2 = 0.86) for SSC levels of 49-533 mg/L. Despite these models being developed with a smaller sample and a lower SSC range in their respective river systems, our results show that the red to green band reflectance ratio from multiple Landsat missions can be used to estimate SSCs of certain reaches in the Upper Mississippi River. This result is consistent with the findings of Markert et al. [48] and Peterson et al. [50], both of which reported SSC-SR models with R 2 ≈ 0.5 using SR data from multiple Landsat missions (5, 7, and 8). Nonetheless, some biases and errors in our SSC estimates were related to the low sediments in the water column captured by the Landsat 5 sensor. As such, other band combinations should be explored when using Landsat 5 surface reflectance data. Note that the errors may also be due to other variables, such as higher blue and red absorption with the presence of chlorophyll and algae during low flows, and the presence of organic materials upstream in the Upper Mississippi River [59,65].
Studies estimating annual sediment discharge using gauge records and composited remotely sensed data (250-m MODIS) have shown good agreement with observed data with a ±2% mean relative difference [17,42]. Compared with these studies, we demonstrated a "per event" or "per image" estimation approach where we approximated the sediment discharge based on gauge records and 30-m Landsat data. Extreme events, such as the "Great Flood of 1993," a 100-year flood event, were also covered in our simulations for the nine gauge sites of the Upper Mississippi River. Given that the approach is highly dependent on the accuracy of the discharge and SSC estimates, the uncertainties in their estimates can either increase the magnitude of errors or improve the final sediment discharge estimates, as illustrated in Figure 13. For instance, study site 5 (Grafton, IL), shows a negative relative bias (−0.19) in Q estimates and a positive relative bias (0.35) in SSC estimates. As a result, the Q and SSC errors offset each other to yield an unexpectedly good match between sediment discharge estimates and measurements. Unlike study site 5, site 7 (St. Louis, IL) has a minimal relative bias (<0.15) in both Q and SSC estimates. Thus, we can reasonably expect that this site will have good overall sediment discharge estimates with a low relative bias and better model fit. For this reason, we regard the St. Louis site as the best site within the Upper Mississippi River for estimating sediment discharge using Landsat data. Future efforts should be devoted to improving the Q and SSC estimation to further advance the utility of Landsat data or other optical remote sensing platforms for sediment discharge estimation in river systems. A river's Q can be estimated based on channel width, but prior discharge information is required to perform the Bayesian-AMHG inference. Readily accessible prior mean discharge of river reaches, e.g., from Lin et al. [66], may be examined to apply discharge estimations to ungauged rivers. Alternative approaches to estimating SSC, such as non-linear regression and machine learning estimation [50], should be explored. Last, even without gauge records for calibrations, a method estimating water and sediment discharge in rivers using remotely sensed data would still be useful in evaluating and screening for those susceptible to extreme flooding events or exhibiting excessive sedimentation.

Conclusions
We explored the potential of using remotely sensed images to estimate the discharge of water and sediment in a river system. In many cases, traditional monitoring operations are hampered by accuracy requirements, cost, and safety concerns regarding their labor-intensive field sampling methods. For this reason, the number of active gauging stations in many areas across the US continues to decline, despite increasing demands to monitor flow, sediment transport, and other aspects of river health. Our alternative approach, using space-based observations to assess the status of river systems, helps to address this need.
Conclusions from this study are: 1. Width outputs from RivWidthCloud can be effectively used for Bayesian-AMHG inference of river discharge. 2. Discharge estimations are influenced by both prior information and morphologic features along the river. 3. Higher biases and errors in SSC estimates tended to result from Landsat 5 sensors capturing Future efforts should be devoted to improving the Q and SSC estimation to further advance the utility of Landsat data or other optical remote sensing platforms for sediment discharge estimation in river systems. A river's Q can be estimated based on channel width, but prior discharge information is required to perform the Bayesian-AMHG inference. Readily accessible prior mean discharge of river reaches, e.g., from Lin et al. [66], may be examined to apply discharge estimations to ungauged rivers. Alternative approaches to estimating SSC, such as non-linear regression and machine learning estimation [50], should be explored. Last, even without gauge records for calibrations, a method estimating water and sediment discharge in rivers using remotely sensed data would still be useful in evaluating and screening for those susceptible to extreme flooding events or exhibiting excessive sedimentation.

Conclusions
We explored the potential of using remotely sensed images to estimate the discharge of water and sediment in a river system. In many cases, traditional monitoring operations are hampered by accuracy requirements, cost, and safety concerns regarding their labor-intensive field sampling methods. For this reason, the number of active gauging stations in many areas across the US continues to decline, despite increasing demands to monitor flow, sediment transport, and other aspects of river health. Our alternative approach, using space-based observations to assess the status of river systems, helps to address this need.
Conclusions from this study are:

1.
Width outputs from RivWidthCloud can be effectively used for Bayesian-AMHG inference of river discharge.

2.
Discharge estimations are influenced by both prior information and morphologic features along the river.

3.
Higher biases and errors in SSC estimates tended to result from Landsat 5 sensors capturing scenes with low sediment levels in the water column. This suggests that for Landsat 5 surface reflectance data, additional band combinations should be explored.

4.
Landsat imagery-based estimates of Q and SSC can yield reasonable sediment discharge estimates. In this study of the Upper Mississippi River, estimates had a relative bias of −25.4, MAE of 6.24 × 10 4 t/d, RRMSE of 1.21, and NSE of 0.49.

5.
Because sediment discharge estimates are the product of two other independent estimates (water discharge and SSC), biases and errors from these component estimates can either increase or decrease the magnitude of errors in the sediment discharge estimates. 6.
Even without gauge records for calibrations, this method can be used to estimate water and sediment discharges of rivers, and to evaluate and screen for rivers susceptible to extreme flooding or exhibiting excessive sedimentation.
This study demonstrates the potential of estimating water and sediment discharges-crucial hydrological information-using remotely sensed data as an alternative to labor-intensive field methods. Future efforts should be devoted to refining the Q and SSC estimation to further advance the utility of Landsat data for sediment discharge estimation in river systems.