Methodology for classifying the ecosystem integrity of forests in Germany using quantified indicators

Atmospheric deposition of nitrogen and climate change can have impacts on ecological structures and functions, and thus on the integrity of ecosystems and their services. Operationalization of ecosystem integrity is still an important desideratum. A methodology for classifying the ecosystem integrity of forests in Germany under the influence of climate change and atmospheric nitrogen deposition is presented. The methodology was based on 14 indicators for six ecosystem functions: habitat function, net primary function, carbon sequestration, nutrient and water flux, resilience. It allows assessments of ecosystem integrity changes by comparing current or prospective ecosystem states with ecosystem-type-specific reference states as described by quantitative indicators for 61 forest ecosystem types based on data before 1990. The method developed enables site-specific classifications of ecosystem integrity as well as classifications with complete coverage and determinations of temporal trends as shown using examples from the Thuringian Forest and the “Kellerwald-Edersee” National Park (Germany).


Background
Climate change and atmospheric nitrogen (N) inputs can alter the integrity of ecosystems, i.e., their dominant structures and functions, and thus limit their bene fits for humans, i.e., the ecosystem services. Therefore, action 5 of the EU Biodiversity Strategy to 2020 foresees that Member States will map and assess the state of ecosystems and their services in their national territory. To this end, an operational guidance to the EU and the Member States on how to assess the condition (or the state) of Europe's ecosystems was developed [13]. Accordingly, ecosystem condition should be measured using indicators and specified for the national level of the EU member states [13]. For Germany, Jenssen et al. [7] and Schröder et al. [23] laid the foundations for a spatially explicit and nationally applicable concept for the classification of changes in ecosystem integrity. This methodology was further deepened and developed by Schröder et al. [22]. It enables an integrative assessment of changes in ecosystem integrity, taking into account the effects of climate change in combination with atmospheric nitrogen (N) deposition. Characteristics of ecosystem integrity concerned are self-organizational capacity, functionality and compliance of abiotic and biotic properties with the natural site potential (identity). The methodology was based on an extensive vegetation database, nationwide available data from digital maps and long-term monitoring programs. It was complemented by dynamic modeling of Open Access *Correspondence: winfried.schroeder@uni-vechta.de 2 Chair of Landscape Ecology, University of Vechta, P.O.B. 1553, 49364 Vechta, Germany Full list of author information is available at the end of the article future climate and soil conditions. The ecosystem condition was assessed on the basis of the criteria of functionality, chemical and biological characteristics, and stress tolerance to anthropogenic nitrogen inputs and climate change. The methodology allows the identification and mapping of potential natural ecosystem types and current near-natural ecosystem types. For certain climate scenarios and atmospheric nitrogen inputs , possible ecosystem developments can be projected and evaluated in the future. The concept complements existing assessment methods for ecosystem conditions by taking abiotic environmental factors and their changes into account as drivers of biological changes and ecosystem functions. At the same time, it should serve to identify the causes of disturbances as early as possible and to derive suitable measures for the preservation and development of certain ecosystem conditions.
For the development of the methodology presented in this paper, the Federal Environment Agency has attached importance to use data from monitoring programs and to cover three spatial levels: the forest stand level as well as the regional and national levels. Thereby, the Germanwide map of hemeroby [25] could not be used since it does not address ecological functions and "is inappropriate for a more accurate calculation of spatial extent and thus the monitoring of local and regional developments" [26]:2.
A fundamental component of the methodology is a classification of Germany's semi-natural ecosystems. Their concordance with other ecosystem classifications for which no spatial concretisation has been carried out nationwide (European Nature Information System EUNIS, [20], habitat types according to Annex I of the Habitats Directive) has been achieved. Thus, the developed ecosystem classification is connectable with other approaches and enables ecologically founded interpretation and spatial differentiation. For 61 selected ecosystem types, a historical reference condition was quantified based on data from the period 1961-1990 [7,23]. The reference condition was defined as a type-specific condition of ecosystems, the characteristics of which are characterized by intervals of historical ecosystem condition variables . These conditions are relatively least affected by substance inputs and climate change, which can be adequately quantified with measurement data.
For selected ecosystem functions (habitat function, net primary function, carbon storage, nutrient flow, water flow and adaptability), indicators were selected with which current and modeled future ecosystem conditions can be compared with the respective ecosystem type-specific reference conditions. The indicators were quantified with data from monitoring programs and from the Waldkunde Institut Eberswalde (W.I.E.) database, whereby the focus was on the effects of changes in the abiotic systemic bases of development.
The reference states were quantified for 40 near-natural forests and 21 cultivated forests. Cultivated forests are distinguished from near-natural forests by a tree species composition that has been changed significantly compared to the potential natural state. The reference states refer to the period up to 1990, mainly from 1960 onwards, but in individual cases to data dating back to the 1920s and 1930s. For each ecosystem type, its reference status was indicated by a data sheet with the following information: 1. Ecosystem code: 1st digit = climate ecological coordinate, 2nd digit = water balance type, 3rd digit = substance cycle type (for description see [22], vol. 3), 2. Name of ecosystem type, 3. EUNIS class, 4. Biotope type BfN [20], 5. Vegetation type according to common plant sociological classifications, 6. Photo, 7. Habitat type according to the Fauna-Flora-Habitat Directive [24], 8. Position in the two-dimensional ecogram with the coordinates soil moisture and base saturation, 9. Location factors: soil shape, soil type, terrain, macroclimate, 10. Habitat function: characteristic species association with continuity and mean quantity development of the soil cover, maximum Kullback distance of the individual records to the mean species quantity distribution, minimum similarity of the individual records with the mean species quantity distribution, 11. Net primary production (NPP): above-ground average annual NPP at the time of culmination in tree wood, leaf/needle mass, ground vegetation and total mass, upper stand height at age 100 as comparative parameter, 12. Carbon storage: carbon stock in humus (C org in humus layer and in soil up to 80 cm depth), 13. Nutrient flow: pH in 1/10 KCl, base saturation V in % and C/N ratio in the uppermost 5 cm from H to Ah horizon (interval of mean value and standard deviation), humus form, nutritional characteristics N%, P%, K%, Ca%, Mg% in the assimilation apparatus of trees in g/100 g of leaf/needle dry matter (August, interval of mean value and standard deviation), 14. Water flow: soil moisture index (interval from mean value and standard deviation) as well as 15. Adaptation to changing environmental conditions: maximum proportions of natural site tree species in self-organized development stages.
The objective of this contribution is to present the methodology of quantifying the ecosystem functions referred to in bullet points 10-15 and to show how ecosystem integrity was classified on this basis.

Methods and results
This paper presents for the first time the following methodical issues for classifying the ecosystem integrity of forests by example of Germany: the operationalisation of ecosystem functions by quantitative indicators (Sect. "Determination of indicator characteristics") and, based on this, the quantification and classification of ecosystem integrity (Sect. "Classification of ecosystem integrity") including its site-specific analysis (Sect. "Site-specific classification"), regionalisation (Sect. "Area-related classification") and temporal trends (Sect. "Determination of temporal trends"). The data used and the results produced as well as software tools developed were published as data and software papers [8,9,16,17].

Indicators of habitat function
The habitat function is simply indicated by the composition of the vegetation according to quality (higher plant species as well as species of soil-dwelling mosses and lichens) and quantity (coverage percentage). For this purpose, the Kullback distance [7,11,23] of the vegetation composition of the study area was calculated from the distribution of the mean species quantities of the reference state (Fig. 1).
For each of the individual vegetation relevés representing the reference state, the Kullback distance to the mean species quantity distribution was calculated and from the sum of the mean value and standard deviation of the totality of these distances, a value was calculated characterizing the reference state, referred to as the maximum Kullback distance of the individual relevés to the mean species quantity distribution, and documented in the data sheets [8] for each ecosystem type. A comparison of the Kullback distance of the vegetation composition of the investigated area with this "limit value" allows a statement on the extent to which the vegetation composition corresponds to the reference condition or not.
In addition to the Kullback distance, an index is calculated which shows the quantitative correspondence of the current species composition of the vegetation with the mean species quantity distribution of the type [4]: This similarity index S is calculated analogously to the Kullback distance [22], Vol. 2: Sect. 2.3). It allows a comparison to be made with the "limit value", also identified in the data sheets [8] as the minimum similarity of the individual relevés representing the reference state with the mean species quantity distribution, which was calculated as the difference between the arithmetic mean of all similarity indices of the individual relevés representing the reference state and their standard deviation. Due to its formal structure as entropy, the Kullback distance emphasizes differences in characteristic combinations of several species, each with medium quantity development, while the similarity index is influenced mainly by agreement of the highly continuous dominating species. This difference may be relevant to the interpretation of habitat function for different groups of plant and animal species in different ecological domains, and therefore both indicators are considered and illustrated in the following by the ICP 1 Forests Level II Location 1605 (Großer Eisenberg, Germany): [19: Table 3 [9] "C4-6d-B1_Vegetationsgesamttabelle.xls") are calculated in an analogous manner, the sum of the mean value and standard deviation of these KD values is obtained as the value for the maximum Kullback distance of the individual relevés to the mean species quantity distribution KD max = 0.53, which is also shown in Jenssen et al. [8]. An analogous calculation of the percentage similarity index S yields values of S(1960) = 69.4%, S(2001) = 50.1% and S(2006) = 60.6% compared to a minimal similarity of the individual relevés representing the reference state with the mean species quantity distribution of S min = 65%. Thus, both calculated indicators for the habitat function show that the typical species composition of the investigated area is clearly disturbed after 2000, with a reversible development being observed between 2001 and 2006.

Indicator net primary production
The indicator net primary production refers to the net NPP of above-ground wood biomass relevant from a forestry point of view, which is recorded in the form of growth. To make it possible to compare this size, which fluctuates greatly with the age and treatment of the stand, with the reference state valid for the respective ecosystem type, the average net primary production of wood biomass at the time of its culmination is calculated from the wood measurement monitoring data using the method outlined in Fig. 2 and described below.
1. Calculation of the mean stand height HG as the height of the circular mean trunk from the individual tree data recorded on the study area for breast height diameter d i and height hi according to the formula 2. Calculation of the relative height class from the stand age and HG according to the routine BON_REL [9]. 3. Calculation of the growth trend depending on the age of the portfolio according to the routine "Growth" [9] and determination of the average total growth DGZ depending on the age of the portfolio. The "Growth" routine in turn accesses the "Stock" routine and, if necessary, the "Diameter" routine [9]. 4. Determination of the culmination point of the average timber growth (DGZ). Multiplication of the DGZ max with the density of the wood species leads to the indicator maximum of the average net primary production DNP max .
The calculation of the indicator for the year 1995 is exemplified here by example of ICP Forests Level II Location 1605 (Großer Eisenberg, Germany): the Level II data set for the area shows individual timber data for 141 trees of the spruce species for 1995. The data sheet on the reference condition of the Rohhumus-Fichten-Hochbergwald (C4-6d-B1; Rawhumus spruce forest on the altimontane level) [9] shows an average stand height of 22 m at the age of 100 years. This results in a relative height class rating of BON_REL (4,100; 22) = 3.97 for the mean reference condition. Using the GROWTH (4; ALT; 3,97) routine, the average annual net primary production of tree wood DNP max = 5.58 t/ha is now calculated using the method described in point 3, assuming that the reference condition would have the stocking density assumed in the yield table model. In fact, however, it is a natural spruce forest in the ecological battle zone between a closed high forest and an open grove vegetation with a stocking density significantly reduced compared to the forest yield table model, which in turn leads to a proportional reduction of stock and hectare-related growth. The data sheet on the reference condition of the Rohhumus-Fichten-Hochbergwald (C4-6d-B1; Raw-humus spruce forest on the altimontane level) [9] shows an average annual net primary production of DNP max = 2.2 t/ha. This results in a reduction factor of 2.2/5.58 = 0.394, by which the average annual net primary production of tree wood DNP max = 4.93 t/ha determined under point 4 must be multiplied, so that an average annual net primary production of tree wood DNP max = 1.94 t/ha is obtained as an indicator of the net primary production of the tree stock of the monitoring area "Großer Eisenberg".

Indicator carbon storage
As an indicator for carbon storage, the carbon stored in the humus of the organic layer and in the mineral soil between 0 and 80 cm deep is calculated (Fig. 3).
The quantities of C org in g/kg given for the individual soil horizons are multiplied by the respective bulk density (kg/m 3 ) and converted into stock values per hectare using the respective horizon thickness data.
If the data for individual horizons do not contain information on bulk density, the volume-related C reserves can alternatively be calculated on the basis of an empirical relationship between C content and litre weight of the fine soil according to [2] The sample calculation was performed using the Level II data for 2009 ( Table 1). The mean value for C org [g/ kg] was calculated for each of several measurements given per layer. For the layers of the mineral soil M01 (0-10 cm), M12 (10-20 cm), M24 (20-40 cm) and M48 (40-80 cm), the bulk density was calculated according to the empirical formula given above. The hectare stocks of organic carbon obtained by multiplying C org [g/kg] by bulk density [kg/m 3 ] were summed across all layers.

Nutrient and water flow indicators
Nutrient and water flow indicators are calculated using indicator value models (Fig. 4).
Indicator value model for calculating the C/N and pH of the topsoil The indicator value model calculates for a given vegetation survey a probability distribution over the C/N ratio and the pH (KCl) of the topsoil (top 5 cm of the  humus layer or mineral soil). The C/N ratio serves as an indicator of nutrient availability similar to the N number according to Ellenberg et al. [1]. From this distribution, an expected value for the C/N ratio and the pH of the topsoil is calculated. A complete documentation of the model is contained in Jenssen et al. [9]. The basis for the modeling of the C/N ratio and the pH are the probability distributions of the most frequent plant species of the Central European forest vegetation, taking into account their stratum affiliation and quantity development ( [9]: Tables 3 and 6, respectively). These distributions are multiplicatively linked to a probability distribution for the ecotope characterized by the vegetation uptake. From the resulting distribution, the characteristic values can be assigned to the ecotope. In the applications performed, the expected value assigned to the ecotope was the arithmetic mean of the class values of the C/N ratio or the pH weighted with the class probabilities of the resulting probability distribution ( [9]: Tables 1 and 4).
The following model algorithm has been implemented: 1. Reading a table (tblVEG) with the vegetation relevé including all occurring species separated by tree layer, lower and upper shrub layer, field layer and the corresponding percentage cover values. 2. Reading the class mean values C/N and pH [9]: Tables 1 and 4. for 20 classes, respectively. 3. Calculation of the probability densities of the occurring species taking into account stratification and cover value class using the function Example: ICP Forests Level II site 1605 (Großer Eisenberg, Germany).
The calculation should be performed using an executable program that implements the algorithm described     above based on the data documented by Jenssen et al. [9]. For a better comprehensibility of the model algorithm, a spreadsheet calculation was carried out in Table 2 for the vegetation survey of the ICP Forests site LII-1605 from 1960, which is documented under the area designation STO 180 in Jenssen et al. [8] ("C4-6d-B1_Vegetationsgesamttabelle.xls"). The first two columns contain the plant species taken from the vegetation survey separately by stratum and the corresponding percentage cover values (step 1 of the model algorithm). The first row contains the class mean values C/N [9]: Table 1 for 20 classes each (step 2). The inner fields of the table contain the probability densities calculated according to steps 2 and 3 for the occurring species over the respective classes, provided that a density function (parameters in [9]) exists for the species. The probability densities were normalized so that their sum over all classes of the C/N ratio (row sum, stored in the last column) results in one. In the row "Column product", the probability densities of the different plant species were multiplied by the respective C/N classes (step 4). In the next row "Column product, weighed with class frequencies" the probability densities are weighed with the case numbers (absolute frequencies) of the individual classes from Jenssen et al. [9]: Table 1 (step 5). This is to ensure that the two extreme classes with a significantly lower number of underlying measured values and correspondingly lower statistical representation are given a lower weighting in the calculation of the expected value (explanations in [9]). In our example, however, these classes have a zero probability, so step 5 has no effect on the result ( Table 2). The modal value of the distribution remains unchanged above the class with the class value C/N = 28.2. The corresponding class probability is 43% after normalization across all classes. In the last line, the expected value for the topsoil C/N ratio of the investigated area is calculated by multiplying the class values with the respective class probabilities and then summing all classes (step 6). The result for the expected value is C/N = 27.4. A completely analogous calculation is carried out to determine the pH (KCl) in the topsoil. pH = 2.8 is obtained.
Indicator value model for calculating the base saturation of the topsoil So far, no probability density functions have been created for base saturation (V). Therefore, with the help of 787 measurements of the V value available in the W.I.E. database and on the basis of a close correlation to the pH values, a total of 838 forest plant species were assigned mean values of the base saturation of the topsoil [9]: Table 9. For the calculation of an area-related base saturation, an average value weighted with the cover values was calculated. This value is subject to corresponding uncertainties compared to the C/N ratios and pHs.
Example: ICP Forests Level II site 1605 (Großer Eisenberg, Thuringian Forest, Germany). Table 3 shows the calculation of the base saturation for the ICP Forests site LII-1605 (Großer Eisenberg, Thuringian Forest, Germany) in 1960 on the basis of the vegetation survey of the site LII-1605 from 1960, which is documented under the area designation STO 180 in Jenssen et al. [9] ("C4-6d-B1_Vegetationsgesamttabelle.xls"). The first two columns contain the plant species taken from the vegetation survey and the corresponding percentage cover values, whereby the values r = 0.01% and + = 0.1% were set. The third column contains the mean V values of the individual plant species taken from Jenssen et al. [9]: Table 9. The last column contains the products from these V values and the coverage percentage, i.e., the products from the two previous columns. In the last row, the column sum of these for each of the 20 classes. You get a vector pd that contains the resulting probability density for each of the 20 class values. Normalization of the resulting probability density across all classes to one. 7. Calculation of a surface-related topsoil moisture as mean value over the class values of the topsoil moisture weighed with the probability densities pd.
Example: ICP Forests Level II Location 1605 (Großer Eisenberg, Thuringian Forest, Germany). Table 4 Table 8. From this, the parameters of a normal distribution according to step 3 were derived. The first row of Table 5 contains the class mean values of 20 equidistantly distributed classes of topsoil moisture between the extremes 0 and 10 (step 2). The inner cells contain the class probabilities (probability densities) calculated for each type of shrub layer and ground vegetation occurring and for each class mean according to the formula given in step 3 (probability densities) normalized to one across all classes (last column, step 4). Table 6 shows the class probabilities multiplied by the respective percentage cover values (second column) (probability densities, step 5). The row "Column product" contains the product of the class probability pages calculated for each moisture class and weighted with the cover values of the species. In the line below, the class probabilities (probability densities) were normalized to one across all classes. The greatest probability of 45% is calculated for the class between the moisture index 5.5 (medium to permanently fresh) and 6 (permanently fresh). These values were multiplied in the lowest line by the respective class values of the moisture index. The sum of the lines is the expected value of the moisture content of the topsoil on the test area and gives the ratio 5.6 (step 7).

Indicator of adaptability to changing environmental conditions
As an indicator of adaptability to changing, unpredictable environmental conditions, the percentage similarity of the current proportions of tree species with the spectrum of natural site tree species [5] is used: The p i denote the percentage amount shares of the tree species indexed with i for the current stocking on the area to be valuated, whereby the quantity shares are summed over several possibly existing tree layers and set to 100%: The P i max describe the maximum percentage areas of natural site-specific tree species that are not exceeded in the course of self-organized development stages under the respective site conditions (explanations in [7], Chapter 4).
The P i max were determined on the basis of the knowledge of the natural distribution and the site requirements as well as the growth and competition behavior of the native tree species for the different ecosystem types and documented in the data sheets of the reference conditions [9]. The p i are determined from the actual cover values of the natural site-specific tree species on the receiving area (Fig. 5).
For a raw humus pine beech forest (Eb-4n-B2), for example, the data sheet for reference conditions [9] gives 80% for the red beech, 40% for the grape oak, 30% each for sand birch and pine and 5% for the mountain ash as the maximum proportions of natural site-specific tree species in self-organized development stages. At present, 70% pine, 20% red beech, 5% grape oak and 5% sand birch are found on a vegetation survey area assigned to this reference condition. The indicator is, thus, calculated as follows.
The current proportions of tree species are, therefore, 60% similar to the spectrum of natural site tree species, i.e., those tree species which are not exceeded in the course of self-organized developmental stages (in contrast to the classical definition of potential natural vegetation (PNV), the spectrum covers not only the main stages but also temporary intermediate forest stages) at this site. Despite excessive amounts of pine, the stock still has sufficient potential to adapt to possibly P = min(70%, 30%) + min(20%, 80%) + min(5%, 40%) + min(5%, 30%) = 60%.
changing environmental conditions through self-organized development.
For ecosystem type C4-6d-B1, the data sheet of reference states [9] gives the maximum proportions of natural site-specific tree species in self-organized development stages as 100% for spruce, 5% for mountain ash, 2% for bog birch and 5% for silver fir as maximum proportions of natural site tree species. In 1960 and in all other monitoring years, 100% spruce was found on the vegetation survey area. Thus, the indicator is calculated to P = min(100%, 100%) = 100%. The current proportions of tree species are, therefore, 100% similar to the spectrum of natural site tree species, i.e., those tree species which are not exceeded in the course of self-organized development stages at this site.

Classification of ecosystem integrity
The analysis and estimation of ecosystem conditions and their development over time is basically carried out by comparison with a functionally and structurally determined historical reference condition [7,23]. The change in ecosystem integrity is to be classified as higher the more the status parameters deviate from those of the respective ecosystem type-specific reference status. To classify ecosystem integrity, a distinction is made between 5 levels in comparison with the reference state, namely on three levels (Fig. 6): 1. Deviations from the reference situation for individual indicators, 2. Deviations from the reference state for individual functions (based on indicators), 3. Changes in ecosystem integrity across functions for the ecosystem type under consideration.
The deviations between the interval limits of the ecosystem-specific reference condition and the lower or upper limit of the ecosystem-specific total span are subdivided into the levels very low (= within the reference span) and low, medium, high and very high. The classifications for a total of 13 indicators assigned to 6 ecosystem functions are based on the principles presented in Table 7.
The variance levels for the indicators of habitat function, net primary production and carbon sequestration are based on the historically determined reference span and a quartering of the residual interval within the total span. For the nutrient balance indicators whose reference range is defined by the respective ecosystem mean value of ± the simple standard deviation (see above): indicator value models for calculating the C/N and pH the topsoil, the base saturation of the topsoil and the moisture index of the topsoil), are assigned low and medium with values into the ranges smaller than 2 or 3 times the standard deviation. The steps high and very high were defined by halving the remaining intervals within the overall spans specific to the material or material cycle type. For ecosystem types where the triple standard deviation is already outside the total span, the deviation steps high and very high are omitted (e.g., for the pH of the Moder-Tannen-Buchen-Bergwald (D2-6d-C2; Moder fir and beech forests of the montane level), or the base saturation of the Magerrohhumus-Sand-Kiefernwald (Ed-2n-A2; Sandy meager-raw-humus pine forests). According to Jenssen and Hofmann [5], near-natural stands are characterized by values of 60-100%. This span is used as a reference (= very small). The deviation interval below 60% is quartered and distributed equally for all ecosystem types to the levels low (= 59-45%), medium (= 44-30%), high (= 15-29%) and very high (= 0-14%).
In the case of ecological functions, which are described by only one indicator (net primary production, carbon storage, water flow, adaptability), the classification is directly linked to the assessment of the individual indicators. For functions with several assigned indicators (habitat, nutrient flow), the individual estimates are aggregated via the modal value (level 2 in Fig. 6. The modal value also provides good orientation for the overall functional classification of changes in ecosystem integrity, with all criteria being considered equally. For the interpretation, the meanings of the 5 levels of deviation and change are described in Table 8. With regard to the interpretation of detected deviations from the reference condition, various threshold values can be used to derive different needs for action on the basis of Mitchell et al. [14] (Fig. 7). Under the simplified Very low The values of the indicators habitat function, net primary production and carbon storage correspond to those of the reference range. For the respective ecosystem type there are also no or only very minor changes in the values for the physical-chemical conditions (nutrient flow, water flow) compared to the values that characterize the historical reference condition . The similarity of the current tree species composition with the spectrum of natural site-specific tree species is very high

Low
The indicators for habitat function, net primary production and carbon storage show small deviations from the historical reference values. The physical-chemical conditions also deviate only slightly from the values of the reference type. There is a high similarity between the current tree species composition and the spectrum of natural site-specific tree species Medium The values of ecosystem function indicators differ moderately from those normally associated with the historically determined reference status. The values give indications of moderate deviations and show significantly stronger interference than was the case under historical conditions. There is a moderate similarity between the current tree species composition and the spectrum of natural site-specific tree species High Ecosystem conditions where indicators point to major changes and differ significantly from those normally associated with the reference status. There is little similarity between the current tree species composition and the spectrum of natural site-specific tree species Very high Ecosystem conditions in which the indicators point to very strong changes and deviate very strongly from the historical reference condition. There is very little similarity between the current tree species composition and the spectrum of natural site-specific tree species assumption of a linear decrease in ecosystem integrity beyond the characteristic range of the reference state, continued measures of environmental monitoring are recommended from the stage low, at the latest from the stage medium, if the temporal trend indicates an increase in deviations from the reference state. Management measures are recommended from the high level onwards, at which conditions outside the ecosystem-specific, natural variability are given (here: mean value of ± three times standard deviation). Appropriate concretizations can be justified with the changes identified at indicator, function and/or ecosystem level.

Site-specific classification
The assessment forms for 61 forest and forest ecosystem types [9] can be used to classify ecosystem integrity at individual sites.
A. Deviation levels for the indicators At the first level, it is assumed that the ecosystem type was determined by means of a determination key or computer-aided comparison of a vegetation survey with the reference conditions documented in Jenssen et al. [9,10,22], vol. 2, chapter 2) and that the characteristics of the indicators for the 6 ecosystem functions (Sect. ""Determination of indicator characteristics) were determined. On this basis, the deviation from the respective reference status is classified using the assessment sheet specified for the ecosystem type identified (level 1 in Fig. 6).
For the Rohhumus-Fichten-Hochbergwald (C4-6d-B1; Raw-humus spruce forest on the altimontane level) at ICP Forests Level II site 1605 in the Thuringian Forest, the interval of the reference condition is 12.9-19.7%. A base saturation value of 25% was estimated on the basis of vegetation cover (Sect. "Indicator value model for calculating the base saturation of the topsoil"). The deviation from the reference condition is classified as 'medium' on the basis of the rules documented in the evaluation sheet. Continuous observation of the examination area is recommended. Fig. 6), deviations from the reference status are determined by evaluating the characteristics of the assigned indicators. Since net primary production, carbon storage, water flow and adaptation to changing environmental conditions are each described by a singular indicator, the level of deviation for the function is identical to that of the indicator. With regard to habitat function and nutrient flux, deviations from the reference status are determined by aggregating the levels of deviation of the associated indicators. Aggregation is done with the modal value providing good orientation. All steps of the deviation from the Example: ICP Forests Level II site 1605 (Großer Eisenberg, Thuringian Forest, Germany)-"nutrient flow" at the Rohhumus-Fichten-Hochbergwald (C4-6d-B1; Rawhumus spruce forest on the altimontane level; Fig. 8).

B. Deviation levels for the functions In the case of ecosystem functions (level 2 in
For the Rohhumus-Fichten-Hochbergwald (C4-6d-B1; Raw-humus spruce forest on the altimontane level) at the ICP Forests Level II 1605 site, the determined values of the nutrient flow indicators result in different levels of deviations from the reference condition: high for the pH, medium for the base saturation and very high for the C/N ratio, furthermore high (N content needles), very low (P content needles), low (K content needles), high (Ca content needles), medium (Mg content needles). The aggregation of these 8 classifications results in an overall deviation from the reference status high for the category nutrient flow at the level of ecosystem functions, which indicates an increased need for countermeasures (Fig. 7).

C. Stages of change for ecosystem types
In a third step, the change in ecosystem integrity is classified across functions (level 3 in Fig. 6.
For the habitat function, there are low or medium deviations from the reference status, low deviations from the net primary production, very low deviations from the C-storage, high deviations from the nutrient flow, low deviations from the water flow and very low deviations from the reference status for the adaptability to changing environmental conditions. This results in 2 denominations of the level very low, 3 denominations of the level low, one denomination of the level medium and one of the level high. The subsequent assessment leads to the classification slight change in ecosystem integrity, to take into account the modal value of the deviation from the reference condition in this example. With regard to ecosystem integrity, continued observation of the study area is recommended.

Area-related classification
For an area-related classification of ecosystem integrity at the regional level, data from vegetation-reception areas (e.g., in nature reserves) are often available or-in contrast to comparable soil data-can be collected with relatively little effort. Deviations of current ecosystem conditions from the respective reference conditions can be spatially generalized either on the basis of mapping of the ecosystem types in the area under investigation or by division of the area into a regular grid (e.g., 2.5 km × 2.5 km) in each case in connection with a representative selection of vegetation survey areas as a sample.
Method 1: classification of ecosystem integrity for mapped ecosystem types In a first step, the quantification of the topsoil parameters from the vegetation structures according to Sect. "Nutrient and water flow indicators" is carried out. Based on the variances of the reference data to the ecosystem types [22], at least 5 vegetation relevés per ecosystem type are recommended for the regional level. Using the wise values of the plant species, a site-specific classification of the deviations from the reference condition (1: very low, 2: low, 3: medium, 4: high, 5: very high) for the corresponding indicators of ecosystem integrity (Sect. "Site-specific classification") is then made for each vegetation survey area. Next, the ecosystem type-specific median of the deviation levels is determined as the central tendency and the maximum deviation. The medians are then aggregated via the modal value to levels of deviation and change at the levels of ecological functions or ecosystems (Sect. "Site-specific classification"). The spatial transfer of the ecosystem-specific classifications determined in this way can finally take place by allocation to For the area of the Kellerwald-Edersee National Park (Hesse, federal state of Germany), an ecosystem type mapping on a scale of 1:5000 and 70 vegetation surveys for five different ecosystem types are available [18] (Table 9). After determining the topsoil parameters, for example for the Braunmull-Buchen-Bergwald (D1-6d-D1, Brownmull beech forests of the montane level), the value range of the moisture indicators (DKF) originating from eight indicator value calculations is between 4.7 and 6.9. This results in eight site-related classifications of the deviations from the reference condition between very low and high Fig. 9. The median value as the central tendency is medium. Together with the median values for the other indicators (here: Kullback distance, similarity of species abundance distribution, pH, base saturation, C/N ratio), at the level of ecosystem functions the modal value results in the levels low for habitat function, medium for nutrient flow and medium for water flow. The modal value of these three deviation steps results in a mean change of ecosystem integrity for the Braunmull-Buchen-Bergwald (D1-6d-D1,Brownmull beech forests of the montane level). This results in the recommendation to further monitor the development of ecosystem integrity of D1-6d-D1 in the basement forest in the future, especially with regard to base saturation, C/N ratio and soil moisture.

Method 2: classification of ecosystem integrity in an area
grid Alternatively, if no information on the spatial distribution of ecosystem types is available for an area, the mean deviation from the reference state can be determined on the basis of an area grid. At a resolution of 2.5 km × 2.5 km per raster element, for example, at least two representative vegetation surveys are recommended for evaluating indicators of soil condition. As in Method 1, a site-specific classification of the deviations from the ecosystem-specific reference condition for the indicators is first carried out for all the sites to be surveyed. For each raster element, the arithmetic mean of the deviation levels is determined, also to show tendencies within the deviation level. On this basis, as in Method 1, the medians are aggregated via the modal value to deviation and change levels at the levels of ecological functions or ecosystems (Sect. "Site-specific classification"). This finally makes it possible to transfer the site-related classifications to the study area grid.
Based on the 70 vegetation surveys available for the Kellerwald-Edersee National Park, the variance levels for five different ecosystem types and the indicators of ecosystem integrity were determined. For the moisture index as an indicator of water flow, there are 70 site-related classifications of deviations from the reference condition between very low and high. For each grid, the mean value of the deviation levels is calculated and used as the basis for the grid-related classification (Fig. 10). Information losses with regard to the variance and the maximum deviation in each grid cell can be avoided by displaying the classifications at the individual locations. The example of the Kellerwald [18,19] shows that a higher spatial differentiation can be achieved by using an area grid compared to the ecosystem type-specific classification, since the area is strongly dominated by the Moder-Buchen-Bergwald (D1-5n-C2, Moder beech forests of the montane level). Together with the variance levels for the other indicators (e.g., Kullback distance, similarity of species quantity distribution, pH, base saturation, C/N ratio), the variance and change levels at the levels of ecological functions or ecosystems are aggregated via the modal value, as in method (1).
To make the application of the classification model more effective for larger amounts of data, the OESI tool can be used as a functional extension of ArcGIS ® Desktop 10.2 [16]. It was implemented using the Python programming language and supports the classification of deviations from the reference state on the basis of 2711 rules for 60 ecosystem types and 6 indicators to date: Kullback distance, similarity of species quantity distribution, pH, base saturation, C/N ratio and moisture index.

Determination of temporal trends
From repeated recordings of the ecosystem condition, temporal trends can be derived and interpreted. This is demonstrated by the example of the ICP Forests Level II program In 1960, the stand matched all ecological parameters of the associated ecosystem type (Rohhumus-Fichten-Hochbergwald). This type corresponds to the natural type of forest that forms in the ridges of the Thuringian Forest under today's climatic conditions in self-organisation.
The data from 1995 onwards could be attributed to the effect of liming, whereby both the vegetation and the topsoil data point to a decreasing effect until 2009 with at least partial reversible development to the original spruce forest type. This becomes particularly clear in the development of base saturation. In vegetation, this effect is reflected in the occurrence of Rubus idaeus and Oxalis acetosella from 2001 and in their decline with simultaneous increases in Trientalis europaea and Deschampsia flexuosa in 2006. The significant narrowing of the C/N ratio of the topsoil and the increased N value in the needles compared to the reference condition observed since 1995 can be attributed essentially to the increased N release due to mineralisation of the topsoil due to liming.
In fact, at the end of the 1980s, the Institute of Forest Sciences Eberswalde (IFE) developed a method to remedy the Mg deficiency symptoms of the spruce in the Thuringian Forest that had appeared on a large scale at the time (personal communication by Prof. Dr  from the air, which, however, showed no effect due to insufficient quantities. As a result, it was decided to fly large areas of the Thuringian Forest, including the ridges around the Großer Eisenberg, with Kamsdorfer Mg marl (high proportions of CaCO 3 , MgCO 3 ).
With regard to the assessment of the longer-term effects of liming, it is remarkable that the proportion of basic cations in the spruce needles between 1995 and 2009 remains well below the reference state. A temporal trend was not detectable in the data. This obviously reflects the effect of nitrous gases, which in the eighties had led to extensive needle yellowing in the Thuringian Forest. The nutrient disharmonies in the needles induced by the N effect were not eliminated by the marl fertilization, but were obviously maintained as a result of the increased N mineralization and N uptake by the roots. The relative increase in DNP at the time of culmination, calculated from the development of the mean level, indicates accelerated growth due to the enhanced N input. These cause-effect relationships indicated here should be further analyzed against the background of the process data collected in the EU measurement program. It is clear that with the proven ecological changes compared to the 1960s, the habitat function protected by the FFH habitat type has been adversely affected. The observed reversible development of vegetation and topsoil condition, on the other hand, is positive.
Since the 1960s, the moisture level indicated by the observed vegetation formation has shifted half a degree towards drier condition and is now at the lower interval limit of the reference condition. A further warming to be expected as a result of climate change could lead to a further decrease in the spread of the natural spruce forest type in the ridges of the Thuringian Forest and a significant expansion of the spectrum of natural site tree species.

Discussion
Over the last couple of years, there has been a broad discussion about sustainable development, ecosystem integrity, ecosystem services and biodiversity in science and the public. However, aside from a bunch of definitory treatises, operationalized approaches which are based on indicators quantifying ecosystem functions and structures with data from monitoring programs of competent authorities are lacking [21].
Ecosystem integrity or related notions are referred to in several national and international biodiversity and ecosystem policies. However, it is still poorly defined and operationalised. Based on a broad literature review Roche and Campagne [21] identified five forms of  The vegetation survey from 1960 does not originate from the Level II data set, but from the W.I.E. database and was recorded by H. SCHLÜTER. Regarding the year in which the stock was established, there were deviations between the information on the website of the Thuringian State Institute for Forest, hunting and fishing (http:// www. thuer ingen. de/ imper ia/ md/ conte nt/ folder/ waldo ekolog/ waldz ustan dsueb erwac hung/ eisen berg09. pdf ) and the Level II dataset, that could not be fully clarified. For the calculation of net primary production, an age of 80 years in 1995 was assumed on the basis of the latter information. This may explain the low absolute values of the calculated NPP ecosystem integrity: 1. ecosystem integrity of wilderness, 2. ecosystem functional and structural integrity, 3. ecosystem stability and resilience, 4. ecosystem condition and 5. ecosystem quality and value. The concept of hemeroby is associated with form 1 proposing that natural state or ahemeroby can be defined by the absence evidence of past and actual human management [12]. The values of Hemeroby index are determined by the degree of occurrences of human pressures, generally indicated by land use, landscape patterns and species assemblages. Walz and Stein (2014) published a Hemeroby map of Germany linking some few surface covering data: the CORINE land cover data, the Base-Landscape Model of the Authoritative Topographic-Cartographic Information System and the Digital Land Cover model for Germany. A seven-point scale was used to classify land use by degree of hemeroby. Ecosystem functions were not taken into account, nor were data from ecological environmental monitoring. Following Walz and Stein [26], their hemeroby mapping approach "is inappropriate for a more accurate calculation of spatial extent and thus the monitoring of local and regional developments. " Therefore, in our investigation, we developed a methodology based on 14 indicators for six ecosystem functions (habitat function, net primary function, carbon sequestration, nutrient and water flux, resilience) by example of Germany. It allows assessments of ecosystem integrity changes by comparing current or prospective ecosystem states with ecosystem-type-specific reference states as described by quantitative indicators for 61 forest ecosystem types based on data before 1990 [10]. Advantages of the rule-based method are the increase in reproducibility and effectiveness in assessments of ecosystem integrity and the coverage of three spatial scale: the forest stand level as well as the regional and national level. As a limitation of the methodology, it should be noted that the use of the mode in the aggregation scheme is associated with a leveling within the range of the ratings. For special issues (e.g., early warning) other linkage algorithms (e.g., maximum) may be more appropriate. Due to the formal classification at the levels of the indicators and ecosystem functions, the evaluation always remains comprehensible in detail. In addition, thresholds derived from the classification can provide orientation for deriving existing needs for action.

Conclusion
A methodological gap in the operationalization of the ecological integrity of forests on the basis of generally available data in Germany was closed by the method presented in this article. Following Roche and Campagne [21] who identified five major forms of ecosystem integrity concepts, the approach presented operationalises the ecosystem functional and structural integrity. Thereby, the functional and structural indicators were quantified by data from monitoring programs of competent authorities. Opportunities for the further research are the extension of the methodology to agrarian ecosystems in Germany and forest ecosystems of Europe. Following meteorology, the method presented refers to a historic reference period. Climatological reference periods usually cover 30 years so that the statistical parameters of the various climatological parameters can be determined with satisfactory accuracy. Longer periods are not used, because then climatic changes influence the series and also in many cases the data basis becomes too scarce. The World Meteorological Organization (WMO) has defined the periods 1931-1960, 1961-1990, 1991-2020 and 2021-2050 1961-1990 as the valid international climatological reference periods. Since meteorological characteristics are the most powerful drivers of ecosystem development, the temporal parallelization of climatic and ecological monitoring periods is technically required. This is all the more true as climate change models are calculated for climate reference periods and a coupling of climate change models with ecosystem change models requires identical time references.
The deviance of a current or potential future ecosystem condition from that of an ecosystem in the reference period does not indicate a drop in ecosystem integrity but a transition of the integrity of an ecosystem type in the reference period to the integrity of a current or future ecosystem condition. This transition is a period with decreasing stability of the ecosystem (integrity) in the reference period which may end up with an upcoming new ecosystem type replacing the former one of the reference period. So, the method ordinates the transition of the integrity of an ecosystem type in the reference period to the integrity of a current or future ecosystem condition. This transition may end up with an upcoming new integer ecosystem type replacing the former one of the reference period.