IntroductionTop

Precision agriculture technologies and variable rate application (VRA) of inputs can play a fundamental role in the farm management, improving production and nutrient use efficiency (Monaghan et al., 2013; Hedley, 2015). Traditionally, agricultural fields have been managed as single units, although it has long been known that soil condition and crop yield are not homogeneous within them (Frogbrook & Oliver, 2007; Vitha­rana et al., 2008; Alletto et al., 2010; López-Lo­zano et al., 2010). Traffic of agricultural machines may vary in terms of intensity and geographical distribution on the field. As a consequence, high variability of soil physical properties and crop yield may occur, even in soils characterized by homogeneous distribution of physical properties (Mouazen et al., 2003; Servadio, 2010). Spatial and temporal variability of soil properties may affect crop growth, yield and quality at the within-field scale (Diacono et al., 2012). The magnitude and structure of such variability may suggest the suitability of site-specific management, with the aim of increasing both profitability of crop production and environment protection (Godwin & Miller, 2003; Mzuku et al., 2005; Vitharana et al., 2006). Site-specific management can improve the energy efficiency of the farm (Servadio & Bergonzoli, 2015) by optimizing the application of inputs. Therefore, it can reduce the negative impacts of pollution due to over-application of chemicals (Di Fonzo et al., 2001; Basso et al., 2011). Moreover, the enlargement of single management units, resulting from the enlargement of arable lands, can encourage the application of non-uniform management techniques (Sylvester-Bradley et al., 1999), including soil tillage (Servadio et al., 2014).

The subdivision of the field in management zones (MZs) is based on the knowledge of the spatial variability of soil parameters that are, generally, stable with respect of time and related to crop yield (Schepers et al., 2000, 2004). Once a data set has been acquired, cluster analysis can be performed to define the management units (Taylor et al., 2003; Fleming et al., 2004) by implementing, for instance, fuzzy k-means or Gustafson-Kessel algorithms (Höppner, 1999; Stafford et al., 1999; Vrindts et al., 2005; Guo et al., 2013).

To delineate these zones, various parameters have been evaluated in literature. For example, MZs were defined considering yield (Vrindts et al., 2005; Xiang et al., 2007; Diacono et al., 2012), soil fertility (Ortega & Santibáñez, 2007; Davatgar et al., 2012; Van Meirvenne et al., 2013), or soil electrical-hydraulic properties (Moral et al., 2010; Keller et al., 2012; Naderi-Boldaji et al., 2013; Doolittle & Brevik, 2014). A combined use of different sets of parameters, such as a combination of physical and chemical soil parameters (Servadio et al., 2017), could lead to an in-depth investigation into spatial heterogeneity and to a more comprehensive knowledge of soil plant system (Beni et al., 2012). Among the physical parameters, soil strength influences many aspects of the cultivation, such as tractors performance during tillage (Servadio & Bergonzoli, 2015; Servadio et al., 2016) and root growth. Furthermore, when compaction occurs, soil permeability and regeneration can be reduced (Manuwa & Olaiya, 2012). Variations of soil texture may also have a significant effect on soil management, as studied in previous investigations (Vitharana et al., 2006; Gooley et al., 2014; Havaee et al., 2014).

The novel contribution of this paper consists of a new combination of limiting-governing factors on three different soil texture (clay loam, sandy clay loam, clay), in order to identify and investigate potential MZs. More specifically, soil texture, organic matter and water content were selected as limiting factors for the yield and as governing factors for the soil mechanical status. The three different soil types were classified as Cambisol (FAO, 2006). The soil mechanical status was described by means of the following parameters: structural stability, shear strength, and penetration resistance. The limiting-governing factors were used as input parameters in a fuzzy clustering algorithm, in order to: i) identify potential management zones; ii) analyze the influence of the soil attributes on soil mechanical properties and on crop yield; iii) analyze and compare to each other the results related to different soil types; iv) provide information for site-specific management within the field. The results of this study provide the basis for the generation of prescription maps, therefore for VRA of inputs and variable soil tillage.

Material and methodsTop

Site and data acquisition

The tests were performed in three different fields belonging to a farm in center of Italy (locality Maccarese, province of Rome, Lazio). As shown in Fig. 1, the three fields, labeled F1-N [41°52’59” Latitude (N), 12°13’37” Longitude (E)], F2-C [41°50’54” Latitude (N), 12°14’48” Longitude (E)], and F3-S [41°49’28” Latitude (N), 12°14’10” Longitude (E)], are located in the northern, central, and southern sector of the farm, respectively. The chosen fields represent three different soil types. For each field, soil texture, soil classification and field size are reported in Table 1.

Figure 1. Fields location in a farm of central Italy.

In this study, three cropping seasons were consi­dered, and meteorological data were acquired for each season by a weather station type Davis Vantage Pro2, located inside the farm [41°53’13” Latitude (N), 12°11’18” Longitude (E)]. More specifically, monthly rainfall and minimum and maximum temperatures were recorded from October 2007 to July 2010. The climatic conditions, typical of the Mediterranean climate, were reported and discussed in a previous study (Servadio et al., 2017). The same management technique was applied for each field. Every year, the soil was ploughed at 0.40 m depth and then harrowed at 0.20 m depth with a rotary harrow. After the preparation of the seedbed, wheat (Triticum durum L.) was sown in each field in the second half of October, and immediately was fertilized with 150 kg/ha of ammonium nitrate. A second fertilization was performed at the mid-season growth stage with 200 kg/ha of urea. Finally, the crop was harvested in mid-July for each growing season, using a combine harvester (New Holland, NH CX860 equipped with PLM System). To obtain grain yield values, the harvester was equipped with GPS sensor (for the acquisition of EGNOS signal) and a grain mass flow sensor. The acquired data were processed with the Precision Land Management Software (Case New Holland, Inc.).

Sampling test and measurements

To determine some soil physical properties, and quantify within field spatial variability, a georeferenced soil sampling was performed at the end of the third growing season. Soil samples were collected on a 30×30 m grid pattern, from 0 m to 0.20 m depth (for composited samples, 3 sub-samples were collected from the vertices of a 3 m side triangle). The defined grid led to the definition of 20 points per field, corresponding to about ten samples per hectare of soil. This number of samples is considered acceptable (Diacono et al., 2012). Each sampling location was georeferenced using a GPS device (mod. GEO XM, Trimble navigation, CA, USA) and differential correction was performed using the post processing software Trimble Pathfinder Office (vers. 4.20, Trimble navigation, CA, USA). For each sampling point, structural stability, water content, soil texture and organic matter were determined. Furthermore, as soil strength indicators, measurements of shear strength and penetration resistance were performed, from 0 m to 0.20 m depth. Structural stability of soil aggregates was determined on the 0.25 mm fraction using a wet sieving apparatus (8 sieves) Eijkelkamp, by means of the Kemper method (Kemper & Chepil, 1965). This method is based on the principle that unstable aggregates break down more easily than stable aggregates when immersed into water. The sieves were filled with a certain amount of soil aggregates and placed in a can filled with water, moving upward and downward for a fixed time. Unstable aggregates fell apart, passed through the sieve and were collected in the water-filled can underneath the sieve. After, the cans were removed and replaced by new water-filled cans. When all aggregates were destroyed, sand and plant roots remained on the sieve, and only aggregates were considered. After drying the cans with the aggregates, the weight of both stable and unstable aggregates were measured, and the weight of stable aggregates was divided by the total weight of the aggregates.

Regarding the water content, the soil samples were extracted with a corer sampling ring, weighed and then dried until they reached a constant weight. The soil texture was determined by sedimentation after separation of clay, sand, and silt, whereas the organic matter content was derived from the total organic carbon (C × 1.72). Considering the soil strength indicators, shear strength was measured using a field inspection vane tester, from 0 to 260 kPa (Eijkelkamp model 14.05). According to the ASAE standard, the cone index was determined using a self-recording electronic penetrometer (Eijkelkamp penetrologger model 06.15.SA), with a 60° cone and base area of 100 mm2, driven into the soil at a constant velocity (5 cm/s).

Identification of management zones

In order to define site-specific management zones, a cluster analysis based on fuzzy c-means algorithm (FCM) was conducted to: (i) identify, for each field, potential management zones; (ii) analyze, for each field, the influence of some soil attributes on soil mechanical properties and on crop yield; (iii) compare the results of each field and identify the most important factors affecting soil mechanical status and crop yield.

Input and output parameters of the cluster analysis

Soil texture, organic matter and water content were selected as limiting factors for the yield and as governing factors for the soil mechanical status, represented by the following set of soil mechanical parameters: structural stability, shear strength and penetration resistance. To avoid the spurious correlations due to the compositional nature of the textural fractions (individual elements sum to 100%), only the clay fraction was considered as input (Vitharana et al., 2008). Therefore, clay percentage, organic matter and water content were selected as input parameters in the cluster analysis, whereas structural stability, shear strength, penetration resistance and yield were selected as test parameters.

The previous parameters were chosen for the following reasons.

1. Input parameters:

a) Clay content was selected as soil strength indicator and as indicator of soil compaction susceptibility (Marsili et al., 1998; Servadio et al., 2005). According to Kumar et al. (2012), for both no-tillage and conventional tillage, cone index values depends on clay fraction values, whereas Eneje & Adanma (2007) observed a trend and a degree of relationship between clay content and aggregate stability.

b) Water content was chosen as yield limiting factor (Melo Damian et al., 2016). This parameter, generally, varies with respect to time. However, despite large variation over time and space of the specific points of the field, Vachaud et al. (1985) showed that spatial patterns of soil moisture changed little with respect to time. Temporal stability of soil moisture was recognized by some authors (da Silva et al., 2001; Starr, 2005; Cosh et al., 2008; Guber et al., 2008; Brocca et al., 2009; Coppola et al., 2011). Since water content is related to soil hydraulic and mechanical properties and it follows quasi-steady spatial patterns with respect to time, it was considered as input parameter in the cluster analysis. In particular, water in soil acts both as a lubricant and as a binding agent among the soil particulate materials, then it influences soil structural stability (Dane & Topp, 2002). It also affects soil shear strength (Bláhová et al., 2013) and cone index values (Kumar et al., 2012; Servadio, 2013).

c) Organic matter was selected as soil quality indicator. The results from the study conducted by Ekwue (1990) showed that the influence of organic matter on soil shear strength depends on whether or not it improves soil aggregate stability. Moreover, Soane (1990) highlighted that: i) strength of aggregates is closely related to the presence of organic matter, roots and fungal hyphae; ii) maintenance of aggregate stability through appropriate crop management practices will enhance soil resistance to compaction loads; iii) soil compaction is sensitive to even quite small change in the amount of organic matter.

2. Test parameters:

Soil structural stability, shear strength and pene­tration resistance were selected as test data set. To evaluate the performance of the analysis, yield was also included in this set, rather than in the input set for delineating MZs. The grain yield is the resultant of the complex interactions among root growth and soil, in term of soil structure, soil strength, water content and soil nutrient availability (Córdoba et al., 2013). Soil compaction reduces the soil permeability to water, reduces regeneration of the soil (compromising metabolic activities of roots), and increases the mechanical strength of the soil (obstructing root growth) (Tracy et al., 2011). The relative importance of these factors depends on the soil-water regime. Shear strength and penetration resistance, considered as soil strength indicators, were widely used to assess trafficability and workability of soil, and to define the status of soil in function of the crop growth (Servadio et al., 2014, 2016; Servadio & Bergonzoli, 2015). A reduction in soil compaction resulted in an increase of yield, attributed to the lowering of the mechanical impedance to root growth (Hamblin, 1985). Therefore, soil penetration resistance and shear strength can help to identify areas where soil mechanical characteristics are negatively affecting yield. Both shear strength and penetration resistance depends on clay fraction (Kumar et al., 2012) and on moisture content (Servadio et al., 2005; Manuwa & Olaiya, 2012; Bláhová et al., 2013).

Cluster analysis algorithm

Fuzzy classification produces a continuous grouping of objects by assigning partial class membership values, according to the properties variability in the soil continuum (Bezdek et al., 1984; Cannon et al., 1986; Pal et al., 2005). In the present investigation, fuzzy c-means algorithm was applied by using the FuzMe software functions (Minasny & McBratney, 2002) in MATLAB (www.mathworks.com).

In order to perform the cluster analysis, the number of zones and the fuzziness exponent have to be defined. In many investigations, the number of classes range from two to eight or nine (Fridgen et al., 2004; Vrindts et al., 2005; Reyniers et al., 2006; Davatgar et al., 2012). Regarding the fuzziness exponent f, its value can be chosen between 1 and infinity. A value close to 1 represent a FCM approaching the hard c-means algorithm, while a value approaching infinity represents a solution approaching the highest degree of fuzziness. Odeh et al. (1992) set f equal to 1.35, while Bezdek suggested [1.5-3] as range of possible values (Bezdek et al., 1984). Many researchers proposed f = 2 (Yu et al., 2004).

In this investigation, calculations were performed by setting equal to 8 the maximum number of clusters, and equal to 2 the fuzziness exponent. The convergence tolerance was set equal to 0.0001. The Mahalanobis distance was used to calculate the distance of data points to cluster center points (Vrindts et al., 2005; Reyniers et al., 2006). Regarding the clustering validation, fuzziness performance index (FPI), normalized classification entropy (NCE) and separate fuzzy validity (S) were used to evaluate the optimum number of classes. More specifically, the FPI function estimates the degree of fuzziness generated by a specified number of classes, while the NCE function estimates the degree of disorganization created by a specified number of classes (Gorsevski et al., 2006). The S function is a measure of the ratio of variance within the clusters to the variance between the clusters (Xie & Beni, 1991). The optimal number of clusters were selected by minimizing the FPI, NCE, and S indexes. In case of more than one minimum, the option with the smallest number of clusters was selected, following a rule of parsimony (Lark & Stafford, 1997; Reyniers et al., 2006).

Statistical methods

Descriptive statistics including mean, standard deviation (SD), coefficient of variation (CV), and standard error (SE) were determined for the input parameters (clay content, organic matter, water content), for the test parameters (structural stability, shear strength, penetration resistance), and for the yield of the three seasons (2007/08, 2008/09, 2009/10).

ResultsTop

Wheat yield

The mean values of wheat yield (referred for each field, to the three consecutive growing seasons 2007–10), together with standard deviation, coefficients of variation and standard error are listed in Table 2.

Considering the temporal variability, the data show that, both in 2007–08 and 2009–10, yield mean values were in the range 5.00–7.03 t/ha, except for F2-C (2.74 t/ha in 2009-10, lower with respect to the others means). In 2008–09, the mean yield values were lower (2.80–3.81 t/ha) with respect to the above reported values. Such temporal variability could be ascribed to the meteorological pattern recorded during the three seasons. In particular, referring to the rain precipitation, in the period from October to July mean values of rainfall were higher in 2009–2010 (824 mm) and 2008–2009 (642 mm) in comparison with 2007–2008 (446 mm). More specifically, during February, March, and April 2009, rainfall values resulted very low (27, 51 and 15 mm, respectively) with respect to the same period of 2008 (91, 34 and 75 mm, respectively) and 2010 (68, 67, and 73 mm, respectively). Furthermore, in May 2009, a very high rainfall value (52 mm) was recorded with respect to the same month of 2008 (28 mm) and 2010 (8 mm). See Servadio et al. (2017) for more details.

Considering the spatial variability, F2-C and F3-S were more productive with respect to F1-N in 2008-09, whereas F1-N and F3-S were more productive with respect to F2-C in 2009-10. Since the weather conditions were the same, such variability could be ascribed to the soil physical-chemical properties.

Average values of the three growing seasons show that the recorded yield was higher in the F3-S (clay) and in F1-N (clay loam). Both fields were characterized by high content of clay and organic matter (as shown in Table 3) with respect to F2-C (sandy clay loam).

Soil parameters

The mean values of soil water content (WC), organic matter content (OMC) and clay content (CC) (used as data sources in cluster analysis), together with standard deviation, coefficient of variability and standard error, are listed in Table 3. Results show that CC was higher in F3-S and in F1-N (583 and 345 g/kg, respectively) with respect to F2-C (307 g/kg). WC in Field F1-N was higher (219 g/kg) in comparison with F2-C and F3-S (198 and 197 g/kg, respectively), whereas OMC was higher in F1-N and F3-S (23.4 and 24.3 g/kg, respectively) with respect to F2-C (13.5 g/kg).

Results of mean values of soil structural stability (STA), shear strength (STR) and cone index (CI), together with standard deviation, coefficient of variability and standard error, are shown in Table 4. STA was higher in F2-C (683 g/kg) with respect to F3-S (496 g/kg) that, in its turn, was higher with respect to F1-N (332 g/kg). STR and CI were higher in F3-S (215 kPa and 2.47 MPa, respectively) with respect to F2-C (58.1 kPa and 1.18 MPa, respectively), which also were higher with respect to F1-N (33.4 kPa and 0.39 MPa, respectively). The analyzed soil parameters and their fields variability can be considered as good indicators of the soil status and then as input parameters for the definition of MZs.

Definition of management zones

Regarding the cluster analysis performed on the F1-N, F2-C and F3-S fields, the optimal clustering results were obtained with a number of classes equal to 4, 5 and 5 respectively. The values of FPI, NCE, and S, with respect to the number of classes c, are shown in Fig. 2.

Figure 2. Values of the indices FPI, NCE and S, with respect to the number of classes c: (a) feld F1-N, (b) feld F2-C and (c) feld F3-S.

Tables 5 and 6 list, for each field and for its corresponding MZs, the mean values and the CV of the input parameters and of the test parameters, respectively. The defined management zones of the fields F1-N, F2-C and F3-S are depicted in Figs. 3a, 3b, and 3c, respectively.

Figure 3. Management zones maps of the three felds: (a) F1-N; (b) F2-C; (c) F3-S.

- Management zones of the Field F1-N

Cluster analysis of F1-N produced four management zones, as shown in Fig. 3a. The clay-loam field was characterized by 4,80 t/ha three years average yield (Table 2), 219 g/kg WC, 23.4 g/kg OMC, and 345 g/kg CC (Table 3). As reported in Table 5, with respect to the mean values of the whole field, zones 1 and 2 achieved high values of CC (368 and 343 g/kg, respectively), OMC (23.6 and 25.0 g/kg, respectively) and WC (228 and 224 g/kg, respectively). Zone 3 achieved a high value of CC (346 g/kg) and the lowest value of OMC (21.7 g/kg). Regarding the data set parameters showed in Table 6, zones 1, 2 and 3 showed the highest values of STA (449, 368 and 358 g/kg, respectively), and zones 2 and 3 achieved the higher values of yield (5.09 and 4.93 t/ha, respectively). A trend of very low soil strength values in term of STR and CI was found among the zones.

- Management zones of the Field F2-C

Cluster analysis of F2-C produced five management zones, as shown in Fig. 3b. The sand clay loam field was characterized by 4.29 t/ha three years average yield (Table 2), 198 g/kg WC, 13.5 g/kg OMC and 307 g/kg CC (Table 3). As reported in Table 5, with respect to the mean values of the whole field, zones 4 and 5 achieved the highest values of CC (371 and 340 g/kg, respectively) and WC (219 and 220 g/kg, respectively ). Also, OMC was equal to 14.2 g/kg in both zones, very close to the maximum value registered in Zone 3 (14.4 g/kg). Regarding data set parameters in Table 6, zone 2 showed the higher values of STR (73.9 kPa), CI (1.34 MPa) and STA (762 g/kg). Zones 3 and 4 showed high STA values (720 and 718 g/kg, respectively) and zones 4 and 5 achieved highest values of yield (4.68 and 4.64 t/ha, respectively). A trend of medium soil strength values in term of STR and CI was found among the zones, and zones 4 and 5 showed the lower values of CI (1.09 and 1.11 MPa, respectively).

- Management zones of the Field F3-S

Cluster analysis of F3-S produced five management zones, as shown in Fig. 3c. Clay Field 3 (F3-S) was characterized by 5.25 t/ha three years average yield (Table 2), 197 g/kg WC, 24.3 g/kg OMC and 583 g/kg CC (Table 3). As reported in Table 5, with respect to the mean values of the whole field, zones 1, 2 and 3 achieved the highest values of CC (588, 594, and 586 g/kg, respectively). Zones 3 and 4 presented high values of OMC (26.5 and 25.2 g/kg, respectively), and zones 2, 3 and 5 presented high values of WC (239, 207 and 207 g/kg, respectively). Regarding the test parameters in Table 6, zones 2, 3 and 4 showed high values of STA (500, 517 and 497 g/kg, respectively), and zones 1, 4 and 5 achieved the highest yield (5.26, 5.34 and 5.34 t/ha, respectively). A trend of high soil strength values in term of STR and CI was found among the zones: zones 4 and 5 showed the highest values of shear strength (233 and 262 kPa, respectively), whereas zones 2 and 5 showed the highest values of CI (2.64 and 2.60 MPa, respectively).

DiscussionTop

The delineation of MZs represents a management tool allowing farmers to improve the environmental sustainability and the energy efficiency of their farm. In a previous study (Servadio et al., 2017), the definition of three management zones within a field showed that shear strength and structural stability were the most significant limiting factors for wheat yield. Furthermore, results showed that, in the field conditions of the tests, nitrogen fertilization could not be considered as a limiting factor for the yield.

In the present study, considering the three different fields, yield mean values (Table 2) showed a temporal variability between the three consecutive growing seasons (2007–10). This could be ascribed to the meteorological conditions (Basso et al., 2009), in particular to the difference in rainfall distribution recorded during the three seasons (Servadio et al., 2017).

According to Melo Damian et al. (2016), moisture could hardly be considered a limiting factor for wheat growth in 2008-09 growing season. Considering source parameters and test parameters, values associated with the whole field are in agreement with the values obtained from average yield of the three growing seasons. In fact, as shown in Table 2, yield was higher in the fields F3-S (clay) and F1-N (clay loam), which are both characterized by a higher content of clay and organic matter (as shown in Table 3) with respect to F2-C (sandy clay loam). According to Kumar et al. (2012), in F3-S, due to the higher content of clay and to the lower content of water, soil strength (in terms of STR and CI) was higher both with respect to the F2-C and to the F1-N (as shown in Table 4). Furthermore, higher percentage of structural stability was found in F2-C (sandy clay loam).

Considering the results of the cluster analysis, from the defined classes within each field emerged that input soil attributes influenced soil physical-mechanical status and crop yield. In particular, in field conditions of F2-C (zones 4 and 5), high values of CC, OMC and WC correspond to high values of yield and to low/mid values of STA. On the contrary, in zone 2, low values of CC, OMC and WC correspond to low values of yield and to high value of STA.

Similar results can be observed for F1-C (zones 2 and 3), where mid/high values of CC and WC correspond to high values of yield and to mid values of STA. However, OMC was high in zone 2 and low in zone 3. Also, in zones 1 and 4, low values of yield correspond to high CC in zone 1 and to low CC in zone 4.

In F3-S field conditions, the five zones were characterized by close means values of CC and yield. In zone 3, high values of CC, OMC and WC correspond to high values of STA. Furthermore, the high soil strength values found in field 3 seems to not affect yield.

The results obtained from cluster analysis can provide the basis for the application of VRA of input and variable tillage (VT). In particular, minimum tillage and/or sod seeding instead of deep plough, the latter being often applied in Central Italy (Servadio et al., 2014). VRA and VT can reduce total energy employed, fossil-fuel energy requirements, carbon dioxide emissions, and total cost of the crop cycle, at the cost of a small reduction of the crop yield (Servadio et al., 2014, 2017; Servadio & Bergonzoli, 2015).

As conclusions, in this paper, cluster analysis based on fuzzy c-means algorithm was applied to delineate site-specific management zones within three fields characterized by different soil texture. Crop yield and various soil parameters, both physical-chemical (soil structural stability, clay fraction, water content, and organic matter) and mechanical (shear strength and penetration resistance) were mapped and processed. Results obtained from each field were analyzed and compared to identify the most important factors affecting the soil mechanical status and the crop yield. These results highlighted: (i) spatial variability of the analyzed soil parameters; (ii) spatial and temporal variability of the grain yield between three consecutive growing seasons (2007–10), ascribed to the difference in rainfall distribution recorded during the three seasons; (iii) agreement between the whole-field values of source and test parameters and values obtained from average yield of the three growing seasons. Generally, high values of clay and organic matter contents corresponded to high values of grain yield and, in some zones, were effective to enhance structural stability. In the clay field, higher clay content related to the lower content of water increased soil strength (in terms of shear strength and cone index). Finally, the high soil strength values found in clay field seems to not affect grain yield.

The results obtained in this investigation stressed the importance of the definition of management zones as a tool to address spatial and temporal variability of soil properties and crop yield. The delineation of these zones lay down the basis for the application of further precision agricultural practices (variable rate and variable tillage applications), which can improve the energy efficiency of the farm and its environmental impact.