Spectral Unmixing of Forest Crown Components at Close Range, Airborne and Simulated Sentinel-2 and EnMAP Spectral Imaging Scale

Forest biochemical and biophysical variables and their spatial and temporal distribution are essential inputs to process-orientated ecosystem models. To provide this information, imaging spectroscopy appears to be a promising tool. In this context, the present study investigates the potential of spectral unmixing to derive sub-pixel crown component fractions in a temperate deciduous forest ecosystem. However, the high proportion of foliage in this complex vegetation structure leads to the problem of saturation effects, when applying broadband vegetation indices. This study illustrates that multiple endmember spectral mixture analysis (MESMA) can contribute to overcoming this challenge. Reference OPEN ACCESS Remote Sens. 2015, 7 15362 fractional abundances, as well as spectral measurements of the canopy components, could be precisely determined from a crane measurement platform situated in a deciduous forest in North-East Germany. In contrast to most other studies, which only use leaf and soil endmembers, this experimental setup allowed for the inclusion of a bark endmember for the unmixing of components within the canopy. This study demonstrates that the inclusion of additional endmembers markedly improves the accuracy. A mean absolute error of 7.9% could be achieved for the fractional occurrence of the leaf endmember and 5.9% for the bark endmember. In order to evaluate the results of this field-based study for airborne and satellite-based remote sensing applications, a transfer to Airborne Imaging Spectrometer for Applications (AISA) and simulated Environmental Mapping and Analysis Program (EnMAP) and Sentinel-2 imagery was carried out. All sensors were capable of unmixing crown components with a mean absolute error ranging between 3% and 21%.


Introduction
The spatial and temporal distribution of forest parameters are indispensable inputs into ecological models, quantifying the exchange of energy and matter between land surface and atmosphere [1].Precise monitoring of the variability in canopy cover also provides important information for forest management practices.In this context, parameters of particular interest include foliage, crown closure, leaf pigmentation and turgescence, as well as vitality, disturbances and species distribution [2].Generally, the relevant biophysical parameters can be measured directly in the field but are not easily extrapolated to large areas.Remote sensing has traditionally been used for these tasks, because they require frequent and relatively inexpensive large-scale monitoring.However, this implies that it is possible to derive the parameters of interest from air-and spaceborne data with sufficient accuracy.
For the last four decades, the normalized difference vegetation index (NDVI) [3,4] has been widely used as an estimate for vegetation density.The NDVI has often been correlated with different biophysical parameters in a variety of land-cover classes, for example with biomass in tropical landscapes [5] or with the leaf area index (LAI) in boreal conifer forests [6].However, as the NDVI saturates at a certain vegetation density, its straightforward application in regression analysis is hampered in dense and structurally complex vegetation assemblages, such as in the tropic regions [7] or temperate broadleaved forests [8].Modifications to the NDVI, such as the enhanced vegetation index (EVI) [9] or the soil adjusted vegetation index (SAVI) [10], do not sufficiently compensate for these tendencies.Thus, the relationship between vegetation indices and biophysical parameters like the LAI is negatively affected by this saturation effect [11].
Using imaging spectroscopy data, it is possible to derive more information than that contained in broadband indices.Studies using imaging spectroscopy data in forest ecosystems have proven this for the estimation of various forest stand variables (e.g., for LAI and crown volume) [12].Airborne and satellite imaging spectroscopy sensors have also been successfully used for assessing canopy biochemistry and physiology in tropical and sub-tropical forests [13,14].However, it is clear that even with hyperspectral sensors a simple relation between the mixed forest canopy signal and biophysical parameters is still difficult to establish.
Various techniques have been used to extract the relevant information from spectroscopy data.One of these techniques is multiple endmember spectral mixture analysis (MESMA).MESMA provides a potential means to circumvent the saturation issue and improve the leaf/foliage fraction estimates.This spectral unmixing technique determines the fractional cover of different materials (endmembers) within one pixel [15].MESMA has been successfully applied in various settings for vegetation and geological analysis.It differs from other spectral umixing approaches in that it allows for a pixel-by-pixel variation of endmembers.Examples of MESMA applications include the mapping of different vegetation species [16,17], the improved estimation of tree and shrub LAI in peatlands [18] and the estimation of defoliation in Eucalyptus plantations [19].
However, despite the potential of MESMA to estimate the bark fraction in addition to the sub-pixel leaf fraction, typical MESMA implementations only rely on soil, green vegetation, non-photosynthetic vegetation (a mixture of senescent plant material and bark), and shade (e.g., [20,21]).Many studies either use MESMA approaches for the mapping of different vegetation species [16,17], or they focus on the quantitative abundance of green vegetation and only unmix the vegetation and soil endmembers [22,23].For these studies, the typical vegetation endmember is a composite of both fractions, leaf and bark.
Depending on species distribution, not only soil but also bark represent a relatively high fraction of the overall reflectance within the image pixels, if viewed at nadir [24].Asner [25] stated that stem material played a small but significant role in determining canopy reflectance in woody plant canopies, especially those with an LAI smaller than 5.0.When considering intra-annual phenological change, the bark-endmember fraction becomes even more pronounced, increasing to 20% in the defoliated state [26].Nevertheless, the contribution of exposed bark to the overall reflectance of trees or canopies is rarely considered in remote sensing studies [27].As MESMA is able to determine the fractional occurrence of several components simultaneously, it allows for an improved insight into the within-crown composition of canopies, compared to conventional vegetation index approaches.
Radiative transfer (RT) modelling is also a particularly well-suited method for the derivation of quantitative information on the composition of forest canopies or on the retrieval of biophysical forest parameters, e.g., [28].Especially the more recent generation of ray tracing models, which simulate the photon interactions within the vegetation layer, can give detailed insight into the distribution of all plant parts, regardless of their visibility from the direction of observation.Due to the complex character of most RT models, analytical solutions for the parameter retrieval are rarely possible.Inversion strategies need to be adopted and often introduce errors [29].The radiation transfer model intercomparison (RAMI) activity aims at assessing the reliability of physics-based RT models under controlled experimental conditions.Recent evaluations of this initiative revealed clear discrepancies among different RT models, especially when complex canopies are examined [30].
Indeed, several applications could benefit from a detailed characterization of the different canopy components, including not only leaf and soil fractions, but also bark fractions.
Especially in hydrological modelling, parameters such as interception and stemflow are major drivers.Tree structure and the distribution of exposed bark within a canopy can dominate the water storage capacity of vegetation because these retain precipitation [31].The area-wide derivation of these parameters requires a separation of the woody biomass fraction inside canopies, which could potentially be provided by unmixing techniques.
The insufficient distinction between bark and leaf might be a further cause for the saturation effect in vegetation indices.The relation of an index to the LAI neglects the fact that the spectral signal and the subsequently-derived index include all parts of the tree (e.g., branches).Thus, the relation is instead established between a field-based LAI measurement and a remote-sensing derived proxy for the plant area index (PAI) [32], which might influence the saturation effects mentioned above.Therefore, a more detailed insight into the composition of canopy components other than leaf (especially bark) could enhance the causal connection of vegetation indices to biophysical parameters.
So far, there has been insufficient research showing the precise spectral composition of forested ecosystems.This is partly due to the fact that reference data for crown component fractions cannot be collected easily.As opposed to other types of vegetation, close range analysis and measurements of crown components are associated with a variety of technical and practical challenges due to the size of the objects and their complex structure [33].In order to get a nadir view of a canopy, devices like zeppelins [34], cherry pickers [35] or towers, situated in forest [36] have been used in experimental setups.The scale-dependent interactions between laboratory measurements of small scale objects like single leaves and of crown-or landscape-scale imagery can be well observed from towers or cranes, because both the leaf and the crown level can be reached by this means.Gong et al. [37] stated that the use of laboratory or field data for canopy level endmember characterization may be significantly impeded by differences between spectral signature determinations at the branch level.The crane-based experimental setups enable the use of the same sensor on both the canopy and canopy component levels.Thus, differences between field data and image data are circumvented.
Nevertheless, for an application over larger areas, the field measurements and unmixing techniques need to be transferable between remote sensing data sets from different sensors [38].This step includes the analytical comparison of reference data, which is presumed to represent the target value for satellite-sensor-derived products [39].This "upscaling" can show the potential of different sensors for applications based on spectral unmixing.However, this process might also be prone to errors [40].A major source of inaccuracies in the upscaling process might be due to a decreased spectral resolution or a decreased wavelength range of the applied sensor.A number of studies have shown that bare soil and non-photosynthetic vegetation cannot be easily separated in visible and near-infrared (NIR) wavelength regions (e.g., [10,21,41]).While the fractional cover of photosynthetic vegetation, non-photosynthetic vegetation and bare soils can be estimated with imaging spectroscopy sensors (e.g., [21]), the accuracy decreases with multispectral sensors (e.g., [42][43][44]).Further sources of errors in the upscaling process may be sensor-specific (e.g., radiometric resolution, signal-to-noise ratio) [45], related to ambient effects caused by illumination and observation angles, or atmospheric effects or artifacts caused by atmospheric corrections [46].More uncertainty is introduced by co-location errors, which are caused by a non-ideal fit of the reference measurement points to image pixels, e.g., [47].
For forest ecosystems in particular, there is a need to derive model-based information.Thus, this study focuses on determining the composition of structural components of the forest canopy.The close range analysis facilitated by a crane measurement platform enables a detailed determination not only of the fractions of leaf and soil, but also of bark within the canopy on four different spectral scales.To our knowledge, a comparison of the ability to quantify the fractions of within-crown components of temperate deciduous forests using spectral mixture models has not yet been evaluated for spectral sampling schemes using different sensors.
The research questions of the present study can be summarized as follows: (a) Does spectral unmixing have the potential to simultaneously estimate the fractional distribution of different canopy components, including exposed bark?(b) Does the inclusion of a bark endmember improve the accuracy of the unmixing process?(c) Can spectral unmixing methods achieve an enhanced estimate for leaf cover when compared to the NDVI in densely vegetated ecosystems?(d) What is the effect of different spectral scales on the unmixing process-Analytical Spectral Device (ASD), Airborne Imaging Spectrometer for Applications (AISA), Environmental Mapping and Analysis Program (EnMAP), Sentinel-2)?

Study Area
A deciduous temperate forest, situated in northeastern Germany at 53°55′N and 12°57′E, served as test area for this study (see Figure 1).The study site is part of the Northeastern German Lowland Observatory belonging to the TERENO (Terrestrial Environmental Observatories) initiative [48].The study area was equipped with a 40 m high Liebherr top-slewing crane to serve as the measurement platform for the close range data collection.Roughly 300 trees grow in the area below the 45 m radius of the crane's extension arm.The lower parts of the area surrounding the crane are influenced by the water regime of the river Trebel.The soil is characterized by a high humus content and is permanently waterlogged or flooded.In parts of the area covered by the crane, Alnus glutinosa (black alder) forms a very uniform alluvial forest, characterized by a homogeneous crown-height of about 25 m.There is no shrub layer present, but a very dense sedge undergrowth of about 1 meter in height is formed by Carex appropinquata (Fibrous tussock sedge).During data collection in 2011 and 2012, the LAI was approximately 4.5, thereby ignoring the LAI of the dense undergrowth.
The average annual temperature varies between 8 and 8.5 °C and average rainfall is between 550 and 600 mm.For a more detailed description of the study area, see [8].

Data
For the current study, spectral signals over four different spectral and spatial scales were used (see Table 1 for more detail): (a) a close range data set with nadir measurements obtained by a portable Full Range spectroradiometer (Analytical Spectral Devices ©, Boulder, CO, USA), (b) an airborne scene taken by the AISA Eagle and Hawk imaging spectroscopy sensors (Specim, Oulu, Finland), (c) a simulated spaceborne imaging spectroscopy scene with the spectral resolution and radiometric characteristics of EnMAP, and (d) a simulated multispectral Sentinel-2 scene.All four sensors generally encompass the 400 nm to 2500 nm spectral domain, and all images cover the area around the crane.These pixels or measurements are referred to as mixed signals in this article.The pure endmember signals were measured with the ASD spectroradiometer from the crane platform or on the ground below the crane.

Mixed Pixel Spectra-Close Range ASD Data Set
For the close range data set, 23 ASD measurements were obtained along a 35 m transect over an alder forest stand (see Figure 1).All instruments were installed on the working platform fixed to the crane hook.The platform was lifted and moved to approximately 10 m above the canopy depending on tree height (see Figure 2a).The ASD optics and a digital camera for the determination of the reference fractions were mounted to the crane platform side by side (see Figure 2b).Both devices simultaneously took spectral measurements and digital photographs.The ASD spectroradiometer was equipped with a fore-optic to reduce the field-of-view from 25° to 8° in order to minimize angular effects.These underlying conditions determined a measurement diameter of 1.4 m at crown level and 4.9 m at ground level, respectively.All measurements were performed at nadir, 15 to 10 minutes before solar noon during cloud-free conditions.A differential GPS permanently logged the current position.All measuring devices and the crane operating unit were controlled by a laptop computer installed on the platform, which performed an automatic data acquisition routine.The main water absorption bands centered at 1349 to 1439 nm, 1779 to 1979 nm and 2369 to 2500 nm were excluded from the measurements, resulting in 1727 remaining bands used in the study.This study is built upon transect measurements and field-based endmember measurements, which were taken with the same ASD spectroradiometer.The situation of having field measurements and mixed measurements from exactly the same sensor is a very convenient case, rarely occurring in practical applications.More common are studies comparing laboratory measurements to simulated imaging spectroscopy scenes or to air-or spaceborne imagery from different sensors, or studies with relatively "pure" endmembers taken directly from the imagery [49,50].The latter typically examine trees or forest stands as a single object, because the spatial resolution is not sufficient for the distinction of canopy components like single leaves and bark.

Mixed Pixel Spectra-Airborne AISA Data Set
In order to assess the suitability of the ASD-based method for remote sensing imagery, a corresponding approach was carried out with mixed pixels based on airborne and simulated spaceborne imagery in this study.On the airborne scale, an AISA scene was acquired over the area on 29 June 2011.The ground-sampling resolution was three meters.Using the GPS positions of the reference measurements, the corresponding AISA pixels were located and used as mixed pixels.Visual assessment insured an accurate fit of field and airborne data.In order to match the AISA ground-sampling resolution during the spatial scaling, three reference points had to be aggregated for each AISA pixel.This reduced the database to only eight mixed pixels.Therefore, the statistical significance of this fact has to be taken with caution.After the removal of the main water absorption features as mentioned earlier, 286 bands were remaining.

Mixed Pixel Spectra-Spectrally Simulated Spaceborne Data Sets
The AISA pixels were used as input for the simulation of the two spaceborne sensors, EnMAP and Sentinel-2.Simulated data for the forthcoming EnMAP sensor (see [51] in this special issue) represent the spaceborne imaging spectroscopy data.In addition, simulated data for the multispectral Sentinel-2 sensor were included in the study.Sentinel-2, which was launched in June 2015, is a sensor with 13 bands, particularly suitable for vegetation studies [52].
Based on the AISA data, EnMAP and Sentinel-2 L2a products were simulated using the sensor end-to-end simulation tools EeteS [53] and S2eteS [54].These tools simulate the entire image data acquisition, calibration and processing chain from spatially and spectrally oversampled data to intermediate Level-1a raw data and to the final EnMAP and Sentinel-2 products, such as Level-1b, Level-1c and Level-2a data.The data acquisition consists of a sequential processing chain represented by four independent modules, namely the atmospheric, spatial, spectral, and radiometric modules.These modules allow flexible customization of a wide range of simulation input parameters.They are coupled with a backward simulation branch consisting of calibration modules, such as non-linearity, dark current, and absolute radiometric calibration, and a series of pre-processing modules, such as radiometric calibration, co-registration, orthorectification, and atmospheric correction.A further spatial scaling to EnMAP and Sentinel-2 ground sampling size would have additionally reduced the number of reference points.Thus, the spatial resolution of the AISA scene was kept throughout the subsequent processing, and the simulation of the two spaceborne sensors is only performed in terms of spectral and radiometric characteristics.No spatial simulation was performed with EeteS and S2eteS, keeping the AISA ground sampling unchanged in the simulated image products for this study.The atmospheric simulation was performed for an image acquisition at the end of June using identical atmospheric parameters for all pixels (aerosol optical thickness: 0.2, columnar water vapor: 1.5-2.5 g•cm −2 , rural aerosol model).
In the EnMAP data set, the main water absorption bands from 1335 to 1418 nm, 1783 to 2007 nm and 2370 to 2439 nm were excluded from the measurements, resulting in 200 remaining bands used in the study.From the Sentinel-2 data, the three 60 m resolution bands, which are mainly dedicated to atmospheric correction and cloud screening, were excluded from the study (bands centered at 443 nm, 940 nm and 1375 nm) [52].See Table 1 for an overview of the characteristics and origins of the data.

Endmember Spectra
The spectral reflectance of all endmembers was measured in the field with the same portable ASD spectroradiometer, which was used for the mixed pixel acquisition.The spectrometer fore-optic was positioned 10 cm above the measurement targets, resulting in a 1.4 cm sensor field-of-view.Soil reflectances as well as a share of the leaf and bark reflectances from the lower parts of the forest stand were collected in the area below the mixed-pixel transect.Most trees in the study area have multiple trunks with an inclined growth in the lower parts.Bark spectra from tree trunks could be collected at nadir from these inclined parts.Most leaf and additional bark reflectances of the exposed branches were collected from above the top of the canopy with the help of the crane working platform.For this purpose, the measurement platform was moved directly above or into the trees.To avoid shading effects, the ASD fore optic was removed from its fixed position on the measurement platform, and the handheld ASD pistol-grip was used for endmember collection.All measurements were acquired at nadir during the years 2011 to 2014.To avoid phenological differences, only measurements from appropriate dates around June were used for this study, resulting in a set of 168 spectral reflectances for the three endmembers.As all endmember reflectances were taken in field, nonlinear mixing effects are inherently included in the spectral curves.The effects of shading and shadowing are included as well and involve brightness difference in the endmember reflectances (see Figure 3).The three endmember classes are characterized by a relatively high interclass spectral similarity.All classes show typical vegetation attributes.Due to the high humus content and undecomposed plant parts, even the soil reflectance partly resembles the vegetation reflectance.The dark humic soil is partially flooded with shallow water, causing a very low spectral reflectance.Laboratory measurements of soil probes ensured that the low overall brightness of the soil spectra is not due to measurement mistakes.
As the crown components never occur as whole pure pixels, there was no spatial aggregation of the endmembers in the scaling process to the airborne and simulated spaceborne data sets.AISA, EnMAP and Sentinel-2 spectral resolutions were simulated from the field-based ASD measurements using the specific spectral response functions.There was no further spatial aggregation of the information contained in the endmember spectra.

Reference Fractions
In order to determine the reference ground cover fraction of the endmembers in each measurement, digital RGB photographs were taken simultaneously with the ASD transect measurements.The camera was a Canon EOS Mark II, a digital single-lens reflex camera with approximately 21.1 effective megapixels.The large single-plate CMOS sensor provided 14 bit RAW images.The camera was equipped with a fast lens (Canon EF 85 mm f/1.8 USM (ultrasonic motor)) to limit exposure times during possible swaying of the crane.For depth-of-field, the aperture was adjusted to 5.6 to ensure sharpness of the image on crown and ground level.Both devices were installed on the crane measurement platform side-by-side and adjusted to point in the same direction, to insure a sound match between ASD measurements and digital nadir images.Both optical axes were orientated parallel to each other at a distance of 7 cm.The camera object lens captures a 24° field-of-view (see Figure 3).
The photographs were reduced to a circular area matching the 8° field-of-view of the ASD measurements.They were segmented and classified semi-automatically to estimate the fractional abundances of leaf, bark, and soil (see Figure 4).A supervised hierarchical approach was chosen for the classification [26].The feature information used for the defined objects included mean layer values, ratios of the layer values, area, asymmetry, compactness, elliptic fit, roundness, and shape index.The components "sunlit leaf", "shaded leaf", "undergrowth", "bark (stem)", "branch", "lichen" and "soil" were first classified separately and then combined to three classes: "leaf", "bark" and "soil".Subsequently, obviously misclassified segments were manually assigned to the appropriate class.As the point-spread function of the ASD can be roughly approximated by a two-dimensional Gaussian model, the sampling intensity is not stable over the whole field-of-view [55].This can introduce errors when comparing the photographic fractions to the fractions modelled from ASD measurements.In order to handle the non-uniform measuring intensity across the field-of-view of the ASD measuring device, all photograph pixels were weighted with a Gaussian filter resembling the ASD measuring intensity (see [35]).
Total leaf fractions are very high in all 23 measured mixed pixels of the transect (see Figure 5a), which vary between 72% and 91%.The lowest fractions correspond with larger crown gaps at the beginning and at the end of the transect.Bark has the second highest share with a relatively stable fraction of 7% to 13% (average 10%).The soil fraction has the lowest share with an average of 6%, but it shows a higher variability (1% to 18%).
In accordance with the rescaling of the mixed pixel spectra to the airborne and simulated spaceborne scale, the reference points had to be aggregated to match the AISA pixels.The reference points located in each AISA pixel were averaged and assumed to represent the pixel.In comparison to the ASD reference data, the spatially aggregated data set still contains a similar range of values (see Figure 5b).

Spectral Unmixing
This study uses a spectral unmixing approach for the estimation of leaf, bark and soil components in a closed forest canopy.A well-established approach is Spectral Mixture Analysis (SMA) [56].It models a mixed spectrum (r) as the sum of the pure spectral signatures of its constituent components (i.e., endmembers), each weighted by their a priori unknown sub-pixel fractional coverage as follows: where M is a matrix containing N columns corresponding to the number of endmembers, ƒ is a N*1 column vector containing the subpixel cover fractions of each of the N endmembers and ϵ is the residual error.With a set of appropriate endmembers and their pure spectral signatures, sub-pixel cover distribution maps can be generated by model inversion [57].
A drawback of the SMA technique, however, is that it fails to account for any endmember variability between the pixels of a scene.To account for this, Multiple Endmember Spectral Mixture Analysis (MESMA) was introduced [15].Instead of one fixed spectral signature for each endmember, endmember libraries are used.This iterative technique evaluates different possible endmember combinations.The model with the lowest root mean square error (RMSE) between modelled and mixed signals is assigned to each pixel individually.
In this study, spectral unmixing of close-range spectral measurements was carried out using MESMA (see Figure 6).
The study aims at evaluating the different impacts of the soil and bark endmembers on the unmixing process in a very densely vegetated forest ecosystem.Thus, the unmixing process was carried out using different combinations of the endmembers "leaf", "bark" and "soil".To account for brightness differences, a shade endmember with zero percent reflectance over all wavelengths was included (e.g., [58]).The resulting combinations had three different levels of complexity: one two-endmember-model (leaf and shade (L-Sh)), two three-endmember-models (leaf, soil, shade (L-S-Sh) and leaf, bark, shade (L-B-Sh)), and one four-endmember-model (leaf, bark, soil and shade (L-B-S-Sh)).All possible combinations without a leaf endmember are not further considered in this study, as "leaf" represents the dominant fraction in all measurements.The reference fractions always sum to 100%.As the shade endmember is not included in the reference data, a shade normalization was carried out.MESMA was run with a constraint of non negativeness and a maximum allowable fraction of 100%.In order to evaluate the different model types, MESMA was first run with fixed levels of complexity and fixed types of models (e.g., L-B-Sh only).In a next step, in order to fully exploit the capabilities of the MESMA technique, MESMA selected the type of models and the level of complexity for each pixel individually.In three attempts, models with only two endmembers (2 EM), models with two or three endmembers (2-3 EM) and models including all endmember combinations (2-3-4 EM) were tested.The selection criterion was the lowest modelling error, calculated as the lowest RMSE between measured and modelled spectral reflectance.
The selection of a more complex model over a simpler model generally leads to a decreased modelling error, because mathematical models composed of more endmembers will generate lower RMSE values.As the modelling error decreases, the fractional error in the MESMA results will often increase at the same time, as endmembers may be included which are not actually present in the pixel [59].To overcome this phenomenon, multilevel fusion thresholds were introduced.These determine the minimal improvement in RMSE for selecting a model of a higher complexity level.The two opposing trends of modeling error and fractional error were not visible in the ASD data set, presumably due to the fact that all endmembers were known and no false endmembers could be included.Therefore, the multilevel fusion threshold was set to zero in this study.
This unmixing process was used to model the endmember fractions of the 23 mixed pixels in the ASD transect over the alder forest.The results were validated using fractions derived from nadir photographs taken simultaneously with the ASD measurements.The relation between modelled and observed fractions of leaf, bark and soil endmembers were examined using the mean absolute errors.Finally, for comparison of the complex unmixing result, the results of a simple application of the NDVI were used as a benchmark.
In the present study, the MESMA routines implemented in the Visualization and Image Processing for Environmental Research (VIPER) Tools 2 were applied.Thesetools were originally developed at the UC Santa Barbara's Geography Department in California and rewritten at the KU Leuven Department of Earth & Environmental Sciences (Belgium).Starting at version 1.5, the VIPER Tools software is licensed under the GNU General Public License (GPLv3).

Close Range Data Set
In Table 2, the mean absolute errors of the three different endmember models are shown as an estimate of accuracy for the ASD-based unmixing.The mean absolute error measures the average magnitude of the errors without considering their direction or the variability of their magnitude.A value of zero indicates a perfect fit between reference and modelled data.The mean absolute error can exceed 100%.The overall error averages the mean absolute errors of all three crown components, presuming that they are present in the related model.The modelling RMSE describes the error between modelled and observed spectra, regardless of the accuracy concerning the reference fractions.
The fractions for all three crown components could be unmixed with a good accuracy.A clear trend is visible and the inclusion of more endmembers into the unmixing process enhances the results.The error for the leaf fraction decreases from 11.9% in the L-Sh model to 8.1% in the L-B-S-Sh model.Bark and soil fractions can be determined with an error of 5.9% and 6.9%, respectively.
The scatterplots of Figure 7 give more detailed insight into the distribution of the modelled and observed fractions.The leaf endmember was sufficiently well-modelled with all four combinations of endmembers, but shows a clear underestimation of some values in the L-Sh model.As soon as the bark or soil endmember was introduced, this effect disappeared.The bark endmember showed very stable results in both models.The soil endmember was generally slightly overestimated but showed an improvement in the more complex model.The shade endmember could not be depicted here as there were no corresponding reference values.The shade fractions showed a high correlation with the soil fractions.Table 3 summarizes the results of the models, including different complexity levels.In addition to the fractional errors for all crown components and the modelling errors, the models which were finally selected by the MESMA algorithm are listed.The 2 EM model incorporates only the leaf and shade endmember and is identical to the L-Sh model presented in Table 2.
Again, the more complex models showed an increased accuracy.Only marginal improvements can be expected by detecting the optimal multilevel fusion threshold.Compared to the straightforward calculation of broadband indices for the estimation of leaf coverage, the application of the MESMA algorithm requires a marked additional effort.In order to evaluate whether this effort is justified, the NDVI was calculated for all ASD transect measurements.There was no further regression established between the NDVI and the leaf cover fraction.The NDVI was directly compared to the results of the MESMA model, including all complexity levels (see Figure 8).In this figure, a good fit of reference and modelled data would be indicated when all values were situated close to the diagonal one-to-one line.This would, in turn, cause the trend line to be orientated similarly to the diagonal, which was only the case for the MESMA results.Although the NDVI results have a better correlation to the trend line, the orientation of the trend line indicates saturation.Even in the dense vegetation structure, characterized by high fractional abundances of the leaf endmember, the unmixed fractions still showed a clear trend, whereas the NDVI was saturated and showed no dynamics.

Airborne and Simulated Spaceborne Data Sets
Although the unmixing results of the data sets based on the AISA overpass show decreased accuracies, the same patterns and trends are visible; the leaf endmember could be unmixed with accuracies ranging between 7.6% and 21.7%, the bark endmember could be unmixed with a mean absolute error between 3% and 8.5%, while the soil endmember resulted in higher errors (6.3% to 20.1%).Here again, the inclusion of additional endmembers enhanced the accuracies.In Table 4, the mean absolute errors of all mixture models for the airborne and simulated spaceborne data sets are listed.Figure 9 summarizes the overall mean absolute errors on the different levels of complexity for all four sensors.As expected, the transfer from the close range data set to the airborne data set and the further simulation of spaceborne data resulted in increased errors.The results achieved with the simulated Sentinel-2 wavebands make an exception to this trend.Generally, the fractional overall error decreased with an increasing complexity.Sentinel-2 data achieved the lowest errors of the whole study for the L-B-Sh model.The multispectral sensor outperformed even the original ASD measurements.
Nevertheless, in the Sentinel-2 data set, the inclusion of a fourth endmember did not change results.No L-B-S-Sh model was selected by the MESMA algorithm.In the EnMAP data set, the inclusion of a four endmember model even decreased the accuracy.

Discussions
The present study investigates the potential of MESMA to derive sub-pixel crown component fractions in a temperate deciduous forest ecosystem with a complex structure and a high proportion of foliation.This setting includes high intra-class variability and a relatively low inter-class variability, compared to other studies differentiating green vegetation and soil in open canopies.

Close Range Data Set
The mean absolute error for the estimation of leaf fractions ranged between 6.3% and 11.9% for the different endmember combinations, which was a reasonably good result for such a challenging setting (see Table 2).The inclusion of bark and soil endmembers enhanced the results, compared to a model considering leaf and shade only.The shade endmember alone was not capable of explaining the variation in the leaf fractions.As depicted in Figure 7, leaf fractions were often underestimated when only leaf and shade endmembers were included.The inclusion of additional endmembers clearly enhanced the results.This trend could also be observed in the models incorporating different levels of complexity.The best results were achieved with the 2-3-4 EM model, which was capable of describing the given setting more realistically, and, therefore, was more accurate than the models of a lower complexity.
The experimental setup allowed for a very good fit between reference data, mixed pixel data and endmember data; thus, only minor errors may have been induced by sensor-specific, pre-processing and environmental conditions.Inaccuracies due to the different spatial extents of the reference nadir photographs and the spectral measurements were kept marginal, as both instruments were mounted on the crane side by side and were automatically triggered simultaneously.Nevertheless, the ASD instrument had an integration time of around three seconds, so minor sway of the crane could have caused differences in the measuring location of both instruments.However, no sensor-specific deviations could have occurred, as the same instrument was used for endmember collection and mixed pixel collection.A major source of errors may have been the fact that the endmembers were collected on similar dates, but for different years.Classification errors of the validation data set also need to be taken into account.
The contribution of exposed bark to the overall reflectance of trees or canopies has rarely been considered in remote sensing studies.For plant canopies with an LAI below 5.0, stem material was shown to have a small but significant influence on the canopy reflectance [25].In the presented study, the spectral reflectance of bark was measured in the field with the help of a crane measurement platform.This insures a very good comparability between endmember and mixed pixels.Furthermore, the utilization of a crane facilitated the accurate estimation of the reference data, which has not been feasible for the majority of studies.
In the forest canopy observed in this study, the fractional share of the bark endmember was around 10% on average, whereas soil mostly had a lower share visible at nadir.With phenological changes over the year, the fraction of the bark endmember becomes even more pronounced, increasing to as high as 20% in the defoliated state [26].
The results demonstrate that the inclusion of a bark endmember is an advantageous step for the unmixing of fractional abundances in dense deciduous forest ecosystems.At the moment, the mixed signal of bark, soil, shade, and leaves is being used for several applications.Firstly, the mixed signal is used to estimate the bio-physically active biomass of forest (e.g., for the estimation of biomass in the Reducing Emissions from Deforestation and Forest Degradation (REDD) process).However, the underlying signal varies.In particular, the amount of bark changes with tree-type, phenological season, or age and structure of the forest [26].An estimation of the share of bark from the signal could be used to correct for these errors.Secondly, the amount of bark within a forest is a very important parameter for different types of hydrological modelling.The share of stems and branches influences relevant factors, such as stemflow or interception.Presently, hydrological models integrate these factors using very rough estimations [60].The first attempts to integrate remote sensing have shown enhanced modelling results, but still fail to integrate the amount of bark [61].Therefore, the inclusion of the endmember class bark into forest unmixing processes, as shown in this study, is a necessary step to provide high quality parameters for further applications.However, the mean absolute error of 5.9%, as shown here, is still high in the context of the relatively low average bark fraction of 10%.
Many studies report a decreasing accuracy for the separation of soil and non-photosynthetic vegetation, when cover fractions were particularly low (e.g., [47]).This effect can also be observed in the distribution of errors for the soil fractions in this study.
Other studies have reported the main confusion to occur between non-photosynthetic vegetation (e.g., dry grass, leaf litter, and woody material) and soil [20].A high similarity between endmembers leads to a high correlation between the endmember spectra, which in turn results in an unstable inverse matrix and a decrease in accuracy [46].As the soil is very humic and rich in dead plant material, it bears some features typical for non-photosynthetic vegetation.Even more striking is the similarity to the shade endmember with zero percent reflection over all wavelengths.During the collection of mixed pixels as well as during the collection of endmember spectral reflectances, the soil was waterlogged or partly flooded with shallow water, which led to heavy photon absorption and corresponding low reflectance.
Modelling results also depend on the location of the crown component within the canopy [41].Hence, soil contributes less to the overall signal as it is mostly shaded, whereas bark in the upper parts of the canopy is more likely to be sunlit and less exposed to nonlinear mixing effects.The similarity of soil to shade in this study may have enhanced the unmixing results, especially for the model only including leaf and soil.At the same time, the similarity might have led to confusion between soil and shade endmembers, which in turn caused the slight decrease in soil accuracy after shade normalization.
It could also be demonstrated that the NDVI, representing vegetation indices in general, could not reproduce the range of fractional occurrences of the leaf endmember.Due to the high leaf density, the index saturated (see Figure 8).In contrast, MESMA was still capable of discriminating different levels of foliation cover even in a dense forest canopy.

Airborne AISA Data Set
In the second part of the study, the same MESMA models were applied to mixed pixels of different airborne and simulated spaceborne sensors.As expected, the unmixing accuracy generally decreases (see Figure 9).The transfer from ASD to AISA is the only step in this study, which includes a spatial scaling.Due to the averaging of reference points located in the same AISA pixels, the reference data was reduced from 23 ASD measurement points to only eight cases during the upscaling process.As a consequence, the validity of the statistical evaluation of the upscaling results is limited and has to be considered with caution.Nevertheless, for future applications, it is still valuable to consider the transferability of the study's results to remote sensing imagery on different scales.
All reference points on the transect were thoroughly checked for their location in the correct AISA pixel to keep co-location errors marginal.However, the ground sampling area of the AISA pixel was still approximately three times greater than the reference sampling area, leading to a possible source of error.Another major source of error might have been the time shift of one year between the collection of reference data and the AISA image acquisition.Nevertheless, in comparison to many other vegetation types, the structure of the alder stand can be assumed to be fairly stable after this time span, as the tree morphology remains relatively unchanged and the total amount of foliage was similar between the two dates.The LAIs of both dates were matching and the distinct positions of the leaves on the branches cannot be considered to be crucial for the estimation of sub-pixel fractions.The soil was waterlogged on both dates; therefore, no soil moisture differences can have altered the reflectance.The dense undergrowth formed solely by Carex appropinquata can be considered as very stable, as the spatial expansion is mainly influenced by the local relief.The endmember library included reference measurements of both years and can be assumed to be equally appropriate for both dates.
Further reasons for decreased accuracy may include radiometric and spectral differences, artifacts caused by the pre-processing of the AISA data and the signal to noise ratio, co-registration errors between the AISA Eagle and Hawk sensors, and further changes in the surrounding conditions.Furthermore, since the transect was not located in the center of the AISA flight strip, bidirectional effects could have altered the proportional visibility of the endmembers.This effect was not considered in the reference data.
In spite of all potential sources of error, which usually accompany the transfer to airborne data, the overall mean absolute error increased only slightly, to values ranging between 8.9% and 14.3% compared to 7.0% to 11.9% for the ASD data set.Additionally, the L-B-Sh model performs slightly better than the L-S-Sh model, possibly due to the similarity between soil and shade in this data set.The model including both soil and shade might have been impaired by their similarity, whereas the L-B-Sh model might even have been enhanced by this similarity, as the missing soil component was partly compensated by the shade endmember.This trend is visible in all three data sets based on the AISA imagery.As the ASD results do not show this trend, it might be caused by differences between the two sensors.In contrast to the ASD dataset, the four endmember model (L-B-S-Sh) showed a slightly decreased overall accuracy when compared to the less complex L-B-Sh model, but it still performed better than the L-S-Sh model.Here, the additional inclusion of bark was beneficial, whereas the inclusion of soil introduced an increased error.If MESMA was run without a fixed model type, the 2-3-4 EM model slightly outperformed the less complex models (see Figure 9).

Spectrally Simulated EnMAP and Sentinel-2 Data Sets
Compared to the AISA results, the accuracy further decreased to values between 11.8% and 21.7% in the EnMAP unmixing results.As the EnMAP data were simulated from the AISA data, no further co-location errors or time shifts could have caused this decrease.Again, the L-B-Sh model performed best, but all models of a high complexity outperformed the simple L-Sh model.Like in all other data sets, the inclusion of bark improved the results more than the inclusion of soil.In the EnMAP data set, the decrease of accuracy in the step from three to four endmembers becomes more pronounced than in the AISA dataset.A more suitable multilevel fusion threshold should at least be able to stabilize the error on the level of the less complex 2-3 EM model.
As the spatial extent of the reference area is very limited, the further transfer to simulated EnMAP and Sentinel-2 data was carried out in terms of spectral and radiometric resolution as well as simulated atmospherical characteristics, but not in terms of ground-sampling size.
Theoretically, a further spatial downscaling to 10 m or 20 m for Sentinel-2 or to 30 m for EnMAP would not be expected to decisively alter the results.Numerous studies, examining the effect of spatial scale on remote sensing applications (e.g., [62,63]), conclude that the impact of spatial scale originates from the relation between the size of the objects under observation and the size of the pixels (e.g., [40]).The pixel size offering the highest information content for the various applications is often determined by means of statistical analysis of the local variance of the pixel neighborhood.When the spatial resolution is considerably finer than the objects of interest, the neighboring pixels become highly correlated.This is most likely due to their depiction of the same object, which causes the local variance to be low.A similar homogeneity is induced when the spatial resolution increases until many ground features become aggregated in a single pixel.The peak of local variance between these two cases occurs when the pixel size is close to the dimension of the objects.This situation was found to correlate with the optimal spatial resolution for classifying the specific objects [40].
In this study, the peak in local variance can be expected to be related to the size of a single leaf or branch, and any coarser scales can be assumed to imply a loss of information.This effect is assumed to be already present in the ASD measurements.Provided that the distribution of crown components remains fairly homogeneous over large areas, the variance can be assumed to stabilize over a range of spatial resolutions, causing the further loss of information to be marginal.With the decrease to Sentinel-2 or EnMAP spatial resolution, spectral detail is progressively aggregated, causing the variance between pixels to be lower than in the ASD mixed pixels.Nevertheless, the mean value should remain stable over all spatial scales [40], and unmixing results based on these pixels should not be impeded noticeably.As the objects are already far smaller than the pixel sizes or the field-of-view of the ASD, the intention of this study is not a classification of the objects but rather a subpixel estimation of fractional abundances.The discussion of the influence of scale on subpixel unmixing approaches has not been comprehensively concluded.Thus, the predominant effects, which occur when the coarser pixels exceed the boundaries of the alder forest, cannot be stated with sufficient reliability.
Although the spatial scaling process might not distinctly decrease the accuracies for the fractional abundances, major errors could be introduced by the fact that the reference data may not sufficiently represent the mixed pixels.In the case of an EnMAP pixel, the validation data set only covers less than one 30th of the mixed pixel.The detailed determination of crown component fractions for validation does not seem feasible for the extent of one or even more EnMAP pixels.A complete match between the extents of the mixed pixels and reference data, as is the case for the ASD data set in this study, is therefore not possible.However, with the EnMAP data set the leaf fractions were underestimated by 10% to 25%.Based on personal field observations, it can be assumed that this underestimation is mainly due to modeling error and not to a validation error.
The multispectral Sentinel-2 data were also simulated from the AISA imaging spectroscopy data and, therefore, bears all errors involved in the AISA data.Nevertheless, the Sentinel-2 results were comparable to those of the ASD data set and, in some cases, even outperformed the ASD models.The three-endmember model, including leaf, bark and shade, achieved the best results of all data sets, with an error of 7.6% for the leaf endmember and 3.0% for the bark endmember.This is remarkable, as only 10 Sentinel-2 bands were used in this study in contrast to 1727 ASD bands.This positive trend cannot be explained by the simulated rescaling process itself, as it bears no benefit in terms of atmospheric conditions, illumination, observation angle, artifacts or radiometric and spectral resolution.The successful results of the Sentinel-2 models could therefore indicate that a reduction in wavebands may be advantageous for the unmixing process.This would be in accordance with other studies that increased the unmixing accuracy through waveband selection techniques, focusing on absorption features [64] or zones indifferent to changes in spectral variability [65].Data reduction has frequently been used in image spectroscopy data processing to facilitate an efficient analysis and to improve feature extraction [66].A considered selection of wavelengths, which is focused on the signal and minimizing the noise, could thus significantly improve SMA accuracy [21].Evidently, the Sentinel-2 bands are well placed for this kind of evaluation.The arrangement of the Sentinel-2 bands highlights the red-edge wavelengths, and exclude some spectral regions less relevant for the discrimination of vegetation.The spectral unmixing results benefitted from this reduction, even though co-location errors and atmospheric noise were added.In a brief test, MESMA was run with only 10 of the 242 EnMAP bands, each positioned in the center of a Sentinel-2 band.The accuracy clearly increased to a level similar to the Sentinel-2 results.This emphasizes the potential or even the necessity of a waveband reduction for this kind of evaluation.In contrast, the reduction of the ASD bands to only ten bands decreased the accuracy.The main reason for these opposing effects may be found in the high precision of fit of the ASD dataset (e.g., same sensor and date), which is not given for the air-and spaceborne imaging spectroscopy data sets.The same effect can be observed when the EnMAP wavelengths are resampled (not simulated) directly from the ASD data.With this setting the EnMAP fractional errors were comparable with the original ASD results, but a reduction to the 10 Sentinel-2 wavebands decreased the overall accuracy from 7% to 10%.
Nevertheless, in contrast to the ASD and AISA data sets, the Sentinel-2 results showed additional errors linked to the increase in complexity from the three-endmember to the four-endmember models.This also led to the stagnation of the positive trend visible in Figure 9. Here, the comprehensive information regarding the imaging spectroscopy data sets was advantageous for the more complex unmixing model, whereas it was not necessary for the less complex models.

Conclusions and Outlook
This paper proposes the utilization of a bark endmember for spectral unmixing of crown components in dense deciduous forest ecosystems.The share of bark within the spectral signal is a valuable contribution to a more accurate retrieval of biophysical parameters.Moreover, the amount of bark is an important factor for hydrological modelling, as it relates to variables such as stemflow or interception.This study evaluates the different impacts of soil and bark endmembers on unmixing models and demonstrates that the inclusion of bark and soil endmembers markedly improve their accuracy.Regardless of the sensor used, not only could the bark and soil endmembers be included without any loss of accuracy compared to leaf and shade models, but the inclusion of these additional endmembers also enhanced the accuracy of the leaf fraction estimates.The overall absolute error of the ASD based MESMA model was approximately 11% for the model including leaf and shade endmembers only.The error could be decreased to 8.7% by the inclusion of a soil endmember or to 7.1% by instead including a bark endmember.A MESMA model with the optional inclusion of all four endmembers performed best with an overall error of 6.9%.
Even though the given vegetation structure of this study was characterized by high fractional abundances of the leaf endmember, the unmixed fractions still followed a clear trend, whereas the NDVI saturated and showed no dynamics.The more complex unmixing method is therefore preferable over the straightforward application of broadband vegetation indices like the NDVI.Moreover, the unmixing approach has the further benefit of being able to simultaneously estimate the sub-pixel fraction distribution of additional canopy components.Still, the MESMA models do not include the full physically-based interactions within a crown canopy structure.Future work could include the advanced parameterization of an RT model such as the Invertible Forest Reflectance Model (INFORM) [67] to incorporate the empirical findings in an enhanced RT model parameterization.
The wet and humic soil present in this study is especially challenging, as it resembles both the shade endmember as well as the leaf endmember.This issue could be less pronounced in settings with brighter and less humic soils.Overall, the utilization of the shade endmember is still an open question that needs to be further examined.To contribute to the transferability of the results, other particularly dense forest ecosystems should be examined.
The transfer to airborne and simulated spaceborne data sets revealed two general trends; for the less complex three-endmember models, the spectral and radiometric reductions led to better results, whereas the four-endmember models were more dependent on higher spectral resolutions.The best models for the air-and spaceborne data sets were the models including leaf, bark and shade.With the LBSh models, the overall error was 8.3% for the AISA data set, 11.8% for the EnMAP data set and only 5.3% for the Sentinel data set.The successful performance of sensors with fewer bands in this study suggests that the unmixing of densely vegetated ecosystems could benefit from a waveband selection technique.This is in good accordance with the conclusions of previous studies.
The models based on airborne and spaceborne data sets benefitted from a waveband reduction, whereas the ASD based models did not.In order to improve the application on remote sensing imagery, the factors influencing these effects should be further identified and examined.
As remote sensing data include more observation angles than nadir, future research should evaluate how the composition of the spectral signal changes, if angular measurements are unmixed instead of nadir measurements.

Figure 1 .
Figure 1.Mixed pixel data at different scales (RGB).(a) Aerial image with crane and alder ASD-transects marked as asterisks and.(b) AISA subset with transect.White circles represent ASD fields of view on upper crown level.

Figure 2 .
Figure 2. Crane measurement setup.(a) Crane with measurement platform and 8° field-of-view of the ASD fore-optic and (b) Camera and ASD fore-optic mounted on movable table attached to the platform.

Figure 3 .
Figure 3. Maximum, mean and minimum spectral signatures of all endmember classes derived from ASD field measurements; bark (left), leaf (middle) and soil (right).

Figure 4 .
Figure 4. (a) Typical digital nadir photograph of an ASD-measurement, (b) associated classified photo-segmentation and (c) derived percentages of the components and the corresponding ASD-measurement.

Figure 5 .
Figure 5. Fractional abundances of all three endmembers (a) along the 23 mixed pixels of the ASD transect and (b) along the AISA mixed pixel transect.

Figure 6 .
Figure 6.Flowchart of the study design.

Figure 7 .
Figure 7. Scatterplot with one-to-one line and modelled and observed fractions for different endmember combinations (in columns) and the different endmembers (in rows) for the ASD data set.No fraction can be estimated for some combinations, because the endmember classes were not part of the associated model.

Figure 8 .
Figure 8. Scatterplot with one-to-one line (solid) and trend line (dashed) of the reference leaf fraction against modelled fractions (a) for the MESMA model including all three complexity levels and (b) for the NDVI.

Figure 9 .
Figure 9. Overall mean absolute error for different levels of complexity for all four sensors.

Table 1 .
Characteristics and origins of data.

Table 2 .
Mean absolute errors and modelling root mean square error (RMSE) of analytical spectral device (ASD) models.

Table 3 .
Mean absolute errors and modelling RMSE of ASD models on different complexity levels.

Table 4 .
Mean absolute errors and modelling errors of airborne imaging spectrometer for applications (AISA), environmental mapping and analysis program (EnMAP) and Sentinel-2 models.