A Method to Analyze the Potential of Optical Remote Sensing for Benthic Habitat Mapping

Quantifying the number and type of benthic classes that are able to be spectrally identified in shallow water remote sensing is important in understanding its potential for habitat mapping. Factors that impact the effectiveness of shallow water habitat mapping include water column turbidity, depth, sensor and environmental noise, spectral resolution of the sensor and spectral variability of the benthic classes. In this paper, we present a simple hierarchical clustering method coupled with a shallow water forward model to generate water-column specific spectral libraries. This technique requires no prior decision on the number of classes to output: the resultant classes are optically separable above the spectral noise introduced by the sensor, image based radiometric corrections, the benthos’ natural spectral variability and the attenuating properties of a variable water column at depth. The modeling reveals the effect reducing the spectral resolution has on the number and type of classes that are optically distinct. We illustrate the potential of this clustering algorithm in an analysis of the conditions, including clustering accuracy, sensor spectral resolution and water column optical properties and depth that enabled the spectral distinction of the seagrass Amphibolis antartica from benthic algae. OPEN ACCESS Remote Sens. 2015, 7 13158


Introduction
A fundamental issue with benthic classification of remotely imaged shallow water environments is determining the appropriate definition and number of benthic classes that: (i) optimizes the classification accuracy and precision [1]; and (ii) standardizes classification maps for ease of inter-comparisons [2].The selection of classes and their descriptive resolution, i.e., biological detail [3], for a spectral library or training dataset, along with the image classification method, have typically been scene and sensor specific [4].For multispectral imagery, the optimum spectral library could consist of image-derived or in situ spectra of pure endmembers (a discrete taxanomic unit i.e., high descriptive resolution) if the benthos in the scene is either homogeneous or was imaged with very high spatial resolution [5][6][7].For a scene with patchy or heterogeneous benthos or that was imaged with moderate to low spatial resolution, the spectral library could contain classes pertaining to merged endmembers (i.e., mixed benthic assemblages) and low descriptive resolution [3,8,9].Due to the spectral characteristics of the sensors used (e.g., SPOT, Landsat 5-7, IKONOS; QuickBird), the above studies mostly utilized supervised classification schemes with or without water column correction to produce benthic habitat maps (additional examples include [10][11][12]).
Hyperspectral sensors, in contrast, have enough spectral resolution and bands to potentially facilitate the spectral unmixing of an image spectrum based on the fractional cover of a subset of pure endmembers [13].This has often led to the implementation of a spectral library of pure endmembers which would be linearly mixed either to pre-defined proportions such as in look-up table methods [14,15] or during spectral optimization in shallow-water inversion methods [16][17][18].These methods have consequently achieved moderate to high benthic classification accuracies [14,[19][20][21].
Previous studies with multispectral imagery have shown an inverse relationship between the number of benthic classes and the classification accuracy [1,3,8].This same relationship has also been shown for hyperspectral imagery, see [21,22].This raises the question of how to optimize the class selection, or equivalently how to merge classes of higher descriptive resolution in order to achieve both accurate and useful maps.Karpouzli et al. [23] observed that agglomerating pure endmembers based on their genera to obtain average spectra of coral, seagrass, macroalgae and sand, reduced the classification accuracy.Specifically high classification accuracy was obtained with more classes at finer descriptive resolution.This was attributed to the fact that the intra-class variability exceeded the inter-class separability-as it has been noted that some species of corals are spectrally similar to macroalgae and vice versa [24].Clearly, averaging pure endmembers based on their genera, which seems a logical and ecologically meaningful step, may not maximize spectral separability between classes, and potentially leads to a higher probability of confusion during classification.
Furthermore increasing the descriptive resolution of the benthic classes may lead to a decrease in the precision of the resultant classification, particularly if sensor and environmental noise is taken into account.Such noise is a component of remotely sensed imagery originating from the sensor and from atmospheric, sunglint and air-water interface corrections that maybe imperfect at times [25,26].The combined impact of sensor and environmental noise and spectral variability of a given taxonomic species observed in local and regional scales [27,28] act to degrade the inter-class spectral separability as described in Hedley et al. [26].
Minimizing confusion or uncertainty arising from spectrally similar classes necessitates a procedure to re-define the spectral library of endmembers into more distinct classes.Clustering is an approach that can merge those spectrally similar classes thereby reducing the spectral confusion in a spectral library.A variety of indices can be used to define the spectral similarity between two classes including root mean square error, spectral angle mapper [29] or spectral information divergence [30].Though clustering is a common procedure in unsupervised classification [31][32][33][34], it is seldom performed on the endmember spectra that constitute the spectral library (for exceptions see [35][36][37]).Particularly on benthic spectra that are further modulated by the physical processes occurring in an optically variable water column.In this study we present a hierarchical clustering algorithm, based on linear discriminant coordinates, tailored for shallow-water inversion models that uses the intra-class variability to merge those classes in the benthic spectral library that overlap, i.e., are optically indistinguishable.Here, the intra-class variability incorporates the individual benthos' natural spectral variability plus image-based sensor and environmental noise.The hierarchical clustering ceases when there is no spectral overlap in the variance between groups; and thus outputs a set of classes that are spectrally distinguishable under actual operational conditions.In combination with a shallow water forward model this hierarchical procedure was used to develop depth and water-column specific spectral libraries.The endmembers of which are spectrally distinguishable above the attenuating properties of an optically variable water column at depth.In addition we investigate the effect spectral resolution has on benthic class separability by analyzing the clustering from simulations using the following three sensors: Hyperspectral Imager for the Coastal Ocean (HICO); HyVista's HyMap, and; Worldview-2 (WV2)-see Table 1 for a list of symbols and acronyms used in the text.Such information on the number and type (i.e., definition) of distinguishable classes at any given depth and water column optical properties helps to understand the potential of benthic classification from shallow water remote sensing imagery.

Methods
The overall goal is the development of a procedure that can quantitate the number and type of classes that are spectrally distinct for any given water column optical properties and depth from a spectral library of representative benthic species.We begin in Section 2.1 with a description of the benthic irradiance reflectance, ρb, library collected from the field and we then describe the hierarchical clustering algorithm in Section 2.2.Using the measured ρb data, the clustering algorithm outputs a library of endmembers that are spectrally separable above the total system's variability but only suitable just below the air-water interface as the attenuating effects of a water column were not modeled.Given an estimated or known range of depths and water column optical properties in a scene, a semi-analytical shallow water forward model coupled with the clustering algorithm can be used to predict the benthic classes that are optically distinct.This procedure is detailed in Section 2.2.3.

Benthic Reflectance Library
A spectral library of pure endmembers was derived from in-air irradiance reflectance measurements of 22 benthic species collected from the Point Peron (32.2715°S, 115.6865°E) study site, Western Australia, on 22 August 2014 (Table 2).Benthic samples were collected with enough material to cover the 25° field of view of the ASD field spectrometer (1 nm resolution, 350-2500 nm coverage).However, the quantities of Posidonia sp. and Metagoniolithon stelliferum were too small to cover the field of view and these were placed over a substrate (sediment and rock/rubble respectively) that would normally be underneath.These samples were stored overnight in a tank that had aerated, flowing filtered seawater to minimize the effect of pigment degradation.The collected samples were then spectrally analyzed the following day on a bench top set up outside on a flat roof.All samples were placed on top of a large flat, non-reflective black tray.The ASD fiber optic was held using a retort stand at a constant 30 cm above the sample and at a small angle off nadir (approx.10°) and a white reference plate was used to compute the reflectance for each spectral measurement at the same angular geometry.Furthermore, each benthic sample was stirred and mixed every three spectral measurements to capture the sub local-scale taxanomic spectral variability.With this setup approximately 30-40 in-air reflectance measurements were recorded for each (wet) benthos collected.At the time of sample collection (winter) the leaves of Posidonia sp. were in the dormant phases, and new leaves were not visible in the collection area.The leaves that were collected, however, were already detached and slightly senescent.The benthic irradiance reflectance spectra, ρb, between 400 and 680 nm were then convolved with the spectral response functions (SRF) of the HICO, WV2 and HyMap sensors to generate three spectral libraries (Figure 1).The clustering algorithm was used on these spectral libraries to assess the impact of reducing spectral resolution.Here, wavelengths past 680 nm were not considered due the high absorption of light in the water column.This reduced the number of bands of HICO and HyMap to 49 and 16 respectively.For WV2 the five spectral bands between 400 and 680 nm were used.

Measure of Interclass Overlap
A measure of the degree of overlap or misclassification between a pair of classes was used in the clustering algorithm to identify which classes to merge per iteration.Specifically, the closest (most similar) pair of classes whose misclassification proportion, τm, exceeded a user-defined amount were combined.We have defined τm between two classes i and j as the proportion of misclassified spectra in Linear Discriminant (LD) coordinates as given by Rencher and Christensen [38], where ni and nj are the number of LD points (i.e., spectra) in classes i and j, whilst nmi and nmj are the number of misclassified LD points from those respective classes.We consider a point misclassified if it is closer or equal in distance to the mean of a different class than to its own class mean.Here the straight line distance between a point and a class mean in LD coordinates is computed using the RMSE, where similar spectra have low RMSE.This approach is equivalent to inserting a separating plane between two classes and counting the number of spectra on their incorrect side as done in Hedley et al. [26].This per-point per-class comparison was performed for all possible pairs of benthic classes in the clustering procedure per iteration.

Hierarchical Clustering Using Linear Discriminant Coordinates (HDC)
The measured ρb spectra (initial spectral library) may contain many benthos that are indistinguishable in terms of their reflectance when sensor and environmental noise NEΔrrs [18] and the species natural spectral variability are taken into account.To assess this, an agglomerative centroid hierarchical clustering algorithm based on linear discriminant coordinates (hereafter referred to as HDC) was developed.Here classes were merged that had a misclassification proportion (τm) greater than a user-specified threshold.This threshold, set at 5% in this study, is inversely proportional to the number of benthic classes generated from the HDC algorithm.We define the class-merging (or clustering) accuracy as 100 − τm, and thus would allow a user to set a threshold for the level of class separability.This is a more quantitative and intuitive approach than predetermining, arbitrarily, the number of classes that the clustering should output.The HDC algorithm can be summarized in the following three steps (flowchart A, Figure 2): (1) the addition of NEΔrrs to the measured ρb spectra; (2) the transformation from spectral space to LD coordinates; and (3) subsequent hierarchical clustering.The advantage of this approach is that no prior decision on the number of resultant clusters to output is needed.This is a common issue faced by conventional hierarchical clustering [39,40].Here the number of output clusters is dependent on the following: (a) the magnitude of the error given by the covariance matrices-pertaining to NEΔrrs and the taxanomic spectral variability; (b) the spectral resolution of the sensor; and (c) the τm threshold.To account for the total noise in the system just below the air-water interface in step (1), we perturb the measured ρb spectra with NEΔrrs using the procedure developed by Hedley et al. [41,42].As such spectral noise caused from the sensor and any spatial noise from atmospheric, sunglint and air-to-water interface corrections were included.Briefly, the spectral covariance matrix, CSE, was extracted from an imaged homogeneous deep-water region of subsurface remote sensing reflectance (rrs) imagery.A set of pseudo-random spectral noise terms, δρSE, were computed from CSE and were added to a single ρb spectrum to generate a number of noise-perturbed spectra, ρb ± δρSE.This number varied for a given class but was such that it would generate 2000 noise-perturbed benthic reflectance spectra.A total of 44,000 spectra from the 22 classes were consequently produced.A detailed description of the computation of δρSE (presented as δrrs) is given by Garcia et al. [43] (see "Optimization and uncertainty propagation" section).
In this developmental study, we estimated the spectral shape and magnitude of NEΔrrs from an imaged homogeneous deep-water region of the HICO-derived rrs image of Shark Bay 14 December 2011.Interpolation to the wavelength centers of HyMap and WV2 were performed to keep the magnitude and spectral shapes of NEΔrrs consistent between sensors.The spectral covariance matrices CSE WV2 and CSE HyMap were then computed and used in the HDC algorithm for the respective sensors.The atmospheric, sunglint and air-to-water interface corrections used on the HICO image of Shark Bay to convert top-of-atmosphere calibrated radiances to rrs are described in Garcia et al. [44].
In step (2) the Linear Discriminant Analysis (LDA) procedure for the case of several classes was implemented [38] to convert the 44,000 ρb spectra into LD coordinates.The optimal number of eigenvectors, s, of the resultant LD points was chosen through an iterative approach, where s was successively increased to maximize the number of classes, k, generated from the subsequent clustering.If k remained constant with increasing s for more than three successive increments of s then the iteration was ceased and the eigenvector that produced the last improvement in k was selected.
For a given eigenvector number an iterative centroid-based hierarchical clustering procedure was implemented to merge overlapping benthic classes.Here a single pair of classes was merged per iteration with τm set to 5%.The iterative clustering begins with the computation of the RMSE between all possible pairs of class means, and τm between those pairs.The pair of classes that had the lowest RMSE and whose τm > 5% were merged to form a single class.In this step the LD coordinates of the two merging classes remain unchanged only that they are labeled as one.The ρb spectra of these two classes were used to compute the weighted average of the newly merged class.This iteration continues until the remaining classes all have a τm ≤ 5%; in other words the clustering ceases when there is a 95% accuracy that all classes are optically separable above the sensor and environmental noise and the benthos' natural spectral variability.This clustering accuracy can be reduced (with a corresponding increase in the number of "separable" classes) according to the desire of the user or application.The output of this entire procedure is a spectral library that contains the optimum set of endmembers for classifying the set of substrates of interest.
We compare the resultant classes from the HDC algorithm with an agglomerative centroid-based hierarchical clustering procedure taken from Everitt et al. [45] using only bottom reflectances.For simplicity in this text we refer to this type of clustering as the standard hierarchical clustering.In this standard hierarchical clustering the pair of classes with the lowest RMSE were merged per iteration.Note, the mean class ρb spectra were used to compute the RMSE between all possible pairs, and the weighted average spectrum was computed for a newly merged group.The merging continues until all groups have clustered into one class.We estimate the appropriate number of clusters by locating the knee from a bi-plot of the linkage distance against the number of clusters [46].Here, the linkage distance is simply the RMSE between the pair of groups that are merged at a given cluster iteration.We utilize the definition of "knee of the plot" as the point that experiences an abrupt change in the RMSE as done by Torrecilla et al. [47].

Depth and Water Column Specific Spectral Libraries
Without incorporating the attenuating properties of a water column to ρb the HDC generated spectral library would only be suitable just below the water's surface.Employing such a spectral library into a physics-based inversion, for instance, may not give a realistic representation of benthic classification at increased depth and/or turbidity (discussed in Section 3.4) as fewer classes would be optically separable [26,48].To include the attenuating properties of a water column and water depth we used the semi-analytical (SA) shallow water forward model given by Lee et al. [16].This model estimates rrs from the following scalars: (a) the absorption coefficient of phytoplankton at 440 nm, P (m −1 ); (b) the absorption coefficient of dissolved and detrital matter at 440 nm, G (m −1 ); (c) the backscattering coefficient of suspended particles at 550 nm, X (m −1 ); (d) depth, H (m); (e) and the bottom reflectance, ρ.This SA model was used to compute rrs over each measured ρb spectrum described Section 2.1 from a given set of P, G, X and depth.The set of rrs spectra were then passed through the HDC (flowchart B in Figure 2) to merge those classes that were optically indistinguishable above the total system's noise and attenuating water column.After the clustering based on rrs, the corresponding benthic reflectances (ρb) were merged to produce the mean bottom endmember spectra.Note that in this SA model the phytoplankton spectral absorption shape is almost fixed.Though in reality this spectral shape is likely to change with the presence of different phytoplankton functional groups (spatially and temporally); benthic classification from remote sensing, which this algorithm is tailored towards, is typically performed in relatively clear shallow waters with relatively low chlorophyll concentrations.In these conditions, the shape of the phytoplankton absorption spectrum would likely have very minor impact to the clustering outputs.
The impact of water column depth and optical properties on the classes produced from the HDC algorithm was modeled by selecting specific depths between 0.5 and 20 m for the following water column optical properties: (1) ). κ and ZSD are the attenuation coefficient at 490 nm and the estimated Secchi depth respectively, κZSD = 1.44 [49].Here, ZSD was used as a qualitative and readily interpretable measure of water clarity.An analysis of the conditions (clustering accuracy, spectral resolution, water column depth and optical properties) where the seagrass species A. antartica is spectrally resolvable against benthic algae was performed.For this analysis the clustering was run iteratively over a range of depths (0.25 to 25 m at 0.25 m increments) at the four water clarities mentioned above for eight different clustering accuracies (60% to 95% at 5% increments).

Hierarchical Clustering of Benthic Irradiance Reflectance Spectra
The HDC-derived dendrogram using the basic ρb spectra at HICO bands describes which benthic classes merge and the RMSE of that union (Figure 3).The respective dendrogram generated from the standard hierarchical clustering is presented in Figure 4.The RMSE refers to the spectral difference (in LD coordinates for the HDC algorithm) between one class mean and another, where the larger the value the larger the spectral difference.Based on the simulated spectral variability, occurring just below the air-water interface, the clustering output 18 classes for HICO (Figure 3) and HyMap (Figure 5) and 14 classes using WV2 bands (Figure 6).The clusters that were formed using HICO and HyMap bands were identical except that the merging of the clusters occurred at lower RMSE for HyMap (Figure 5).For example the merging of C. sinuosa with the cluster S. linearifolium/S.spinuligerum occurred at an RMSE of ~0.035 for HICO and ~0.032 for HyMap.This slight decrease in RMSE is likely due to the lower spectral resolution of the HyMap sensor or due to the lack of spectral bands below 450 nm (i.e., less spectral bands to facilitate separation).The dendrogram for HICO (Figure 3) and HyMap (Figure 5 These two clustering algorithms generate different classes because the HDC first identifies pairs of classes with τm > 5% (the user defined threshold) and then merges the pair with the lowest RMSE (Mahalanobis distance in spectral space).The centroid based hierarchical clustering (Figure 4) instead merges the most similar class pairs regardless of whether their variance overlap or not.In other words a pair of classes may be grouped even if they could be considered optically distinct.An added advantage of the HDC is that the number of classes selected is imbedded inside the clustering, where the merging of clusters ceases when the variance overlap between groups is below a user defined threshold.Thus the criterion for stopping is statistically meaningful and aligns with the sensor and environmentally limited system of remote sensing as described and illustrated in Hedley et al. [26] (see Figure 1 in [26]).The user defined threshold describes the proportion of misclassification that a pair of classes should have before they are considered optically indistinguishable and hence merged.For example using HICO bands, the two classes C. germinata and C. flexis had a τm of 20%, i.e., 20% of the total spectra of these two classes were located in the overlap region and hence inseparable.Given that τm > 5% and proximity of the class means the HDC clustered these C. germinata and C. flexis after the second iteration.
The much lower spectral resolution and fewer wavebands of the WV2 sensor only facilitated the separation of 14 classes with the HDC.The dendrogram (Figure 6) shows the formation of the following clusters: (1) all the brown algae classes; (2) the green alga cluster Entermorpha sp./ U. australis; (3) the mixed seagrass-green algae class of Posidonia sp./C.germinata/C.flexis/ B. vestita; and (4) the red alga cluster E. articulata/M.stelliferum.Thus out of the seagrass species, only A. antartica was spectrally separable using WV2 bands.At the time of collection the seagrass Posidonia sp.samples were senescent and less spectrally distinct (using WV2 bands at least), which likely lead to its inclusion into the green algae cluster.Fyfe [50] compared the in-air reflectance of several unfouled and fouled seagrass species and showed that the latter had broader and less reflective peaks between 520-580 nm compared to the unfouled case.It is likely that this seagrass species can be optically separated if their reflectance spectra are collected during their growing season.

Water-Column Specific Benthic Spectral Libraries
Two spectrally distinct benthic classes remain at a modeled water depth and optical properties of 15.0 m, P = 0.05 m −1 , G = 0.10 m −1 and X = 0.010 m −1 using HICO bands (Figure 7).Here the HDC has merged all benthic vegetation species into a mixed vegetation class and merged the remaining two sediment classes.Thus at this water column optical property and depth it is not possible to distinguish between any of the initial benthic vegetation species above the total system noise, and only a bright (i.e., sediment) and dark (i.e., vegetation) substrate can be distinguished, assuming completely filled pixels of each type.Figure 7 shows the possibility of determining a priori which benthic classes are separable in a single image pixel at a fixed depth and optical property.Thus for an image with varying depths and water column optical properties the HDC clustering will need to be applied per-pixel.Physics-based inversion methods can be used to derive the depth and water column optical properties as will be discussed in Section 3.4.
The decrease in the number of optically separable benthic classes for increasing depth and water column turbidity at HICO, HyMap and WV2 bands was quantified with the HDC clustering (Figure8).Clear to turbid water columns (Figure 8, water types 1 and 4 respectively) were modeled to assess their impact on benthic class optical separation.For a given water clarity the number of optically distinguishable benthic classes decreased in a near exponential manner with increasing depth, as has been described in the literature [26,48].Increasing the water turbidity also reduces the number of separable classes dramatically with increasing depth.With HICO spectral bands, 13 classes can be optically distinguished at 5 m depth for water type 1 (48 m Secchi depth, Figure 8a), whereas only 11, 6 and 2 classes are separable for water types 2, 3 and 4 (Secchi depths of 13, 6.2 and 2 m) respectively at that depth.The water column saturation point, that is, where the water column contributes to nearly all of the water leaving radiance, is located when none of the benthic classes can be distinguished (i.e., number of classes = 1).For water types 2, 3 and 4 these occur at depths 20 m, 12 m and 5.5 m respectively using HICO bands.
Decreasing the number of spectral bands also increases the slope of the exponential curve as illustrated when comparing Figures 8a (HICO) and 8c (WV2).For water type 3 at 3 m depth the HICO, HyMap and WV2 wavebands, according to the HDC, can distinguish 12, 11 and 5 classes respectively.Increasing the depth to 5 m reduced the number of classes to 6, 6 and 3 respectively.For depths less than the water column saturation point, the five WV2 bands between 400 and 680 nm can generally separate half as many benthic classes as HICO.Based on the modeling shown in Figure 8 and a comparison between the exponential slopes, there appears to be minor differences in the number of classes distinguishable at depth between HICO and HyMap wavebands.In this analysis, wavelengths between 400 and 680 nm were used which reduced the number of bands to 49 and 16 for HICO and HyMap respectively.This may seem like a large reduction of wavebands without greatly affecting the benthic separability.We should note that unlike HICO, HyMap does not have bands below 450 nm [51], hence both sensors have wavebands that emphasize the spectral differences between the benthic classes analyzed.The literature has shown that specific wavelength ranges predominantly within 520-680 nm are useful for identifying substrate types, such as healthy coral, bleached coral, seagrass, sand and algae [34,39,50,52,53].Therefore, it is likely that the difference in bandwidth is the cause where the spectral resolution of HICO and HyMap are 5.7 and 15 nm [51] respectively.This is in line with what Hochberg and Atkinson [24] showed where the separation of coral, algae and sand using LDA were nearly identical for spectra at 10, 5.5 and 1 nm spectral resolutions.The number of classes within a spectral library can be stratified with depth, for instance using water type 1, the spectral library for WV2 (Figure 8c) has the same number of classes from 5 m to 15 m.Stratification in the number of classes is expected as they are integer values and hence more emphasized in slowly decreasing curves and less observed in areas of rapid change (e.g., 0-3 m depth of Figure 8b for water type 3).Analyses of the spectral libraries from all three sensors within the stratified regions of a given water type indicate that at >5 m the output classes produced are identical though subtle differences exist in the RMSE of when classes merge.At shallower depths, however, the spectral libraries in the stratified regions have minor differences in the classes that merge and thus in the output classes.For HICO spectral bands 13 classes were generated in the stratified region between 3 and 6 m depth with water type 1.The spectral libraries at 3 and 4 m are identical as are the spectral libraries at 5 and 6 m (Table 3).Between these two sets of spectral libraries the following differences occur: (1) the brown alga E. radiata and C. sinusoa from the mixed vegetation class (at 3 m) split to form an individual mixed brown algae class; and (2) the union of E. articulata and M. stelliferum (Table 3).Table 3.The classes from the HICO spectral libraries at depths 3-6 m using water column optical properties of P = 0.01 m −1 , G = 0.01 m −1 , and X = 0.001 m −1 .The superscripts represent the number of species that form that cluster, for example mixed brown algae 4 means that four brown algae species were merged to form one mixed brown algae class.The spectral libraries at 5 and 6 m depth are identical.The classes from the HyMap spectral libraries for depths 3 to 6 m with water type 1 are given in Table 4.A stratified region (Figure 8b) is formed in this depth range and water type, where 12 classes are optically separable at 3, 5 and 6 m depth.At 3 m depth four species of brown algae merge to form one cluster, however at 4, 5 and 6 m two of those brown algae species (E.radiata and C. sinusoa) form a separate mixed class and the other two merge into a mixed brown algae/green algae/seagrass class (Table 4).Other differences include: (1) the separation of the Sediment-Rubble/Sediment class (3 and 4 m) to their individual class at 5 and 6 m; and (2) the union of the green alga Entermorpha sp. and U. australis at >3 m depth.Some clusters are consistently generated within this depth range which include Rock-red coralline algae and six red algae classes E. articulata/M.stelliferum, Ballia sp., H. ramentacea; A. armata; B. callitrichia, and A. anceps.Indeed each spectral library between 3 and 6 m has six red algae classes, >1 sediment class, >1 green algae class and a mixed vegetation class.At 4 m depth the HDC procedure generated 11 classes using three eigenvectors, analysis has shown that using six eigenvectors facilitates the separation of 12 classes identical to that produced at 5 and 6 m (Table 4).Here six eigenvectors were not used because the fourth and fifth eigenvectors did not increase the number of clusters generated and thus the iterative eigenvector selection procedure chose the first three.Like for HICO (Table 3) there are minor differences between spectral libraries in these stratified regions, as such we propose the use of a single spectral library for such regions.

Resolving Seagrass Species from Algae
An important ecological application of optical remote sensing is mapping seagrass meadows [54].Using the HDC algorithm, we consider the clustering accuracy, sensor spectral resolution, water column depth and optical properties that enable the optical distinction of the seagrass species A. antartica from green and brown alga.Figure 9 shows the modeled depth at which we can no longer resolve A. antartica from the green algae species (Table 2) for clear to turbid water clarities (water types 1 to 4 respectively) using the different sensors.Thus for HICO bands A. antartica is optically distinguishable from green algae to a depth of approximately 4.25 m at 90% clustering accuracy for a clear water column (water type 1, Figure 9a).Beyond this depth at 90% accuracy A. antartica cannot be spectrally distinguished.Note that the results for WV2 in Figure 9c relates to the ability of distinguishing A. antartica from the Posidonia sp./green algae cluster, as the latter seagrass species could not be resolved using bottom reflectances (Figure 6).
For a given water type, decreasing the clustering accuracy enables the ability to spectrally resolve seagrass from green algae to greater depths.The reason being that the HDC allows for greater misclassification proportion between pairs of classes before they are considered "indistinguishable".For example, using HyMap bands A. antartica can be distinguished at 2.5 and 3.75 m for water type 2 at 85% and 75% clustering accuracies respectively (Figure 9b).Therefore if the user desires the ability to map A. antartica up to a depth of 4 m with HyMap, then 75% clustering accuracy should be used for water columns with Secchi depths greater than 13 m.A reduction of the water clarity also decreases the ability to resolve species at depth, as observed by the ability of WV2 to distinguish A. antartica to 2 m and 0.75 m for water types with Secchi depths 48 and 2 m (Figure 9c) respectively at 80% accuracy.Therefore as the modeled water depth or turbidity is increased the distances between classes decreases such that neighboring classes experience more overlap [26].Hence a greater τm accounts for the decrease in the inter-class distance, where classes are considered optically distinct despite their increase in overlap.We have low confidence for results that used clustering accuracies less than 70% as, for example, HyMap bands can distinguish A. antartica to 2.25 m for water type 4 (Figure 9b) that has an estimated Secchi depth of 2 m.The advantage of using hyperspectral sensors is clearly demonstrated in Figure 9, where HICO and HyMap have the ability to resolve seagrass from green algae to much greater depths for a given clustering accuracy and water clarity compared to WV2.Indeed at 85% clustering accuracy, HICO and HyMap spectral bands can distinguish A. antartica to an estimated depth of 6.5 and 3.75 m respectively for water type 1 (48 m Secchi Depth), whilst WV2 can only resolve it to a depth of 1.5 m.On average HICO bands can detect A. antartica to depths 2.7 times greater than WV2 bands across the clustering accuracies and water types.In most conditions analyzed there is close similarity in the depth limits at which A. antartica can be resolved using HICO and HyMap bands, where on average HICO can detect this species to depths 1.2 times greater than HyMap.Table 4.The classes from the HyMap spectral libraries at depths 3-6 m using water column optical properties of P = 0.01 m −1 , G = 0.01 m −1 , and X = 0.001 m −1 .The superscripts represent the number of species that form that cluster, for example mixed brown algae 4 means that four brown algae species were merged to form one mixed brown algae class.The spectral libraries at 5 and 6 m depth are identical.

Clusters
The ability to spectrally resolve mixed seagrass/green algae from brown algae was analyzed for the same conditions presented in Figure 9.However, a direct comparison between seagrass and brown algae was not possible for HICO, HyMap and WV2 at any of the water column conditions presented.The reason being that Posidonia sp. and A. antartica both merge with a mixed green algae class.This seagrass/green algae mix then clusters with a mixed brown algae class at deeper depths.The distinction between mixed seagrass/green algae and brown algae (effectively between green and brown colored benthic species) follows the same trend as that presented in Figure 9 expect that they can be distinguished to much deeper depths.

Implications to Shallow Water Habitat Mapping
Modeling for a variety of water turbidities and depth as displayed in Figure 8 enables predictions on the number and type of benthic classes that can be expected from a hyperspectral or multispectral image at a constant clustering accuracy.The level of accuracy can be changed according to the needs of the user or application with corresponding changes to the number and type of benthic classes defined by the HDC. Figure 10 illustrates that the number of optically separable classes increases as the clustering accuracy decreases for a constant set of water column optical properties at HICO wavebands.These results are similar to those produced by Andrefouet et al. [1] where the thematic accuracy decreased with increasing number of habitat classes used in the benthic classification scheme.Although thematic and clustering accuracy are applied in different contexts, Figure 8 does imply that if a benthic habitat map at a constant classification accuracy is required, then the number of classes used in the classification scheme must change with the range of depths and water column optical properties within an image.Prior knowledge of the depth and optical properties in a given scene therefore becomes a requisite in applying the HDC algorithm.Note the depth, P, G, and X can be determined by applying a shallow water inversion model for sensors with enough spectral bands in the visible [16,17].
Here the full benthic spectral library could be used solely to optimize these parameters.Benthic classifications of remotely sensed imagery using physics-based inversion models to date have used a spectral library of representative benthic species (fine descriptive resolution, see Dekker et al. [20]).Depending on the shallow water model used, the inversion process iterates through unique combinations of benthic mixtures and optimizes for their fractional coverage, water column optical properties and depth.The benthic combination and fractional coverage that generated the best fit to the input reflectance are assigned.Therefore distinct benthic species would be assigned to image pixels at any given depth and water column optical property, no matter if it were unlikely to distinguish one benthic species from another.
The research presented here indicates that this type of benthic classification may not be appropriate, particularly when mapping regions of varying water depths or when propagating sensor and environmental noise in the inversion scheme [41,42,44,55].This can be illustrated by the dendrogram presented in Figure 7, where the at-surface spectral differences between the benthic vegetation species do not make observable differences to rrs above the total system noise at the depth and P, G and X modeled.Clearly if a spectral library consisting of the initial 22 benthic classes were used for an inversion of say a pixel having 100% coverage of E. radiata with P, G, X and depth used in Figure 7; then accordingly a shallow water inversion model would likely not be able to distinguish the forward modeled rrs derived from E. radiata from any other benthic vegetation species when NEΔrrs is incorporated.In this situation it is likely that equal probability of assignment across the benthic vegetation classes will result.This was illustrated in Figure 5 by Hedley et al. [42] where the uncertainty of the bottom type increased with depth and increasing number of endmembers in the spectral library.Furthermore, the assignment of a specific benthic species (even correctly identified E. radiata) cannot be made with sufficient certainty, and a classification assignment of "mixed vegetation" would be more appropriate for these optical properties and depth.In fact such a class name would convey that level of uncertainty and lack of spectral separability.Indeed, analysis of HDC-derived output classes has showed that as the water column becomes more turbid or deeper, the ability to distinguish distinct benthic species is lost, where resultant classes consist of mixed clusters.Again in these situations, assigning the appropriate mixed cluster would help convey the level of certainty and optical separability into the classification scheme.
Presently the HDC algorithm assumes 100% substrate coverage of a pixel and as a post-inversion tool would be suited for high spatial resolution image data.In addition, its application would be more appropriate towards constraining the benthic classification of the HOPE inversion model [16], if the endmembers of a high descriptive resolution spectral library were cycled through during the optimization process.Note that the HOPE model considers one substrate endmember rather than a linear mixture of two [18] or three [17].The requirement of moderate to high spatial resolution image data does limit the applicability of the HDC algorithm as a post-inversion tool to mostly airborne sensors such as AISA Eagle, HyVista's HyMap, CASI-2, AVIRIS and Ocean PHILLS.From the current and planned satellite sensors, WV2, Sentinel-2 and VENμS SSC [56] possess high spatial resolution (≤10 m) and enough spectral bands in the visible domain to facilitate optimization using current shallow water inversion models (e.g., [42,57]).
Analysis has showed that the HDC-derived water column specific libraries are sensor and scene specific.Sensor specific because the number and type of HDC-derived benthic classes depend on the sensor's spectral resolution and signal to noise ratio.Scene specific because: (1) the environmental noise in the imagery, such as sunglint which even if corrected [58] leaves residual spatial noise; and (2) The number and type of output classes are dependent on the representative benthic species present.The latter is implied by the clustering of the seagrass Posidonia sp. and the green algae C. germinata, C. flexis and B. vestita at WV2 bands (Figure 6).If none of the green algae species were present in the scene then Posidonia sp. would be spectrally distinct, and would lead to different benthic class mixtures.
As the clustering is based on LD coordinates the ability to distinguish individual classes potentially decreases when the number of initial classes increase.This is illustrated in Figure 11, which shows the LD coordinates of the four brown algae E. radiata, S. linearifolium, S. spinuligerum and C. sinuosa.Three clusters resulted when the ρb spectra (at HICO bands) of these four brown algae species were passed through the HDC algorithm (see Figure 11).Recall that two brown algae clusters were considered distinct (S. linearifolium/S.spinuligerum/C.sinuosa and E. radiata) when the initial spectral library of 22 classes where passed through the HDC (see dendrogram in Figure 3).The reason for this is that the LDA seeks a projection that maximizes the distance between classes whilst minimizing the within-class variance [59,60].Thus in the case of Figure 11, an optimal projection was found that enabled the distinction of three brown algae clusters.These results imply that only a minimal set of initial benthic classes should be used in the clustering.Although a sensor's spectral characteristics, total system noise and water column attenuation affect substrate separability; spatial resolution, not accounted by the HDC algorithm, also plays a significant role.Particularly since coastal and coral reefs can be spatially heterogeneous, where a diverse array of spectrally distinct substrates can occur at sub-meter scales [55].Previous research involving thematic mapping and Object-Based Image Analysis of coral reefs has shown improved benthic classification using imaging platforms with high spatial resolution (<10 m) compared to those with moderate resolution (30 m) [61].Thus although a sensor, such as HICO, with high spectral resolution can spectrally discriminate between many benthos to deeper depths, its spatial resolution (~90 m) will significantly affect the ability to resolve individual substrates compared to HyMap or even WV2 with resolutions less than 10 m.For moderate spatial resolution sensors broad classification of pixels into key substrate components would be a more feasible approach (e.g., [17,62]).
As a final note, a sensor's radiometric calibration and the atmospheric correction can significantly affect the accuracy of the retrievals from inversion models [20], and hence on subsequent results from the HDC-algorithm as a post-inversion tool.Image pre-processing steps have been used, prior to shallow water inversion models, to minimize the effect of residual under-or over-atmospheric correction [17,63].In Klonowski et al. [17] the Rrs of an image were vertically shifted to make Rrs(750) = 0, this approach was later modified where the Rrs were subtracted by the median value between 650 and 800 nm [20].In Lee et al. [63] the Rrs were adjusted by setting the Rrs(750) to a value computed from an empirical relationship.These approaches are therefore useful if the atmospheric parameters were assumed constant for an image and the atmospheric correction produced a systematic offset.Atmospheric correction that deduce atmospheric parameters per-pixel can potentially use different aerosol models in an image and thereby cause variation in the spectral shape of the derived Rrs.If imprecise the pre-processing steps will not be able to account for an incorrect spectral shape of the Rrs.Radiometric sensor calibration can also significantly affect both the shape and offset of sensor-derived Rrs.Indeed, Lewis et al. [64] showed substantial improvement in match ups between in situ and HICO-derived water leaving radiance after vicarious calibration.Minimizing errors from calibration and atmospheric correction is particularly crucial over dark targets such as deep water or submerged aquatic vegetation.This is due to a reduced signal-to-noise ratio and where a higher proportion of the at-sensor radiance is contaminated by the atmospheric path radiance compared to bright shallow targets.Hence small errors in the calibration or deviations of the actual atmospheric parameters compared to what is used during atmospheric correction can produce negative reflectances [65][66][67].Improvements to a sensor's SNR, accurate calibration and atmospheric correction will lead to improved benthic classification from shallow water inversion models.We note that the use of the HDC algorithm to examine the ability of spectrally resolving a set of representative substrates for various environmental conditions takes into account noise introduced from inaccurate atmospheric or sunglint correction.

Conclusions
We have presented an agglomerative hierarchical clustering using linear discriminant coordinates (HDC), which when used on a spectral library of benthic endmembers outputs those classes that are optically distinguishable above the total system's variability.The HDC clustering differs from other hierarchical clustering schemes by the inclusion of the system's total variability onto each data point and the clustering based on pairwise distance and misclassification proportions.The iterative merging of classes in the HDC ceases when the remaining classes have a misclassification rate below a user specified threshold.Thus the user does not need to decide beforehand or post-clustering the number of classes to output-a common problem faced in the application of clustering procedures.Rather the total system's variability controls this outcome where a noisier system would lead to fewer output classes for a given misclassification threshold and vice-versa.In the context of benthic classification from shallow water remote sensing, variability arises from the sensor and environment (imperfect atmospheric, sunglint and air-water interface corrections) and the spectral variability observed within a given benthic species.Both sources degrade the spectral separability between benthic classes.
The HDC was applied to a set of measured reflectance spectra of 22 benthic classes collected from the Point Peron study site, Western Australia.This initial spectral library contained the sub-local scale spectral variability of each benthic class and with the addition of sensor and environmental noise taken from a HICO image, accounted for the total system's variability.Analysis of the clustering on datasets convolved to the wavebands of WV2, HyVista's HyMap and HICO sensors showed a reduction in spectral resolution reduced the spectral separability between classes and hence reduced the number of distinct classes.For instance, at WV2 and HICO wavebands 14 and 18 classes were spectrally separable just below the water's surface respectively, from the initial 22.Dendrograms from the HDC have showed that benthic species from a given genera preferentially merge first before clustering with species of other genera.
Depth and water column turbidity have been demonstrated in the literature as the most influential factor inhibiting benthic classification, where fewer classes are optically separable with increasing depth and/or turbidity.Using the measured benthic reflectance spectra, a shallow water forward model was used to simulate the water-leaving subsurface reflectance at any specified depth and water column optical property.Combined with the HDC, depth and water column specific spectral libraries were generated and subsequent analysis have quantified that the number of benthic classes decreases with increasing depth and water column turbidity at a constant clustering accuracy.Furthermore, the clustering analysis identifies those benthic classes that merge as a consequence of water column attenuation and the total system's noise.
If the irradiance reflectance spectra of the representative benthos in a scene of interest are known, then the HDC algorithm offers the ability to explore the potential of benthic classification from shallow water remote sensing, where the number and type of classes can be quantified for any given sensor, clustering accuracy, water column optical property and depth.From a management perspective, the HDC can be used to pre-determine which sensor to use so that the mapping requirements match the expected outcomes from the resultant benthic classification map.This was illustrated through the analysis of conditions that include clustering accuracy, water clarity and depth and spectral resolution that enable A. antartica to be spectrally resolved from benthic algae.Such an analysis can aid in an assessment of the feasibility of detecting phase shifts in seagrass meadows, for example, from a variety of sensors; but does not incorporate the effect of spatial resolution, which can decrease spectral distinctions.In situations where local knowledge of the spectral library is absent then the algorithm presented here cannot be applied, as the derived class names are scene specific.
Further research includes the validation of the HDC algorithm with vertical profiling optical measurements and above water radiometry.Such studies would include a comparative analysis of the clustering algorithm to optical measurement based cluster dendrograms using the rand index [68] or cophenetic correlation [69], in an analogous manner to Torrecilla et al. [47].Further development includes the application of the HDC-derived spectral libraries into physics-based shallow water inversion models.The purpose of which would be to obtain a benthic classification map with more appropriate class names that help convey the level of certainty and optical separability.It is likely that the resultant benthic classification map will be depth dependent (and hence layered), as the class names (and numbers) change with depth.Change detection is also an avenue of future research, particularly in assessing the amount of temporal change of a benthic species when it is not spectrally distinct at depth.In other words, assessing which benthic species/clusters we can detect temporal changes in.Mixed benthic assemblages typically occur at spatial resolutions greater than 2 m, how this affects benthic classification using HDC-derived spectral libraries should be investigated.

Figure 1 .
Figure 1.The average irradiance reflectance spectra of the benthos collected, convolved to the spectral resolutions of (a) HICO; (b) HyMap; and (c) WV2.

Figure 2 .
Figure 2. The HDC algorithm using the basic bottom reflectances (flowchart A) and coupled with a shallow water forward model to include a water column (flowchart B).
) show the formation of two green algae clusters (C.germinata/C.flexis and B. vestita/C.duthieae) and a brown algae cluster (C.sinuosa/S.linearifolium/S.spinuligerum) with the other benthic classes remaining optically separable.The dendrogram produced from the standard hierarchical clustering [45] (Figure 3) illustrates the potential indecision of what RMSE to cut the dendrogram and extract the relevant clusters.Here clustering continues to a single class with no indication of the optimum number of classes.At ~0.003 RMSE, 18 classes were extracted with the following clusters: (a) seagrass-green algae A. antartica/C.germinata/C.flexis; (b) mixed brown algae E. radiata/S.spinuligerum; and (c) brown-green algae S. linearifolium/B.vestita.A post clustering approach was utilized to estimate the appropriate number of classes from the standard hierarchical clustering (Figure 4).Here, nine classes was chosen as the optimum number, where the following general clusters remained: (a) two separate green algae species (Entromorpha sp. and U. australis); (b) two separate sediment classes (sediment/rubble and sediment); (c) mixed red/brown algae; (d) mixed brown, red, green algae and seagrass; and (e) a mixed red algae class.The number of classes selected by this post-clustering approach is much lower than that selected by the HDC.Moreover the HDC preferentially merges classes of the same genera first.For example the centroid hierarchical clustering merges the seagrass Posidonia sp. with the mixed C. germinata/C.flexis algae green class at the second iteration (Figure 4, RMSE ~0.0135), whereas the HDC considers these two classes optically distinct.

Figure 3 .
Figure 3. Dendrogram from the HDC of ρb spectra using HICO's spectral bands.Here the iterative selection of eigenvectors chose the first six, which produced 18 classes.

Figure 4 .
Figure 4. Dendrogram from the standard hierarchical clustering of ρb spectra using HICO's spectral bands.The top right panel shows the linkage distance vs. number of clusters, where nine classes were selected as the optimal based on the location of the knee of the curve.

Figure 5 .
Figure 5. Dendrogram from the HDC of ρb spectra using HyMap spectral bands.The iterative selection of eigenvectors chose the first seven, which produced 18 classes.

Figure 6 .
Figure 6.Dendrogram from the HDC of ρb spectra using WV2's spectral bands.The iterative selection of eigenvectors chose the first four, which produced 14 classes.

Figure 7 .
Figure 7. Dendrogram from the HDC of modeled rrs spectra using HICO's spectral bands.The water column was modeled with depth 15 m and P = 0.05, G = 0.1 and X = 0.010 m −1 .Here the iterative eigenvector selection chose only one eigenvector.

Figure 8 .
Figure 8. Number of optically separable classes vs. depth using the spectral bands and resolutions for (a) HICO; (b) HyMap; and (c) WV2.Four water-columns were modeled with varying water turbidity representing very clear to turbid waters.The exponential curves were derived through a least squares fit.

Figure 9 .
Figure 9.The estimated modeled depth (m) at which A. antartica can no longer be distinguished from green algae at various clustering accuracies for (a) HICO; (b) HyMap and (c) WV2 wavebands.Clear (water type 1) to turbid (water type 4) water clarities were modeled.The depths given here are estimates as the clustering at every 0.25 m depth increment were analyzed.

Figure 10 .
Figure 10.Number of optically separable classes vs. depth at three different clustering accuracies using the spectral bands and resolutions for HICO.The water column was modeled with P = 0.05, G = 0.1 and X = 0.010 m −1 .

Figure 11 .
Figure 11.LD coordinates of the four brown algae species using HICO wavebands.The HDC-derived dendrogram is shown on the top left corner.Here one eigenvector (z1) allowed the optimal distinction of the following three clusters: (1) S. linearifolium/S.spinuligerum; (2) E. radiata, and (3) C. sinuosa.

Table 1 .
List of symbols and acronyms used through the text.

Table 2 .
List of the benthic species collected from the Point Peron study site.