Mapping spatio-temporal dynamics of the cover and management factor (C-factor) for grasslands in Switzerland

The decrease in vegetation cover is one of the main triggering factors for soil erosion of grasslands. Within the Revised Universal Soil Loss Equation (RUSLE), a model commonly used to describe soil erosion, the vegetation cover for grassland is expressed in the cover and management factor (C-factor). The site-specific C-factor is a combination of the relative erosion susceptibility of a particular plant development stage (here expressed as soil loss ratio SLR) and the corresponding rainfall pattern (here expressed as R-factor ratio). Thus, for grasslands the fraction of green vegetation cover (FGVC) determines the SLRs. Although Switzerland is a country dominated by grassland with high percentages of mountainous regions and evidence for high erosion rates of grassland exists, soil erosion risk modeling of grasslands and especially of mountainous grasslands in Switzerland is restricted to a few studies. Here, we present a spatio-temporal approach to assess the dynamics of the C-factor for Swiss grasslands and to identify erosion prone regions and seasons simultaneously. We combine different satellite data, aerial data, and derivative products like Climate Change Initiative (CCI) Land Cover, Swissimage false-color infrared (Swissimage FCIR), PROBA-V Fraction of green Vegetation Cover (FCover300m), and MODIS Vegetation Indices 16-Day L3 Global (MOD13Q1) for the FGVC mapping of grasslands. In the spatial mapping, the FGVC is extracted from Swissimage FCIR (spat. res. 2 m) by linear spectral unmixing (LSU). The spatially derived results are then fused with the 10-day deviations of temporal FGVC derived by FCover300m. Following the original RUSLE approach, the combined FGVC are transformed to SLRs and weighted with high spatiotemporal resolved ratios of R-factors to result in spatio-temporal C-factors for Swiss grasslands. The annual average C-factor of all Swiss grasslands is 0.012. Seasonal and regional patterns (low C in winter, high C in summer, dependency on elevation) are recognizable in the spatio-temporal mapping approach. They are mainly explicable by the R-factor distribution within a year. Knowledge about the spatio-temporal dynamic of erosion triggering factors is of high interest for agronomists who can introduce areal and time specific selective erosion control measures and thereby reduce the direct costs of mitigation as well as erosion measures.


Introduction
Among all soil erosion risk factors in USLE-type (Universal Soil Loss Equation) and USLE based soil erosion models (e.g., RUSLE Revised Universal Soil Loss Equation), the cover and management factor namely C-factor is the one most sensitive as it follows plant growth and rainfall dynamics (Wischmeier and Smith, 1978;Nearing et al., 2005). The Cfactor represents the effect of cropping and management practices on soil erosion rates by water (Renard et al., 1997). The factor can be expressed as a combination of crop and plant systems, management, and rainfall pattern (Wischmeier and Smith, 1978). Following the USLE-original approach (Wischmeier and Smith, 1978;Schwertmann et al., 1987), a site-and time-specific C-factor is derived by the ratio of T temporal variation among all the RUSLE factors (Zhang et al., 2011;Estrada-Carmona et al., 2016). Thus, the factor can easily alter by a change of policy and farming strategies (McCool et al., 1995;Panagos et al., 2015a). An alteration of the support practice factor (P) (e.g., introducing of stone walls, grass margins, contour farming, terracing) often requires higher financial investments and soil conservation subsidies (Panagos et al., 2015b(Panagos et al., , 2015c. Other important soil erosion risk factors such as rainfall erosivity (R), soil erodibility (K) and topography (LS) are mainly determined by natural conditions and are relatively more independent from anthropogenic interventions.
Commonly, remote sensing approaches to determine the C-factors (Vrieling, 2006;Zhang et al., 2011;Panagos et al., 2014a) are not calculating SLRs but frequently assess the C-factor directly without weighting SLRs with the intra-annual distribution of rainfall erosivity to assess C-factors in the original sense of (R)USLE. Remote sensing methods for C-factor determination are often based on vegetation indices like the Normalized Difference Vegetation Index (NDVI). NDVIs are directly transformed to C-factors by a linear (de Jong et al., 1998) or exponential regression (van der Knijff et al., 2000) or related to field observations (Karaburun, 2010;Vatandaşlar and Yavuz, 2017). NDVI based C-factor modeling also exists for determining the C-factor for mountainous grasslands (regions of Korea, Lee and Won, 2012;China, Zhang and Li, 2015;Kyrgyzstan, Kulikov et al., 2016;Turkey, Vatandaşlar and Yavuz, 2017). However, drawbacks of that technique are its uncertainty due to the poor correlation with vegetation attributions, the soil reflectance, and the changing vitality of plants (de Jong, 1994;Vrieling, 2006;Asis and Omasa, 2007;Montandon and Small, 2008;Meusburger et al., 2010a;Grauso et al., 2015;Panagos et al., 2015b). As an alternative to NDVI-based approaches, spectral unmixing can estimate the fractional abundance of green vegetation (here called the fraction of green vegetation cover FGVC) and bare soils/bedrock simultaneously (Paringit and Nadaoka, 2003;Guerschman et al., 2009) which are related to C-factors after including rainfall erosivity (Yang, 2014). Spectral unmixing techniques (e.g., linear spectral unmixing LSU) are used in many erosion studies to determine C-factors over the last years (Hill et al., 1995;Ma et al., 2003;Lu et al., 2004;de Asis and Omasa, 2007;de Asis et al., 2008;de Jong and Epema, 2010;Meusburger et al., 2010aMeusburger et al., , 2010b). An advantage of spectral unmixing compared to traditional hard classification methods is the decomposition of mixed pixels in its corresponding component fractions rather than assigning them to a unique single class (Foody, 2006). Under consideration of the NDVI-related disadvantages, de Asis and Omasa (2007), de Asis et al. (2008) and Yang (2014) perform a comparative analysis of C-factors, derived from NDVI-and LSU-approaches, which result in better results for LSU. A relationship between C-factor and canopy cover fraction can be observed in various studies. Zhang et al. (2003) and Gao et al. (2012) determine an exponential decrease of the C-factor with an increase in canopy cover of grasslands. Wischmeier and Smith (1978) also observed a negatively exponential relationship of decreasing C-factors with increasing coverage in their empirical experiments on the USLE plots.
The (R)USLE factors C and R are highly dynamic with a clear annual cycle (Wischmeier and Smith, 1978;Renard and Freimund, 1994;Vrieling, 2006;Vrieling et al., 2014;Möller et al., 2017) in contrast to the rather constant RUSLE-factors K and LS Alexandridis et al., 2015). The status of grasslands is diversified within a year owing to the natural growth cycle, periodical cutting of hay, or pasture farming (Wiegand et al., 2008). Despite, this spatio-temporal variability of the C-factor for grasslands, it is often parameterized without accounting for the spatial variability within a land cover unit (Ozcan et al., 2008;Bosco et al., 2009;Efthimiou, 2016;Mancino et al., 2016) nor for the temporal variations (Wang et al., 2002). Hawkins (1985) stated already that "the complications of time and spatial variations in site properties are usually not considered" by applying the USLE. Alexandridis et al. (2015) conclude that a dynamic approach focusing on C-factors for the four seasons or 12 months of a year might help to reduce errors in the annual soil loss compared to a single annual C-factor. Vrieling et al. (2008Vrieling et al. ( , 2014 follow a multi-temporal and spatial approach to assess the riskiest erosion periods of the year for Brazil and Africa. López-Vicente et al. (2008) capture erosive periods among a year for a study area in the mountains of the Central Spanish Pyrenees by a dynamic approach on a monthly scale. Such time-dependent assessments of soil loss are relevant to support policy makers and farmers to protect the soil more targeted like it was done by López-Vicente et al. (2008). Panagos et al. (2012Panagos et al. ( , 2016 and Panagos (2016, 2017) propose a monthly time-step to be appropriate for soil erosion modeling. The same resolution was already proposed by Wischmeier and Smith (1965). Grazhdani and Shumka (2007) modeled the soil erosion rate for Albania on a monthly scale. A combination of both spatially and temporally varying R-and C-factors lead to a more dynamic soil erosion risk assessment and simultaneously allows the identification of susceptible seasons and regions (Panagos et al., 2014a;Ballabio et al., 2017;Möller et al., 2017). As it is shown in Meusburger et al. (2012), Schmidt et al. (2016) and Ballabio et al. (2017), the Rfactor of Switzerland also has a high intra-annual variability with clear regional patterns.
So far, most of the existing national C-Factor maps either do not include grassland areas (Friedli, 2006;Alexandridis et al., 2015), do not consider the temporal variations of vegetation cover and management (Friedli, 2006;Bosco et al., 2009;Panagos et al., 2015b), nor taking rainfall erosivity for C-factor calculation into account. An assessment following the original approach by Wischmeier and Smith (1978) to derive C-factor maps with a high spatio-temporal resolution based on SLRs and spatio-temporal R-factor ratios does not yet exist on a national scale. We aim to (i) determine the fractional vegetation cover with a linear spectral unmixing of orthophotos (2 m spatial resolution), and (ii) quantify the temporal change of vegetation fraction (10 days temporal resolution) to (iii) assess the spatial and temporal patterns of the C-factor based on SLRs and high-spatio-temporal R-factor ratios.

Swiss grassland areas
Switzerland is a country with a high heterogeneity of climatic, topographic and edaphic conditions. Hills and mountains cover more than one-third of the state. The Swiss elevation ranges can be clustered in elevation zones (in m a.s.l. modified after Ellenberg et al., 2010: Colline zone < 800; Montane > 800-1800; Subalpine > 1800-2300; Alpine > 2300-2700; Subnival > 2700-3100; Nival > 3100), which are typical for the plant development in the Swiss Alps. Owed to these natural conditions, permanent grassland is the predominant land use in about 28% of the territory of Switzerland with a share of 72% of the total agricultural area (Bötsch, 2004;Jeangros and Thomet, 2004;Schmidt et al., submitted). Grassland is the prevailing land use type at elevations above 1500 m a.s.l. (Hotz and Weibel, 2005). Almost half (46%) of the grassland area is therefore designated as alpine grassland (Hotz and Weibel, 2005). Alpine soils have been managed by humans for about 500 years already, but an intensification of the usage and management of grasslands can be observed since the last 50 years (Jeangros and Thomet, 2004;Bätzing, 2015;Alewell et al., 2008). Changes in grassland cover are expected due to land use and climate change.

Datasets for C-factor mapping
We subdivided the datasets of the C-factor mapping approach into data for the spatial and for the temporal assessment. In the spatial modeling approach, we used a high spatial resolution false-color infrared orthophoto (0.25 to 0.50 m; G R NIR) mosaicked of a set of 3432 tiles. This orthophoto mosaic called Swissimage FCIR (Swisstopo, 2010) is recorded with a Leica ADS80 airborne digital sensor, containing the channels green (533-587 nm), red (604-664 nm) and near-infrared (833-920 nm). The production process of Swissimage FCIR is based on an along-track scanning from east to west that generates stripes of aerial photos during each flight. The scheduling of the flights of the used version of Swissimage FCIR was in the years 2012, 2013, 2014 and 2015 between the months March and September. In the preprocessing step, the aerial photos have undergone a georeferencing, orthorectification, mosaicking, and clipping to tiles of 4375 m × 3000 m by Swisstopo. We reduced the file size (original file size 1.17 Gigabytes per tile) and the spatial resolution by resampling to 2 m for a more straightforward data handling.
The temporal variations of grassland cover in Switzerland are derived from time series of 10-day fractions of the green vegetation cover (FCover300m, spatial resolution 300 m; Smets et al., 2017) as a product from PROBA-V. The FGVC is expressed in percentages from 0% (no vegetation cover) to 100% (full vegetation cover). PROBA-V is a satellite with an assembled vegetation (V) instrument to image the global land surface vegetation regularly (Blair, 2013 (Leilei et al., 2014). This information is relevant for normalizing different recording periods of the Swissimage to the date of the peak growing period. A data accuracy modification was applied for MOD13Q1. Not processed or filled data, marginal data, and cloudy grid cells were substituted either by the preceding or succeeding good data or snow/ice data. With this routine, unreliable pixels were adjusted by the temporally closest reliable values.
We used the Swiss National Grassland Map of the year 2015 (Schmidt et al., submitted) for clipping the previously mentioned datasets to the grassland extent. Further, the dynamics of the long-term snow occurrence in Switzerland (Fig. S1) are derived from the Climate Change Initiate (CCI) Land Cover provided by the European Space Agency (ESA) (Arino and Ramoino, 2017). Elevation zones are extracted from the Swiss digital elevation model (SwissAlti3D, Swisstopo, 2017a). An overview of all used datasets is provided in Table 1. Data processing was done in ENVI 5.3., ESRI ArcGIS 10.3.1., and GDAL 2.1.3.

Concept of C-factor mapping for Swiss grasslands
Firstly, we derived the spatial pattern of Fraction of Green Vegetation Cover (FGVC spatial ) by LSU from the high spatial resolution Swissimage FCIR (Section 2.3.1). Secondly, we used FCover300m to estimate the temporal changes in the FGVC (FGVC temporal ; Section 2.3.2). Both approaches, the high spatial and the high temporal one are combined (Chen et al., 2015;Zhang and Li, 2015) via a normalizing procedure to result in a set of monthly FGVC maps for Switzerland (Section 2.3.3). This procedure involves the normalization of the orthophoto mosaic Swissimage FCIR with the temporal variations in vegetation cover of FCover300m to a given base date. The normalized high spatial and temporal FGVC spatio-temp maps of Swiss grasslands were then converted to SLR maps. The relationship of SLR and the fraction of vegetation cover (FVC) is based on measured data in alpine grasslands by Martin et al. (2010) andSchindler Wildhaber et al. (2012). SLRs were derived from the measured sediment yield for the given FVC classes proportional to an uncovered soil surface (SLR 100%; Schwertmann et al., 1987). SLR and FVC describe an exponential relationship (Eq. (1), Fig. 1). The SLRs are multiplied by the corresponding proportion of rainfall erosivity (R ratio ) to result in the C-factor according to the original approach by Wischmeier and Smith (1978) and Schwertmann et al. (1987). Monthly R ratio for Swiss grasslands with a spatial resolution of 100 m can be obtained from Schmidt et al. (2016). The processing workflow and manipulation of data is visualized in Fig. 2 (1) 2.3.1. Spatial modeling of fraction of green vegetation cover (FGVC spatial ) by linear spectral unmixing Spectral unmixing assumes that the spectrum measured by a sensor and represented as a mixed pixel is a combination of the spectra of components within the instantaneous field of view. As such, the reflectance of a mixed pixel is a mixture of distinct spectra (Roberts et al., 1993;Gilabert, 2000;Heidari Mozaffar et al., 2008). In spectral unmixing techniques, the mixed pixel is decomposed into a collection of endmembers and a set of fractional abundances according to the endmembers (Keshava and Mustard, 2002). The image endmembers, also called pure pixels, are at the vertices of the image simplex in an ndimensional space (Smith et al., 1985). Pixels defined as endmembers are relatively unmixed with other endmember signals (Rogge et al., 2007). Among the spectral mixture methods, the LSU is by far the most common type (de Asis and Omasa, 2007). LSU assumes that the incoming radiation only interacts with a single component of surface and is represented in a mixed pixel without multiple scattering between different components (van der Meer, 2010). Although this is a crucial assumption, the effects of intimate association between the components have been found to be relatively minor (Kerdiles and Grondona, 1995). LSU is expressed as the spectral reflectance (R i ) of the mixed pixel in band i as followed (Smith et al., 1990;Hill et al., 1995;de Asis et al., 2008):  et al., 2010a). However, QuickBird data is too cost intensive and heterogeneous for a national assessment and therefore rather applicable for catchment studies like it was done by Meusburger et al. (2010aMeusburger et al. ( , 2010b. Guerschman et al. (2009) use the hyperspectral EO-1 Hyperion in combination with MODIS data to result in a higher variety of endmembers with a spatial resolution of 1000 m. However, that spatial resolution of fractional cover is relatively coarse to explain the spatial patterns of the FGVC, SLRs and C-factors.

SLR e FVC FGVC and
In the present study, orthophotos (Swissimage FCIR) with a national coverage and resampled resolution of 2 m (resampled from 0.25 m to 0.5 m) were used. The spatial assessment for deriving FGVC spatial (see Fig. 2) is based on all three bands of the Swissimage FCIR. ENVI 5.2 provides a Pixel Purity Index tool (PPI) to automatically identify the most spectrally pure pixels of the image, designated to be the mixing endmembers (Pal et al., 2011;RSI Research Systems, 2004). PPI works with an iterative process by counting the number of times a pixel is registered as extreme pixel for each run. Pixels that appear to be extreme most often are then endmembers (González et al., 2010). We performed 10.000 iterations with a threshold value of 2.5 and identified a maximum of 100.000 pure pixels. The application of LSU can result in n + 1 endmembers where n is the number of bands (Phillips et al., 2005). PPI based on the three bands (G, R, and NIR) of Swissimage FCIR and determined the following endmembers namely i) vegetation, ii) bedrock, bare soil, asphalt, and iii) shade. These endmembers are the typical groups of endmembers which are distributed all over the grassland areas in the country (Roberts et al., 1993;Adams et al., 1995;Theseira et al., 2003;Meusburger et al., 2010a). Although the spectrum of water is relatively pure, water was not selected as an endmember since it is occurring only locally (Adams et al., 1995).
Swissimage FCIR has undergone a Minimum Noise Fraction (MNF) rotation before the selection of purest pixel and unmixing (Green et al., 1988). The MNF rotation is a two-step principle component analysis and used to determine the inherent dimensionality of the image data, to improve the signal-to-noise ratio and reduce the processing time (Boardman and Kruse, 1994;Nascimento and Dias, 2005). MNF can improve the quality of the resulting abundance maps by a decorrelation of the bands (van der Meer and de Jong, 2000). Furthermore, since the spectra are neither purposed to be linked to laboratory and field reflectance spectra nor to be meant for temporal approaches, a transformation of encoded-radiances in digital numbers (DN) was not required in this study (Adams et al., 1995;van der Meer, 2002).
A well-known problem of FGVC mapping is its underestimation due to the presence of dry vegetation (Meusburger et al., 2010a(Meusburger et al., , 2010b. This problem can either be addressed by long-wave spectral bands in hyperspectral sensors at the expense of spatial resolution (Guerschman et al., 2009) or by a calibration of the approach. As we aim to explain the spatio-temporal dynamics in soil erosion for Switzerland, we decided to preserve the high spatial resolution of our dataset (Swissimage FCIR) and followed the second option by using 1000 calibration points (FGVC cal ) to calibrate the FGVC spatial (based on the LSU) and to identify potential biases in the automated assignment of vegetation abundances. These points are randomly set for grassland areas. The FGVC cal is estimated user-driven for each point based on the 0.25 m resolved Swissimage FCIR and RGB. Besides that, the types of vegetation (photosynthetic and non-photosynthetic grassland, clipped grass, forest) or non-vegetation (shade, asphalt), slope degree and exposition are recorded. Although the calibration procedure assesses dry vegetation, it is not to be differentiated from bare soil in the LSU approach. Thus, the endmember of bare soil includes e.g. non-photosynthetic grassland. Thereby, the unmixed vegetation cover can be calibrated by the biases of dry vegetation. The density of optimization points is  Table 1) to result in spatiotemporal C-factors of Swiss grasslands. 37 km 2 , corresponding to one optimization point for each 6 to 6 km on average. An acquisition of ground truth data with a representative distribution in the field is hardly feasible on a national scale.

Temporal mapping of fraction of green vegetation cover (FGVC temporal and FGVC deviation )
Temporal variations of the fraction of green vegetation cover (FGVC temporal ) are provided within the FCover300m dataset. We averaged three files of the same date by the years 2014 to 2016 to a shortterm mean fraction of green vegetation (FGVC temporal ; see Fig. 2; Smets et al., 2017). Each of the three years of FCover300m is represented by a set of 36 files (108 files in total) in a 10-day resolution from 10th of January to 31st of December. The deviation of FGVC temporal to a base date is determined on a per pixel scale (FGVC deviation ) to be used for normalizing the FGVC spatial in the following Section 2.3.3. The processing of the FCover300m data is done within the Copernicus program where FCover is derived from the leaf area index and further canopy structural variables (Smets et al., 2017). Concerning its computation, FCover300m is more robust than classical vegetation indexes like NDVI which has stronger dependencies on geometry and illumination of surface cover (Weiss et al., 2000;Fontana et al., 2008).
A series of 253 NDVI datasets from 2005 to 2015 of the MOD13Q1 (Didan et al., 2015) were used for determining this respective base date as mean peak growing season indicated by the maximum NDVI within a year (Leilei et al., 2014). Fontana et al. (2008) demonstrate that the relationship between plant growth records in alpine grasslands and NDVI is quite remarkable. Busetto et al. (2010) use a time series from 2005 to 2007 of MOD13Q1 to determine the start and the end of the growing season of larches in the alpine region. For more robust results we averaged all ten years by each specific recording date to derive a mean NDVI per recording date for Switzerland. A correction of snow cover like it was done by Busetto et al. (2010) was neglected in the study as we are only focusing on the assessment of the peak growing season and not on minimum NDVI. The maximum NDVI of all the averaged datasets was selected for each cell and the corresponding DOY assigned to the associated cell. If a cell contained a no data value, it was skipped and the averaging done over the cells of the remaining year(s).

Merging of spatial and temporal fraction of green vegetation cover (FGVC spatio-temp )
As Swissimage is a mosaic of tiles recorded at heterogeneous dates, the vegetation cover can be assumed to be different between tiles according to its recording date. We used a normalizing process to make all tiles comparable. Therefore, the FGVC spatial are normalized to a base date. The spatial results, as well as the temporal results, are meant for being combined to spatio-temporal FGVC spatio-temp of grasslands (see Fig. 2). First of all, we extracted the recording dates of each along-track scanning stripe, and spatial joined the dates with the 3432 image tiles. In cases of multiple recording dates, we used the mode to extract the most common date. Tiles with same recording dates were aggregated to a multiple tile mask (Fig. S2) and later used to clip the FGVC spatial according to their recording dates.
Each FGVC spatial tileset of a specific DOY i can be normalized to that base date by weighting it with the the relative change of the FGVC deviation to the same base date as expressed in Eq. (4): Thus, tiles recorded early in the season where the plant growth can be assumed to be low are weighed by a greater FGVC deviation factor compared to an image tile recorded close to the base date.
All FGVC norm are merged to a new raster which represents a national map of FGVC at the defined base date. The normalized composite raster of the base date can then be recalculated to other dates.

Spatio-temporal mapping of grassland C-factors by considering soil loss ratios (SLRs) and rainfall erosivity (R-factor)
Originating from the FGVC spatio-temp , the SLR can be calculated with the relationship proposed in Eq. (1). SLRs express the ratio of soil loss of an area with a certain plant development relative to an uncovered surface (Renard et al., 1997). The SLRs are weighted with the ratio of the total annual rainfall erosivity (R ratio ) of the same period to result in the C-factor. The R ratio can be derived from monthly R-factor maps which exist with a high spatial resolution (100 m) for Switzerland (Schmidt et al., 2016). Monthly rainfall erosivity maps (100 m spatial resolution) for Switzerland are generated by regression-kriging of 10min rainfall records at 87 automated gauging stations (19.5 yrs. measuring sequences) and with the use of up to five spatial covariates. The 12 maps have a mean R 2 of 0.51 and a mean RMSE of 93.27 MJ mm ha −1 h −1 month −1 with highest uncertainties in winter due to generally low rainfall erosivity. The authors have discussed the variability of monthly R-factors for Switzerland in detail. R ratio can be assessed by calculating the monthly fraction of R-factor of the sum of all 12 maps. For the present purpose of Swiss grasslands, the monthly national maps of the R-factor are clipped to the extent of the improved Swiss National Grassland Map (Schmidt et al., submitted). The R ratio maps are multiplied with the SLR maps for grassland to result in monthly C-factor maps with a high spatial resolution. For each month we averaged the three corresponding FGVC spatio-temp maps to monthly FGVC maps to comply with the temporal resolution of the R-factor maps.

Spatial pattern of the fraction of green vegetation cover of Swiss grasslands
The optimized LSU of the Swissimage FCIR enables the differentiation of the FGVC spatial as well as the fractions for bare soil and bedrock. Spatial patterns of FGVC spatial are visualized on a national scale as well on a local level (Fig. 3). Such an analysis of the degree of fractional vegetation cover is of high relevance when categorizing land use for potential hot spots of erosion since it is more likely that an erosion process starts from the uncovered or bare soil.
The dimensionality of the Swissimage FCIR stays unchanged after noise segregation by MNF. The estimated ranges of FGVC spatial had 0.56% outliers outside the LSU constrained range of 0 to 1 (100%), which indicates that one or more of the endmembers chosen for the analysis is probably not well-characterized or that additional endmembers might be missing (RSI Research Systems, 2004). These outliers were omitted. They predominantly consisted of constructed environments (buildings, streets) that could not be masked in the grassland areas (Schmidt et al., submitted). The RMSE of the LSU for Switzerland is 22.6%. Higher uncertainties generally occur in the valleys of the Alpine foothill (Fig. 4). One reason for the high RMSE is the incorrect separation of grassland from arable land due to the coarse resolution (300 m) of the grassland map based on CCI Land Cover.
The mean FGVC of the 1000 calibration points (FGVC cal ; 61%) identifies a systematic underestimation of the mean FGVC spatial (39%) by 22% which is close to the mean RMSE. The highest discrepancy between FGVC cal and FGVC spatial mainly arises by an erroneous classification of non-photosynthetic vegetation (33% deviation), shades and artifacts (42% deviation), and forested areas (46% deviation). The segregation of non-photosynthetic vegetation and bare soil is impeded due to the very similar spectral characteristics. Shaded areas and artifacts disrupt the spectral signal of vegetation cover which is visually detectable but automatically assigned with a very low degree of coverage. The pattern of discrepancy between FGVC cal and FGVC spatial show a strong dependency to slope exposition. Highest deviations up to 34% are present at northern exposed slopes. All FGVC spatial were calibrated by adding the amount of mean underestimation to each grid cell. Subsequently, we used the calibrated FGVC spatial for all further calculations. The accuracy of the LSU approach could be increased with a more accurate grassland map and a higher number of spectral bands as it was already discussed in Meusburger et al. (2010a). A new orthophoto of Switzerland (Swissimage RS; Swisstopo, 2017b) with four spectral bands (NIR, R, G, B) is about to be released in 2020. Such an increase in bands could result in an additional endmember and might improve the LSU.

Temporal variation in the green vegetation cover of Swiss grasslands
The annual distribution of the mean FGVC temporal for Swiss grasslands visualizes the seasonal dynamic of grasslands with periods of dormancy and growing (Fig. 5). Higher FGVC temporal lasts until the end of October (approx. DOY 304) in lower elevations (Colline and Montane zone) of northern Switzerland. According to FCover300m, an FGVC temporal below 40% is present for most of the Swiss grasslands from December to February. The annual distribution of the FGVC temporal is comprehensive and complies with the typical expectable grassland plant growth cycle (Fontana et al., 2008;Filippa et al., 2015;Inoue et al., 2015). The lack of FCover300m data mainly covers the northern latitudes of Switzerland. According to the high solar altitude in summer, missing values are relatively rare during that season. Winter records are comprised of a higher number of no data values due to snow cover (especially in the Nival zone), sun path and cloudiness (Camacho, 2016). Thus, erosion in winter continues to be a blank spot, because we can neither observe changes in FGVC below the snow cover (which will affect the SLR and C-factor) nor assess the erosivity induced by snow movement and snowmelt (which will affect the R-factor) (Ceaglio et al., 2012;Meusburger et al., 2014;Stanchi et al., 2014). We excluded no data pixels (indicating snow) from the dataset if they are presented in all the three averaged years. The FCover300m still is in demonstration mode and has only undergone a validation over Europe yet (Camacho, 2016). Therefore, uncertainty could be introduced in the absolute fraction of green vegetation cover. Nevertheless, as all the 10-day data are assessed identically, the relative deviation of the values can be deemed correctly.
Based on the MOD13Q1 data, the long-term (2005 to 2015) maximum NDVI of the most considerable proportion of pixels is DOY 177 (26th of June, Fig. S3). We used the 30th of June (DOY 181) as the base date as this date has a high temporal proximity to the maximum NDVI of our analysis. This is in agreement with Jonas et al. (2008) who proposed the 6th of July as the mean date of the maximum height of grassland cover for elevations between 1560 and 2545 m a.s.l.. According to model results by Garonna et al. (2014), the growing season in the alpine zone starts at DOY 118 and lasts until DOY 266. The FGVC deviation in relation to DOY 181 marks a positive trend from DOY 181 to DOY 232 which determines the peak growing season for the national grassland area ( Table 2). The minimal FGVC in relation to DOY 181 is met on DOY 20 with a reduction of 58% in green vegetation cover.

Spatio-temporal patterns of the fraction of green vegetation cover of Swiss grasslands
The mean FGVC spatio-temp of Swiss grasslands on DOY 181 (30th of June; Fig. 6) is 60%. Grasslands next to the border of Austria (Cantons Appenzell and St. Gallen) have the lowest FGVC spatio-temp . These Cantons (see a map of Swiss cantons in Fig. S4) are fully dominated by meadows and alpine pastures (Table 3; Federal Statistical Office Switzerland, 2017). As the management of these grasslands is very intense (grazing, fodder), the FGVC spatio-temp is comparatively low. Intense grazing causes a significant limitation in grass growth (Bilotta et al., 2007;Mayer et al., 2009) which results in a degradation of vegetation cover (Yong-Zhong et al., 2005). These regions have one of the highest mean livestock unit (1.7 per ha; Table 3) and mean share of grazing livestock farming (78.8%). Hence, most of the areas in the region are already mowed at the 30th of June (typical mowing period for St. Gallen is DOY 166 to DOY 196;Zwingli, 2017). The whole Switzerland experienced a land use intensification of grassland over the last decades. It is apparent by an increase in stocking rates (~50% increase of sheep numbers during 40 years) and an alteration in grazing systems (permanent shepherding replaced by uncontrolled grazing, Troxler et al., 2004).

Spatial and temporal hot-spots of C-factors on Swiss grasslands
The monthly maps (Fig. S5) are averaged to seasonal maps of Cfactors for grasslands (Fig. 7). They represent the high temporal and spatial variability of the C-factors for grasslands throughout a year. According to the modeling results, relative high C-factors in winter can only be observed in the Jura mountain at the border to France and the western Alps. These patterns are mainly controlled by the ratio of the annual rainfall erosivity (R ratio ; Fig. 8). The whole alpine range experiences increased values in spring. The distribution of C-factors in summer for Swiss grasslands is relative diffuse with a spatial cluster in the north-eastern region of Switzerland (Cantons Appenzell and St. Gallen) which is a result by the low FGVC due to intense grassland land use (see Section 3.3) and the high rainfall erosivity. Absolute C-factors are decreasing in fall but with regional pattern of high C-factors at the southern and eastern Alps. The minimum C-factors within a year are covering the lowland areas of Switzerland in winter. Maximum C-factors are observable in the previously mentioned region of the Cantons Appenzell and St. Gallen (close to the border of Austria) in summer.
The mean annual C-factor for Switzerland is 0.012 (Table 4). Lowest mean C-factors of Swiss grasslands can be observed in January (0.003), highest in the summer months July (0.024) and August (0.025) (Fig. 9). The maximum C-factor in August is about 8 times higher than the minimum C-factor in January. The trend marks an abrupt increase of Cfactors from April to August with a decrease in its low winter values. The natural plant growth cycle determines the annual trend of FGVC. As the C-factor is not solely related to FGVC but further a product of SLR and weighted R-factor ratios, the trend of the C-factor is influenced by the regional and temporal rainfall erosivity pattern.
The rainfall erosivity, as well as the FGVC, is controlled by elevation level (Fig. 10). The C-factors per month and elevation zone follow typical patterns. Highest C-factors can be observed in the Alpine zone. The Subalpine, Alpine and Subnival zone show more than one peak with highest C-factors. The Colline and Montane zone have only one maximum in August. The C-factors in all elevation zones are lowest in the winter months January and February. FGVC in winter is low due to the reduced plant growth. The here excluded presence of snow cover in winter results in a delay of increasing FGVC with elevation after meltout. The typical melt-out at elevations between 1560 and 2545 m a.s.l.  is recorded by Jonas et al. (2008) and Fontana et al. (2008) at DOY 147. Large areas of Switzerland show a snow occurrence in winter (Fig. S1). Protection of grassland soils by plant cover is relatively low in winter but simultaneously affected by only very low rainfall erosivity. However, the tremendous impact of snow gliding on exposed soil surfaces during winter might be a crucial impact . Although the fraction of vegetation cover is increasing in summer for all the grasslands, the weighting with the R ratio causes a high C-factor. As discussed in Schmidt et al. (2016), a significant fraction of the annual rainfall erosivity is within the time window between June and September. The predominantly glaciated Nival zone (> 3100 m a.s.l.) could not be considered due to a small proportion of grassland areas (0.6% of the zone). Cantons in the east of Switzerland (Fig. S6) have slightly higher C-  factors in the month May to December which is also related to the differences elevation level (mean elevation of eastern cantons 1122 m a.s.l., western cantons 865 m a.s.l.) and different ratios of Rfactors. The elevation patterns become also visible by comparing the northern and southern cantons (mean elevation 928 m a.s.l. and 1795 m a.s.l., respectively). The capturing of the relationship between C-factor and elevation zone meets our expectations and confirms the plausibility of the input parameters and modeling approach. Bosco et al. (2009) already observed a relationship of C-factors and elevation level based on literature values. Kulikov et al. (2016) studied the temporal variations of C-factors of Kyrgyz mountain grasslands. They observed decreasing C-factors from April (immediately after snowmelt) to June in both of their study areas. They assess the months April and May with the highest potential soil loss owing to high C-factors with simultaneous high rainfall erosivity. A soil erosion assessment for a watershed in Brazil (de Carvalho et al., 2014) reveals highest soil loss in the rainy season where rainfall erosivity is high and the C-factor low. Another combination of dynamic Rand C-factors, done by Panagos et al. (2014b) for Crete in Greece, assesses March as a month with high rainfall erosivity and low fractional vegetation cover. Thus, it is important for C-factor assessment to consider the relative timing of peak C-and peak R-Factor. Panagos et al. (2015b) derived C-factors for grasslands for the 28 European Union member states from FCover300m and ranges of literature values. Their results present a mean European grassland Cfactor of 0.0435 which is about 3.5 times higher than the one for Switzerland. However, C-factors in Mediterranean regions, which are included in the mean European C-factor, are substantially higher than ones in Central Europe. The surrounding countries of Switzerland have mean national values between 0.0345 (Austria) to 0.0421 (Germany). Switzerland's nationwide C-factor for grasslands (0.012) is 70% lower than the mean of the four neighboring countries (0.0396). A different seasonal trend and lower values compared to Panagos et al. (2015b) and Kulikov et al. (2016) can be explained by the different methods to compute C-factors and the neglecting of the rainfall erosivity.
Extensive pasture systems might have a positive effect on a dense vegetation cover. Furthermore, rotation grazing systems or reduced stocking rates supports the development of a better-closed vegetation cover (Troxler et al., 2004). The exclusion (e.g., by fencing) of susceptible soils or spots with a reduced growth period due to a late meltout could effectively prevent soils from being mobilized. The regeneration time of a degraded sward will take many years, and as long as then the soil surface remains uncovered, it will be fragile and highly prone to an expansion of soil degradation in the form of erosion.
The study of the dynamic soil erosion is of high importance as growing seasons in the European Alps are about to be extended under futures changing climates and shortened snow-cover periods (Defila and Clot, 2001;Studer et al., 2005;Bänninger et al., 2006;Fontana et al., 2008;Frei et al., 2017). A long-term effect of the prolonged growing season for alpine plants would be the favoring of higher and faster-growing plants with enhanced biomass production. More biomass production increases the vegetation cover and lowers the C-factor in summer (Rammig et al., 2010). Simultaneously, the warmer climate and heavy precipitation events during fall and winter will result in higher R-factors (after snowmelt; Fuhrer et al., 2006;Rajczak et al., 2013;Rajczak and Schär, 2017). Sparsely covered soils in late fall (before snow cover) and early spring are then more susceptible to erosion by water. A significant increase and intensification in the coldseason precipitation is already observable for Switzerland (Widmann and Schär, 1997;Schmidli et al., 2002;Schmidli and Frei, 2005).

Conclusion and outlook
We derived Swiss C-factor maps of grasslands from soil loss ratios weighted with R-factor ratios in using the most state-of-the-art remote sensing products for Switzerland (e.g., national orthophoto with an original spatial resolution of 0.25 m (Swissimage FCIR) and a 10-day time series of fractional green vegetation cover (FGVC, FCover300m)). The assessment enables the nationwide quantification of the C-factor of grasslands in its dynamic throughout a year. C-factors are much higher in winter than in summer due to the relation to rainfall erosivity ratio and show the expected dependency on elevation gradient. The mean annual C-factor of Swiss grasslands is 0.012 which complies with the Cfactor of October. An improved spectral resolution will be available with the future Swissimage RS product which might increase the accuracy and quality of the linear spectral unmixing results. However, the present results can help to implement soil conservation strategies of an adopted land use management. The identification of regions in Switzerland and periods of the year with high C-factors in combination with the dynamic R-factors might help agronomists to introduce selective mitigation strategies for erosion control of Swiss grasslands. The mitigation potential of soil erosion particularly relies on the C-factor since the R-factor is climate driven and not directly to be altered by human interventions. The utilized grassland areas of Switzerland are of particular interest since grazing might degrade soil functions and stability and has an impact on soil cover. Grazing in alpine environments usually takes place during the most susceptible season. As sediment yield is reduced to a minimum under closed vegetation cover, priority should be on keeping the vegetation coverage of grassland high. The FGVC can be increased, and thus the C-factor lowered by avoiding grazing on highly susceptible grassland or at least by paying more attention to the choice of the grazing animal species and stocking numbers/ diversity. To capture the spread of degraded surfaces, the automated identification and classification of bare soil spots with a higher spectral resolution is envisaged for future studies. Beyond the current state of C-factors, the models can be linked to land use and climate scenarios to get an idea of future impacts of soil erosion. As we demonstrated the usefulness and applicability of the C-factor and its relation to the R-factor, this study also highlights the advantages of USLEtype modeling. Individual computation and assessment of every single factor result in a high transparency and verifiability of USLE-based erosion models. Each individual factor does not only have the advantage to be adjusted and evaluated on its own but also deliver valuable conclusions for other environmental issues. Fig. 8. Monthly ratio maps of the annual rainfall erosivity (R-factor) of Swiss grasslands. Monthly R-factor ratios are the fraction of R-factor related to the total annual R-factor sum. Rainfall erosivity maps of Switzerland are based on Schmidt et al. (2016). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.) Table 4 Mean C-factors of Swiss grasslands per month.

Acknowledgement
This work was supported by the Swiss Federal Office for the Environment (FOEN) (grant numbers N°N222-0350 and N°P182-1535). The authors would like to thank all data providers, namely Swisstopo, Copernicus Global Land Services, European Space Agency (ESA), and National Aeronautics and Space Administration (NASA) for making their data available. Furthermore, we would like to thank the anonymous reviewers for their valuable comments and suggestions to improve the quality of the paper.

Funding
This work was supported by the Swiss Federal Office for the Environment (FOEN) (grant numbers N°N222-0350 and N°P 182-1535).

Declaration of interest
None.