Mapping cropland-use intensity across Europe using MODIS NDVI time series

Global agricultural production will likely need to increase in the future due to population growth, changing diets, and the rising importance of bioenergy. Intensifying already existing cropland is often considered more sustainable than converting more natural areas. Unfortunately, our understanding of cropping patterns and intensity is weak, especially at broad geographic scales. We characterized and mapped cropping systems in Europe, a region containing diverse cropping systems, using four indicators: (a) cropping frequency (number of cropped years), (b) multi-cropping (number of harvests per year), (c) fallow cycles, and (d) crop duration ratio (actual time under crops) based on the MODIS Normalized Difference Vegetation Index (NDVI) time series from 2000 to 2012. Second, we used these cropping indicators and self-organizing maps to identify typical cropping systems. The resulting six clusters correspond well with other indicators of agricultural intensity (e.g., nitrogen input, yields) and reveal substantial differences in cropping intensity across Europe. Cropping intensity was highest in Germany, Poland, and the eastern European Black Earth regions, characterized by high cropping frequency, multi-cropping and a high crop duration ratio. Contrarily, we found lowest cropping intensity in eastern Europe outside the Black Earth region, characterized by longer fallow cycles. Our approach highlights how satellite image time series can help to characterize spatial patterns in cropping intensity—information that is rarely surveyed on the ground and commonly not included in agricultural statistics: our clustering approach also shows a way forward to reduce complexity when measuring multiple indicators. The four cropping indicators we used could become part of continental-scale agricultural monitoring in order to identify target regions for sustainable intensification, where trade-offs between intensification and the environmental should be explored.


Introduction
Agricultural expansion and intensification have led to marked increases in agricultural production since World War II (Rudel et al 2009, Tilman et al 2002, albeit at substantial environmental costs (Tscharntke et al 2012). Arguably, demand for agricultural products will have to increase in the future as the world's population grows, diets change, and bioenergy becomes more important (Beringer et al 2011, Krausmann et al 2013. Demand-side strategies such as reducing meat consumption and curbing food wastage (Bajzelj et al 2014, Erb et al 2009, Tilman et al 2011 are promising, but will likely not be enough. How production increases could be achieved while curbing the environmental impacts of agriculture remains unclear (Butchart et al 2010, West et al 2014. Expanding agriculture further into the last remaining undeveloped fertile lands in South America and Africa would entail drastic environmental costs (e.g., substantial carbon emissions and biodiversity loss) (Laurance et al 2014, Licker et al 2010, Ramankutty et al 2002. Crop production can also be increased by intensifying agriculture on existing cropland, yet underperforming cropland (Godfray et al 2010, Johnston et al 2011. Such intensification could, for example, entail an optimization of crop rotations (e.g., less fallow) and higher resource efficiency (e.g., nutrient, water use) (Mueller et al 2012, Ray and Foley 2013, Siebert et al 2010b. Since the environmental impacts of intensification can be substantial (Licker et al 2010, Matson et al 1997, careful, context-specific assessments of the risks and opportunities of intensification are required (Garnett et al 2013).
Information on spatial and temporal patterns of cropland use at multiple geographic scales is required to better understand the potential for intensification. Unfortunately, existing data on cropland-use intensity are mostly coarse in scale, heavily rely on uncertain cropland maps (Fritz et al 2011, Fritz et al 2013, or are based on national statistics, which themselves may contain uncertainties (Verburg et al 2011, Zaks and. Many existing datasets represent snapshots in time and cannot reflect the often highly dynamic management intensity of agricultural land , Siebert et al 2010b. Substantial progress in mapping indicators of croplanduse intensity has been made recently, including yield gaps (Johnston et al 2011, Monfreda et al 2008, West et al 2014, fertilizer use (Potter et al 2010), human appropriation of net primary production (HANPP) (Haberl et al 2012), field size (Fritz et al 2015) or the extent of irrigated agriculture (Siebert et al 2010a, Thenkabail et al 2009 or tillage (Johnson 2013). Yet, our understanding of cropping patterns (i.e., temporal dynamics of cropland use), such as crop rotations, multi-cropping (i.e., number of harvests per year), crop duration (i.e., fraction of the year in which the cropland is covered with crops), cropping frequency (i.e., the number of cropped years ), or the fallow land extent , Li et al 2014, Portmann et al 2010, Siebert et al 2010b is limited.
Dense time series of medium-resolution satellite images, such as from the Moderate Resolution Imaging Spectroradiometer (MODIS), can help identifying cropland dynamics across broad geographic extents , Ganguly et al 2010, Rogan and Chen 2004, Siebert et al 2010b (Li et al 2014), or the mid-western United States (Wardlow and Egbert 2008). These studies highlight the potential of medium-resolution sensors to map cropping indicators, but used only single indicators over relatively short time periods to describe often highly dynamic and heterogeneous agricultural systems.
Our main goal here was to characterize European cropping systems by mapping four MODIS-based cropping indicators and by identifying typical cropping clusters. Europe is interesting for assessing patterns of cropland-use intensity for at least three main reasons. First, Europe is characterized by a wide range of agricultural systems, spanning strong environmental and socio-economic gradients, resulting in systems with different degrees of cropland-use intensity (Herzog et al 2006, Rounsevell et al 2012, Rudel et al 2009. Second, Europe has a long land-use history, with most land-use change nowadays happening along gradients of intensity change. Third, Eastern Europe experienced a dramatic declined in croplanduse intensity after the breakdown of socialism (Kuemmerle et al 2013, Prishchepov et al 2012) and is consistently highlighted as a candidate region for sustainable intensification . Therefore Europe is a prime example to develop methods that allow capturing and mapping cropland-use intensity and changes therein.
Specifically, we assessed the following research questions: 1. What were the spatial patterns of cropping intensity in Europe from 2001 to 2012, as measured by cropping frequency, multi-cropping, fallow cycles, and crop duration ratio?
2. What are regions of similar cropping systems across Europe?

Data and methods
Our study area included the entire European continent and Turkey to define the cropland extent within this region, we used the GlobCorine land-cover map from 2005 (Defourny et al 2010). To overcome scale differences between GlobCorine (300 m) and the MODIS time series (231.6 m), we downscaled the prior to the latter using a nearest neighbour algorithm. Based on this map we focussed on all classes containing cropland, specifically rainfed and irrigated cropland, complex cropland, and mosaic cropland/natural vegetation, and we masked all other classes (figure S1 in the supplementary material).
To map cropping indicators, we used a NDVI time series pre-processed in a previous study following three steps (Estel et al 2015). First, using satellite images from both Terra and Aqua satellites from 2000 to 2012 we reduced effects from clouds, water, snow, and ice by excluding poor-quality observations based on the MODIS quality information, land surface temperature, the land-water mask and interpolating missing values. Second, due to the strong climate gradient across Europe and the resulting varying phenology (e.g., earlier green-up and shifted vegetation peak in the Mediterranean, higher seasonality in the north), we normalized the NDVI time series to make them more comparable. The normalization was twofold and included an accounting for the shifted vegetation peak in Mediterranean environments in Europe, and harmonization the NDVI time series in regards to vegetation maxima and amplitudes across Europe. The normalization procedure is described in detail in Estel et al (2015). Third, we gathered an extensive training dataset on active and fallow cropland by interpreting the NDVI time series and high-resolution images from GoogleEarth, and classified each cropland pixel into active (i.e., managed) and fallow farmland (i.e., unmanaged) for each year between 2001 and 2012. The annual active/fallow maps had an average overall accuracy of >90% based on independent validation data (Estel et al 2015).

Mapping cropping indicators
We used the pre-processed NDVI time series and the annual fallow/active maps from Estel et al (2015) to map four cropland-use intensity indicators: (1) cropping frequency, (2) multi-cropping, (3) fallow cycles, (4) crop duration ratio at a spatial resolution of 231.6 m. To provide an intensity measurement and to better describe and compare the mapped patterns of intensity we divided each indicator into a high, medium and low intensity class (table 1). We furthermore excluded all pixels labelled as abandoned or permanently fallow based on Estel et al (2015). Our final cropland mask included an area of 400 Mha. For the indicators crop duration ratio and multi-cropping, we considered only non-fallow years. We calculated the cropping frequency as the number of years a cropland pixel was cropped over the observation period (see supplementary material). Higher cropping frequencies thus signifies higher cropland-use intensity (table 1). Multi-cropping refers to the number of harvests within a single year (i.e., growing season, Spera et al 2014). In Europe, either single or doublecropping occurs. To identify double-cropping, we counted the number vegetation peaks per growing season using TIMESAT (Jönsson andEklundh 2004, Li et al 2014), which detects double peaks based on the amplitude ratio between the primary and the secondary peak (Jönsson and Eklundh 2004). We derived annual single vs. double-cropping maps and summarized the number of double-cropped years (table 1).
Fallow cycles refer to recurring periods of fallow cropland. Longer and frequent fallow periods thus signify less intense land management. We defined fallow as cropland without management (i.e., not sown, cropped, or ploughed) (Estel et al 2015). 'Active fallow' (e.g., cultivation of legumes for nitrogen-fixation) was not considered. To identify fallow cycles, we screened the active/fallow time series for 'chain segments', i.e., a certain number (1, 2 or 3) of consecutive fallow years. We identified chain segments consisting of one (FC1), two (FC2), and three (FC3) fallow years between active years. We counted the occurrence of these chains across the entire time-series per pixel (see figure S3), and summarized all chain segments using a weighting scheme. Weights were calculated as the ratio of the total number of years in the time series (12) and the number of maximally possible chain segments of a particular cycle type (see supplementary material). The resulting index provides information about the level of cyclicity, and thus management intensity. We considered only time series with at least two chain segments (figure 3).
The crop duration ratio is the time a field is cropped in relation to the total length of the growing season (Siebert et al 2010b). We derived the total length of the growing season as the number of days with a land surface temperature above 5°C, i.e., the time between the earliest and the latest MODIS acquisition date when plants are assumed to actively grow (Hickler et al 2012, Zhang et al 2004. The time period a pixel was cropped was defined as a vegetation signal of at least half the peak of the phenological curve (see supplementary material). The half-maximum is frequently used as a phenological marker for leaf unfolding and the loss of canopy structure of natural vegetation (Bradley et al 2007, Fisher et al 2006 and we used it as a proxy for crop green-up and harvesting. We then computed the crop duration ratio for each year and calculated the average crop duration ratio from 2001 to 2012. Since cropping cycles are more dynamic than natural vegetation (e.g., varying timing of ploughing, sowing and harvesting), we carried out a sensitivity analyses to test the robustness of our results in relation to the choice of threshold, and derived crop duration ratios for thresholds of 40%, 45%, 50% (=half-maximum), 55%, and 60%. We then calculated the standard deviation crop duration (see supplementary material).

Mapping typical clusters of cropland-use intensity
To identify similar cropping systems, we used selforganizing maps (SOMs). SOMs are an unsupervised clustering technique based on competitive learning that reduce a high-dimensional dataset to a two dimensional map by grouping observations according to their similarity (Skupin and Agarwal 2008). To identify the optimal number of clusters, we applied a sensitivity analyses with SOM clusters varying from 2×2 to 4×4 clusters (Maulik and Bandyopadhyay 2002) and used the Davies-Bouldin index that compares intra-and inter-cluster variability (Davies and Bouldin 1979) to pick the optimal cluster number. We z-transformed our indicators prior to the clustering, and calculated average values across a 1×1 km 2 grid (see supplementary material).
Multi-cropping was widespread in the study area. About 56% (table 1) of all croplands in Europe were double-cropped at least once during 2001-2012. Areas where multi-cropping was high (i.e., 9-12 doublecropped years) accounted for only 2% of the total croplands, and were most widespread in central Europe (i.e., north-eastern Germany, central Poland) and Russia (i.e., Black Earth regions). Medium multi-cropping (i.e., 5-8 double-cropped years) accounted for about 9% of all croplands, mainly in central Europe, central Spain, southern Ukraine, Russia (figure 2).
Our fallow cycle mapping showed that about 27% of all European croplands had one of the three fallow cycles (figure S3). Out of all croplands about 1% had a high fallow cyclicity (upper tercile of index values) and these areas occurred predominately in the southern Iberian Peninsula, north-eastern Turkey, and eastern Europe (table 1). About 26% had a medium (i.e., mid tercile) or low fallow cyclicity (i.e., lower tercile) and occurred all over Europe with concentrations on the Iberian Peninsula, eastern Europe (i.e., northern Ukraine, Russia) and some Mediterranean areas.
Mapping crop duration ratio (figure 4) revealed that about 20% of all European cropland was characterized by high crop duration ratios, mainly in central Europe (e.g., Germany, Poland, eastern France, and southern Hungary). About 69% of all croplands showed medium crop duration ratios (i.e., 0.44-0.70), occurring mainly in European Russia, Ukraine, UK and eastern France (table 1). A few regions (∼12%) showed lower crop duration ratios (i.e., <0.43), mainly in Spain, Turkey, Italy, Greece and southern    European Russia. Our sensitivity analyses showed the robustness towards alternative definitions, with standard deviations in crop duration ratio <20%. For some Mediterranean areas, where crop duration ratio is lower than elsewhere in Europe, we found a higher, but still moderate sensitivity (see supplementary material).
We identified six clusters of similar cropping systems (figure 5) using SOMs and a sensitivity analysis across varying cluster numbers (see supplementary material). To describe the magnitude and direction of the different cropping indicators in each cluster (C1-C6), we provided here the deviation (±) from mean z-score (=0, table S1). Positive and negative numbers thus signify above and below average values respectively, whereas values close to zero mean that a specific indicator is close to the overall mean of the study area (figure 6). Cluster 1 was determined by high cropping frequencies (+0.82%), crop duration ratios slightly above average (+0.22%), and a very low fallow cyclicity (−1.25%). This cluster occurred in northern France, England, Italy, and around the Black Sea in Romania, Ukraine, Russia, and Turkey. Cluster 2 was mainly determined by very high multi cropping (+2.76%), high crop duration ratios (+0.80%) and cropping frequencies (+0.72%), but a low fallow cyclicity (−0.67%). This cluster 2 was mainly located in Germany, central Poland, southern Russia, and western France. Cluster 3 was characterized by very low cropping frequencies (−1.78%), infrequent multicropping (−0.80%), a marked fallow cyclicity (+0.73) and occurred mainly in eastern Europe (e.g., European Russia, Baltics, Belarus, and northern Ukraine), the Mediterranean, and in mountain areas (e.g., Alps, Pyrenees, Caucasus). Cluster 4 had very low crop duration ratios (−1.62%) and a very low fallow cyclicity (−1.23%), and the highest cropping frequencies (+0.85%) of all cluster. This cluster occurred mainly in the Mediterranean, southern Ukraine and southern Russia. Cluster 5 was characterized by a high fallow cyclicity (+0.73%), low multi-cropping (−0.53%) and low crop duration ratios (−0.35%). This cluster occurred mainly in Ukraine and Russia, on the Iberian Peninsula, western France, and Turkey. Cluster 6 had the highest fallow cyclicity (+0.74%) and the highest crop duration ratios (+1.00%) of all cluster. This cluster occurred mainly in central Europe and southern Ukraine.

Discussion
Understanding spatial patterns in cropland-use intensity is important for identifying target regions for intensification, and for assessing its potential environmental trade-offs. Better information on cropping systems, including cropping frequency, fallow cycles, multi-cropping, and crop duration, are important indicators in this context; yet agricultural census do not cover them at fine scale. Using a 12-year MODIS NDVI time series, we mapped these cropping indicators at the continental scale for Europe. The cropping patterns we find correspond well with other indicators of agricultural land management in Europe, highlighting the potential for satellite-based cropping measures to support agricultural monitoring. Moreover, our indicators capture important aspects of agricultural intensity, highlighting how management intensity varies in space and time. The identified cropping systems can be explained by agro-environmental conditions (e.g., soil quality, water availability), socio-economic conditions (e.g., rural depopulation), management practices (e.g., irrigation), and cropspecific management (e.g., rice growing in northern Italy). Our satellite-based indicators of cropland-use intensity and the mapping of cropping systems may help to identify baselines and candidate region for intensifying croplands sustainably.
The cropping patterns captured by our indicators correspond well with those from other indicators of agricultural management in Europe. For example, existing maps of cropland net primary production, yields, and yield gaps consistently show the highest output intensity where we identified intensive used multi-cropping systems and a long crop duration (Monfreda et al 2008, Mueller et al 2012, Neumann et al 2010. Maps of fertilizer usage also show highest fertilizer use in these regions (Potter et al 2010, Temme and. In contrast, these maps suggest low input and output cropland intensity in eastern Europe, especially in European Russia, Ukraine, and the Baltic States, congruent with our cropping indicators. A more quantitative comparison among our maps and other land-use intensity indicators is not feasible, given different resolutions and time periods covered, and considering that cropland-use indicators related to inputs and outputs usually represent downscaling agricultural statistics (usually national scale). Although multi-cropping has increased due the biofuel expansion (in Europe mostly based on wheat) in the last decades it is not considered in agricultural statistics (Langeveld et al 2014).
That our satellite-based indicators identify, on a general level, the similar spatial patterns of high and low intensity is encouraging, and highlights the potential of satellite-based to observe and monitor agricultural systems more directly and with finer spatial detail (Kuemmerle et al 2013, Zaks and. The spatial congruence between our cropping indices and alternative measures of agricultural management, such as fertilizer and yields, also underlines the value of satellite-based cropping indices to more directly measure cropland-use intensity. Yet some of our indicators, for example cropping frequency, (i.e., more cropped years indicates higher intensity) are more closely related to intensity than others (e.g., crop duration ratio). A few points of caution need to be mentioned when interpreting our individual indicators. First, our multi-cropping indicator captured whether there are one or two vegetation peaks in a given year, but not in all cases a second peak refers to a second crop, nor does it attest to whether the second crop was actually harvested or left in the field as green manure. Second, our fallow cycle indicator captured cropping cycles over multiple years, but rests on a reliable identification of active vs. fallow cropland. Although our accuracy assessment of the active/fallow maps suggests these maps are reliable (>90% overall accuracy, Estel et al 2015), uncertainty for some regions (e.g., with small fields) may be higher than for others. Third, we measured crop duration ratio as the length of the vegetation signal on cropland, but lower crop duration ratios may be due to management or climate. To interpret crop duration ratio as a measure of cropland intensity, agro-environmental conditions and irrigation structure should be considered in future work. Finally, validation of our cropping indicators is challenging because retrospective field-level data on cropping do not exist for larger areas, and in situ data are currently not feasible to gather at the continental scale. However, we note that the annual fallow/active maps derived in Estel et al (2015) and used as input data for the cropping intensity indicators cropping frequency and fallow cycle has been validated extensively (Estel et al 2015), and sensitivity analyses for the crop duration indicator attest to the robustness of these indicator.
SOMs were a useful tool to identify regions with similar cropping and thus to help interpret and reduce complexity in our multi-dimension indicator dataset. The six clusters of cropping systems appear to be related to distinctly different agro-environmental and socio-economic conditions across Europe. Cluster 1 characterized moderately intensive rain-fed cropping and included some irrigated areas in southern Europe along rivers and reservoirs (e.g., Ebro-basin in northern Spain, Po Valley in northern Italy, Black Sea area) (Siebert et al 2006, Salmon et al 2015. Irrigation here led to an uncoupling from climate constraints, allowing for central/western European cropping systems in these areas.
Cluster 2 was clearly linked to the most intensified rain-fed cropping in highly favorable agro-environmental conditions (e.g. in Germany and Denmark, Neumann et al 2010). Interestingly, this system was also found in Europe's east, particularly, in southern Ukraine, Romania and southern Russia, where some of the world's most fertile soils are found (Fischer et al 2000). High double cropping rates that characterized this cluster (e.g., southern Russia, northern Germany) are linked to the cultivation of winter wheat as main crop. The specific cropping time of winter wheat allows preceding or subsequent crops (e.g., rape, summer wheat) (Gienapp et al 2012, Schierhorn et al 2014. Cluster 3, characterized by frequent fallow years, occurred mainly in water-limited regions (e.g., Extremadura and northern Andalusia in Spain, southern Portugal). Fallowing is characteristic for semi-arid regions to maintain soil moisture and fertility (Boellstorff and Benito 2005). Interestingly, this cropping system also occurred in eastern Europe, although agro-environmental conditions are more favorable there. The higher fallow frequencies there are possibly a legacy of the breakdown of socialism and the subsequent restructuring of agricultural sectors (Fischer et al 2000, EU 2005, leading to widespread farmland dis-intensification and abandonment (Rey Cluster 4 was clearly related to rain-fed cropping under water limitation (e.g., the Mediterranean), allowing only for shorter growing seasons (Fischer et al 2000). This does not necessarily indicate lower cropland-use intensity since yields in this regions can be high (Mueller et al 2012). For example, rice growing areas in the Po Valley (northern Italy) had very short crop durations but are considered among the most intensively used croplands in Europe (Blengini andBusto 2009, Erb et al 2013).
Finally, cluster 5 contained much cropland without major constraints (e.g., accessibility, soil fertility), yet that is currently not used to its full potential. Cluster 6, characterized by high fallow rates, occurred mainly in central and Eastern Europe in favorable agro-environmental conditions (e.g., southern Romania, southern Ukraine). Both cluster could entail particularly candidate regions for sustainable intensification.
In sum, our four satellite-based cropping indicators, as well as our cluster analyses to identify similar cropping systems across Europe, appear valuable for broad-scale agricultural monitoring. Satellite-based indicators could therefore complement ground-data that must always remain sample-based with wall-towall observations of agricultural management. Our study showed that time series of satellite images can also help to better characterize cropland-use intensity, and thus to assess baselines and potentials for intensifying croplands sustainably. Such information is critical to scrutinize the possible socio-economic and environmental trade-offs, as well as synergies, of intensifying croplands.
Commission (VOLANTE, No. 265104 and HER-CULES, No. 603447), and the DFG (LUCC-BIO, No. 3210 3109 99). We also would like to thank the ESA GlobCORINE Project for making the GlobCORINE map available.