Recent and Future Changes in Rainfall Erosivity and Implications for the Soil Erosion Risk in Brandenburg, NE Germany

The universal soil loss equation (USLE) is widely used to identify areas of erosion risk at regional scales. In Brandenburg, USLE R factors are usually estimated from summer rainfall, based on a relationship from the 1990s. We compared estimated and calculated factors of 22 stations with 10-minutes rainfall data. To obtain more realistic estimations, we regressed the latter to three rainfall indices (total and heavy-rainfall sums). These models were applied to estimate future R factors of 188 climate stations. To assess uncertainties, we derived eight scenarios from 15 climate models and two representative concentration pathways (RCP), and compared the effects of index choice to the choices of climate model, RCP, and bias correction. The existing regression model underestimated the calculated R factors by 40%. Moreover, using heavy-rainfall sums instead of total sums explained the variability of current R factors better, increased their future changes, and reduced the model uncertainty. The impact of index choice on future R factors was similar to the other choices. Despite all uncertainties, the results indicate that average R factors will remain above past values. Instead, the extent of arable land experiencing excessive soil loss might double until the mid-century with RCP 8.5 and unchanged land management.


Introduction
Soils are a fundamental resource for life on Earth and provide numerous goods and services for the human society [1]. The degradation of soils poses a global threat to our well-being, mainly due to soil erosion by water [2,3] (henceforth soil erosion) exceeding the 'tolerable' natural formation rate of around 1 t ha 1 a 1 [4]. Among other consequences, soil erosion can hamper the sustainable agricultural production [5], impair water quality and habitats [6,7], and reduce the lifetime of reservoirs [8].
Unsustainable agriculture is a key driver of soil erosion, not only in Europe [3,9,10], especially when heavy rainfall meets inappropriate management. To mitigate agricultural impacts on soils, the Common Agricultural Policy (CAP) of the European Union introduced the Good Agricultural and Environmental Conditions (GAEC) in 2003 as a set of environmental standards and rules on cross compliance for financially supported farmers. Accordingly, German laws and regulations on both federal and state levels address soil erosion, not only to implement the GAEC, but also within the European Water Framework Directive. To estimate potential soil erosion, derivatives of the universal soil loss equation (USLE) are commonly used. For instance, German federal states apply an adapted version of the USLE [11,12] to identify areas prone to soil loss and to impose countermeasures on farmers. The USLE estimates the soil erosion rate from five factors, namely rainfall erosivity (R factor), soil erodibility (K factor), slope length and steepness (LS factor), soil cover (C factor), and soil conservation (P factor).
The R factor of the USLE is the long-term average annual sum of the rainfall erosivity. The erosivity expresses the capacity of a rainfall to induce soil erosion and integrates its duration, amount, and intensity (cf. Appendix A, also for the units of the USLE factors). The calculation of R factors requires long time series of precipitation data at high temporal resolution which are often unavailable. Many studies thus rely on daily to annual data for extrapolation in space and time [13][14][15]. Accordingly, the German industrial norm (DIN) on soil erosion [11] lists linear regression models to estimate R factors from summer or annual rainfall. While they were derived at the state level, the explained variability of the calculated R was unsatisfying in Brandenburg and other federal states in the North German Plain, with Pearson's r being lower than in other parts of Germany [11]. The DIN models were established during the early 1990s, so they also reflect the climate from the 1960s to the 1980s (cf. [16]). The lack of accuracy and the age raise the question whether the DIN model can be recommended to estimate current and future R factors for land and water management in Brandenburg.
Climate change can exhibit direct e↵ects (changes in the amount, intensity, and distribution of rainfall) and indirect e↵ects on soil erosion (rainfall and temperature changes a↵ect biomass production, soil moisture, and the growing season) [17]. Numerous trend studies have assessed past changes in precipitation, especially changes in extreme precipitation across Europe [18][19][20][21], in (Northeastern, NE) Germany [22][23][24][25] and neighboring countries [23,26,27]. These studies consistently found increasing winter rainfall in Central Europe, while the trends for summer precipitation are less coherent. The latter studies indicate less rainfall and longer dry spells during summer, although this holds not necessarily true for extreme rainfall [23] and when recent data is included [26]. The few regional to national studies on rainfall erosivity in Central/Western Europe show a multi-decadal variability with an increase in erosion risks since the end of the 20th century [28,29], although the strength and direction of regional trends are not necessarily valid everywhere [30].
Additionally, studies on future R factors show contradictory changes for NE Germany. While a European assessment proposed that R factors might approximately double from 2010 to 2050 even under moderate climate change scenarios [31], a German-wide analysis reported an average increase of 10% for 2011 to 2041 and a decline in the same order of magnitude for 2041 to 2071 compared to the reference period 1971-2000 [32]. Both studies did not assess how the choices of their respective climate model and rainfall indices to estimate the R factors a↵ected their scenario results. According to more recent German ensemble studies, there might be a tendency of decreasing summer totals, but increasing summer extremes [33,34]. Thus, di↵erent choices of aggregated rainfall indices might change the direction of future R changes, with implications for discussions of a sustainable regional land and water management under climate change.
The main objective of this study is to assess the impact of climate change on rainfall erosivity and, subsequently, on the potential soil erosion risk in Brandenburg, without considering further impacts on land use and vegetation cover. Using an ensemble of climate scenarios, our study addresses the following main questions: • Which aggregated rainfall index is best to estimate current R factors and their recent change in NE Germany? • How does climate change a↵ect regional R factors and the risk of soil erosion? • How does the rainfall index a↵ect future trends and how does the impact compare to other sources of uncertainty such as the choice of climate model, bias correction, and RCP scenario? Three aggregated rainfall indices were derived from daily data and tested to explain the variability of the calculated R factors with linear regression models. We compared the total sum from May to September (P sum ) to the total sum of rainfall on heavy rainfall days. As no common criteria for "heavy rainfall day" exists, we considered exemplarily the 10 highest daily values (P max10 ) and daily values above 11.8 mm (the 33th percentile of all calculated erosive rainfall events, which is above the threshold of 10 mm d 1 used for a "heavy rainfall day" e.g., by Reference [37]). In contrast to the DIN and similar empirical models, we excluded the October because rainfalls were less erosive than in summer months  Figure 2 in the Results section). For each climate station, P sum , P max10 , and P 11.8 were calculated as multi-annual means for the years in Table 1.
We applied the new regression models to the German-wide dataset of regionalized daily station data (REGNIE, [36]) to create new state-wide R maps for the years 2001-2015. These maps replaced the fixed R factor of 50 kJ m 2 mm h 1 of the current erosion map [38], e.g., used by the State O ce of Environment to identify risk areas, i.e., arable land with erosion rates above 1 t ha 1 a 1 . Likewise, we obtained the current extent of risk areas. For this reason, we used the REGNIE data to establish the regression models. As REGNIE is derived from station data, di↵erences between both only occur if stations are not considered for the regionalization. and daily values above 11.8 mm (the 33th percentile of all calculated erosive rainfall events, which is above the threshold of 10 mm d 1 used for a "heavy rainfall day" e.g., by Reference [37]). In contrast to the DIN and similar empirical models, we excluded the October because rainfalls were less erosive than in summer months (cf. Figure 2 in the Results section). For each climate station, Psum, Pmax10, and P11.8 were calculated as multi-annual means for the years in Table 1.
We applied the new regression models to the German-wide dataset of regionalized daily station data (REGNIE, [36]) to create new state-wide R maps for the years 2001-2015. These maps replaced the fixed R factor of 50 kJ m 2 mm h 1 of the current erosion map [38], e.g., used by the State Office of Environment to identify risk areas, i.e., arable land with erosion rates above 1 t ha 1 a 1 . Likewise, we obtained the current extent of risk areas. For this reason, we used the REGNIE data to establish the regression models. As REGNIE is derived from station data, differences between both only occur if stations are not considered for the regionalization.

Climate Scenarios
We used climate model data from the Coordinated Regional Downscaling Experiment which provides the most recent climate data for Europe (EURO-CORDEX domain, Reference [39]). Harmonized daily datasets of climate parameters are available in high spatial resolution from various GCM (general circulation models), RCM (regional circulation models), and RCP  (Table 1), all stations with daily data to estimate R factors, bias correction of the climate models, and scenario analyses, black: stations with data for 1971-2015, hollow: stations with shorter time series. To harmonize the time series and to fill gaps we extracted the time series from regionalized station data.

Climate Scenarios
We used climate model data from the Coordinated Regional Downscaling Experiment which provides the most recent climate data for Europe (EURO-CORDEX domain, Reference [39]). Harmonized daily datasets of climate parameters are available in high spatial resolution from various GCM (general circulation models), RCM (regional circulation models), and RCP (representative concentration pathways). This combination makes the CORDEX data most suitable for assessing uncertainties related to the selection of climate models and RCP scenarios in regional studies. The data availability in the CORDEX database di↵ers among the RCP. To obtain a su ciently large, yet homogeneous ensemble for the uncertainty assessment, we used 15 combinations of GCM and dynamical RCM in the highest possible spatial resolution of 0.11 (approx. 12 km, CORDEX domain EUR-11) available for RCP 4.5 (i.e., moderate climate change) and RCP 8.5 (i.e., extreme climate change) ( Table 3). We downloaded the CORDEX data for the years 1971 to 2100 and extracted the time-series for 188 climate stations in and around Brandenburg (black circles in Figure 1). To reduce the deviations to station data, we applied a bias correction (distribution mapping with linear scaling) to adjust the monthly frequency distributions and means for the period 1971-2015. This period includes hindcasted data where both RCP are similar  and climate change projections where the RCP deviate (2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015). We used an adapted version of the software CMHyd (References [40,41], see Reference [42] for the conceptual and methodical background) to perform the bias correction. The software also transformed the 360-days calendar of the GCM MOHC-HadGEM2-ES to a 365-days calendar.
While most stations had su ciently long time series for the bias correction ( Figure 1), we considered additional stations with shorter time series until 2015 to increase the density for the spatial interpolation (Chapter 2.4). Thus, we used again the REGNIE data instead of the station data. P sum , P max10 , and P 11.8 were calculated from the climate model output and the REGNIE data. We used a simple ranking scheme based on the average Kling-Gupta e ciency (KGE, [43]), the absolute percent bias (APB, in %), and the root-mean square error (RMSE) to assess the performance of the bias-corrected climate models for the period 1971-2015. The ranking considered the • Spatial rank: For each station, the long-term average rainfall indices were calculated. KGE, APB and RMSE were obtained for each climate model.

•
Trend rank: For each station, we calculated the annual indices from REGNIE and the climate models and determined KGE, APB and RMSE from the linear trend. For each climate model, we averaged KGE, APB, and RMSE for the ranking.
The rank of each climate models was determined with the function rank.avg in MS Excel which returns the average rank for equal values, e.g., eight if all 15 models would perform equally well for an index. To rank almost similar values equally, we rounded the KGE to the nearest 0.05 and APB and RMSE to the nearest integer. The five best-performing climate models were selected (subset) and compared to the whole ensemble of 15 models to evaluate the e↵ect of model choice on R trends.

Impact of Climate Change on R Factors, Uncertainties, Consequences for the Extent of Erosion-Risk Areas
For all 188 climate stations and the 8 combination of climate models (whole ensemble or subset), bias correction (with or without), RCP (4.5 or 8.5), and index (heavy-rainfall or total sum), we derived the change in R relative to 2001-2015 for other periods of 15 years overlapping by five years (i.e., 2011-2025, 2021-2035, etc.) to assess how modelling choices influenced future trends. Each choice resulted in two di↵erent R factors (R and R') for all the stations and periods. We considered this di↵erence as uncertainty. We averaged the ratios of max(R, R') and min(R, R') over four periods in the first and the second half of this century  to compare the impact of each choice.
Finally, the relative station values were interpolated using thin splines [44] (a) to assess the sub-regional pattern of R changes, and (b) to quantify the change in the extent of arable land [45] at risk of erosion rates exceeding the threshold of 1 t ha 1 a 1 . For the latter, we applied the interpolated data as correction factor to the revised state erosion map (Chapter 2.2).

Rainfall Indices and the Spatial Variability of Current R Factors
The calculated R factors considerably exceeded the values estimated with the DIN regression model (Figure 2a), especially in the south-eastern part where the highest R factors were calculated. While the average calculated value was 76.8 kJ m 2 mm h 1 , the DIN model estimated only 42.7 kJ m 2 mm h 1 . The fixed value used for the existing erosion map corresponds to the maximum of the estimated values but to the minimum of the calculated values. In consequence, the current state-wide erosion map underestimates the potential erosion rate and the extent of risk areas.
About 90% of the annual rainfall erosivity occurred between May and September, 50% in July and August ( Figure 2b). During the rest of the year, rainfall events were far less erosive. As aggregated rainfall indices cannot capture the lower peak intensities in these months and would overestimate their contribution to R factors, we did not consider these months in our rainfall indices for the regression analyses.
Water 2019, 11, 904 6 of 18 an index. To rank almost similar values equally, we rounded the KGE to the nearest 0.05 and APB and RMSE to the nearest integer. The five best-performing climate models were selected (subset) and compared to the whole ensemble of 15 models to evaluate the effect of model choice on R trends.

Impact of Climate Change on R Factors, Uncertainties, Consequences for the Extent of Erosion-Risk Areas
For all 188 climate stations and the 8 combination of climate models (whole ensemble or subset), bias correction (with or without), RCP (4.5 or 8.5), and index (heavy-rainfall or total sum), we derived the change in R relative to 2001-2015 for other periods of 15 years overlapping by five years (i.e., 2011-2025, 2021-2035, etc.) to assess how modelling choices influenced future trends. Each choice resulted in two different R factors (R and R') for all the stations and periods. We considered this difference as uncertainty. We averaged the ratios of max(R, R') and min(R, R') over four periods in the first and the second half of this century  to compare the impact of each choice.
Finally, the relative station values were interpolated using thin splines [44] (a) to assess the sub-regional pattern of R changes, and (b) to quantify the change in the extent of arable land [45] at risk of erosion rates exceeding the threshold of 1 t ha 1 a 1 . For the latter, we applied the interpolated data as correction factor to the revised state erosion map (Chapter 2.2).

Rainfall Indices and the Spatial Variability of Current R Factors
The calculated R factors considerably exceeded the values estimated with the DIN regression model (Figure 2a), especially in the south-eastern part where the highest R factors were calculated. While the average calculated value was 76.8 kJ m 2 mm h 1 , the DIN model estimated only 42.7 kJ m 2 mm h 1 . The fixed value used for the existing erosion map corresponds to the maximum of the estimated values but to the minimum of the calculated values. In consequence, the current state-wide erosion map underestimates the potential erosion rate and the extent of risk areas.
About 90% of the annual rainfall erosivity occurred between May and September, 50% in July and August (Figure 2b). During the rest of the year, rainfall events were far less erosive. As aggregated rainfall indices cannot capture the lower peak intensities in these months and would overestimate their contribution to R factors, we did not consider these months in our rainfall indices for the regression analyses.  The heavy-rainfall indices P max10 and P 11.8 explained the variability of the 22 calculated R factors better (r 2 = 0.46-0.50) than the total sum P sum (r 2 = 0.21, and 0.08 if calculated for May to October as used in the DIN). However, the performance was strongly a↵ected by the most extreme event, its erosivity (EI 30 value) being more than twice the second highest value (494 kJ m 2 mm h 1 , Appendix C). Without this event, the corresponding R factor (of station 430, Table 1) decreased by 33% resulting in r 2 values of around 0.6 (Equations (2) and (3)) and 0.3 (Equation (1)).
Using heavy-rainfall indices instead of total sums a↵ected the spatial pattern of the estimated R factor (Figure 3). While the outcomes of Equations (1) (horizontal line). Both, calculation and estimation of R factors according to the German norm DIN; (b) average intra-annual distribution of the rainfall erosivity. Rainfall from October to April was far less erosive than summer rainfall and therefore not considered in the revised regression models for R factors in Brandenburg.
The heavy-rainfall indices Pmax10 and P11.8 explained the variability of the 22 calculated R factors better (r² = 0.46-0.50) than the total sum Psum (r² = 0.21, and 0.08 if calculated for May to October as used in the DIN). However, the performance was strongly affected by the most extreme event, its erosivity (EI30 value) being more than twice the second highest value (494 kJ m 2 mm h 1 , Appendix C). Without this event, the corresponding R factor (of station 430, Table 1) decreased by 33% resulting in r² values of around 0.6 (Equations (2) and (3)) and 0.3 (Equation (1)).
Using heavy-rainfall indices instead of total sums affected the spatial pattern of the estimated R factor (Figure 3). While the outcomes of Equations (1)

Ranking of Climate Models
No single GCM or RCM was unanimously ranked as "best matching" (Figure 4, Appendix D), since 4 combinations of GCM and RCM were equally best performing (ID 2, 5, 10, and 11, Table 3). The 5th rank was ambivalent with different results for both RCP. However, with the chosen model 7, the ensemble subset comprised 4 out of the 5 GCM and 3 out of the 6 RCM.

Ranking of Climate Models
No single GCM or RCM was unanimously ranked as "best matching" (Figure 4, Appendix D), since 4 combinations of GCM and RCM were equally best performing (ID 2, 5, 10, and 11, Table 3). The 5th rank was ambivalent with di↵erent results for both RCP. However, with the chosen model 7, the ensemble subset comprised 4 out of the 5 GCM and 3 out of the 6 RCM.   Table 3. The lower the rank, the better the model performance. Selected climate models for the ensemble subset in black.

Climate Scenarios and Soil Erosion Risk
Only the subset of the 5 top-ranked climate models captured, albeit underestimated, the general increase of the R factors during the past four decades (R values relative to 2001-15 below 1.0 in Figure 5). In contrast, the whole ensemble did not reveal this increase (past values 1.0). On average, the RCP 4.5 values (relative R of 0.81) were closer to the REGNIE data (relative R of 0.71) than the RCP 8.5 values (0.90, RCP not differentiated in Figure 5).
The variability was largest for the whole uncorrected ensemble, partly due to negative regional R factors obtained with the climate model ID 9 and RCP 4.5. For instance, the coefficient of variation (CV) was here above 22,500% for Psum (Equation (1), but 22.9% without the negative values) and above 30% for Pmax10, P11.8, and RCP 8.5. The ranking helped to reduce the CV to 12%-19%, which was close to the reference data (11%).   Table 3. The lower the rank, the better the model performance. Selected climate models for the ensemble subset in black.

Climate Scenarios and Soil Erosion Risk
Only the subset of the 5 top-ranked climate models captured, albeit underestimated, the general increase of the R factors during the past four decades (R values relative to 2001-15 below 1.0 in Figure 5). In contrast, the whole ensemble did not reveal this increase (past values ⇡ 1.0). On average, the RCP 4.5 values (relative R of 0.81) were closer to the REGNIE data (relative R of 0.71) than the RCP 8.5 values (0.90, RCP not di↵erentiated in Figure 5).
The variability was largest for the whole uncorrected ensemble, partly due to negative regional R factors obtained with the climate model ID 9 and RCP 4.5. For instance, the coe cient of variation (CV) was here above 22,500% for P sum (Equation (1), but 22.9% without the negative values) and above 30% for P max10 , P 11.8 , and RCP 8.5. The ranking helped to reduce the CV to 12-19%, which was close to the reference data (11%).   Table 3. The lower the rank, the better the model performance. Selected climate models for the ensemble subset in black.

Climate Scenarios and Soil Erosion Risk
Only the subset of the 5 top-ranked climate models captured, albeit underestimated, the general increase of the R factors during the past four decades (R values relative to 2001-15 below 1.0 in Figure 5). In contrast, the whole ensemble did not reveal this increase (past values 1.0). On average, the RCP 4.5 values (relative R of 0.81) were closer to the REGNIE data (relative R of 0.71) than the RCP 8.5 values (0.90, RCP not differentiated in Figure 5).
The variability was largest for the whole uncorrected ensemble, partly due to negative regional R factors obtained with the climate model ID 9 and RCP 4.5. For instance, the coefficient of variation (CV) was here above 22,500% for Psum (Equation (1), but 22.9% without the negative values) and above 30% for Pmax10, P11.8, and RCP 8.5. The ranking helped to reduce the CV to 12%-19%, which was close to the reference data (11%).  Regarding future trends, all ensembles and regression models predicted increasing R factors towards the end of the 21st century ( Figure 6), except for Equation (1) in combination with the whole uncorrected ensemble and RCP 4.5 (Figure 6a). However, this combination was again hampered by negative R factors.
As expected, the model ensemble and the RCP influenced the decadal variability. The ensemble subset amplified the di↵erences between RCP 4.5 and 8.5. Until mid-century, the estimated regional R factors would generally remain above the reference period if the whole ensemble was used (solid lines in Figure 6). With the subset (broken line), the overall changes might be negligible and even negative for RCP 4.5, but increase more strongly for RCP 8.5. Regarding future trends, all ensembles and regression models predicted increasing R factors towards the end of the 21st century ( Figure 6), except for Equation (1) in combination with the whole uncorrected ensemble and RCP 4.5 (Figure 6a). However, this combination was again hampered by negative R factors.
As expected, the model ensemble and the RCP influenced the decadal variability. The ensemble subset amplified the differences between RCP 4.5 and 8.5. Until mid-century, the estimated regional R factors would generally remain above the reference period if the whole ensemble was used (solid lines in Figure 6). With the subset (broken line), the overall changes might be negligible and even negative for RCP 4.5, but increase more strongly for RCP 8.5. However, switching from Psum to the heavy-rainfall indices also resulted in an increase of future R factors. The effect was larger for the uncorrected than for the bias-corrected climate models because the bias correction raised the ratio of heavy and total rainfall (Figure 7). The impacts of index choice on trends were comparable to the other modelling decisions (Figure 8).  However, switching from P sum to the heavy-rainfall indices also resulted in an increase of future R factors. The e↵ect was larger for the uncorrected than for the bias-corrected climate models because the bias correction raised the ratio of heavy and total rainfall (Figure 7). The impacts of index choice on trends were comparable to the other modelling decisions (Figure 8). Regarding future trends, all ensembles and regression models predicted increasing R factors towards the end of the 21st century ( Figure 6), except for Equation (1) in combination with the whole uncorrected ensemble and RCP 4.5 (Figure 6a). However, this combination was again hampered by negative R factors.
As expected, the model ensemble and the RCP influenced the decadal variability. The ensemble subset amplified the differences between RCP 4.5 and 8.5. Until mid-century, the estimated regional R factors would generally remain above the reference period if the whole ensemble was used (solid lines in Figure 6). With the subset (broken line), the overall changes might be negligible and even negative for RCP 4.5, but increase more strongly for RCP 8.5. However, switching from Psum to the heavy-rainfall indices also resulted in an increase of future R factors. The effect was larger for the uncorrected than for the bias-corrected climate models because the bias correction raised the ratio of heavy and total rainfall (Figure 7). The impacts of index choice on trends were comparable to the other modelling decisions (Figure 8).  The combination of heavy-rainfall index, RCP 4.5, and ensemble subsetting suggested very low overall changes in R factors around mid-century (2041-2065) compared to 2001-2015 ( Figure 6). However, the spatial aggregation masked distinctively opposite trends in different sub-regions. While R factors in some parts of Northern Brandenburg might increase by more than 15%, they can decrease by more than 10% in the South-Eastern part (Figure 9). Due to the spatial variability of R changes, of the other USLE factors, and of the distribution of arable land, changes of R and the extent of risk areas on arable land can differ. In Brandenburg, the risk areas might in fact increase more than the R factor, for instance from currently 3% to 6% by The combination of heavy-rainfall index, RCP 4.5, and ensemble subsetting suggested very low overall changes in R factors around mid-century (2041-2065) compared to 2001-2015 ( Figure 6). However, the spatial aggregation masked distinctively opposite trends in di↵erent sub-regions. While R factors in some parts of Northern Brandenburg might increase by more than 15%, they can decrease by more than 10% in the South-Eastern part (Figure 9). The combination of heavy-rainfall index, RCP 4.5, and ensemble subsetting suggested very low overall changes in R factors around mid-century (2041-2065) compared to 2001-2015 ( Figure 6). However, the spatial aggregation masked distinctively opposite trends in different sub-regions. While R factors in some parts of Northern Brandenburg might increase by more than 15%, they can decrease by more than 10% in the South-Eastern part (Figure 9).  (2) and (3). Values interpolated with thin splines. The overall change for this combination is close to 0%.
Due to the spatial variability of R changes, of the other USLE factors, and of the distribution of arable land, changes of R and the extent of risk areas on arable land can differ. In Brandenburg, the risk areas might in fact increase more than the R factor, for instance from currently 3% to 6% by  (2) and (3). Values interpolated with thin splines. The overall change for this combination is close to 0%.
Due to the spatial variability of R changes, of the other USLE factors, and of the distribution of arable land, changes of R and the extent of risk areas on arable land can di↵er. In Brandenburg, the risk areas might in fact increase more than the R factor, for instance from currently 3% to 6% by mid-century (2041-2065) with RCP 8.5 under status-quo land management, if the bias-corrected ensemble subset was applied (Figure 10), while the average R factor may raise by about 40% (Figure 6, left column). In accordance to the R factor, the heavy-rainfall indices favored a more pronounced change in risk areas than the total sum.  (Figure 10), while the average R factor may raise by about 40% ( Figure  6, left column). In accordance to the R factor, the heavy-rainfall indices favored a more pronounced change in risk areas than the total sum.

Indices to Estimate Current R Factors
The variability of total summer rainfall does not match the variability of R factors in Brandenburg. This explains the lack of accuracy of the existing DIN model, which cannot be resolved by a re-calibration. Instead, heavy-rainfall sums are more suitable for regional assessments of soil erosion risks. By changing the rainfall index, the estimated R factors and potential soil loss in the Southern part increase relatively to the northern part which corresponds to its higher calculated R factors (the 4 stations with R > 90 kJ m 2 mm h 1 are located in the South-Eastern part).
Given the low performance of the existing DIN models across the Central Plains ecoregion [11], this adjustment can help to improve the delineation of risk areas in other regions as well. However, the application of such simple empirical relationships is inherently impaired by the aggregated indices, and much of the unexplained variability of 40% has to be attributed to our use of daily data. Soil erosion often occurs during short peak intensities and the extreme erosivity of single events can considerably affect local R factors as multi-year averages. Likewise, interpolated station data cannot capture the high variability of rainfall erosivity even over short distances [46,47]. Better regional and national datasets are unavailable. Nonetheless, the data quality is continuously improving, including radar data [48,49].

Climate Change Impacts on R Factors and the Soil Erosion Risk
The significant underestimation of the calculated recent R factors by the existing DIN regression model indicates that the climate-based erosion risk increased since the 1990s in NE Germany. This corresponds well to rising total and heavy rainfalls in summer and is confirmed by an independent multi-decadal dataset with steadily increasing maximum and minimum annual values ( Figure 11) and lower historical values for NE Germany [16]. Similar increases during the last decades were also reported for Western Germany [28] and recently deducted for the whole country [49].
Despite all model uncertainty, the climate scenarios consistently suggest that the rainfall erosivity and the erosion risk for the whole of Brandenburg might not return to levels of 2-3 decades

Indices to Estimate Current R Factors
The variability of total summer rainfall does not match the variability of R factors in Brandenburg. This explains the lack of accuracy of the existing DIN model, which cannot be resolved by a re-calibration. Instead, heavy-rainfall sums are more suitable for regional assessments of soil erosion risks. By changing the rainfall index, the estimated R factors and potential soil loss in the Southern part increase relatively to the northern part which corresponds to its higher calculated R factors (the 4 stations with R > 90 kJ m 2 mm h 1 are located in the South-Eastern part).
Given the low performance of the existing DIN models across the Central Plains ecoregion [11], this adjustment can help to improve the delineation of risk areas in other regions as well. However, the application of such simple empirical relationships is inherently impaired by the aggregated indices, and much of the unexplained variability of 40% has to be attributed to our use of daily data. Soil erosion often occurs during short peak intensities and the extreme erosivity of single events can considerably a↵ect local R factors as multi-year averages. Likewise, interpolated station data cannot capture the high variability of rainfall erosivity even over short distances [46,47]. Better regional and national datasets are unavailable. Nonetheless, the data quality is continuously improving, including radar data [48,49].

Climate Change Impacts on R Factors and the Soil Erosion Risk
The significant underestimation of the calculated recent R factors by the existing DIN regression model indicates that the climate-based erosion risk increased since the 1990s in NE Germany. This corresponds well to rising total and heavy rainfalls in summer and is confirmed by an independent multi-decadal dataset with steadily increasing maximum and minimum annual values ( Figure 11) and lower historical values for NE Germany [16]. Similar increases during the last decades were also reported for Western Germany [28] and recently deducted for the whole country [49].
Despite all model uncertainty, the climate scenarios consistently suggest that the rainfall erosivity and the erosion risk for the whole of Brandenburg might not return to levels of 2-3 decades ago-as observed during the 20th century in Western Germany [28] and suggested before to occur mid-century in NE Germany [32]. While the GCM HadGEM2 appeared twice in our ensemble subset, the rather low mid-century impact of a moderate climate change (RCP 4.5 in Figure 9) is nonetheless in contrast to the significant increase as derived in a European analysis [31]. However, the latter two studies relied only on a single and statistically downscaled GCM rather than an ensemble of dynamical RCM as used by CORDEX. The RCM type was recently shown to a↵ect the direction of rainfall trends in Germany [33] (and also noted in Reference [32]).
A moderate climate change might, thus, open a window of opportunity to discuss and establish climate change adaptation, especially in southern Brandenburg. More soil conservation and better soil coverage on arable land are e↵ective agricultural measures to mitigate the on-site and o↵-site impacts of intensified soil erosion under climate change.
Water 2019, 11, 904 12 of 18 ago-as observed during the 20th century in Western Germany [28] and suggested before to occur mid-century in NE Germany [32]. While the GCM HadGEM2 appeared twice in our ensemble subset, the rather low mid-century impact of a moderate climate change (RCP 4.5 in Figure 9) is nonetheless in contrast to the significant increase as derived in a European analysis [31]. However, the latter two studies relied only on a single and statistically downscaled GCM rather than an ensemble of dynamical RCM as used by CORDEX. The RCM type was recently shown to affect the direction of rainfall trends in Germany [33] (and also noted in Reference [32]). A moderate climate change might, thus, open a window of opportunity to discuss and establish climate change adaptation, especially in southern Brandenburg. More soil conservation and better soil coverage on arable land are effective agricultural measures to mitigate the on-site and off-site impacts of intensified soil erosion under climate change.

Impact of Index Choice on Future R Factors and other Sources of Uncertainty
The choice of the R estimator, i.e., the replacement of the currently used total sum by heavy-rainfall indices, can affect the trend similarly to other modelling decisions as sources of uncertainty in scenario analyses. The impact on trend strength-the increase of R and the erosion risk was stronger for heavy-rainfall than for total sums-supports findings, e.g., by References [33] and [34]: future summers might be characterized by more concentrated rainfall and longer dry spells, especially with RCP 8.5 ( Figure 12).

Impact of Index Choice on Future R Factors and other Sources of Uncertainty
The choice of the R estimator, i.e., the replacement of the currently used total sum by heavy-rainfall indices, can a↵ect the trend similarly to other modelling decisions as sources of uncertainty in scenario analyses. The impact on trend strength-the increase of R and the erosion risk was stronger for heavy-rainfall than for total sums-supports findings, e.g., by References [33] and [34]: future summers might be characterized by more concentrated rainfall and longer dry spells, especially with RCP 8.5 ( Figure 12). ago-as observed during the 20th century in Western Germany [28] and suggested before to occur mid-century in NE Germany [32]. While the GCM HadGEM2 appeared twice in our ensemble subset, the rather low mid-century impact of a moderate climate change (RCP 4.5 in Figure 9) is nonetheless in contrast to the significant increase as derived in a European analysis [31]. However, the latter two studies relied only on a single and statistically downscaled GCM rather than an ensemble of dynamical RCM as used by CORDEX. The RCM type was recently shown to affect the direction of rainfall trends in Germany [33] (and also noted in Reference [32]). A moderate climate change might, thus, open a window of opportunity to discuss and establish climate change adaptation, especially in southern Brandenburg. More soil conservation and better soil coverage on arable land are effective agricultural measures to mitigate the on-site and off-site impacts of intensified soil erosion under climate change.

Impact of Index Choice on Future R Factors and other Sources of Uncertainty
The choice of the R estimator, i.e., the replacement of the currently used total sum by heavy-rainfall indices, can affect the trend similarly to other modelling decisions as sources of uncertainty in scenario analyses. The impact on trend strength-the increase of R and the erosion risk was stronger for heavy-rainfall than for total sums-supports findings, e.g., by References [33] and [34]: future summers might be characterized by more concentrated rainfall and longer dry spells, especially with RCP 8.5 ( Figure 12).  Additionally, past and future R factors were more prone to model uncertainty if estimated from total instead of the heavy-rainfall sums. While a careful selection of the climate models and a bias correction can also help to reduce this uncertainty, both implicitly assume that model deviations are stationary which might not necessarily hold true under climate change conditions [50]. While certain combinations of GCM and RCM depicted the spatio-temporal variability of rainfall indices particularly well in the past, none can unanimously be recommended as "the best performing". Firstly, the CORDEX ensemble does not comprise all combinations of the available GCM and RCM. Secondly, the performance and ranking of climate models (and bias corrections) vary among ecoregions and depend on the hydro-climatic variables under consideration [41].
The similar results obtained with the two conceptually di↵erent heavy-rainfall indices show that the exact definition of "heavy rainfall days" is of minor relevance for the estimation of current and future R factors. The underestimation of the recent increase and the trend towards more and more heavy rainfall in summer point to the unknown e↵ect of climate change on peak intensities of storm events. Should extreme events, such as the maximum value in our time-series (Appendix C) or the reported values by Reference [51], occur more frequently, even our revised empirical relationships would underestimate the future erosion risk and would thus become obsolete. Acknowledgments: We thank the German Meteorological Service (DWD) for providing the precipitation data and for the revision of suspicious values. We very much acknowledge the World Climate Research Programme's Working Group on Regional Climate and the Working Group on Coupled Modelling, former coordinating body of CORDEX and responsible panel for CMIP5 as well as the climate modelling groups for producing and making available their model outputs (cf. Table 3). We also appreciate the fast review process and the helpful comments of two anonymous reviewers. J.K. and A.G. express their gratitude to H. Rathjens for providing the CMHyd software and the support with the bias correction.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Appendix A
The R factor of the USLE is the long-term average annual sum of EI 30 , the rainfall erosivity of a storm event. EI 30 is the product of the total kinetic energy E (in kJ m 2 mm h 1 ) and the peak 30-min intensity (I 30 , in mm h 1 ). Erosive storm events exceed the thresholds of 10 mm total rainfall or of 10 mm h 1 peak intensity within 30 min. Each event is separated by gaps of more than 6 h. Following [11,12], E was derived from the intensity I (mm h 1 ) and precipitation P (mm) for each time step i (Equation (A1)). Note, that the unit kJ m 2 mm h 1 (or N h 1 ) for our R factors equals 0.1 MJ mm ha 1 h 1 , e.g., used by [15]. Accordingly, the calculation of the soil loss in tons ha -1 requires USLE K factors in tons ha 1

Appendix B
The 10-min data was aggregated and compared to daily and hourly sums from [36]. As the control data is not a priori error-free, we asked the DWD to revise deviations (gaps and spurious values) with the strongest impact on the R factors ( [52,53], Table A1). The missing data at station 880 was estimated from rainfall radar data (RADOLAN RW) [36] adjusted to the conventionally measured value of 85.7 mm/6 h by assuming linear changes during each hour ( Figure A1). The EI 30 value increased from 5.6 to 163.7 kJ m 2 mm h 1 .

Appendix B
The 10-min data was aggregated and compared to daily and hourly sums from [36]. As the control data is not a priori error-free, we asked the DWD to revise deviations (gaps and spurious values) with the strongest impact on the R factors ( [52,53], Table A1). The missing data at station 880 was estimated from rainfall radar data (RADOLAN RW) [36] adjusted to the conventionally measured value of 85.7 mm/6 h by assuming linear changes during each hour ( Figure A1). The EI30 value increased from 5.6 to 163.7 kJ m 2 mm h 1 .

Appendix C
The highest and most intense rainfall event in the dataset occurred on 25 August 2006 at Berlin Tegel (station ID 430 in Table 1). During 70 minutes, between 17:10 and 18:19, more than 122 mm was registered which corresponds to more than 20% of the total precipitation in this year (597 mm, Figure A2). The peak intensity of I30 = 135 mm/h resulted in an erosivity of 494 kJ m 2 , or 86% of the annual EI30. The R factor for the years 2000-15 was 94 kJ m 2 mm h 1 with and 62 kJ m 2 mm h 1 without this event. Due to the much lower amount (cf. Figure A3) and intensity of rainfall, the daily erosivity at other climate stations was below 10 kJ m 2 mm h 1 . Information on the weather conditions in Central Europe at this time are available online [54] (in German).

Appendix C
The highest and most intense rainfall event in the dataset occurred on 25 August 2006 at Berlin Tegel (station ID 430 in Table 1). During 70 minutes, between 17:10 and 18:19, more than 122 mm was registered which corresponds to more than 20% of the total precipitation in this year (597 mm, Figure A2). The peak intensity of I 30 = 135 mm/h resulted in an erosivity of 494 kJ m 2 , or 86% of the annual EI 30 . The R factor for the years 2000-15 was 94 kJ m 2 mm h 1 with and 62 kJ m 2 mm h 1 without this event. Due to the much lower amount (cf. Figure A3) and intensity of rainfall, the daily erosivity at other climate stations was below 10 kJ m 2 mm h 1 . Information on the weather conditions in Central Europe at this time are available online [54] (in German).  The rainfall amount on this day was significantly lower elsewhere, as well as the rainfall erosivity (not shown). Data source: RADOLAN RW [36]. The linear features are artifacts due to radar shadows. Table A2 lists the aggregated rankings as shown in Figure 4. The spatial variability of Psum was, as expected, almost identical between the bias-corrected (KGE 1) climate models and the REGNIE data, and in general very good for the two other variables (KGE > 0.7). Although the mostly negative KGE for the linear trends are not satisfactory, the bias correction could significantly reduce RMSE and APB (Table A3).   The rainfall amount on this day was significantly lower elsewhere, as well as the rainfall erosivity (not shown). Data source: RADOLAN RW [36]. The linear features are artifacts due to radar shadows. Table A2 lists the aggregated rankings as shown in Figure 4. The spatial variability of Psum was, as expected, almost identical between the bias-corrected (KGE 1) climate models and the REGNIE data, and in general very good for the two other variables (KGE > 0.7). Although the mostly negative KGE for the linear trends are not satisfactory, the bias correction could significantly reduce RMSE and APB (Table A3).  The rainfall amount on this day was significantly lower elsewhere, as well as the rainfall erosivity (not shown). Data source: RADOLAN RW [36]. The linear features are artifacts due to radar shadows. Table A2 lists the aggregated rankings as shown in Figure 4. The spatial variability of P sum was, as expected, almost identical between the bias-corrected (KGE ⇡ 1) climate models and the REGNIE data, and in general very good for the two other variables (KGE > 0.7). Although the mostly negative KGE for the linear trends are not satisfactory, the bias correction could significantly reduce RMSE and APB (Table A3).