Remote Sensing of Turbidity in the Tennessee River Using Landsat 8 Satellite

: The Tennessee River in the United States is one of the most ecologically distinct rivers in the world and serves as a great resource for local residents. However, it is also one of the most polluted rivers in the world, and a leading cause of this pollution is storm water runoff. Satellite remote sensing technology, which has been used successfully to study surface water quality parameters for many years, could be very useful to study and monitor the quality of water in the Tennessee River. This study developed a numerical turbidity estimation model for the Tennessee River and its tributaries in Southeast Tennessee using Landsat 8 satellite imagery coupled with near real-time in situ measurements. The obtained results suggest that a nonlinear regression-based numerical model can be developed using Band 4 (red) surface reﬂectance values of the Landsat 8 OLI sensor to estimate turbidity in these water bodies with the potential of high accuracy. The accuracy assessment of the estimated turbidity achieved a coefﬁcient of determination (R 2 ) value and root mean square error (RMSE) as high as 0.97 and 1.41 NTU, respectively. The model was also tested on imagery acquired on a different date to assess its potential for routine remote estimation of turbidity and produced encouraging results with R 2 value of 0.94 and relatively high RMSE. to study surface water quality in the Tennessee River and its surrounding water bodies. This study aims to investigate the potential of remote sensing technology to study surface water quality in the watersheds of Southeast Tennessee using satellite observations coupled with ﬁeld measurements.


Introduction
Satellite remote sensing technology has been providing an opportunity for synoptic and multitemporal viewing of water quality for more than 25 years [1][2][3][4]. Sensors aboard satellites can measure the amount of solar radiation at various wavelengths reflected by surface water which can be correlated to different water quality parameters [5][6][7]. This constitutes an alternative means of estimating water quality, offering three significant advantages over ground sampling. First, the near-continuous spatial coverage of satellite imagery allows for synoptic estimates over large areas. Second, the global coverage of satellites allows for the estimation of water quality in remote and inaccessible areas. Third, the relatively long record of archived imagery (e.g., Landsat since the early 1970s) allows estimation of historical water quality when no ground measurements can possibly be performed. Although there are challenges with low spatial resolution imagery when compared with in situ measurements, the availability of spatially and temporally distributed datasets can be very useful to understand the trend of water quality distribution at basin or regional scales.
Currently no remote sensing-based algorithm is available to study surface water quality in the Tennessee River and its surrounding water bodies. This study aims to investigate the potential of remote sensing technology to study surface water quality in the watersheds of Southeast Tennessee using satellite observations coupled with field measurements. Path radiance (L p ) is the radiance that is recorded by a sensor resulting from downwelling solar and sky radiation. This is usually unwanted as it never reaches the water. The free-surface layer or boundary layer radiance (L s ) is the radiance that reaches the air-water interface. It only penetrates about 1 mm and is then reflected from the water surface. The spectral information about the near-surface characteristics of the water is contained in this reflected energy. The subsurface volumetric radiance (L v ) is the radiance that penetrates the air-water interface and interacts with the organic/inorganic constituents in the water, and then exits the water column without encountering the bottom. It provides information about the internal bulk characteristics of the water column. The bottom radiance (L b ) is the radiance that reaches the bottom of the waterbody, reflects from it, propagates back through the water column, and then exits the water column. This is useful to infer information about the bottom (e.g., depth, color) [8].
The water-quality studies using remotely sensed data are usually interested in measuring the subsurface volumetric radiance (L v ) exiting the water column toward the sensor. The characteristics of this radiant energy are a function of the concentration of pure water (W), inorganic suspended minerals (SM), organic chlorophyll a (Chl-a), dissolved organic material (DOM), and the total amount of absorption and scattering attenuation that takes place in the water column due to each of these constituents, c( λ ), as shown in Equation (2) [8]: It is important to look at the effect that each of these constituents has on the spectral reflectance characteristics of the water column. To estimate water quality parameters (WQPs) from different satellite data, three different approaches are generally used [3,4,9]: Empirical approach, which is based on the development of bivariate or multivariate regressions between remote sensing data and measured in situ WQPs. Digital numbers (DN) or radiance values at the sensor, as well as their band combinations, are correlated with in situ measurements of different WQPs, usually acquired at near real-time of the sensor overpass. A summary of empirical approaches for lakes can be found in Lindell et al. [10].

2.
Semi-empirical approach, which may be used when spectral characteristics of the parameters of interest are known. This knowledge of spectral characteristics is included in the statistical analysis by focusing on well-chosen spectral areas and appropriate spectral bands used as correlates. An example of a semi-empirical approach with different sensors is reported by Härmä et al. [11], who investigated Finnish lakes.

3.
Analytical approach, where WQPs are related to the bulk Inherent Optical Properties (IOPs) via the Specific Inherent Optical Properties (SIOPs). The IOPs of the water column are then related to the Apparent Optical Properties (AOPs) and hence to the Top of Atmosphere (TOA) radiance, such as described by the radiative transfer theory [12,13]. The analytical method involves inverting the above-described relations to determine the WQPs from remote sensing data. An example of such an approach, using Landsat over lakes, can be found in Dekker et al. [14] for the total suspended matter retrieval.
A recent comprehensive review on water quality remote sensing can be found in Gholizadeh et al. [15], while a brief discussion on the current state of the art of the estimation of common WQPs using remote sensing is provided below.

Turbidity
Turbidity is a measure of water clarity derived from how much light is scattered by material in a body of water. Turbidity is a commonly assessed water quality parameter due to its optically active properties and was predicted as having high potential for correlation before satellite imagery was readily available [16]. In the last four decades, quantitative predictions of turbidity, and other water transparency parameters, have been successfully carried out using spaceborne satellites in lacustrine, estuarine, reservoir, and coastal environments and validated against in situ measurements [17][18][19][20][21][22][23][24][25][26].
Regression analyses with single or multi-band algorithms are common methods for assessing correlation between spectral reflectance and in situ turbidity measurements [17,18,20,21,27,30]. These methods are simple, efficient, and inexpensive for predicting water quality parameters with near-continuous spatial coverage. Other methods have found significant potential for correlation as well. Gene-expression programming (GEP) and geographically and temporally weighted regression (GTWR) have been proven to predict turbidity more accurately than that by multiple linear regression (MLR) among different reservoirs [23,26].
Beyond quantitative studies, qualitative studies are also common for prediction of turbidity when there are little to no in situ data available [19,29,31]. These qualitative predictions can be useful for general assessment of water quality conditions when no quantitative predictions are necessary. For example, one study developed an index for qualitative analysis of turbidity, the Normalized Difference Turbidity Index (NDTI). This index used the green and red band reflectance values instead of in situ measurements and was successfully utilized during the COVID-19 pandemic when field sampling could not take place [29,31].
However, few studies have used spaceborne satellites to study turbidity in rivers. River waters, which contain inherent absorption and scattering properties, remain poorly understood. They consistently contribute to markedly greater error when they are correlated with spectral reflectance [28]. Broad predictions of turbidity levels in rivers can still be made with statistically sound results when using specialized algorithms and regression models [27,28]. A qualitative assessment of fluvial turbidity was conducted by Greg et al. [29].
Chl-a is an optically active form of chlorophyll used in photosynthesis. It is one of the most studied remotely sensed water quality parameters. The majority of studies have found strongly correlated values detecting remotely sensed Chl-a in lakes, reservoirs and coastal areas using data from Landsat 1-3, Landsat 5, Landsat 8, RapidEye, Sentinel-2, and MODIS satellites [1,5,6,20,21,[32][33][34][35][36][37][38]. Successful correlations have also been found in past studies that utilized airborne remote sensing [30,39]. Limited success has been found in estuarine environments using imagery from Landsat 5 satellite and may need fine resolution imagery, as shown in [19,40] where IKONOS imagery was used.
An alternative way of assessing Chl-a is by using the maximum chlorophyll index (MCI) which was developed to detect algal blooms specifically with MEdium Resolution Imaging Spectrometer (MERIS) imagery. Planet's RapidEye satellite constellation was successfully used on lakes to create an equivalent product to that of MERIS detecting MCI [1].

Suspended Sediment Concentration (SSC)
SSC values are derived by measuring the dry weight of sediment from a known volume of water. SSC has optically active properties that can be easily correlated against in situ measurements using regression analysis techniques. SSC has been successfully estimated in lacustrine, fluvial, and estuarine environments by several modern spaceborne sensors such as IKONOS, Landsat 8, Sentinel-2, RapidEye, and MODIS [34,[40][41][42][43][44].
SSC has also been estimated using turbidity as a surrogate. There is a need to provide continuous records of SSC among ambient water gauge stations, particularly in polluted fluvial environments, as collecting in situ data is inefficient and costly [45]. In situ turbidimeters are more practical and ubiquitous in water gauge stations and have proven to be a successful surrogate for indirectly measured SSC [46,47]. Many techniques have been published to reduce turbidity-derived SSC errors, including using artificial neural networks, calibrating turbidity estimates using remotely sensed imagery for seasonal variation differences, considering physical properties of SSC, and understanding the effects of hysteresis [48][49][50][51].

Total Suspended Solids (TSS)
TSS are optically active and they differ from SSC only by the differing laboratory analyses used to assess each. When compared to SSC, quantifying TSS has proven to be unreliable for samples measured from natural waters as methods were originally designed for analyzing wastewater [52]. Regardless, TSS have been successfully estimated in lacustrine and estuarine environments using regression analysis techniques based on in situ measurements and reflectance data obtained from spaceborne optical sensors such as Landsat 5 and Landsat 8 [5,[19][20][21]. Airborne remote sensing has been successful in assessing TSS of reservoirs and rivers, noting specific advantages over spaceborne satellites including less atmospheric interference, adaptable spatial resolution, and fewer temporal restraints [30].

Salinity
Salinity is the total dissolved salt in water measured in parts per thousand. Salinity is known to be optically active, but there are also closely associated parameters, such as suspended solids, that can be used as a spectral surrogate in remote sensing studies [19]. Lagerloef et al. [53] presented a remote sensing salinity investigation prior to 1995 with early efforts using the dielectric constant of saline solutions as the physical basis for microwave remote sensing of ocean salinity. Airborne low frequency microwave radiometer experiments have shown success in measuring coastal surface salinity [54]. Linear regression models are commonly used to estimate salinity in estuarine and coastal waters using imagery acquired by Landsat 5 and Landsat 8 satellites [6,19,33,55,56].
SDD is a measure of water transparency, specifically the depth at which the opaque Secchi disk itself loses visibility from the surface. SDD has been successfully measured for decades via satellite and airborne remote sensing methods in a variety of environments including estuaries, rivers, reservoirs, and lakes. The most common method is to use varying techniques of regression analysis for predicting SDD as tested against in situ measurements [20,30,32,39,40].

Total Dissolved Solids (TDS)
TDS are a measure of dissolved inorganic and organic substances in water and is related to electrical conductivity. TDS have been successfully estimated within lacustrine environments using regression techniques based on in situ measurements and reflectance data obtained from spaceborne optical sensors such as Landsat 5 and Landsat 8 [21,57].
1.1.8. Dissolved Oxygen (DO) DO in a system is affected by many factors such as water temperature, respiring and decaying organisms, photosynthesizing plants, stream flow, and aeration. DO is optically inactive, but water temperature is a major influence of DO as oxygen more readily dissolves in cold water than in warm water. No sensor has been established to perform consistently accurate measurements of remotely sensed DO [15]. A Landsat 5 TM-based study conducted by Qiu et al. achieved the best coefficient of determination (R 2 ) value of 0.68 by using a DO inversion exponential model [58].

Dissolved Organic Carbon (DOC)
DOC is an optically inactive water quality parameter and few papers have been published attempting to predict DOC with satellite sensors against in situ data [15]. Mapping dissolved organic carbon for inland waters through remote sensing was shown to be successful when utilizing common band ratio algorithms using Sentinel-2 imagery [36]. DOC transport estimates along a fluvial environment were also found successful using a simple ratio algorithm using SeaWiFS imagery derived ocean color data [59].

pH
It is not feasible to attain a direct correlation between surface reflectance and pH values, because pH is optically inactive [15]. Using RapidEye imagery, Yigit Avdan et al. [60] attempted to develop a model to predict pH by comparing in situ measurements and sensor reflectance values obtained over a lake but were unable to find any correlation among the data. Huang et al. [61] developed a method to indirectly assess pH in water bodies through correlation with other water quality parameters such as salinity, temperature, and dissolved oxygen.

Total Nitrogen (TN) and Total Phosphorus (TP)
TN and TP concentrations in water are optically inactive. However, they have been estimated indirectly by correlating them with optically active chlorophyll a, which was easily determined using Landsat 8 imagery [62]. Older studies attempted to directly assess TP in lacustrine and estuarine environments using remotely sensed imagery (Landsat 1-3 and 5) correlated with in situ data to very limited success [19,32]. An anomalous paper found very high correlation values for both TN and TP using multiple regression analysis with in situ measurements and reflectance data obtained from Landsat 5 sensor [57].

Colored Dissolved Organic Matter (CDOM)
CDOM is a large optically measurable portion of DOC. Inland bodies of water have high correlation values when assessing CDOM via remote sensing using a variety of sensors such as Sentinel-2, Landsat 8, and MODIS [34,36].

Objectives of the Current Study
Since no remote sensing-based algorithm is available to study surface water quality in the Tennessee River and its surrounding water bodies, the Geological and Environmental Remote Sensing Laboratory (GERS-Lab) at the University of Tennessee at Chattanooga recently started developing research programs to study surface water quality in the watersheds of Southeast Tennessee using satellite observations coupled with field measurements. As part of this initiative, this study's objectives were twofold: first, conducting a comprehensive literature review on water quality remote sensing to develop a knowledge base about the state of the art of different water quality parameter estimation techniques using remotely sensed imagery, and second, developing models to estimate turbidity for the Tennessee River and its surrounding tributaries, taking the first step towards satellite observation-based water quality estimation model development for this area.
The literature review conducted in this study shows that turbidity is one of the most widely studied water quality parameters by remote sensing technology. Therefore, the objective of the current study is focused on the development of remote sensing-based numerical models to estimate turbidity for the Tennessee River and adjacent creeks. Landsat satellites are the most frequently used data sources for this purpose; thus, this study uses Landsat 8 Operational Land Imager (OLI) imagery coupled with concurrently acquired field measurements to develop the turbidity estimation models.

Study Site
The Tennessee River is the largest tributary of the Ohio River. It is approximately 1049 km long and is located in the southeastern United States within the Tennessee Valley [63]. Its watershed is the fifth-largest in the nation, encompassing parts of six bordering states as well as Tennessee [64]. The Tennessee River is formed at the confluence of the Holston River and the French Broad River in present-day Knoxville, Tennessee. From Knoxville, it flows southwest through East Tennessee into Chattanooga before crossing into Alabama. The Tennessee River defines the boundary between two of Tennessee's grand divisions: Middle and West Tennessee [63].
The Tennessee River is one of the most ecologically diverse rivers in the world and serves as a great resource for the residents living around it. However, it is also one of the most polluted rivers in the world [64]. In fact, it is ranked as the fourth most polluted waterway in the United States [65]. A leading cause of this pollution is storm water runoff. Runoff from rain and snow that consequently picks up fertilizers, insecticides, and other chemicals, such as phosphates (found in soap used to wash cars), carries these substances into wetlands and waterways, thereby polluting neighboring bodies of water such as the Tennessee River.
Water pollution from runoff subsequently causes concerns of eutrophication in the Tennessee River due to the nutrients that are loaded from urban areas. These excessive nutrients lead to algal blooms and ultimately result in a reduction of dissolved oxygen [66]. This means that nutrients, such as nitrogen and phosphorous compounds needed for plant growth, mix with water in the Tennessee River by surface runoff, resulting in blooms of algal species [67]. Algal blooms are concerning because their subsequent decay can produce an excessive oxygen demand resulting in bad odors, taste, and reduced oxygen levels, thereby harming other aquatic life [67,68]. Therefore, it is vital to explore the main source of pollution, storm water runoff, and its consequences such as eutrophication, the pollution's effects on drinking water quality, and seek solutions to maintain a common resource, i.e., the Tennessee River.
The city of Chattanooga, TN, located along the Tennessee River, has grown substantially during the last several decades and has become the center of a series of urbanized sub-watersheds [69]. Environmental impacts of this growth, especially the quality of surface waters, have become a major concern for the sustainable development of the greater Chattanooga areas [69]. The successful development of a turbidity estimation algorithm would provide a powerful tool to aid in the study of the impacts of land use and land cover change on the surface water quality in the watersheds of Southeast Tennessee.
This study selected the part of the Tennessee River that flows through the city of Chattanooga as the study site. The selected site also includes a section of the South Chickamauga Creek adjacent to the Tennessee River. Figure 2 shows the location and extent of the study site.

Data Collection
This study uses a time series of multispectral satellite imagery and concurrently acquired in situ water quality data. The description of the obtained data and the data acquisition techniques are briefly described below.

Satellite Imagery
The satellite images were obtained from the Landsat 8 mission. The images were acquired over the study site by the Operational Land Imager (OLI) sensor onboard the Landsat 8 satellite on 2 December 2018, 13 February 2019, and 15 August 2019. The imagery products were ordered and obtained through the U.S. Geological Survey (USGS) Earth Explorer. The Landsat 8 OLI Level-2 products, which provide the calculation of surface reflectance by the USGS, were used for spectral analysis. The characteristics of Landsat 8 imagery are provided in Table 1. Figures 3 and S1-S3 show the Landsat 8 OLI images acquired and used for this study.

Water Quality Data
In situ turbidity data were collected in different parts of the study site using a research vessel (a 16-foot-long Jon Boat) and a Hydrolab HL7 multiparameter sonde. This research specifically used the data acquired by the turbidity sensor on the sonde. The HL7 sensors were calibrated in the Geological and Environmental Remote Sensing Laboratory (GERS-Lab) at the University of Tennessee at Chattanooga (UTC) on the morning of field sampling.
The calibration was done using the standard methods described in the Hydrolab operating software manual.
The data were collected concurrently with the Landsat imagery acquisition, which occurred on three dates: 2 December 2018, 13 February 2019, and 15 August 2019. Measurements were acquired at a total of 108 sample points across these three dates, with 35 samples on 2 December, 41 samples on 13 February, and 32 samples on 15 August.
Each sample was collected at similar depths within the water, extending about 20 cm from the surface. An individual sample was the average of 10 recordings for each parameter, collected at 1-second intervals. Parameter readings were stored in the Hydrolab Surveyor. The coordinates for each sample were collected concurrently using a Garmin eTrex 30x GPS unit. Hand-written records of each sample number, coordinates, and parameter values were also recorded in a field book.
Logged data files for the sample parameters were downloaded from the Hydrolab Surveyor using the supplied Hydrolab operating software. Downloaded data were imported into Microsoft Excel and associated with the GPS coordinates for further processing. Figure 4 shows the distribution of sampling locations. Table 2 provides the statistics of the acquired in situ turbidity measurements.

Turbidity Estimation Model Development
The turbidity estimation models developed in this study used an empirical approach based upon simple regressions between in situ measurements and surface reflectance data acquired by the Landsat 8 satellite mission. The statistical analyses were carried out using Microsoft Excel. The geospatial analyses were carried out using ESRI's ArcGIS Pro software (versions 2.1-2.4). The simplicity of using single band regression models could create an avenue for developing routine turbidity estimations for the study site using Landsat 8 OLI satellite imagery. The flowchart in Figure 5 shows the steps taken to build the turbidity model.

Data Pre-Processing
The acquired in situ turbidity measurements for all samples were separated by date, organized in a spreadsheet, and assigned unique Sample ID numbers beginning at and increasing sequentially by 1, with associated latitude and longitude coordinates in decimal degrees utilizing the WGS 1984 geographic coordinate system.
A point GIS data layer (shape file) was created for the locations of each date's in situ turbidity measurements using the 'XY Table to Point' geoprocessing tool. The points are shown in Figure 4.
All seven multispectral bands shown in Table 1 were stacked together to prepare a seven-band multispectral composite image for each image acquired (Figures 3 and S1-S3). The 'Composite Bands' geoprocessing tool was used to do this. Two scenes of Landsat 8 OLI images were available to cover the study site for 2 December 2018 and 15 August 2019. Both scenes were mosaicked using the 'Mosaic to New Raster' geoprocessing tool (Figures 3, S1 and S3). There was only one scene available to cover the study site for 13 February 2019 (Figures 3 and S2); therefore, mosaicking was not required for this date.
The scenes were then assessed and geo-rectified using ground control points. This was done to ensure positional accuracy, particularly along shorelines. Once each scene was correctly rectified, the pixels for land were then clipped to the USGS National Hydrography dataset of waterbodies, using the 'Clip Raster' tool. This removed extraneous pixels from land that would extend processing demands and improperly influence the water quality analysis.
The produced GIS point files for each set of in situ turbidity measurements were displayed on top of the processed multispectral images in true color corresponding to each date, i.e., 2 December 2018; sample points were overlayed upon imagery of the same date. This was done to identify the presence of any samples obscured by cloud coverage and proximity of stream banks. This way, only the measurements available without any influence of cloud and land could be selected for the model development. Due to this reason, three measurements from the December dataset and two measurements from the February dataset were removed from the final analysis. Seven measurements from the August dataset had to be removed since part of the image was covered by scattered clouds, as seen in Figures 3 and S3. Finally, the total number of in situ turbidity measurements used for this study is 96, which includes 32, 39, and 25 measurements for December, February, and August image acquisitions, respectively.

Preparing Model Inputs
This study used the December and February datasets for both model development and testing. The August dataset was used exclusively for testing the model as a separate dataset without any influence in model development. The measurements for the December and February image acquisitions were divided into two sets: training and testing. The dataset used for training included 80% of the measurements whereas the dataset used for testing included 20% of the measurements. This provided 25 and 33 measurements for training the model using December and February image acquisitions, respectively. The dataset used for testing the model included 7 and 6 measurements for December and February image acquisitions, respectively. The selection of the measurements was done randomly using the 'Subset Features' geoprocessing tool available in the Geostatistical Analyst module in ArcGIS Pro software. As we can see in Table 2, the turbidity measurements acquired in December did not have much statistical variation; however, the measurements acquired in February had more variation. Due to this, the datasets used for training the model for December and February image acquisitions were combined for model development. This made the total number of sample points used for training 58. All 25 measurements acquired during August image acquisition were used for further testing the model in addition to the datasets used for testing for December and February image acquisitions. Figure 6 shows the distribution of the sample points used for training and testing for each date. Spectral reflectance values were extracted using the 'Extract Multi-Values to Points' geoprocessing tool at the associated pixel for each sample's coordinate location from each dataset used for training the model. This tool stores spectral reflectance for each spectral band, at each sample site, in the attribute table associated with each respective feature class. In order to conduct statistical analysis, the table for each respective feature class was converted to an Excel spreadsheet using the ' Table to Excel' tool. Tables S1 and S2 show the extracted surface reflectance values for December and February image acquisitions.

Model Development
The extracted surface reflectance values (Tables S1 and S2) from the sample points used for model training were exported to Microsoft Excel (Version 2108, Build 14326.20404) spreadsheets for conducting regression analysis. The 'Table to Excel' geoprocessing tool was used for this task. Single band exponential, power, and linear regressions for the dataset used for training were processed by generating scatter plots with spectral reflectance for each band on the x-axis plotted against corresponding in situ turbidity measurements on the y-axis. These regressions determined the spectral bands used as variables and corresponding equations for turbidity estimation. This also provided a value of coefficient of determination (R 2 ) associated with each best fit trend line used to determine the fit between the calculated regression line and the data. Table 3 shows all the obtained regression equations and corresponding R 2 values. Figure 7 shows the regression scatterplots for all bands (1-7).  The regression equations shown in Table 3, scatterplots shown in Figure 7, and the corresponding R 2 values were carefully assessed and the equation with the best R 2 value was selected to develop the model to estimate turbidity for the study site. The model was developed using the 'ModelBuilder' module available in ArcGIS Pro software. The 'ModelBuilder' module used several geoprocessing tools including 'Make Raster Layer' and 'Raster Calculator'. The 'Raster Calculator' tool used the selected regression equation as algebraic expressions in Python, using the specific spectral bands as input variables (e.g., ρBand value). The output of this model was one raster, for each image acquisition date, which contained predicted turbidity values at each pixel of water within the study site.

Accuracy Assessment
The accuracy of the selected algorithm for turbidity estimation was assessed by using R 2 values and root mean square error (RMSE) values for each date from which the datasets used for model training were obtained. The corresponding test datasets, used for assessment of performance, are shown in Table 4. The model's capability to estimate turbidity during other Landsat 8 image acquisition was evaluated using the data collected on 15 August 2019 ( Table 5). The assessment was done using the obtained R 2 values and RMSE values. The R 2 values were calculated using the best fit line drawn on the scatterplots generated using the observed and predicted values. RMSE was calculated using the following equation: where N = number of samples used for accuracy assessment, x i = observed values, x i = predicted values.

Results
Careful observations of all the regression equations (Table 3), scatterplots (Figure 7), and the corresponding R 2 values suggest that surface reflectance values recorded by Band 1 (Coastal), Band 2 (Blue), Band 3 (Green), and Band 4 (Red) are moderately to strongly correlated with the turbidity in the water of the study site as the majority of the obtained R 2 values range from 0.66 to 0.95. On the other hand, the surface reflectance values recorded by Band 5 (NIR), Band 6 (SWIR1), and Band 7 (SWIR2) are not well correlated with the turbidity in the water of the study site as the obtained R 2 values range from 0 to 0.44. For developing the turbidity estimation model, the nonlinear power regression (Equation (4)) that used Band 4 (Red) was selected since it achieved the highest value of R 2 , 0.95: where T = Turbidity and ρ red (Landsat 8) = Landsat 8 satellite red band surface reflectance values.
The Landsat scenes come packaged with Pixel Quality Assessment (QA) rasters that identify, among other classifications, cloud coverage and shadow, water, and land within a scene through an assigned value. The model-estimated turbidity pixels were extracted using a water mask prepared with the provided Landsat Pixel QA rasters. This water mask removed the influence of either clouds or both clouds and land. For visualization purposes, the masks that removed influence of both clouds and land were further clipped to the Tennessee River study site extent using the Tennessee River shapefile available in the GERS-Lab at UTC. Figures 8 and 9 show the turbidity estimated on 2 December 2018 and 13 February 2019, respectively, by the model.  The model-derived turbidity estimation shows higher values of turbidity concentrated in the two major tributaries in the study site: South Chickamauga Creek and the Hiwassee River. This pattern was observed in both dates. However, the overall turbidity values were higher in February compared to those observed in December. Due to cloud coverage, it was not possible to estimate turbidity for the entire study site for the February image.
Using the test points shown in Figure 6 and Table 4, the predicted turbidity values were extracted from the newly created turbidity raster data for both February and December image acquisitions in order to perform accuracy assessments of the model. This was done by calculating the RMSE and R 2 of predicted values against the observed values. The RMSE values were calculated by Equation 3. The R 2 values were calculated by generating a best fit line on the scatterplots between the in situ turbidity measurements and corresponding estimated turbidity values. The scatterplots are shown in Figure 10. The results are shown in Table 6.  As shown in Figure 10 and Table 6, the turbidity estimation model performed better on the February image acquisition with a R 2 value of 0.97. The model only achieved a R 2 value of 0.31 for the December image acquisition. During image acquisition in February, a higher variation of turbidity values was observed; however, during image acquisition in December, a lower variation of turbidity values was observed ( Table 2). This difference appears to be what is causing inconsistency in model performance. Although the RMSE value for December was lower (1.414 NTU) compared to that of February (11.47 NTU), it can also be attributed to the nature of variation inherent to the acquired in situ turbidity measurements for these dates.
The turbidity estimation model's performance was also assessed on a different date, 15 August 2019. This was done to see how the model performs on a Landsat 8 image for which no data are available for training the model. Equation 4 was applied on the processed surface reflectance values recorded by Band 4 (red) of the image acquired on 15 August 2019. The estimated turbidity image is shown in Figure 11. Using the test points shown in Figure 6 and Table 5, the predicted turbidity values were extracted from the newly created turbidity raster data for August image acquisition to perform accuracy assessments of the model. This was done by calculating RMSE and R 2 of predicted values against the observed values. The RMSE values were calculated by Equation 3. The R 2 values were calculated by generating a best fit line on the scatterplots between the observed test points and corresponding predicted turbidity values. The scatterplot is shown in Figure 12. The results are shown in Table 7.  Estimation of August turbidity provided results with generally high accuracy as evidenced by the obtained R 2 value of 0.94. The obtained RMSE value of 18.08 NTU was relatively higher than expected; however, this could be due to a few outliers present within the test point dataset as seen in Figure 12.

Discussions
Ideally, in situ measurements would be obtained for the entire study site to improve the model accuracy, but this was not possible due to cloud coverage and limited resources. Consequently, this caused the data to be more skewed with less variation. Accuracy was improved by integrating the training data acquired during two field campaigns. However, since the test data were used separately, the difference of data variation influenced the results of accuracy assessments significantly. The test data with higher variation among observed values produced a higher R 2 value and a higher RMSE value. Conversely, the test data with limited variation produced a lower R 2 value and a lower RMSE value. This indicates that the availability of more in situ measurements and concurrent cloud free images would increase the consistency of model performance as it would then account for more variability. Therefore, it is strongly recommended to expand this research by incorporating more in situ measurements and concurrent image acquisitions.
Usually, it is desirable to acquire data for training and testing a model during both base flow and after a storm event to capture the variability of turbidity in stream water. However, due to cloud coverage issues and the image acquisition revisit calendar (Landsat 8 acquires images every 16 days) it is not always possible to achieve this. In this study it was possible to acquire training data after two similar precipitation events. The study site received rainfall in the amount of 20. The study site received a very low amount of rainfall (1.02 mm) on the day before the image acquired (on 15 August 2019) for testing the model's performance for routing turbidity estimation capacity. The recorded value of discharge during the image acquisition time was also very low, at 7.11 m 3 /s. This could be the reason for having higher RMSE values in the test results. So, it is recommended to perform this test again for the images acquired during similar hydrologic conditions comparing to the model's training data. This would also provide us with the opportunity to see how different hydrologic conditions correlate with the remote-sensing-based turbidity estimations.
The majority of the remote-sensing-based inland turbidity studies are conducted in lacustrine environments. Spectral bands within the red and NIR wavelengths are most often chosen, either independently or combined with other bands, for developing turbidity estimation models [17][18][19][20][21][24][25][26][27][28][29][30][31]. This study found that the red band performed the best for developing the model to estimate turbidity in the Tennessee River. The model achieved a R 2 value of 0.95. The models developed using the coastal, blue, and green bands also performed well for turbidity estimation. However, the models developed using the NIR band did not perform well as it did not show a strong relationship with turbidity. The obtained R 2 values were 0.16 to 0.44. This finding aligns with the study conducted by Lobo et al. [27] in the Tapozos River Basin using Landsat MS/TM/OLI images. According to Lobo et al. [27], when water turbidity is high, suspended particles are easily detected by the NIR spectrum. However, higher content of clay particles in water may increase the surface reflectance of the red spectrum compared to sand-rich waters. On the other hand, surface reflectance values should be close to zero in NIR channels for perfectly clear waters. However, in narrow water channels, adjacency of land can erroneously increase NIR surface reflectance, especially in the case for satellites with a low spatial resolution. The visible bands are largely unaffected by this adjacency influence. Therefore, it is suggested to avoid using the NIR band to derive turbidity in narrow fluvial environments for this reason. Since the study site of this research includes narrow fluvial environments, it is possible to have similar issues described by Lobo et al. [27] with NIR bands. More investigations are needed to draw conclusions in this regard.
This research aimed to develop a turbidity estimation model using a single band to simplify the regression equation and to make it suitable for easier implementation of the model for routing operation. However, since other studies [17,18,20,21,27,30] showed success in using multiple bands for turbidity estimation, those methods should also be explored in the future for this area.
Strong winds can generate waves and cause variation in surface water turbidity. However, this is expected not to have any significant influence on the overall surface reflectance values since one pixel of a Landsat 8 image covers about 900 m 2 on the ground. Surface reflectance values recorded by one pixel is the average of the reflectance values recorded over an area of 900 m 2 .

Conclusions
Remote sensing technology has been used successfully for water quality estimation for many years. Turbidity is one of the most frequently used water quality parameters that can be estimated from remote sensing, with Landsat satellites being the most widely used spaceborne technology in this regard.
Currently, no remote-sensing-based water quality estimation model is available for the Tennessee River, albeit it being one of the most ecologically diverse and polluted rivers in the world. This research developed a turbidity estimation model using imagery from the Landsat 8 satellite mission, providing us with the opportunity to estimate turbidity in the Tennessee River every 16 days at a 30 m spatial resolution.
At present, there are only a few scattered water quality measuring stations in the Tennessee River and its surrounding tributaries. With this limited water quality measuring capability, it is very challenging or even impossible for some areas to adequately study the water quality in this area. The model developed by this study using Landsat 8 satellite imagery would provide with the turbidity estimation capability in the Tennessee River and its surrounding tributaries at every 900 m 2 surface area. This can be acquired whenever a cloud-free Landsat 8 OLI image is available.
As the city of Chattanooga, TN is being continuously urbanized, the environmental impacts, especially the quality of surface waters due to this growth, have become a major concern for the sustainable development of the greater Chattanooga area [69]. This remotesensing-based turbidity estimation model would be very useful to monitor turbidity in the Tennessee River and its surrounding tributaries throughout the year.
The results obtained in this study indicate that a simple nonlinear regression-based model could be developed using the surface reflectance values recorded in Band 4 (red) of the Landsat 8 OLI sensor to estimate turbidity in the Tennessee River and its tributaries. This is strongly supported by the obtained R 2 values for model development (Table 3 and Figure 7). The accuracy assessments of the estimated turbidity images (Figures 8 and 9) also indicate promising prospects in this regard, as the highest R 2 value was 0.97 and the lowest RMSE was 1.414 NTU ( Table 6).
The model was also tested to evaluate its potential for routine operational use for turbidity estimation in the Tennessee River and its tributaries. The test results were promising as the obtained R 2 and RMSE values were 0.94 and 18.01 NTU, respectively. However, since the model was only tested on one image acquisition date, 15 August 2019, it would be beneficial to test this model on data acquired from other image acquisition dates with near real-time in situ turbidity measurements.

Supplementary Materials:
The following are available online at https://www.mdpi.com/article/ 10.3390/rs13183785/s1, Figure S1: One of three Landsat 8 OLI images acquired over the study site on 2 December 2018. The image is subset to the study site extent displayed in true color. The inset regions represent areas of high turbidity variability due to the convergence of lakes and minor streams with the Tennessee River, Figure S2: One of three Landsat 8 OLI images acquired over the study site on 13 February 2019. The image is subset to the study site extent displayed in true color. The inset regions represent areas of high turbidity variability due to the convergence of lakes and minor streams with the Tennessee River, Figure S3: One of three Landsat 8 OLI images acquired over the study site on 15 August 2019. The image is subset to study site extent displayed in true color. The inset regions represent areas of high turbidity variability due to the convergence of lakes and minor streams with the Tennessee River, Table S1: Model training data for 2 December 2018 image acquisition, Table S2  Informed Consent Statement: This is not relevant to this study. Data Availability Statement: Some of the data are already provided in the supplementary materials section of this publication. Other data presented in this study are available on request from the corresponding author.