Confidence Levels, Sensitivity, and the Role of Bathymetry in Coral Reef Remote Sensing

Remote sensing is playing an increasingly important role in the monitoring and management of coastal regions, coral reefs, inland lakes, waterways, and other shallow aquatic environments. Ongoing advances in algorithm development, sensor technology, computing capabilities, and data availability are continuing to improve our ability to accurately derive information on water properties, water depth, benthic habitat composition, and ecosystem health. However, given the physical complexity and inherent variability of the aquatic environment, most of the remote sensing models used to address these challenges require localized input parameters to be effective and are thereby limited in geographic scope. Additionally, since the parameters in these models are interconnected, particularly with respect to bathymetry, errors in deriving one parameter can significantly impact the accuracy of other derived parameters and products. This study utilizes hyperspectral data acquired in Hawaii in 2000–2001 and 2017–2018 using NASA’s Classic Airborne Visible/Infrared Imaging Spectrometer to evaluate performance and sensitivity of a well-established semi-analytical inversion model used in the assessment of coral reefs. Analysis is performed at several modeled spatial resolutions to emulate characteristics of different feasible moderate resolution hyperspectral satellites, and data processing is approached with the objective of developing a generalized, scalable, automated workflow. Accuracy of derived water depth is evaluated using bathymetric lidar data, which serves to both validate model performance and underscore the importance of image quality on achieving optimal model output. Data are then used to perform a sensitivity analysis and develop confidence levels for model validity and accuracy. Analysis indicates that derived benthic reflectance is most sensitive to errors in bathymetry at shallower depths, yet remains significant at all depths. The confidence levels provide a first-order method for internal quality assessment to determine the physical extent of where and to what degree model output is considered valid. Consistent results were found across different study sites and different spatial resolutions, confirming the suitability of the model for deriving water depth in complex coral reef environments, and expanding our ability to achieve automated widespread mapping and monitoring of global coral reefs.


Introduction
Deriving bathymetry from remote sensing data has long played an important role in the study of coastal regions, coral reefs, inland lakes, waterways, and other aquatic areas. Bathymetry has of multiple data sets. Additionally, the great effort added by expert analysis unfortunately often injects user bias into the analysis, and there are also practical difficulties associated with collecting the required numbers of statistically significant validation points that are independent of the training points used for developing the map [33]. In a similar fashion, many studies deriving bathymetry from multispectral data still rely on two-band ratio methods [4,5], even when incorporating new methods (e.g., machine learning) [34] and new satellites [35]. While these ratio methods can be highly accurate, they are based on constant water column properties throughout the image and output is thus highly localized. As such, even with all the recent progress, automating the classification of coral reefs and extracting bathymetry at regional or global scales remains difficult, and commonly requires significant ancillary data and/or expert knowledge to be incorporated into the image analysis workflow [20].
Hyperspectral remote sensing is being increasingly used in efforts to develop more broadly applicable analysis techniques that can be robustly applied across diverse areas with little to no ancillary information. Studies using hyperspectral data to develop classification and bathymetric maps for optically shallow coral reefs predominantly use semi-analytical inversion methods, employing either numerical optimization modeling or look-up tables to derive bottom reflectance and depth [15,[36][37][38][39][40]. Whereas optimization models can have high computational demands, and whereas look-up tables are reliant on predetermined surface reflectance spectra generated from a select set of benthic reflectance spectra, water depths and water properties, both analysis approaches are sensitive to errors in derived bathymetry. More specifically, if the assumed water depths are incorrect, these errors can have significant impact on the accuracy of the subsequent bottom reflectance and benthic classification steps. Therefore, it is important to validate derived water depth, as well as understand how bathymetry impacts overall model capabilities and limitations.
To pursue these application objectives, this study utilizes hyperspectral data from Hawaii, acquired using NASA's Classic Airborne Visible/Infrared Imaging Spectrometer (AVIRIS-C), to evaluate performance and sensitivity of a representative semi-analytical inversion model that simultaneously derives water depth, water properties, and benthic surface reflectance given only remote sensing-based water surface reflectance as input. Previous studies of different inversion methods applied to hyperspectral imagery have proved to be relatively accurate, with higher accuracy gained when local parameters are incorporated in the modeling process [38]. However, uncertainties regarding how a given algorithm applies to a specific geographic area, and how the algorithm might apply to other areas, are rarely discussed. As such, in the absence of prior knowledge and ancillary data, decision tools are needed to automatically determine where a model is most valid and where it is not (e.g., identifying areas that are too deep or too turbid to obtain viable results). This study uses lidar bathymetry to validate the accuracy of hyperspectral derived water depth for the selected inversion model, and is unique in using the results to perform a sensitivity analysis and to develop confidence levels for model output. Furthermore, the methods described were developed in the interest of widespread application, and discussion is focused on implications for automated large-scale coral reef analysis.

Materials and Methods
This study employs a semi-analytical hyperspectral inversion model to gain insight into and assess impact of derived bathymetry on the retrieval of coral reef habitat information in shallow tropical waters. Airborne data from the AVIRIS-C instrument is used as a proxy for hyperspectral satellite data at varying spatial resolutions, and bathymetric lidar data is utilized as the measured ground truth for evaluating model accuracy. The methodology workflow includes a series of preprocessing steps, model implementation, bathymetry validation, development of confidence levels, and a sensitivity analysis to assess how errors in derived bathymetry translate to uncertainty in derived benthic reflectance.

Study Areas
Three separate coral reef locations in the Hawaiian Islands were utilized as study areas in this investigation (Figure 1), each representing differing reef types, environmental conditions, and bathymetric variability: Kaneohe Bay, a partially enclosed~50 km 2 embayment on the northeast coast of Oahu, contains a central barrier reef, fringing reefs, an extensive lagoon/bay area with numerous patch reefs, and a seaward fore reef. Depths in this area range from <1 m on the reef flats to an average of 10-15 m in the lagoon, with rapidly increasing depth offshore. Water clarity, which impacts the ability to derive bathymetry and benthic information using remote sensing, is typically highest in the northern and more exposed seaward portions of the bay, and lowest along the shore, in the enclosed southern portion of the bay, and following significant rainfall-runoff events. 2.
The south coast of Molokai, specifically a~20 km 2 section of coastal area immediately adjacent to Kaunakakai Harbor, consists of a 1-2 km wide fringing reef, with a well-defined reef crest, shallow shoreward reef flat, and seaward fore reef. Depths are <1-2 m across the reef flat out to the reef crest, transitioning to 30-50 m within 0.5 km seaward of the reef crest, and then deepening quickly offshore. With the exceptions of portions of the harbor itself, where water residence time is higher, and following rainfall-runoff events, the majority of this area is regulated by adjacent oceanic water conditions and water clarity is generally high.

3.
French Frigate Shoals, specifically a~100 km 2 subset in the southeast corner of this rather sizeable 900 km 2 coral atoll in the remote Northwestern Hawaiian Islands, consists of a lengthy~35 km crescent shaped reef and reef flats, a semi-enclosed lagoon, a scattering of small islands, and an extensive seaward fore reef. Depths are <1 m at the reef crest, 5-10 m in the eastern and northern portions of the lagoon, transitioning to 15-20 m in the center of the lagoon and 25-30 m in the open western area of the lagoon, and rapidly deepening along the fore reef and offshore. While water clarity along the reef flat and seaward fore reef is regulated by relatively clear oceanic water, areas inside the lagoon, particularly along its protected eastern end, where sediment suspension is generally greater, have higher turbidity and reduced water clarity.
Remote Sens. 2020, 12, 496 4 of 28 Three separate coral reef locations in the Hawaiian Islands were utilized as study areas in this investigation (Figure 1), each representing differing reef types, environmental conditions, and bathymetric variability: 1. Kaneohe Bay, a partially enclosed ~50 km 2 embayment on the northeast coast of Oahu, contains a central barrier reef, fringing reefs, an extensive lagoon/bay area with numerous patch reefs, and a seaward fore reef. Depths in this area range from <1 m on the reef flats to an average of 10-15 m in the lagoon, with rapidly increasing depth offshore. Water clarity, which impacts the ability to derive bathymetry and benthic information using remote sensing, is typically highest in the northern and more exposed seaward portions of the bay, and lowest along the shore, in the enclosed southern portion of the bay, and following significant rainfall-runoff events. 2. The south coast of Molokai, specifically a ~20 km 2 section of coastal area immediately adjacent to Kaunakakai Harbor, consists of a 1-2 km wide fringing reef, with a well-defined reef crest, shallow shoreward reef flat, and seaward fore reef. Depths are <1-2 m across the reef flat out to the reef crest, transitioning to 30-50 m within 0.5 km seaward of the reef crest, and then deepening quickly offshore. With the exceptions of portions of the harbor itself, where water residence time is higher, and following rainfall-runoff events, the majority of this area is regulated by adjacent oceanic water conditions and water clarity is generally high. 3. French Frigate Shoals, specifically a ~100 km 2 subset in the southeast corner of this rather sizeable 900 km 2 coral atoll in the remote Northwestern Hawaiian Islands, consists of a lengthy ~35 km crescent shaped reef and reef flats, a semi-enclosed lagoon, a scattering of small islands, and an extensive seaward fore reef. Depths are <1 m at the reef crest, 5-10 m in the eastern and northern portions of the lagoon, transitioning to 15-20 m in the center of the lagoon and 25-30 m in the open western area of the lagoon, and rapidly deepening along the fore reef and offshore. While water clarity along the reef flat and seaward fore reef is regulated by relatively clear oceanic water, areas inside the lagoon, particularly along its protected eastern end, where sediment suspension is generally greater, have higher turbidity and reduced water clarity.   This project utilizes a unique time-series of airborne hyperspectral data collected over coral reefs in the Hawaiian Islands over 15 years apart between 2000-2001 and 2017-2018. While hyperspectral remote sensing has been broadly utilized for more than three decades, and is becoming increasingly prevalent across a host of application domains, this time-series is one of the only coral reef hyperspectral datasets spanning such a lengthy time period, thus representing an invaluable resource for exploring hyperspectral remote sensing of coral reefs through space and time. Data was acquired using the AVIRIS-C instrument, which collects 224 contiguous wavelength bands from~400-2500 nm, and is typically flown onboard NASA's DHC-6 Twin Otter or Lockheed ER-2 aircraft, providing both low-and high-altitude deployment options, respectively. As was the case for the data in this study, acquisitions were performed using the ER-2, providing data at either~17 or~7.5 m spatial resolution as a function of aircraft altitude. The specific datasets used in the analysis are listed in Table 1.

Lidar Bathymetry Data
Bathymetric lidar data was used as ground truth for validating hyperspectral derived bathymetry in Kaneohe Bay and Molokai. Two different datasets from 2000 and 2013 were utilized in this process. Lidar bathymetry was collected over Kaneohe Bay and the south coast of Molokai in 2000 by the US Army Corps of Engineers (USACE) using the Scanning Hydrographic Operational Airborne Lidar Survey (SHOALS) [41][42][43][44]. SHOALS point data (horizontal accuracy ± 3 m; vertical accuracy ± 0.15 m [45]) was downloaded from the University of Hawaii's Coastal Geography Group website (http://www. soest.hawaii.edu/coasts/data/) for Oahu (Zones 5-7) and Molokai (Zone 12) as ASCII text files. Another lidar data collection, which had greater geographic coverage in each study area, was performed in the main Hawaiian Islands in 2013 by the USACE using the Coastal Zone Mapping and Imaging Lidar (CZMIL) [46,47]. CZMIL point data (horizontal accuracy ± 3.5 + 0.05 * depth m; vertical accuracy ± sqrt(0.25ˆ2 + (0.0075 * depth)ˆ2) m; which, for reference, equates to ±4 m horizontal and ±0.26 m vertical at 10 m depth [48,49]) were downloaded for both Kaneohe Bay and the south coast of Molokai from the NOAA Data Access Viewer (https://coast.noaa.gov/dataviewer/) in .las format.

Hyperspectral Data Preprocessing
To prepare data for analysis, a series of preprocessing steps were utilized to transform AVIRIS-C radiance data to georeferenced surface reflectance at three simulated spatial resolutions. As discussed below, included in these steps were geocorrection, atmospheric correction, spatial resampling, masking, Remote Sens. 2020, 12, 496 6 of 29 and sunglint correction. Where feasible, in the interest of maintaining a standardized workflow, preprocessing included output from NASA's standard operational AVIRIS-C algorithms. This included atmospheric correction for all data, as well as geocorrection and spatial resampling for the 2017 and 2018 data. For the historic 2000 and 2001 data, however, geocorrection and spatial resampling were performed independently as part of our analysis, as were masking and sunglint correction for all 2000, 2001, 2017, and 2018 datasets.

Geocorrection
The modern AVIRIS-C data collection protocol includes acquisition of high accuracy position and attitude information from an onboard GPS/INS system, which enables automated geocorrection using a standardized geocoding algorithm [50]. By integrating the GPS/INS data with a detailed camera model, the algorithm uses ray tracing to explicitly derive locations for each AVIRIS-C pixel. For the 2017 and 2018 data we utilized the geocorrection output, as derived using this model and provided by NASA, with no further modifications. While NASA does not provide specifications for accuracy of the geocorrection process, indications are that geocoding errors are typically on the subpixel scale [51,52]. For the historic 2000 and 2001 data, insufficient ancillary data is available to implement this same algorithm. Hence, geocorrection of this data was performed in ENVI (version 5.5) by co-registering the 2000 and 2001 data to the associated 2017 data for each study area using manually selected ground control points, triangulation warping (to accommodate for the irregular spatial variations of airborne data), and nearest neighbor resampling (to maintain spectral fidelity) (root mean square error: Kaneohe Bay, rmse = 3.3; Molokai, rmse = 0.3; French Frigate Shoals, rmse = 2.3). Following geocorrection, all data for each study area were regridded to a common projection and spatial extent, again using nearest neighbor resampling, to facilitate direct one-to-one spatial comparison and analysis.

Atmospheric Correction
Atmospheric correction is an important step in most remote sensing applications to minimize influence of the atmosphere on the reflectance from the observed surface. For aquatic remote sensing, where the signal strength is typically much lower than terrestrial applications, atmospheric correction is particularly important since even minor errors or improper corrections can result in significant impacts on results. While a number of different atmospheric correction routines have been utilized in coral reef remote sensing [38,39,[53][54][55], in the interest of maintaining a standardized preprocessing approach, this study utilized atmospheric correction output generated by NASA using the most recent AVIRIS-C data processing workflow. This process [56] derives scaled surface reflectance from calibrated AVIRIS-C radiance data using a modified version of ATREM [57]. For those familiar with previous versions of this algorithm, which were developed primarily for terrestrial applications, the current version has been optimized for application to both terrestrial and aquatic target areas. It is also the same algorithm as being used for processing data from NASA's Coral Reef Airborne Laboratory mission (https://coral.jpl.nasa.gov/). To validate the algorithm, derived surface reflectance for the 2000 Kaneohe Bay data was compared against output from an earlier study [53] in order to confirm that results are consistent with another successful atmospheric correction routine using the same data.

Spatial Resampling
To emulate the response of different feasible moderate resolution hyperspectral satellites, spatial resampling was used to generate output for all datasets at 18, 30 and 60 m spatial resolution. As with the geocorrection, spatial resampling was implemented differently for the historic 2000 and 2001 datasets as compared with the 2017 and 2018 data. For the recent data, resampling was performed as an integrated process within NASA's geocorrection software, where the 18 m data was generated using nearest neighbor resampling and the 30 and 60 m data were generated using a Gaussian function to average all pixels contributing to the coarser resolution output pixel [58]. For the 2000 and 2001 data, where resampling was performed independently, the 18 m data was similarly generated using nearest Remote Sens. 2020, 12, 496 7 of 29 neighbor resampling, whereas the 30 and 60 m data were generated using a linear average to aggregate the contributing pixels. While it is acknowledged that differences in the resampling approaches may contribute uncertainties in subsequent analysis, it is expected that the impacts are minor as compared with uncertainty introduced by variabilities in the water column and water surface.

Masking
Since this project is focused on only the water portions of each study area, steps were taken to mask out areas of each image containing land, clouds, cloud shadows, wave breaks, and even boats. To a large extent masking could be achieved using a combination of standard band indices and band thresholds (e.g., NDVI, NDWI, NIR); however, it was generally determined that no single algorithm, or suite of algorithms, could adequately handle all masking needs, particularly for thin clouds and their associated shadows over water. Hence, a certain amount of manual editing was required to generate comprehensive masks for each image at each spatial resolution.

Sunglint Correction
Specular reflection from the water surface (i.e., sunglint) occurs when sunlight is reflected from the surface at the same angle at which it is being viewed. Unless avoided during data acquisition through control of view angles (e.g., flight direction, time of day, pointing), sunglint can occur at varying levels of magnitude depending on proximity of view angle to the reflected angle. When present in an image, sunglint can obscure our ability to derive information about the water column and benthic surface, and must be minimized where possible during both the image acquisition and analysis process. For the AVIRIS-C data used in this project, the 2017 and 2018 data were specifically acquired to avoid sunglint; however, this was not the case for the 2000 and 2001 data, hence requiring that a sunglint correction scheme be included in the preprocessing workflow. As such, to maintain consistency in the overall workflow regardless of the presence/absence of sunglint, and because sunglint correction can also serve to normalize surface reflectance within and between images, sunglint correction was applied to all data.
The sunglint correction algorithm utilized in this study calculates a per-pixel spectral offset as a function of the pixel's reflectance at 640 and 750 nm, where reflectance is set to zero in optically deep water, but allowed to be greater than zero in optically shallow areas (e.g., shallow sand) [53,59]. Using reflectance output from the atmospheric correction, R*(λ), which is transformed to remote sensing reflectance, R rs *(λ), as R*(λ)/π, the sunglint corrected surface remote sensing reflectance, R rs (λ), is calculated as follows: where λ is the wavelength in nm.

Lidar Data Preprocessing
Lidar preprocessing was used to transform the native lidar point cloud data into rasterized grids consistent with the various hyperspectral datasets for Kaneohe Bay and Molokai (i.e., same geographic extent, pixel locations, and spatial resolution). Subsequent analysis was used to cross-compare the bathymetric rasters from 2000 and 2013 to establish that there were not any significant changes in bathymetry at the study areas over this time period.

Rasterization
To rasterize the bathymetric point data, SHOALS data files were imported into R (version 3.2.0) and then converted to raster images matching the same spatial extent and resolution of each of the different hyperspectral data files for Kaneohe Bay and Molokai. Raster image pixel locations and resolutions were set to the same geographic specifications as the hyperspectral images, and output Remote Sens. 2020, 12, 496 8 of 29 pixel values were calculated as the mean of all depth values for the lidar returns located within each given pixel extent. To ensure only measured data values were used in the output, no extrapolations nor interpolations were used to fill any gaps in the coverage. CZMIL data were also similarly imported into R and rasterized in the same manner as the SHOALS data.

Temporal Comparison
The bathymetric lidar raster images from 2000 and 2013 were statistically compared using a linear correlation analysis to assess consistency in water depth over the given measurement period. This analysis revealed strong correlations between the two datasets for the depth range of 0-30 m (Kaneohe Bay: r = 0.998, mean difference = 0.28 m; Molokai: r = 0.997, mean difference = 0.17 m), and hence indicated that bathymetry remained constant in these study areas over this period. By assuming that bathymetry would similarly continue to remain constant through 2018, this preprocessing step was used to establish that, of the two lidar data sets, the 2013 data dataset, which has greater geographic coverage, could be used alone in the subsequent validation analysis for all years of hyperspectral data acquisition.

Hyperspectral Data Analysis
Data processing and analysis employed a semi-analytical inversion model to retrieve water depth, water properties and benthic reflectance, followed by validation of the bathymetry output, development of output confidence levels, and implementation of a sensitivity analysis to evaluate how variations in bathymetry impact the derivation of benthic reflectance. Considering that such models are often utilized as the foundation, or fundamental component, of different methodologies for habitat assessment and monitoring, understanding both model accuracy and sensitivity represents an important aspect in establishing both capabilities and limitations in the subsequent habitat analysis.

Inversion Model
The semi-analytical hyperspectral inversion model used in this study [53,59], as modified from its original formulation [50,51], and with various other adaptations reported elsewhere [15,37,39,40], simultaneously derives estimates of water depth, water properties, and benthic reflectance given only surface reflectance as input. Since the model does not require any local input parameters or other a priori information as input, it is broadly applicable without customization to most shallow coral reef environments.
The model is governed by the generalized relationship that surface remote sensing reflectance, R rs (λ), can be approximated as a function of five parameters: where P (m −1 ) is the phytoplankton absorption coefficient at 440 nm, G (m −1 ) is the absorption coefficient for gelbstoff and detritus at 440 nm, BP (m −1 ) is the backscattering coefficient of suspended particles at 440 nm, B is the bottom reflectance at 550 nm, and H (m) is the water depth. The two primary equations that then define the semi-analytical model based on these parameters are as follows (with the wavelength term, λ, dropped for notational convenience): r rs ≈ r rs dp 1 − exp −κH D u C + 1 cos(θ w )

Water Column Contribution
where r rs (sr −1 ) is the subsurface remote sensing reflectance, r rs dp (sr −1 ) is the subsurface remote sensing reflectance for optically deep water, θ w (rad) is the subsurface solar zenith angle at the time of image acquisition, D u C and D u B are optical path-elongation factors for scattered photons from the water column and bottom, respectively, κ (m −1 ) is an attenuation coefficient, and r b is a generic bottom spectrum normalized to 1.0 at 550 nm. The variables r rs dp , D u C , D u B , and κ are in turn determined as a function of wavelength, λ, coefficients for backscattering, b b (m −1 ), and absorption, a (m −1 ), and the parameters P, G, and BP (refer to [53,[59][60][61] for more details). Output is derived on a pixel-by-pixel basis using a nonlinear optimization routine to iteratively adjust the five model parameters (P, G, BP, B, and H) such that the root mean square error (rmse) between estimated and measured surface reflectance is minimized. The end result is a per-pixel estimation of water depth, water properties, and benthic reflectance at 550 nm throughout a given image.

Lidar Validation
Using a raster-to-raster image analysis, linear correlation coefficients and mean differences of derived bathymetry values versus the measured lidar data for depths 0-15 m were used to statistically assess the inversion model bathymetry output. This was performed for both Kaneohe Bay and Molokai at each of the three different spatial resolutions (18,30, and 60 m). As the first step in the process, since the lidar data values were reported with respect to the Mean Lower Low Water (MLLW) tidal datum, derived hyperspectral bathymetry values were tide corrected prior to performing the comparison (https://tidesandcurrents.noaa.gov/; Kaneohe Bay: Mokuoloe, HI-Station #1612480, verified, accuracy ± 0.02 m; Molokai: Kaunakakai Harbor, HI-Station #1613198, predicted, approximate accuracy ± 0.1 m). Additionally, the hyperspectral output was masked to remove pixels where the optimization routine did not sufficiently converge (i.e., eliminate any additional areas of land, cloud, and wave breaks not already captured by the mask process). The comparison was then performed by evaluating the linear correlation coefficient and mean difference between the derived and measured bathymetry values.

Temporal Validation
As a secondary metric for bathymetry validation, linear correlation coefficients and mean differences were again used as statistical comparisons between the different years of image acquisition for each study area. Based on the comparison of the 2000 and 2013 lidar data, which indicated bathymetry remained constant in these study areas for this time period, it was assumed that bathymetry would remain similarly unchanged through 2018, thereby implying that year-to-year comparisons of derived bathymetry should also remain constant. Accordingly, linear correlation analysis was performed between each of the different image dates to evaluate whether derived bathymetry remained consistent between data sets in different acquisition years.

Confidence Levels
A common challenge when scaling remote sensing analysis from local to regional scales, and beyond, is maintaining model integrity across varying environmental conditions and application domains [21,62,63]. With airborne and satellite derived bathymetry this includes the potential to encounter water conditions and depths that exceed the capabilities of a given modeling approach. While localized knowledge could be used to constrain models to optimal areas, in the interest of broad-scale application it is beneficial to automate this process. With this objective in mind, this study utilized maximum surface reflectance at 600 nm as a threshold to establish confidence levels in model output. Image analysis and exploration found that reflectance magnitude at 600 nm best correlated to model accuracy and hence provided a reliable first-order indicator of the model's capacity for deriving bathymetry. In contrast, reflectance at shorter visible wavelengths were less reliable at predicting accuracy, especially in optically deep water, and reflectance at longer wavelengths tended toward zero in both optically deep and shallow waters. As such, remote sensing reflectance, sr −1 , for each pixel at 600 nm was used to establish five equally distributed confidence levels based on the average reflectance values observed in the imagery used in this study, listed here from highest confidence to lowest: (i) Rrs > 0.01; (ii) Rrs > 0.0075; (iii) Rrs > 0.005; (iv) Rrs > 0.0025; and (v) Rrs < 0.0025.

Sensitivity Analysis
As a means to assess model sensitivity, the measured bathymetric lidar data was used as input to the inversion model, rather than deriving bathymetry as an output. This served to more narrowly constrain the inversion optimization by transforming water depth from an unknown to a known, while still allowing water properties and benthic reflectance to vary in the inversion process. The modified model was then iteratively run using discrete variations in tide-corrected bathymetry input as a means to assess the resulting impacts on derived benthic reflectance values. It is noted that the model could be further constrained using a fixed set of assumed water properties, but since the objective was to examine the overall impacts of varying water depth, water properties were retained as an unknown variable. The resulting analysis process was repeated in discrete 5% increments such that input water depth varied from 50% to 150% of the measured water depth, resulting in a total of 21 different model runs at each spatial resolution for each hyperspectral dataset in Kaneohe Bay and Molokai. The mean standard deviations of derived benthic reflectance as a function of varying input depths were used to assess sensitivity of the model, where high standard deviation results indicated derived benthic reflectance is sensitive to errors in water depth, and low standard deviation results indicated derived benthic reflectance has limited sensitivity to errors in water depth.

Inversion Model Output
The two primary inversion model outputs of concern in this analysis were the derived parameters for water depth and benthic reflectance. When considering each of the different images and each of the three different spatial resolutions, this equated to results from 27 different model runs, 12 for Kaneohe Bay, nine for Molokai, and six for French Frigate Shoals. Example output from the 2000 Kaneohe Bay data is shown in Figure 2.

Sensitivity Analysis
As a means to assess model sensitivity, the measured bathymetric lidar data was used as input to the inversion model, rather than deriving bathymetry as an output. This served to more narrowly constrain the inversion optimization by transforming water depth from an unknown to a known, while still allowing water properties and benthic reflectance to vary in the inversion process. The modified model was then iteratively run using discrete variations in tide-corrected bathymetry input as a means to assess the resulting impacts on derived benthic reflectance values. It is noted that the model could be further constrained using a fixed set of assumed water properties, but since the objective was to examine the overall impacts of varying water depth, water properties were retained as an unknown variable. The resulting analysis process was repeated in discrete 5% increments such that input water depth varied from 50% to 150% of the measured water depth, resulting in a total of 21 different model runs at each spatial resolution for each hyperspectral dataset in Kaneohe Bay and Molokai. The mean standard deviations of derived benthic reflectance as a function of varying input depths were used to assess sensitivity of the model, where high standard deviation results indicated derived benthic reflectance is sensitive to errors in water depth, and low standard deviation results indicated derived benthic reflectance has limited sensitivity to errors in water depth.

Inversion Model Output
The two primary inversion model outputs of concern in this analysis were the derived parameters for water depth and benthic reflectance. When considering each of the different images and each of the three different spatial resolutions, this equated to results from 27 different model runs, 12 for Kaneohe Bay, nine for Molokai, and six for French Frigate Shoals. Example output from the 2000 Kaneohe Bay data is shown in Figure 2. Output data was first visualized for illustrative purposes and as an important basis for preliminary evaluations of overall data quality and model performance. To do so, each of the different Output data was first visualized for illustrative purposes and as an important basis for preliminary evaluations of overall data quality and model performance. To do so, each of the different input images and associated model outputs were visually compared with existing maps and field knowledge prior to performing any statistical analytics. While the majority of output was consistent with known field characteristics, it was observed that the 2018 output (and to some degree the 2017 output) for Kaneohe Bay and Molokai exhibited lower than expected performance at moderate to deep depths. More specifically, when visually compared with data from the other years of acquisition, this data indicated reduced accuracy in areas deeper than~5 m. Given that the inversion model has been successfully utilized in multiple other studies at depths exceeding this limit [15,[37][38][39][40]53,59,[64][65][66], it was hypothesized that these results were a function of suboptimal illumination conditions (i.e., lower total solar radiation due to time of year and lower solar angle due to time of acquisition; Table 1) and/or environmental conditions (i.e., reduced water clarity and increased surface roughness) at the time of acquisition. This reduced accuracy is explored further throughout the analysis.

Lidar Validation
The 2013 rasterized lidar data was used to assess accuracy of water depths derived from the Kaneohe Bay and Molokai hyperspectral data. Figure 3 provides examples of how the model derived output for Kaneohe Bay compares with the 2013 lidar data at 18 m spatial resolution. Based on an initial visual assessment, the model derived water depth compares favorably with the lidar bathymetry, exhibiting similar geographic extents amongst different features, such as the shallow fringing reefs, patch reefs and reef flat (shown in red), the north and south channels (shown in green), and the offshore transition from shallow (shown in orange/yellow) to deep water (shown in blue/purple). Notable exceptions to this consistency occur predominantly in the lagoon (light blue in lidar data versus green/yellow/purple in the 2000 output, green/yellow/orange in the 2017/2018a output, and orange/red in the 2018b output), which is presumed to be a result of persistently poor water clarity in these areas (based on field observations). For the lidar data, this reduced water clarity impedes the ability to collect accurate lidar readings, hence the reason there are large areas in the lagoon and in the channel at the northern end of the bay with no available lidar data (white areas within the extent of measured areas in Figure 3a). In contrast, the hyperspectral model will generate output for these same areas regardless of water conditions, which illustrates the need to develop confidence levels for improved interpretation of output and more robust model applicability. It is also apparent from these results that there is a significant decline in performance for the 2018 data, particularly 2018b, in the deeper areas of the study area (lagoon, channels, and offshore regions). As mentioned above, this is attributed to illumination and/or environmental conditions at the time of acquisition and discussed in more detail in the statistical analysis that follows.
Scatterplots of lidar data compared to the hyperspectral tide-corrected Kaneohe Bay and Molokai bathymetric data from 0-15 m water depth at 18 m spatial resolution are presented in Figure 4, and statistical comparisons are presented in Table 2. Data deeper than 15 m are not included due to the lower amount of available lidar data at deeper depths, which results from the steep offshore drop-off in both study areas, particularly for Molokai, and not as an indication of any model limits. In terms of interpreting the linear correlation analysis, and thereby evaluating the capability for deriving absolute water depth, stronger relationships are indicated by correlation coefficients, r, closer to 1 and mean differences, m, closer to 0 m. However, results are not expected to be perfect, in part because of uncertainties in the model process itself, but also due to errors introduced by geocorrection offsets, lidar inaccuracies, and tidal corrections. Additionally, while a high correlation coefficient indicates a strong predictive relationship for water depth alone, the capacity for predicting actual water depths is equally important, as this can impact the derivation of other interdependent variables for water properties and benthic reflectance. When assessing these scatterplots it is therefore equally important to consider the statistical correlations, as it is difficult to visually assess point density, and hence interpret dominant data trends, in the scatterplots alone. Further, while Kaneohe Bay has varied bathymetric characteristics, Molokai is dominated by a shallow fringing reef flat with a steep drop-off that exceeds modeling limits. The Molokai data thus need only be accurate for depths down to 3-5 m for the majority analysis to be effective, whereas Kaneohe Bay should be accurate to at least 15 m depth. With that in mind, from these plots and the associated statistics it is evident that the 2000, 2001, and 2017 data provide the best overall results, while the 2018 data has the lowest performance, particularly at depths exceeding 5 m. As surmised in the initial visual assessment of the output, and attributed to suboptimal illumination and/or environmental conditions, it is not unexpected that reported accuracy for the 2018 data would be lower than the data from other years.
fringing reefs, patch reefs and reef flat (shown in red), the north and south channels (shown in green), and the offshore transition from shallow (shown in orange/yellow) to deep water (shown in blue/purple). Notable exceptions to this consistency occur predominantly in the lagoon (light blue in lidar data versus green/yellow/purple in the 2000 output, green/yellow/orange in the 2017/2018a output, and orange/red in the 2018b output), which is presumed to be a result of persistently poor water clarity in these areas (based on field observations). For the lidar data, this reduced water clarity impedes the ability to collect accurate lidar readings, hence the reason there are large areas in the lagoon and in the channel at the northern end of the bay with no available lidar data (white areas within the extent of measured areas in Figure 3a). In contrast, the hyperspectral model will generate output for these same areas regardless of water conditions, which illustrates the need to develop confidence levels for improved interpretation of output and more robust model applicability. It is also apparent from these results that there is a significant decline in performance for the 2018 data, particularly 2018b, in the deeper areas of the study area (lagoon, channels, and offshore regions). As mentioned above, this is attributed to illumination and/or environmental conditions at the time of acquisition and discussed in more detail in the statistical analysis that follows.  Scatterplots of lidar data compared to the hyperspectral tide-corrected Kaneohe Bay and Molokai bathymetric data from 0-15 m water depth at 18 m spatial resolution are presented in Figure  4, and statistical comparisons are presented in Table 2. Data deeper than 15 m are not included due to the lower amount of available lidar data at deeper depths, which results from the steep offshore drop-off in both study areas, particularly for Molokai, and not as an indication of any model limits. In terms of interpreting the linear correlation analysis, and thereby evaluating the capability for deriving absolute water depth, stronger relationships are indicated by correlation coefficients, r, closer to 1 and mean differences, m, closer to 0 m. However, results are not expected to be perfect, in part because of uncertainties in the model process itself, but also due to errors introduced by geocorrection offsets, lidar inaccuracies, and tidal corrections. Additionally, while a high correlation coefficient indicates a strong predictive relationship for water depth alone, the capacity for predicting  Focusing closer on the 2000 and 2017 Kaneohe Bay results, whereas the majority of points demonstrate a very strong correlation, there are clusters of points in both datasets that significantly underestimate water depth. Interestingly, these clusters primarily correspond to the aforementioned lagoon areas with known lower water clarity. Previous research comparing the same AVIRIS-C and lidar bathymetric data, both from 2000, did not include these lagoon areas (i.e., did not use the 2013 lidar data as shown in the current analysis) and reported even higher correlations (r = 0.91 compared to r = 0.84 for this study) [53,59]. This again suggests the need for a dynamic methodology to automatically assign confidence levels as a function of input data characteristics.    Focusing closer on the 2000 and 2017 Kaneohe Bay results, whereas the majority of points demonstrate a very strong correlation, there are clusters of points in both datasets that significantly underestimate water depth. Interestingly, these clusters primarily correspond to the aforementioned lagoon areas with known lower water clarity. Previous research comparing the same AVIRIS-C and lidar bathymetric data, both from 2000, did not include these lagoon areas (i.e., did not use the 2013 lidar data as shown in the current analysis) and reported even higher correlations (r = 0.91 compared to r = 0.84 for this study) [53,59]. This again suggests the need for a dynamic methodology to automatically assign confidence levels as a function of input data characteristics.
Additional analysis was also used to compare the lidar data against derived bathymetry from Kaneohe Bay and Molokai hyperspectral data at each of the different spatial resolutions (Table 2). While output indicates a minor increase in the mean difference as spatial resolution is increased, overall results show that model performance remains consistent at different spatial resolutions. Two other trends that can be interpreted from Table 2 and Figure 4, respectively, are: The model appears to consistently, on average, underestimate water depth (i.e., mean differences are negative); and accuracy decreases with increasing water depth (i.e., there is greater scatter with increasing water depth). Given the accuracies achieved in other studies using the same model, it is surmised that the underestimations shown here are in part a function of using single tide stations as the correction datums for each entire study area (this is a common limitation since tide stations are typically sparsely located). Further, the underestimations are also more significantly a function of including lower confidence data in the analysis (e.g., lagoon areas of Kaneohe Bay), where derived depths are biased shallower when reduced water clarity obscures the signal from the bottom, which skews the overall correlation lower. For the other trend, the decreasing accuracy with increasing depth is to be expected, representing typical physical limitations imposed on all bathymetry algorithms derived from passive optical remote sensing data, with both decreasing signal strength from the bottom and increasing uncertainty in water column properties as depth increases.

Temporal Validation
As a measure of relative model validation, results from the different acquisition years were compared against one another for each study area. The assumption is that if bathymetry remains stable over the given study period (2000-2018), then so too should the derived bathymetry, with any differences attributed to variations in model performance, as well as differences in input data characteristics and/or environmental conditions. As two examples of this comparison, Figure 3 above illustrates the similarity between the 2000 and 2017 Kaneohe Bay bathymetry output, and Figure 5 below illustrates the similarity between the 2001 and 2017 bathymetry along Molokai's shallow fringing reef flat. While both examples demonstrate consistency between acquisitions, including detection of subtle spatial variations in water depth, it is also apparent that erroneous output can be generated in deeper water (e.g., as shown in the speckling within the offshore regions of the 2017 Molokai output, which result from differences in sensor view and solar illumination geometries and variations in factors such as surfaces waves and water properties).
From the linear correlation analysis (Table 3) below illustrates the similarity between the 2001 and 2017 bathymetry along Molokai's shallow fringing reef flat. While both examples demonstrate consistency between acquisitions, including detection of subtle spatial variations in water depth, it is also apparent that erroneous output can be generated in deeper water (e.g., as shown in the speckling within the offshore regions of the 2017 Molokai output, which result from differences in sensor view and solar illumination geometries and variations in factors such as surfaces waves and water properties).

Confidence Levels
When preexisting knowledge is available for water properties and bathymetry, it is feasible to manually subset model output to only those areas with an expected high degree of accuracy (i.e., relatively shallow depths and relatively clear water). However, when such preexisting data is not available, an alternative means is needed to appropriately subset the output. Ideally, for widespread model application, this process should be automated. The hyperspectral data used in this study contains multiple situations where this process is needed: Areas of deeper water (offshore Kaneohe Bay, Molokai, and French Frigate Shoals with depths > 30 m); areas with lower water clarity (Kaneohe Bay lagoon); and areas with suboptimal acquisition conditions (2018 Kaneohe Bay and Molokai). For example, as illustrated in Figure 6, a comparison of the 2000 and 2017 bathymetry output for French Frigate Shoals reveals the importance of being able to subset deeper water areas that are beyond the physical limits for extracting depth. Here, whereas the 2000 data correctly delineates deep water from shallow, the 2017 data generates significant speckling in these same deep-water areas. This speaks to the importance of optimizing image acquisition specifications, as well as the need for an automated methodology for assigning confidence levels to output.
Remote Sens. 2020, 12, 496 16 of 28 To address this need for differentiating reliable results from those that are less reliable, a series of five thresholds were applied to the sunglint-corrected remote sensing reflectance at 600 nm to establish five different confidence levels (from highest confidence to lowest: (i) Rrs > 0.01; (ii) Rrs > 0.0075; (iii) Rrs > 0.005; (iv) Rrs > 0.0025; and (v) Rrs < 0.0025). As mentioned previously, preliminary analysis found that reflectance at 600 nm best correlated to model accuracy and hence provided a first-order indicator of model performance. Figure 7 illustrates the application of these thresholds to the 2000 Kaneohe Bay data, where the lowest confidence level (red) corresponds to the offshore deep water, the next two lowest levels (orange and yellow) correspond to the lagoon and moderately deep offshore areas, and the two highest levels (light and dark green) correspond predominantly to the shallow and moderately deep inshore areas. Based on results from the lidar validation analysis, these confidence levels correlate well with expected delineations in accuracy and model performance for this study area. It is further noted that while confidence levels correspond to certain relative depth intervals for this particular study area and this particular image, the relationship between confidence levels and water depth will vary for different locations and different images as a function of water properties, surface conditions, and image acquisition characteristics. Thus, confidence is not simply a function of depth, but also influenced by other factors as well.  To address this need for differentiating reliable results from those that are less reliable, a series of five thresholds were applied to the sunglint-corrected remote sensing reflectance at 600 nm to establish five different confidence levels (from highest confidence to lowest: (i) Rrs > 0.01; (ii) Rrs > 0.0075; (iii) Rrs > 0.005; (iv) Rrs > 0.0025; and (v) Rrs < 0.0025). As mentioned previously, preliminary analysis found that reflectance at 600 nm best correlated to model accuracy and hence provided a first-order indicator of model performance. Figure 7 illustrates the application of these thresholds to the 2000 Kaneohe Bay data, where the lowest confidence level (red) corresponds to the offshore deep water, the next two lowest levels (orange and yellow) correspond to the lagoon and moderately deep offshore areas, and the two highest levels (light and dark green) correspond predominantly to the shallow and moderately deep inshore areas. Based on results from the lidar validation analysis, these confidence levels correlate well with expected delineations in accuracy and model performance for this study area. It is further noted that while confidence levels correspond to certain relative depth intervals for this particular study area and this particular image, the relationship between confidence levels and water depth will vary for different locations and different images as a function of water properties, surface conditions, and image acquisition characteristics. Thus, confidence is not simply a function of depth, but also influenced by other factors as well.  Figure 8 depicts how confidence levels can be utilized to automatically subset output using the high-confidence results for the 2017 and 2018 Molokai data, in this case by retaining the top two confidence levels and masking out the remaining lower confidence levels. These maps show that the deeper water areas with less model confidence are effectively masked and the shallower areas with higher confidence are retained, which indicates that even a straightforward reflectance thresholding approach can serve as a first-order methodology to establish relative constraints on the reliability of hyperspectral derived inversion model output.   Figure 8 depicts how confidence levels can be utilized to automatically subset output using the high-confidence results for the 2017 and 2018 Molokai data, in this case by retaining the top two confidence levels and masking out the remaining lower confidence levels. These maps show that the deeper water areas with less model confidence are effectively masked and the shallower areas with higher confidence are retained, which indicates that even a straightforward reflectance thresholding approach can serve as a first-order methodology to establish relative constraints on the reliability of hyperspectral derived inversion model output.
high-confidence results for the 2017 and 2018 Molokai data, in this case by retaining the top two confidence levels and masking out the remaining lower confidence levels. These maps show that the deeper water areas with less model confidence are effectively masked and the shallower areas with higher confidence are retained, which indicates that even a straightforward reflectance thresholding approach can serve as a first-order methodology to establish relative constraints on the reliability of hyperspectral derived inversion model output.

Sensitivity Analysis
Given the complex interdependent relationship between surface reflectance, as measured by an overhead remote sensing instrument, and the confounding influences of water depth, water properties, and benthic reflectance, output errors in one derived variable are expected to contribute errors to other derived variables. A sensitivity analysis is used here to investigate specifically how variations in bathymetry can impact derived benthic reflectance and thereby impact the subsequent ability to accurately assess benthic composition.
As a foundation for understanding and interpreting the sensitivity analysis, Table 4 provides a summary of the average benthic reflectance for each image at each different spatial resolution for four different depth intervals. These results illustrate that on average benthic reflectance is relatively low in these study areas. While there are areas of bright sand throughout each area (e.g., Figures 2 and 5), the majority of the benthic substrate is dominated by corals, algae, rubble, and other darker reef components; hence, the lower average reflectance. Results also show consistency across spatial resolutions, with a slight decline in average reflectance as spatial resolution is increased, which is again a function of the substrate being dominated by darker components. Interestingly, these results also reveal that reflectance in the deeper depth intervals for some of the imagery (Kaneohe Bay 2018a and Molokai 2017 and 2018) trend towards 0.01, which represents the lower constraint of potential reflectance values used in the optimization routine. Since the combination of deep water and very dark benthic reflectance represents a physical incongruity, this provides further evidence as to concerns with acquisition characteristics for this particular imagery. Results of the sensitivity analysis are displayed as the average standard deviation in benthic reflectance as a function of different unique combinations of spatial resolution, depth interval, and image acquisition (Figures 9 and 10). These findings suggest that the highest model sensitivity is typically found at shallower depths of 0-5 m, with up to 0.05 standard deviation in reflectance occurring as a result of the first ±5% variation in water depth, and a more gradual increase in standard deviation thereafter. Considering that many coral and algae areas have reflectance values of 0.10-0.25, or lower, at 550 nm, a standard deviation of 0.05 can be significant. The deeper depth intervals as a whole exhibit a lesser degree of sensitivity, but still follow a general trend of increasing standard deviation as a function of increasing variation in water depth. The proportionally lower apparent sensitivity with increasing depth is not unexpected given the logarithmic absorption of light in water, whereby less of an absolute change in reflectance is needed to affect the same impact on surface reflectance as would be required at a shallower depth. Put differently, other things being equal, an equivalent change in benthic reflectance at two different depths results in a greater change in surface reflectance for the deeper of the two. However, in terms of inversion model behavior, the relationship between benthic reflectance, water depth, and surface reflectance is not as straightforward since all of these parameters are allowed to vary independently, which adds an additional degree of uncertainty to this relationship. Hence, it is not surprising that results for the different depth intervals in this analysis do not simply show decreasing sensitivity as the depth range increases. This highlights the complex interdependence of the many different environmental variables contributing to coral reef remote sensing.
between benthic reflectance, water depth, and surface reflectance is not as straightforward since all of these parameters are allowed to vary independently, which adds an additional degree of uncertainty to this relationship. Hence, it is not surprising that results for the different depth intervals in this analysis do not simply show decreasing sensitivity as the depth range increases. This highlights the complex interdependence of the many different environmental variables contributing to coral reef remote sensing.  of positive compared to negative errors in depth, and mainly uniform at different spatial resolutions. This indicates that there is no apparent bias associated with over-or underestimating water depth, and that inversion model performance is consistent across different spatial scales. Such consistency implies a high degree of mathematical stability in the inversion model itself, which is not surprising given the physically-based approach used for formulating the model, and is also advantageous when considering the uncertainties associated with modeling this complex environment.

Workflow Analysis
Physically-based analysis methods, such as the semi-analytical model used in this study, can be used to effectively derive and differentiate information regarding water depth, water properties, and benthic reflectance. Illumination characteristics such as incident solar angle, and environmental Notable exceptions to the general trends exhibited in these results are the relatively lower variations exhibited in the deeper depth intervals of the 2018 Molokai and Kaneohe Bay data. This is in part attributed to the previously identified issues with the 2018 acquisitions. Additionally, for the 2018 Molokai data in particular, this is further attributed to the fact that the Molokai study area is dominated by shallow water with a steep drop-off to deeper water beyond modeling limits, whereas in comparison Kaneohe Bay has much more varied bathymetry. Therefore, while lower variability in the sensitivity analysis may at first appear to always be a positive outcome, in this case lower variability correlates only to the model maintaining consistency even when output accuracy is limited.
These results also show that variations in benthic reflectance are mostly symmetric as a function of positive compared to negative errors in depth, and mainly uniform at different spatial resolutions. This indicates that there is no apparent bias associated with over-or underestimating water depth, and that inversion model performance is consistent across different spatial scales. Such consistency implies a high degree of mathematical stability in the inversion model itself, which is not surprising given the physically-based approach used for formulating the model, and is also advantageous when considering the uncertainties associated with modeling this complex environment.

Workflow Analysis
Physically-based analysis methods, such as the semi-analytical model used in this study, can be used to effectively derive and differentiate information regarding water depth, water properties, and benthic reflectance. Illumination characteristics such as incident solar angle, and environmental conditions such as water clarity and surface roughness, play important roles in overall accuracy of derived output products. Specifically, results suggested that lower solar angles (i.e., lower solar altitude) at the time of image acquisition correlate with reduced output quality, whereas higher solar angles correlate with improved output quality. This is attributed to higher solar angles providing greater signal (i.e., incident radiance) at both the water surface (resulting from lower atmospheric transmittance distance and attenuation) and at the benthic surface (resulting from lower aquatic transmittance distance and attenuation). The tradeoff, however, is that higher solar angles can produce greater amounts of sunglint, which in turn can obscure the signal from the benthic surface. Nonetheless, as indicated in previous studies [53,59,67,68], and shown in this study for the 2000 Kaneohe Bay and 2001 Molokai data, sunglint correction algorithms can effectively mitigate these effects. This implies that future satellite and airborne coral reef acquisitions should consider higher solar angles as an important component in optimizing data acquisition specifications, while both still limiting potential sunglint and incorporating sunglint correction into the workflow.
Robust masking routines for differentiating water from land, clouds, cloud shadows, and wave breaks are also important for widespread application of remote sensing methodologies. While numerous algorithms exist for performing these tasks, there is still need for further improvements, particularly with respect to delineating cloud shadows. The relative differences in brightness between land, clouds, and wave breaks compared to water, especially in the near-infrared, enable a variety of different masking methods to be used for these features. In contrast, the relative similarities in brightness between cloud shadows, including shadows on both water and land, compared to areas of darker aquatic surface reflectance (e.g., areas with dark benthic substrate and/or optically deep water) tend to be a challenging issue for most masking algorithms to fully address. Further research should be conducted in this domain to develop cloud masking routines for future hyperspectral satellites similar to those specifically designed for Landsat, MODIS, and SeaWiFS [69][70][71][72].
This research also introduced considerations regarding the relative tradeoffs of retaining bathymetry as a derived parameter compared to utilizing measured bathymetric data from other sources as an input to a given modeling effort. Even in areas where bathymetric data is carefully calibrated to exacting hydrographic survey standards, tide stations are sparsely located at best, and data availability is insufficient to fully account for the nonlinear spatial variations in water levels at the time of image acquisition that result from localized hydrodynamic complexity (e.g., due to topographic and bathymetric variability, wind fetch, river inflows, and waves). Often a single tide station is all that is available for correcting an entire study area, as was the case in this analysis. The question therefore arises as to the tradeoffs associated with using this data, which may contribute both advantages and disadvantages, compared to deriving bathymetry as part of the overall modeling effort. Given the apparent sensitivity of derived benthic reflectance to errors in bathymetry, this is a research area that requires further investigation.
While not directly addressed in this analysis, image distortions resulting from the underwater terrain and variations in view/illumination angles were clearly evident during image geocorrection. For example, a steep south-facing reef slope appears different when viewed (or illuminated) from the south as opposed to the north. Orthorectification is commonly used in terrestrial remote sensing to remove the effects of geometric, sensor, and terrain distortions, but is largely overlooked in aquatic remote sensing applications. This is not unexpected given the lack of readily available algorithms, as well as the physical complexity of performing such corrections involving the overlying water column and water surface. Therefore, as with terrestrial remote sensing, greater positional accuracy and improved image coregistration would result from incorporating orthorectification into the standard coral reef remote sensing workflow.

Bathymetry Analysis
This study was performed with the intent of developing a framework for processing hyperspectral satellite data to derive bathymetry and benthic reflectance of coral reefs for widespread (ideally global) Remote Sens. 2020, 12, 496 23 of 29 application [21][22][23]. The results are consistent across multiple sites in the Hawaiian Islands, suggesting applicability across geographic space and time. Additionally, the main inversion algorithm [60,61] has been applied in coral reefs across the world with high accuracy, including the Caribbean and western Atlantic [73], Indian Ocean [74], and the Red Sea [75], among others. Thus, the method presented here should be applicable to other coral reefs, and future research can further validate the model using data at coral reefs in other diverse locations.
An important aspect of operationalizing any given analysis method is automating processes for quality-assurance (developing a robust standardized workflow) and quality-control (identifying quality of output data). In coral reef remote sensing, the quality-control process relates to differentiating areas where model outputs are considered valid compared to invalid (e.g., optically shallow compared to optically deep water), as well as assigning confidence levels to those outputs which are considered valid. This is particularly important in most coral reef areas where the availability of ancillary data is limited, making it difficult to perform quality-control using in situ information. Recent methods developed for multispectral data have improved the geographic scale of applications and reduced the need for auxiliary data [76], while others have incorporated multiple methods to identify geographical errors and adjust classifications based on these errors [77]. Nonetheless, there remains a need for further research and investigation into establishing confidence levels. In this analysis, a straightforward thresholding scheme was applied to the sunglint-corrected surface reflectance as a first-order approach for performing this quality-control process. This can help not just with interpreting the bathymetry output but also any higher-level derived products that use this data. Note here that while the confidence levels used in this study were formulated as discrete intervals, similar to data quality flags currently found in Landsat and MODIS imagery, they could also be applied as a layer of continuous values that are indicative of relative overall data quality. To further increase the effectiveness of this methodology as a widely applicable analysis tool, further research is recommended to develop more robust confidence levels that are a function of both water surface reflectance and derived water properties. Ideally, this methodology will allow for widespread regional and global application without the need for large amounts of in situ validation or model input data [36,38], ultimately enabling the automatic derivation of benthic maps with quality-assurance and quality-control built into the processing.
Output from any given inversion model is expected to be sensitive to errors in one or more of the input parameters and/or derived variables. This is especially true in deriving water properties and coral reef habitat information from remote sensing imagery, particularly when considering the interdependence of the variables involved. Thus, it is important to explore this sensitivity to develop a better understanding of the accuracy of the derived output products, as well as any subsequently developed products (e.g., habitat composition). For example, the analysis of the inversion model used in this study examined how variations in bathymetry might impact the derived benthic reflectance. Results indicated that significant differences in reflectance are possible in shallow water with only minor changes in depth. Differences in deeper water were less pronounced, but still nonetheless significant (with the additional caveat that overall limits on model accuracy are still present, independent of sensitivity). This confirms the importance of employing a model that accurately derives bathymetry, and also highlights how the interdependence and overall accuracy of other derived parameters can vary nonlinearly. Thus, understanding model sensitivity is an important component in identifying both model limitations and capabilities.
A final observation from this analysis was that the overall accuracy and performance measures of the inversion model were consistent across different spatial scales. Hence, while the level of detail regarding coral reef habitat composition can vary as a function of spatial scale (e.g., higher spatial resolution can equate to greater habitat information), considering that hyperspectral imagery can be unmixed to provide information at sub-pixel spatial scales [38,59,78], this provides support for the value of future moderate resolution hyperspectral satellites for mapping and monitoring coral reefs. This also indicates that such robust physically-based models can serve as a foundation for deriving other fundamental output products, which can then be used as the basis for generating additional higher-level products.

Conclusions
Coral reef remote sensing presents many challenges, not the least of which is resolving the confounding effects of the water column in order to obtain information about the benthic surface. Given the strong absorption of light in water, bathymetry plays an important role in this process, imposing both limits on the depths to which coral reef habitat information can be assessed and constraints on the effectiveness of remote sensing models in different water conditions. Accurately estimating bathymetry and understanding these limitations are thus key components in the coral reef modeling process. Fortunately, the science and technology of remote sensing in general, and coral reef remote sensing in particular, continues to improve both from the data acquisition viewpoint and the data analytics perspective. For example, with the recent advances in the availability of remote sensing data, including Sentinel-2, Landsat, Planet, and UAVs, along with the rapid development of machine learning classification techniques, it is exciting to see these methods exploring new approaches for coral reef habitat assessment and benthic classification [19,20,28,[31][32][33][34][35]. Nonetheless, as these new methods continue to improve output accuracy and advance the science, most still require large amounts of localized validation data, and work yet remains to develop robust, scalable, automated coral reef data processing workflows.
With the objective of advancing global coral reef remote sensing capabilities, this study utilized hyperspectral data acquired using NASA's AVIRIS-C sensor in Hawaii in 2000-2001 and 2017-2018 to evaluate performance and sensitivity of a representative semi-analytical inversion model for deriving water depth and benthic surface reflectance. The study highlights the interdependent role bathymetry plays in coral reef remote sensing analysis, and implements a multi-step workflow for processing hyperspectral imagery at different spatial scales (18,30, and 60 m), thereby representing different potential specifications of future moderate resolution hyperspectral satellite sensors. Analysis evaluated the accuracy of model-derived bathymetry, established confidence limits on the physical extent of model applicability, investigated the sensitivity of derived benthic reflectance to errors in bathymetry, and assessed model performance with respect to large-scale implementation across differing spatial scales and geographic extents. Results confirmed the ability of the given inversion model to accurately derive bathymetry in multiple complex and geographically separate coral reef environments (r = 0.74-0.97), while also underscoring the importance that acquisition characteristics ultimately have on the accuracy of the derived output (reduced correlation for suboptimal illumination and/or environmental conditions; r = 0.54-0.86). Indications were that higher solar angles, which provide greater signal, result in higher accuracy of derived image products, and that future coral reef image acquisitions should consider higher solar angles in conjunction with limiting sunglint as important components in optimizing data acquisition specifications. The processing workflow revealed the importance of further developing methods for masking land, clouds, and particularly cloud shadows, while also provoking questions regarding the use of ancillary bathymetric data as a model input versus only as output validation, as well as the importance of adapting terrestrial orthorectification methods for use in geocorrection of coral reef imagery. The applied methodology was valid across Kaneohe Bay and Molokai, proving versatility in different shallow water environments, and provided comparable results for French Frigate Shoals, further suggesting its potential for widespread application. Sensitivity analysis underscored the nonlinear interdependence of the different modeling parameters, and indicated that derived benthic reflectance is sensitive to errors in bathymetry, particularly at shallower depths (0.03-0.05 standard deviation in benthic reflectance for ±5% variation in depth at depths <5 m). This study also utilized surface reflectance at 600 nm as a threshold to establish confidence levels in model output, thereby providing a first-order, automated quality-control process for evaluating model validity and confidence. While the thresholding for quality-control needs to be developed further, the research suggests that future hyperspectral missions can use this workflow at different spatial resolutions and different locations. Overall, the consistent accuracy across the differing spatial resolutions supports model application to data from moderate resolution hyperspectral satellites, further confirms the suitability of the model for deriving water depth in complex coral reef environments, and expands on our ability to achieve automated widespread mapping and monitoring of global coral reefs.