The Kielstau catchment (Figure 1) is a small rural catchment in northern Germany covering an area of 50 km² (Fohrer et al., 2014; Wagner et al., 2018). It is a sub-basin of the Treene river. The mean annual precipitation and temperature are 893 mm and 8.3°C (DWD, 2009), respectively. The topography is comparatively even, with an elevation ranging from approx. 27 m to 78 m a.s.l. (Figure 1). Soils of Gleysol, Podzol, and Luvisol (Figure 1) dominating this catchment are mainly used as arable land, and grassland or pasture. Kielstau is the main river with a total length of approximately 17 km, flowing through Lake Winderatt about 5 km downstream from the river source. The Lake Winderatt is surrounded by protected areas that are mainly used for moderate grazing (Fohrer and Schmalz, 2012). The Kielstau re- ceives discharge from two main tributaries: the Moorau and the Hennebach. To secure agricultural productivity, subsurface tile drainages have been installed during land reallocation (Riedel and Polensky, 1987), which caused wide, drained areas that were estimated to cover 38% of the entire catchment (Fohrer et al., 2007).

Croplands constitute 63% of the study area and the majority is found in larger land patches as compared to the other land use classes. Among them, winter wheat, winter barley, winter rape, or corn take up a larger proportion, while winter rye, beans, oats, or row crops (sugar beet, potatoes, etc.) constitute a relatively smaller area. Crop rotations are commonly applied, resulting in constant changes within the croplands. Grassland or pasture (20%-21%) is mainly found in proximity to rivers and lakes or in protected areas. Urban, forest, water, and garden or orchard each occupy a smaller and persistent proportion.

Land use maps of the years 2013 and 2017 (Figure 1) are used for spatio-temporal pattern analysis. The maps were acquired through field surveys in spring and early summer of 2013 and 2017. The main base data for this mapping process was extracted from the automated landed property map (ALK-Automatisierte Liegenschaftskarte, 2004) released by the state survey office of Schleswig-Holstein (LVERMA-SH). This map provides the outlines of all land patches. The land use class of every patch was mapped in the field. Furthermore, the consistency of the derived land use data was checked with the help of 0.2 m-resolution orthophotos (LVA, 2013, 2016) that were taken on 01/03/2013 and 26/08/2016.

Croplands comprise 62.8% in 2013 and 63.1% in 2017 of the catchment (Table 1). They are mainly classified into (1) cereals, including winter wheat, winter barley, winter rye, oats, and a smaller fraction of summer wheat and summer barley, (2) energy crops, like winter rape and corn silage, (3) row crops including potatoes and sugar beet, as well as (4) field bean, strawberries, and vegetables. With regard to semi-natural land use, grassland, meadow, mowing meadow, and pasture account for 20.8% in 2013 and 20.3% in 2017 (Table 1), respectively, and are primarily found along the rivers. The forest areas are nearly stable over time. Fallow areas are croplands that are not cultivated in the current year. It is mainly located in the southwestern region. Water areas cover approximately 1.8% including lakes, ponds, rivers and open creeks. Residential sites are located near main road intersections, with a slight increase from 10.5% in 2013 to 10.6% in 2017. Garden or orchard plots are scattered mostly near villages. To avoid very small land use classes and samples, the land use was divided into the following 14 main classes: settlement areas (residential, commercial, and industrial lands), fallow, grassland or pasture (field grass, meadow and pasture), corn, other croplands (strawberry, potato, sugar beet, etc.), forest, winter rye, winter rape, garden/orchard, winter wheat, other cereals (oats, summer wheat, summer barley), field beans, winter barley, and water.

Twenty-six spatially distributed variables have been used to explain the location of land use classes in space. Variables are depicted in Table 2, including topography variables, soil properties, distance and socioeconomic variables, and landscape metrics. All datasets are processed to a 10 m grid resolution. Spatial patterns of slope, silt content, drained soil share, distance to protected areas, population density, and fractal dimension are given as examples in Figure 2. Topography, soil properties, distance and socioeconomic variables are widely used to explain spatial patterns of land use and land use change (Verburg et al., 2004; Qasim et al., 2011). Landscape metrics have been previously used to quantify landscape structure and to assess land use and land cover change.

A 5 m digital elevation model (DEM) derived from topographic map of Schleswig-Holstein (LVA, 1992-2004) has been used to generate elevation, slope, and aspect data at a 10-m resolution. The distance to roads or villages has been calculated from polygon roads and villages shapefiles that were extracted from land use maps in 2013 and 2017. The distance to rivers has been calculated from the main rivers shapefile of the Kielstau catchment (LANU, 2003). The distance to protected areas was derived based on the distribution of protected areas surrounding the Lake Winderatt and the Kielstau river (StiftungNaturschutz, 2016). By combining the land patch map with tile drained areas estimated (Fohrer et al., 2007), the percentage of drained area per patch was calculated. With regard to population density, (i) every 50 m² settlement area is reclassified as one residential site based on a village map that was extracted from the land use map and from the digital basis landscape model (ATKIS-Basis-DLM) (LVA, 2016); (ii) residential sites are interpolated into a residence density raster using a Kernel algorithm (Silverman, 1981; Yan et al., 2011); (iii) the mean value of residents per site in each town is derived by dividing the number of inhabitants per town in June 2013 and December 2015 (Statistik Nord, 2013, 2015) by the amount of residential sites; (iv) the spatial population density is the product of the residence density raster with the mean value of residents per site. Soil properties are acquired from a combined application of soil type distribution derived from a Digital Soil Map (BGR, 1999) and soil attributes generated for a modeling study (Fohrer et al., 2014). As German biogas plants take in cereals, weeds, corn, and sunflowers as feedstocks (Golon, 2009), their distribution potentially affects crop distribution. Distance to biogas plants has therefore been considered and calculated according to the locations of two biogas plants within the catchment and another one near the catchment. Explicit outlines of land patches were extracted from land use maps in 2013 and 2017. Landscape metrics of all patches have been calculated with the Patch Analyst 3.1 extension for ArcGIS 10.3 and all other calculations have been carried out in ArcGIS 10.3.

Logistic regression models are used to analyze spatial patterns of all land use classes in the catchment with 26 explanatory variables (Table 2). Water areas are excluded from the analysis study as they occupy only a small proportion of the area. For each class C, a binary coding is employed, where 1 indicates the presence of class C; 0 indicates the presence of another class. The probability P_i for each pixel i (10m×10m) for the appearance of this land use class is calculated using a function of explanatory variables X_n,i as follows:

where β_n is the regression coefficient for the variable X_n,i. To seek the most important explanatory variables and avoid over-fitting, the number of explanatory variables n is limited to five. Previous studies also showed that five explanatory variables are sufficient to acquire reasonable results (Baumann et al., 2011; Wagner and Waske, 2016).

The relative operating characteristic (ROC) is used to assess model performance. The ROC statistic is the area under the curve of the rate of true positives versus the rate of false positives for a range of threshold values applied to the probabilities to achieve a binary classification. It ranges from 0.5 (random separation) to 1 (perfect discrimination) (Pearce and Ferrier, 2000). An ROC≥0.7 indicates that the independent variables have a strong capacity to model the dependent variable (Pontius and Schneider, 2001). In this study, the ROC value is used to select the model that is more suitable to explain land use patterns. To ensure that the selected variables are not collinear, Pearson’s correlation (r) is calculated for each variable pair. When r exceeds 0.7, the variable that better explains the land use appearance is retained and the other one removed (Baumann et al., 2011; Wagner and Waske, 2016). The logistic regression models of all possible combinations of five explanatory variables from non-collinear variable datasets are set up, by removing non-significant variables that are unable to optimize model performance according to a stepwise approach. Among these combinations, the 13 best logistic regression models with the highest ROC, one for each land use class, are selected and used for the analysis. To avoid spatial autocorrelation, for each specific land use class a stratified random sample is extracted by taking 20% of all pixels from its binary land use raster in 2017, and by excluding adjacent pixels in this sample (Wagner and Waske, 2016). The derived sample is randomly divided into two equal parts: one for calibration and one for validation of the models.

To assess the relative importance of the explanatory variables, all possible models derived from the calibration process by the ROC values are sorted and the first 50 best models for each class are selected. Then, how often each variable is included into these 50 models is counted. The percentage of the inclusion of a variable into best regression models is regarded as variable importance. All calculations, model executions, and analyses are performed in R and with the help of the R packages ROCR (Sing et al., 2005) and raster (Hijmans et al., 2016).

The derived logistic regression models are quantitatively and visually evaluated in space and time (Pontius Jr et al., 2004). First, the validation half of the stratified random samples from the 2017 binary raster are used for spatial validation using the ROC statistic. Second, the model applicability is tested in time. To this end, some of the explanatory variables are updated: landscape metrics in 2013 are updated due to reshaped land parcels. Distance to villages is updated using slightly changed settlement areas in the 2013 land use map. Population density is recalculated by using village residence map and population data for the year 2013. The probabilities for each specific land use class in the entire catchment in 2013 are calculated, by applying the best logistic regression models to the 2013 explanatory variables. The derived probabilities are compared against the observed 2013 land use map with the help of the ROC statistic.

To further assess the plausibility of modeled results, R-G-B composites are used to visualize the competition for land in space between three different land use classes. An R-G-B composite is built by overlaying any three probability maps. In R-G-B composite maps, dark areas indicate low probabilities for these three classes, whereas light colors (near white) correspond to high probabilities for all of them, i.e. strong competition for land. Water areas are masked in white.

All explanatory variables that are included in the best logistic regression models are significant (p<0.05, Table 3). The best models are well above random discrimination as ROCs are greater than 0.5 (Table 3). The ROC for calibration (ROC_cal) ranges from 0.68 for corn to 0.97 for settlement areas. With regard to the spatial validation ROC values (ROC_val_2017) between 0.73 and 0.97 (Table 3) except for a slightly lower ROC of 0.68 for corn, are within the range of reliable precision 0.7-0.9 (Wu et al., 2009) or nearly 100% accurate with ROC value close to 1 (Pontius Jr and Schneider, 2001), highlighting reasonable performances of the derived models. Slightly lower ROC values from temporal validation (ROC_val_2013) indicate that the prediction of the 2013 land use pattern is a little worse (Table 3). Nevertheless, ROC_val_2013 values are greater than 0.71 with the exception of corn. The ROC for calibration differs from that for spatial validation and temporal validation by a maximum of 0.03 and 0.08, respectively. Hence, the derived regression models are reliable and robust both in space and in time. They can be used to accurately predict the spatial distribution of land use classes. Only corn is predicted with a slightly lower accuracy. Spatial patterns of croplands are often less accurately explained than non-croplands. For instance, ROCs for corn, winter wheat, and winter barley are equal to or smaller than 0.75, whereas ROCs for settlement areas, orchard/garden, fallow, and forest are equal to or larger than 0.84.

The majority of explanatory variables are included at least once in the best logistic regression models (Table 3). However, only few soil properties are used due to their higher collinearity, e.g. silt content has a positive correlation with soil organic carbon (r= 0.83) and clay content is negatively correlated with moist bulk density (r= -0.75). As only the variable that better explains land use appearance is included, soil organic carbon, available water capacity, and some other properties were not included in the best models (section 2.4).

The odds ratio, depicted as the value of exponent (β_n) for the n-th variable, represent the impact of this variable on the predictor. The probability for land use presence will increase upon an increase in the n-th variable with odds ratio being greater than 1, whereas probability decreases with an increase in the variable when odds ratio is below 1. Table 3 shows that topography variables tend to affect natural or semi-natural land uses. Specifically, with an increase in slope by 1 degree the probability for fallow, forest, or orchard/garden appear to increase by 35.7%, 34.7%, and 50.3% (odds ratios=1.35698, 1.34681, 1.50291). When terrain ascends by 1 m, the probabilities for fallow and field beans decrease by 17.4% and 15.2% (odds ratios= 0.82567, and 0.84845, respectively).

Soil content variables (clay, silt, sand, or rock content) are mainly used to explain the occurrence of agricultural lands, e.g., with an increase of silt by 1%, winter wheat and winter barley become more probable by 8.9%, and 10.6%, respectively, and winter rye becomes less probable by 2.6%. Grassland or pasture is more probable on areas with a deeper first soil layer (0.2% increase per mm). The drained soil share affects seven land use classes. With an increase in the proportion of drained area by 1%, the probabilities for settlement, garden or orchard, and forest presence decrease by 3.4%, 2.4%, and 2.6%, respectively, while the probability for pasture increases 1.2%. Other crops (3.3%) and winter rye (1.3%) are less probable if the drained area share increases, whereas field beans are more probable (2.0%). The contrary effect on field beans is in agreement with the higher probability for field beans at low elevations, indicating that field beans are primarily found at lowland areas with a high drainage percentage in the catchment.

Distance variables and population density explain all land use classes. Particularly, distance to protected areas and population density are included in the best models for eight and five land use classes, respectively. As further away from protected areas, grassland or pasture becomes less probable (0.03% per meter, Table 3). In contrast, settlement areas and orchard or garden seem more probable (0.03%, 0.04% per meter, respectively) with increasing distance from protected area. There is no clear impact direction on croplands (Table 3). Population density has a clear effect on land use patterns; with an increase of one more person per square kilometer, settlement areas (0.4%) are more likely to be found, whereas fallow (0.9%), forest (0.8%), corn (0.3%), and other crops (0.3%) are less probable.

The influence of patch structure in particular of the fractal dimension is evidently different for croplands and non-croplands. An increase of fractal dimension indicating more complex patch shape (nearer 1 - simpler shape; 2 - irregular and complicated shape (Forman, 2014; Forman and Godron, 1981)) suggests that settlement area and forest (odds ratios > 1, Table 3) are more likely to be found on irregular patches, whereas croplands corn, winter rye, winter rape, winter wheat, other cereals, field beans, and winter barley all with odds ratios < 1 (Table 3) are more likely to be on simpler patches (e.g. large rectangles).

In addition to the best logistic regression model, the 50 best logistic regression models are derived for each land use class to evaluate the importance of explanatory variables. All models achieve a reliable performance with ROCs ranging from 0.72 to 0.97 except for corn (Table 4). Moreover, ROCs of the 50 best models only differ by a maximum of 0.03, indicating similar performance. The variables that are also included in the best model mostly have the greatest importance (percentages are marked in bold in Table 4). Some are even included in all of the 50 best models (100% inclusion). Overall, distance variables and population density are particularly important as they account for 31.7% of the variables that are used in these models, followed by soil properties (29.0%) and landscape metrics (28.7%), whereas topography (10.7%) is less important for rural land use patterns in this lowland region.

Distance to roads, distance to villages, distance to biogas plants, population density, and drained soil share are important for all land use classes, as they are at least included in two of the 50 best models (i.e. ≥4% inclusion, Table 4). Elevation, distance to protected areas, distance to villages, population density, drained soil share, perimeter-area ratio, and fractal dimension are of utmost importance as they are included at least one time in each 50 best models. Among them, patch fractal dimension (Pfd) is quite important for settlement areas (100% inclusion), mixed croplands (100% inclusion), field beans (96% inclusion), winter rye (88% inclusion), and forest (72% inclusion); distance to protected areas affects other croplands, winter rye, winter rape, and other cereals (inclusion ≥98%); drained soil share greatly explains the appearance of settlement areas, other croplands, forest, grassland, and field beans (inclusion ≥90%). Further important influences indicated by 100% inclusion are: elevation for fallow and field beans, distance to villages for orchard/garden, population density for forest, and patch perimeter-area ratio for other croplands.

In general, the variables distance to protected areas (mean inclusion 48%), population density (32%), drained soil share (51%), and patch shape complexity (fractal dimension, 66%) are largely important as they are more frequently included and used for 10 to 13 land use classes. Topography and soil properties contribute to fewer land use classes with a relatively low importance. The lower importance of topography is largely related to smaller terrain variations and similar potentials for crops throughout the catchment. The low impact of soil properties mainly results from the exclusion due to collinearity between soil variables, e.g. organic carbon content and moist bulk density are highly correlated with available water capacity (r= 0.84 and r= -0.85, respectively) and are, therefore, excluded as available water capacity is preferred during variables selection (Table 4). Nevertheless, soil properties indeed provide crucial information for the spatial distribution of cropland, e.g. silt content is an important explanatory variable for the spatial patterns of major cereals (winter rye 34%, winter wheat 98%, and winter barley 58%, respectively, Table 4) and part of the best model for these classes.

The probability maps for each land use class in 2013 are calculated using the best logistic regression model and the respective explanatory variables. For most of the land use classes, the spatial patterns of one or two defining explanatory variables show up clearly in their probability maps. For settlement areas, patterns of population density and drained area share are visible in the probability map: higher probability of settlement in areas of higher population density and lower percentages of drained area (Figures 2 and 3). Soil depth is part of the best model for grassland or pasture, and its outline is prevalent along the Kielstau river and to the east of Winderatt lake in the probability map for grassland. As a component of the best models for fallow and forest, the slope pattern is particularly visible in the southwestern part of the probability maps. For winter rye and winter barley, silt content is a determinative variable and results in clearly silt-analogous pattern in the probability maps. The probability map for orchard/garden exhibits a pattern similar to the variable distance to villages with higher likelihood nearer to villages. Other cereals are clustered near biogas plants, which is generally consistent with the surveyed distribution of other cereals in 2013 (Figure 1). The probability for field beans depends to a great degree on drained soil share. Field borders are explicitly visible in probability maps for all croplands, underlining the great importance of incorporating patch parameters into the cropland models.

R-G-B composite maps have been produced by overlaying probability maps to visualize competition between land use classes. From the composite map of settlement (Red, R), grassland/pasture (Green, G), and forest (Blue, B) as shown in Figure 4a, a majority of areas are found in these colors without mixing, i.e. settlement areas, grassland/pasture, and forest are well separated and do not compete for similar lands. In proximity to water areas, especially near the Kielstau River, green colors dominate indicating the high suitability for grassland or pasture. If the blue channel forest is replaced with winter wheat in Figure 4b, blue colors indicate the suitability for winter wheat. They are mostly found in areas with low probabilities for settlement, grassland/pasture and forest that are depicted in black in Figure 4a. The non-mixing colors red, green, and blue in Figure 4b can imply that no competition exists between dominant agricultural (e.g. winter wheat) and non-agricultural land use classes. In the R-G-B composites of the main croplands in Figures 4c and 4d, the probability maps of corn (R) or winter barley (R), winter rape (G), and winter wheat (B) are compared. An obvious feature is that mixed and light colors are prevailing in the two maps, suggesting that similar areas are simultaneously suitable for two or three crop types. Moreover, mixed colors present in the same patches in both composite maps represent the suitability of similar fields for these four main crops. Due to relatively high simulated probabilities for corn (depicted in red), more purple and magenta colors (mix of red and blue) are dominant in Figure 4c than in Figure 4d indicating higher probabilities for corn as compared to the competition of winter barley with winter wheat and winter rape. The cyan-blue colors dominating Figure 4d represent areas that are suitable for winter rape and winter wheat. The suitability of fields for different crop types is in agreement with crop rotations in the Kielstau catchment. The most common rotation is a mix of cereals and rapeseed, followed by a combination of corn and cereals, and of corn, cereals, and rapeseed (Kandziora et al., 2014).

Our findings indicate which variables are important for the spatial distribution of land use classes in a small agricultural catchment in Germany. Most of these results are in agreement with our general understanding of the catchment. For instance, grassland or pasture have the highest probability near protected areas, the Kielstau river, and Lake Winderatt, which agrees with the principle that grassland or pasture instead of cropland is in direct vicinity to rivers to sustain water resource quality (Gerrish et al., 1995; Qureshi et al., 2013). Obviously, population density is a defining variable for settlement areas (Yue et al., 2013). Simpler patch shapes result in higher probability for croplands while more complex shapes are found in non-cropland areas. This is a reasonable result as crops are usually grown on larger and simpler fields as compared to other land use classes. The visual evaluation of the probability maps further confirms the plausibility of using the logistic regression models to simulate land use competition. The R-G-B composites of probability maps (Figure 4) indicate no strong competition between settlement areas, grassland/pasture, and croplands, which is reasonable in a rural environment that is not exposed to strong population pressure. However, the different crops (winter wheat, winter rape, corn, and winter barley) compete for the same locations, as indicated by the mixing and light colors in the both R-G-B composites of the probability maps Figures 4c and 4d. This is in agreement with the fact that these crops can be grown on the same fields. Crop variations are related to farmer decisions as well as to crop rotation practices of around three years that occur in the study area (Kandziora et al., 2014).

Our analysis shows that the most important variables to explain land use patterns in the Kielstau catchment are distance to protected areas, fractal dimension, drained soil share, and population density. These variables underline the agricultural character of the rural catchment, as fractal dimension, drainage density, and (low) population density are linked to agriculture. Moreover, the course of the Kielstau River affects the land use pattern as the variables distance to protected areas (Figure 2) and distance to rivers are linked to the Kielstau River, and the soil map includes properties of the flood plain (Figure 2). However, this may in part be explained by the fact that our analysis is carried out at the catchment scale. The logistic regression models are based on a dataset with a spatial resolution of 10 m. On a coarser spatial scale the spatial structure of the explanatory variables will be different and also the land use pattern will change, when smaller fields are merged. Consequently, other variables may be more important on the large scale. Nevertheless, we are confident that the 10 m resolution is appropriate for our comparatively small study area.

The quantitative model performance as well as the visual evaluation of the probability maps (Figure 3) indicates that the derived logistic regression models can reasonably explain land use patterns in the Kielstau catchment. The patterns are overall consistent with the land use maps for 2013 and 2017 (Figure 1) and the ROC values for both spatial validation and temporal validation are mostly greater than 0.7 indicating a reasonable performance. The lowest ROC value for corn (0.65) may be attributed to the fact that corn cultivation is possibly more affected than other crops by non-spatial variables like market prices and policies as it may be used for biogas production. The assumption that corn fields can be found near biogas plants was not verified. This might also be explained by the small size of the study area and the small number of biogas plants.

The range of ROC values for the different land use classes can be explained by characteristics of rural land cover, e.g., settlement areas (0.97), orchard/garden (0.89), and forest (0.84) are well distinguished and thus have strong explanatory variables (e.g. population density, distance to villages, drained soil share, respectively), however, the explanatory variables are not similarly defining for one crop or the other and, therefore, yield lower accuracies with ROCs for most agricultural land between 0.7 and 0.8 in Table 3. In many other studies, agriculture is, therefore, only one lump class (Mottet et al., 2006; Yang et al., 2012), whereas this study differentiates eleven agricultural land use classes including eight cropland classes. A possible improvement could be achieved if the land tenure system was included as a spatial variable. This variable could be used to better represent farmers’ decisions as well as the influence of the market and policies.

A comparison of the model performances in 2013 and 2017 shows that the agricultural classes vary in time and that the performances differ, while they are more or less constant for non-agricultural classes. Usually, a slightly worse performance can be expected for temporal validation (Shu et al., 2014). This applies to winter rye, winter rape, and winter barley. But it has to be noted that some slight improvements for e.g. winter wheat and other cereals can be observed, indicating that the regression models fit the 2013 data moderately better. Since models are limited out at five explanatory variables it should not be neglected that if all significant variables are included better results can be achieved.

Landscape metrics have been proven valuable in the context of our agricultural catchment. Hence, land patch structure and shape indicators are suitable to characterize and differentiate land use patterns, which is in agreement with the use of landscape metrics in the context of urban land use change or land cover classification studies (Seto and Fragkias, 2005; Fichera et al., 2012). Our results indicate that the landscape metrics provide important complementary information to the more commonly used biophysical and socioeconomic variables. Due to their explanatory power, they may also be useful in other study areas. The derived regression models are suitable to predict land use patterns and the derived probability maps can be used to visualize competition among land use classes in space, by incorporating the R-G-B composites as a simple model of land use competition. To simulate future land use patterns, the regression models can be incorporated into an intergrated application of CLUE-S (Liu et al., 2017) and Markov chain or cellular automata models (Arsanjani et al., 2013). These models also account for non-spatial variables like policy and market changes to alter the shares of different crops and derive a corresponding land use pattern.