Characterization of Ecological Exergy Based on Benthic

The evaluation of ecosystem health is a fundamental process for conducting effective ecosystem management. Ecological exergy is used primarily to summarize the complex dynamics of lotic ecosystems. In this study, we characterized the functional aspects of lotic ecosystems based on the exergy and specific exergy from headwaters to downstream regions in the river’s dimensions (i.e., river width and depth) and in parallel with the nutrient gradient. Data were extracted from the Ecologische Karakterisering van Oppervlaktewateren in Overijssel (EKOO) database, consisting of 249 lotic study sites (including springs, upper, middle and lower courses) and 690 species. Exergy values were calculated based on trophic groups (carnivores, detritivores, detriti-herbivores, herbivores and omnivores) of benthic macroinvertebrate communities. A Self-Organizing Map (SOM) was applied to characterize the different benthic macroinvertebrate communities in the lotic ecosystem, and the Random Forest model was used to predict the exergy and specific exergy based on environmental variables. The SOM classified the sampling sites into four clusters representing differences in the longitudinal distribution along the river, as well as along nutrient gradients. Exergy tended to increase with stream size, and specific exergy was lowest at sites with a high nutrient load. The Random Forest model results indicated that river width was the most important predictor of exergy followed by dissolved oxygen, ammonium and river depth. Orthophosphate was the most significant predictor for estimating specific exergy followed by nitrate and total phosphate. Exergy and specific exergy exhibited different responses to various environmental conditions. This result suggests that the combination of exergy and specific exergy, as complementary indicators, can be used reliably to evaluate the health condition of a lotic ecosystem.

Exergy, which was introduced into ecology at the end of the 1970s [13,14], represents the difference in the concentrations of specific components between the studied system and a reference state. However, we cannot precisely calculate the total exergy of an ecosystem because it is impossible to measure all of the possible contributions to exergy in a particular ecosystem [15]. Thus, exergy only provides a relative value in comparison with a predefined reference system [16,17]. However, changes in exergy (e.g., comparisons between two different structures) can be indicative of alterations in ecosystem function because exergy, as a holistic indicator, expresses the degree of development and complexity in the ecosystem of interest [17]. For this reason, exergy can be used to assess ecosystem health [18,19]. Exergy represents the total information carried by the biomass, and specific exergy is the relative exergy carried by the units of matter [20]. More developed organisms have more biological information content per unit, and they usually inhabit higher trophic levels and more complicated food webs than do less developed organisms. Therefore, a higher specific exergy is expected in a more complicated ecosystem with more efficient niche utilization [21][22][23]. If the total biomass is constant in a certain system, then variations in exergy largely depend on its structural complexity.
Among the various groups of freshwater organisms, benthic macroinvertebrate assemblages are highly diverse in most streams and play important roles in both bottom-up and top-down processes, such as detritus processing, animal-microbial interactions, herbivory and energy transfer to consumers at higher trophic levels [31,32]. Benthic macroinvertebrates employ a wide variety of structural and functional responses to various factors (e.g., longitudinal changes in energy sources, food webs, nutrient input and pollution gradients) [33] because macroinvertebrate species possess diverse life cycles and generation times and represent several trophic levels and different feeding behaviors. Exergy has typically been used in benthic studies [10] because it has been a useful functional indicator for the evaluation of the impacts of human disturbance on benthic communities [34].
Structural indicators based on community composition (e.g., species richness, Shannon diversity, Simpson diversity and evenness) have been used as surrogates for functional indicators (e.g., metabolism, productivity, emergy and exergy). The patterns (structure) determine the processes (function) and the processes in turn influence the patterns [35]. Structural indicators alone, however, cannot detect all types of impairment [36]. Thus, true functional studies should be developed and included to assess the health and integrity of running-water ecosystems. Previous studies of general system patterns, in terms of their functional and structural aspects, produced contradictory results due to differences in environmental conditions, such as forest composition [37][38][39], pH gradient [38], water abstraction [40,41], temperature [42], nutrient gradients [42] and land-use types [43][44][45]. Therefore, improving our understanding of the relationships between structural and functional indicators will support the monitoring and assessment of lotic ecosystems [46].
The aim of this study is to characterize changes in the functional aspects of lotic ecosystems in terms of exergy and specific exergy using benthic macroinvertebrate communities from the headwaters to downstream. Accordingly, we assume that exergy and specific exergy respond differently to the effects of various environmental characteristics and hypothesize that: (1) the patterns of exergy and specific exergy differ along an upstream to downstream hydrological gradient (i.e., river width and depth gradient), and (2) high nutrient loads can cause exergy to increase whereas specific exergy to decrease.

Ecological Data
Macroinvertebrate data were extracted from the "Ecologische Karakterisering van Oppervlaktewateren in Overijssel (Ecological characterization of surface waters in the province of Overijssel: EKOO)" database in The Netherlands [47]. The EKOO database was compiled from 1981 through 1985 to improve aquatic ecosystem management through the use of water typology to establish a baseline course, including a total of 650 sampling sites classified into 42 different water body types representing both lotic and lentic conditions (e.g., acid springs, springs, temporary streams, upper courses, stream pools, middle courses, lower courses, canals, ditches, lakes, ponds and moorland pools). From the 42 water-body types, we selected 249 sampling sites from four river types (springs: 52 sites, upper courses: 46 sites, middle courses: 76 sites and lower courses: 75 sites) to characterize the longitudinal gradient ( Table 1). The details of the sampling protocols for the EKOO database were described in Verdonschot and Nijboer [47].
In total, 690 species were identified from 249 sampling sites. Community diversity indices (species richness, Shannon index, Simpson index and evenness) were calculated to compare the structural characteristics of the four different river course types.
Feeding strategies of macroinvertebrates reflect the adaptations of species to environmental conditions [48]. Accordingly, a functional feeding classification for macroinvertebrates would improve the knowledge of trophic dynamics in streams [49]. Changes in the distribution patterns of the functional feeding groups (FFGs: scrapers, shredders, gatherers, filterers, predators and piercers), as well as the trophic guilds (carnivores, detritivores, detriti-herbivores, herbivores and omnivores), would reflect the environmental gradient of a lotic ecosystem. Although trophic guilds differ from FFGs, FFGs have been used widely as an alternative to trophic guilds in macroinvertebrate studies [50]. Thus, FFGs have been used to characterize the spatial changes in FFG composition in response to environmental changes in lotic ecosystems, as reviewed in Mihuc [51].
Nineteen environmental variables were selected to reflect the different habitat conditions in the four river course types (Table 1). There were distinct longitudinal gradients in geological and hydrological variables; the slope tended to decrease, and the river width and depth tended to increase along rivers from the springs toward the lower courses. The substrates were classified into two categories (micro substrates: the sum of the relative proportions of peat, clay, silt and sand; and macro substrates: the sum of the relative proportions of gravel and stone). The relative proportion of macro substrate was the highest in the springs, and micro substrate proportions were highest in the middle courses. The sampling sites in springs showed the lowest relative proportion of all macrophyte growth forms, a category that comprises emergent, floating and submerged macrophytes. The relative proportion of floating macrophytes was the highest in the lower courses. The concentrations of ammonium (NH 4 + ), orthophosphate (o-P), total phosphate (t-P), chloride (Cl − ) and calcium (Ca 2+ ) were relatively low in the springs and lower courses, and the concentration of nitrate (NO 3 − ) was highest in the springs. All other water quality variables, except for dissolved oxygen (DO) and NO 3 − , were highest in the middle courses, indicating a high nutrient load or eutrophic conditions [52].

Exergy
Ecological exergy is defined as the amount of work a system can perform if it is brought to thermodynamic equilibrium with its environment or a reference state [13,21]. The exergy of an ecosystem cannot be measured exactly because of the high complexity of almost any ecosystem. However, exergy may be calculated for each system component [53] and thus provide an exergy index as a model of an ecosystem [16]. According to Jørgensen et al. [24], exergy (Ex) can be calculated as follows in Equation (1): where C i is the concentration of component i (i.e., the biomass of species i) in the ecosystem and β i is the weight coefficient that expresses the "quantity of information" embedded in the biomass [54]. Detritus is used as a reference level (β i = 1), and exergy is calculated using this equation in detritus energy equivalents. We multiplied the exergy value by 18.7 to express exergy in units of kJ [55]. In this study, the weighting coefficients were based on trophic groups (carnivores, 47; detritivores, 30; detriti-herbivores, 32.5; herbivores, 35; and omnivores, 41) [56].
Specific exergy (SpEx) is defined as the exergy divided by the concentration of biological content and is calculated using Equation (2): (2)

Self-Organizing Map
To characterize macroinvertebrate communities based on their abundance in the four river-course types, we applied Self-Organizing Maps (SOMs) [57,58]. SOM is an unsupervised learning algorithm based on artificial neural networks that approximates the probability density function of the input data [58]. A SOM consists of input and output layers connected with computational weights (connection intensities). The array of input neurons (computational units) operates as a flow-through layer for the input vectors, whereas the output layer consists of a two-dimensional network of neurons arranged in a hexagonal lattice.
In the learning process of a SOM, the input data (i.e., the abundances of 243 species) were initially subjected to the network. The number of output neurons was set to 77 (= 11 × 7) in a 2D hexagonal lattice derived from previous experience and a preliminary study. The weights of the network were trained for a given dataset. Each node of the output layer computes the summed distance between the weight and input vectors. The output nodes are considered to be virtual units that represent typical patterns of the input dataset assigned to their units after the learning process [59]. Among all of the virtual units, the best matching unit, which has the minimum distance between weight and input vectors, becomes the winner. For the best matching unit and its neighborhood units, new weight vectors are calculated by the SOM learning rule. This procedure trains the network to classify the input vectors according to the weight vectors that are closest to the inputs. For the training SOM, we used the functions provided in the SOM toolbox [60] in Matlab for Windows, version 6.1 [61].
We used two different methods to cluster the trained SOM units into several groups. First, the unified distance matrix algorithm (U-matrix) [62] was applied. The U-matrix calculates distances between neighboring map units, and these distances can be visualized to represent clusters using a grey scale display on the map. Second, a hierarchical cluster analysis with Ward's linkage method based on the Euclidean distance measure [63] was applied to the weights of the SOM output units [59,64]. After the clusters in the SOM analysis were defined, Multi-Response Permutation Procedures (MRPP) were conducted to test whether there was a significant difference among the clusters. The analyses were performed using PC-ORD for Windows version 5.0 [65].
Prior to the analysis of the SOM, we selected 243 species having an occurrence frequency of > 5% in the study sites. Significant correlation was identified between the selected species and all other species through a Spearman correlation analysis (species richness: R 2 = 0.988, abundance: R 2 = 0.992, exergy: R 2 = 0.985 and specific exergy: R 2 = 0.926). The abundance data were natural log-transformed, with a value of one added to the abundance data to avoid taking logarithms of zero.
Indicator species analysis (IndVal) [66] was used to determine the significant indicator species in the clusters based on SOM analysis. IndVal determines indicator species by combining information on the relative abundance of species with their degree of occurrence in a particular group. Indicator values range from 0 to a maximum value of 100. The maximum value indicates perfect indication for a particular group. Species having an IndVal greater than 25 were considered indicator species , individuals of a certain species were present in more than 50% of the samples in one cluster, and the relative abundance of the indicator species in the clusters was more than 50%) [67,68]. The significance of the indicator values for each species was tested with Monte Carlo permutation (9,999 random permutations) using PC-ORD for Windows version 5.0 [65].

Random Forest
The Random Forest method [69] was used to evaluate the relative explanations of environmental variables for exergy and specific exergy. The Random Forest is a non-parametric method for the prediction and assessment of the relationship between a group of potential predictor variables and a response variable [69]. To measure feature importance, we used the mean decrease Gini (MDG), which quantifies the importance of a variable by summing all decreases in Gini impurity due to a given variable when the values of the variable are randomly permuted. This approach is based on a comparison with the original observations. The values of the MDG for each environmental variable were rescaled for exergy and specific exergy within a range of 0 to 100 to represent the relative importance of each environmental variable. Random Forest analysis was repeated 10 times, and the extracted MDG values were averaged for the smallest generalized error. This analysis was performed using the 'randomForest' package in R (http://cran.r-project.org).

Statistical Analysis
An analysis of variance (ANOVA) [70] test was used to evaluate the differences in exergy, specific exergy, community indices and environmental variables among the four river course types and among different clusters defined by the SOM results. A Tukey's multiple comparison test was applied to detect significant differences between factors in the ANOVA test (p < 0.05). A Spearman rank correlation was applied to determine the relationships between community diversity indices, exergy and specific exergy. All statistical analyses were performed with STATISTICA for Windows version 7 [71].

Differences in Biological Indices
Exergy and the community indices tended to increase from springs to the lower courses (Table 1). However, specific exergy was lowest in the middle courses and highest in the lower courses. The relative abundance (%) of shredders decreased from the springs (16.7%) to the lower courses (5.5%), whereas collector-filterers increased from the springs (8.3%) to the lower courses (17.9%) (Figure 1). The proportion of collector-gatherers was relatively high in the springs (47.6%) and the middle courses (45.2%), whereas scrapers (springs: 4.2% and middle courses: 8.8%) and piercers (springs: 2.0% and middle courses: 7.1%) were lower in these two river types.

Patterning Benthic Macroinvertebrate Communities
The study sites were patterned based on the similarities of benthic macroinvertebrate communities in the SOM (final quantization error: 3.989 and final topographic error: 0.012) (Figure 2). The SOM output units were classified into four clusters based on the U-matrix and the dendrogram obtained with hierarchical cluster analysis ( Figure 2). An MRPP showed significant differences between benthic macroinvertebrate communities among the four clusters (A = 0.12, p < 0.01). The size of the letters in the figure is proportional to the number of sampling sites in the SOM map and ranges from 1 to 17 sampling units (Figure 2). The classification of the study sites reflects the characteristics of different river types (e.g., topographical variation, spatial variation and nutrient gradients). The study sites in the springs were classified primarily in cluster 1, and all study sites in the upper courses were classified in cluster 2. Most sites in the middle courses were classified in cluster 3. The sampling sites in the lower courses were primarily classified in cluster 4. Cluster 1 showed low values for community diversity indices, conductivity, NH 4 + , o-P, t-P, Cl − and Ca 2+ (Table 2), whereas clusters 2 and 3 displayed intermediate ranges for community diversity indices and high values for the water quality variables. Cluster 4 displayed the highest community diversity indices and low values for the water quality variables, except for conductivity and Ca 2+ . The relative proportion of collector-filterers increased from cluster 1 to cluster 4, whereas the proportion of shredders decreased ( Table 2). The proportion of piercers was lowest in cluster 1, and collectorgatherers were higher in clusters 1 and 3. The differences in exergy values among the four clusters were similar to the differences in the community diversity indices ( Table 2). Exergy values tended to increase along with the stream size. In contrast, the changes in the specific exergy among the four clusters were different from the changes in the community diversity indices and exergy. Specific exergy was low in cluster 1, representing the study sites in springs, and cluster 3, representing the study sites with high nutrient loads. The highest specific exergy was observed in cluster 4, representing study sites in the lower courses with relatively low nutrient levels. Based on IndVal of the SOM clusters, 113 species were selected as indicator species (p < 0.05) (Appendix 1). The abundance and frequency of rheophilic species inhabiting oligosaprobic rivers were relatively high in cluster 1 and included Nemoura cinerea, Nemurella pictetii, Plectrocnemia conspersa, Hexatominae, Brillia modesta, Micropsectra sp. and Elodes minuta. The decomposition of fallen leaves in forest areas resulted in a relatively high level of detritus and a high abundance of Gammarus pulex (Table 2) in the sites of cluster 1. In cluster 2, the species usually inhabiting alpha to beta mesosaprobic rivers were dominant and included Hydrobius fuscipes, Hydroporus palustris, Anisus leucostoma and Psectrotanypus varius. Piercers were also common in cluster 2, including Haliplus lineatocollis, Hydroporus palustris and Agabus sp. The abundance and the occurrence frequency of detritivores were high in cluster 3, representing the middle courses. Cluster 3 was characterized with Limnodrilus hoffmeisteri, Tubificidae and Procladius sp., and cluster 4, representing the lower courses, was characterized with the species (primarily predators and piercers) inhabiting emergent macrophytes and typically swimming at the water surface. These species were represented by the Coleoptera, including Haliplus fluviatilis, Hygrotos versicolor, Noterus clavicornis, Laccophilus hyalinus and Haliplus sp. and the Hemiptera, including Sigara falleni, S. strigata, S. distincta/falleni/longipalis and Ilyocoris cimicoides.

Prediction of Exergy
Exergy and specific exergy were predicted based on differences in environmental variables through the learning process in the Random Forest model, and the relative importance of environmental variables was evaluated (Figure 3). The values of predictive power for exergy and specific exergy were specified by coefficients of determination (R 2 ) of 0.51 and 0.23, respectively. Exergy and specific exergy responded differently to environmental variables. For example, the longitudinal gradient of stream size (i.e., river width) was the most influential environmental factor for exergy prediction, followed by DO, ammonium and the river width. Orthophosphate was the most important predictor of specific exergy, followed by nitrate and total phosphate.
The Spearman rank correlation coefficients indicated that exergy, specific exergy and the community indices (species richness, Shannon diversity index and Simpson diversity index) were positively correlated, with correlations ranging from 0.267 to 0.998 (p < 0.01). The highest correlation was between the Simpson diversity and the Shannon diversity indices (0.998), while the lowest was between exergy and specific exergy (0.267). In addition, the exergy and community indices showed high correlation coefficients (0.589 to 0.919) when compared with the correlation coefficients between specific exergy and the community indices (0.267 to 0.320) ( Table 3).  Table 1.

Discussion
Lotic ecosystems are complicated and dynamic, and their physical, chemical and biological components interact and influence each other through non-equilibrium, dissipative processes [72]. Changes in environmental factors often cause alterations in ecosystem structure. These alterations include changes in species composition or biodiversity as well as ecosystem function. A high load of organic substances and nutrients (e.g., nitrogen or phosphorus) usually causes changes in species composition and trophic structure. Accordingly, from the water management point of view, there is a need to find indices that express the alteration of the physicochemical environment in the lotic ecosystem as well as the changes in energy flow between the various trophic levels in the food web [73].
In this study, four naturally separated groups are clearly classified by the SOM clusters. A comparison between natural separated groups (Table 1) and SOM clustered groups (Table 2) indicates the effectiveness of the SOM model. Such similar results suggest that SOM is a useful tool for identifying differences between naturally separated groups [59,74].
The exergy and biodiversity indices tended to increase from cluster 1 to cluster 4 in the SOM analysis. However, specific exergy was lowest in cluster 3 with a high proportion of detritivores (i.e., a low β value), whereas specific exergy was highest in cluster 4 with a high proportion of carnivores (i.e., a high β value). Similar results were observed in the lotic and lentic ecosystems. According to Marques et al. [28] and de Wit et al. [75], eutrophication (i.e., a high nutrient load) results in higher exergy and lower specific exergy because the biomass bulk has a low weighting factor. Unlike exergy, specific exergy is not dependent on biomass alone and expresses the ability of the ecosystem to accept and utilize external fluxes of energy while simultaneously serving as an indicator of ecosystem development, reflecting the complexity and level of evolutionary development of species in an ecosystem [10]. Accordingly, specific exergy reflects diversity (i.e., a more complex ecosystem), and a higher specific exergy represents more highly developed organisms (i.e., higher β values that represent more information) [34].
Middle river courses (cluster 3) were directly influenced by nutrients transported in runoff from the adjacent pastures and agricultural areas along the streams together with nutrient inputs from upstream. Furthermore, the relative ratio of micro substrate was higher than in the lower courses. The middle courses receive large amounts of sediment as well as organic matter, resulting in the accumulation of detritus on the stream bottom and thus producing a high abundance of detritivores. The lower courses (cluster 4) showed relatively stable habitat conditions, and such conditions typically support high species complexity and diversity. The nutrient sources for the lower courses originated primarily from upstream nutrient input and represented a dilution of the nutrients coming from the middle courses. Discharge and current velocity were low but stable over all seasons, as is also the case for lentic habitats [52,76].
Even though any environmental factor could potentially influence the exergy and specific exergy in our study, the Random Forest models demonstrated that exergy was more heavily influenced by longterm conditions such as the longitudinal gradient (e.g., river width); whereas short-term disturbances, such as those associated with the nutrient concentrations (e.g., o-P, t-P, NO 3 − and NH 4 + ), were more important in determining the specific exergy. This was congruent with a study by Marques et al. [27], which suggested that exergy reflected slower changes, whereas specific exergy may shift drastically if the species composition of a particular body of water changed abruptly. Those findings indicated that the complementary use of exergy and specific exergy provides information on long-term gradients such as the longitudinal gradient (i.e., river continuum concept), as well as the ecosystem response, including the suppression of the food web and diversity due to environmental stress [77].
Most biodiversity indicators interpret the community composition from structural characteristics, while the exergy and specific exergy focus on the functional characteristics of an ecosystem. In recent ecological studies, exergy values have been used as ecosystem health indicators, reflecting such processes as human activities and eutrophication [10]. However, exergy and specific exergy are promising indicators for shifts in species composition and trophic structures along longitudinal and nutrient gradients. In this study, exergy and the community indices showed the same trends in the different river courses and increased from the springs to the lower courses. In contrast, specific exergy was the highest in the lower courses, with relatively low values in other river types. Such a relationship indicates that the complicated results from different bioassessment methods and that the combined use of structural and functional indicators could therefore provide comprehensive information on ecosystem development and complexity along both longitudinal and nutrient gradients in lotic ecosystems.

Conclusions
It is clearly crucial to obtain an integrative biological indicator to interpret the structure and function of lotic ecosystems. Exergy and specific exergy are potentially promising target functions for this purpose, indicating shifts in species composition, as well as trophic structure along longitudinal and nutrient gradients. Exergy increased from the springs to the lower watercourses, but specific exergy was highest in the lower watercourses. Specific exergy represented ecological complexity based on nutrient concentrations. In this sense, it is advisable to use both indicators to evaluate ecosystem health assessment because exergy and specific exergy relate different and complementary information regarding the functioning of a lotic ecosystem.