Effect of Management on Spatial and Temporal Distribution of Soil Physical Properties

The purpose of this study was to characterize spatial and temporal variations in soil physical properties and to examine possible differences that might occur when the crop rotation shifted to rice cultivation. Soil samples from 0-20 cm and 20-40 cm depth were collected on a 100×100 m grid basis over a 45 ha field, in May 2005 and May 2011. The geostatistical methods were used to model the variance structure of gravimetric water content (GWC), soil bulk density (BD), penetration resistance (PR) and aggregate stability (AgS). The coefficient of variations were 17.5% for GWC, 7.1% for BD, and 25.0% for PR at surface and 12.8% for GWC, 6.7% for BD, and 17.1% for PR at subsurface. PR at subsurface (>2MPa) and BD at surface and subsurface soils (1.47 g cm) exceeded the threshold values in some fields of the study area. Temporal variations of GWC and PR at 0-20 cm and 20-40 cm depths were significant (P<0.01) while BD was significantly varied only at 20-40 cm. Spatial correlation ranges of GWC, BD and PR varied from 72 m to 433 m in 2005 and 40 to 1524 m in 2011. The spatial dependences of GWC, BD and PR at surface were moderate, and AgS had strong dependency at this depth. The spatial distribution maps of soil physical attributes are useful in identifying the limiting factors and take necessary precaution to prevent further loss of soil quality.


Introduction
The appearance of distinct soil physical properties and soil fertility zones are indications of considerable spatial variability inherent to the soil (Wendroth et al 2003).Geostatsitics, a branch of statistics focusing on spatial or spatiotemporal datasets, constitute an important tool in determining spatial effects of soil management practices (Maiorana et al 2001).Identifying the zones with soil compaction problems and developing management options to minimize crop production risks and the harmful impact of traffic on the environment can be achieved by assessing the spatial distribution of the penetration resistance (Özgöz et al 2007).
Considering the spatial variability of soil properties within a given area is important to evaluate the plant production potential of that area, and to optimize profitability, sustainability and protection of the environment.The spatialtemporal variability of important soil properties and associated plant biomass production are also useful developing effective sampling schemes for future site management (Mulla et al 1992).The extent of soil spatial variability depends on the variations of soil forming factors and the management practices applied for a particular crop growth (McGraw 1994;Mulla & McBratney 2000).Soil management practices are important sources of temporal variability than soil type.Therefore, temporal variability of soil physical properties can be even greater than spatial variability in agriculturally managed soil.Many of soil physical attributes particularly the hydraulic properties significantly vary even in a short time period, such as during a crop cycle, especially immediately after tillage (van Es et al 1999).
Soil physical properties such as bulk density, penetration resistance and water content temporally and spatially change as a result of natural phenomena and human activities which are dependent on seasonal climatic conditions, management practices, crop development and biological activity (Reynolds et al 2007).Changes in soil physical properties due to the wetting and drying cycles (Dörner et al 2009), different tillage management practices (Özgöz et al 2007), grazing (Krümmelbein et al 2009) and land use changes (Özgöz et el 2011) have been reported.Dynamic of soil structure has significant effect on many of soil physical properties, thus management systems are main agents for alteration of soil environmental conditions.The effects of several management systems have been clarified by identifying the temporal and spatial changes of soil physical properties (Alletto & Coquet 2009).
Due to the simplicity, low cost of analysis and efficiency in a short period of time, penetration resistance which describes soil strength (Draghi et al 2005) or resistance to root growth (Dec et al 2011) is used to evaluate temporal and spatial changes in soil structure, in addition to soil compaction.
However, since penetration resistance changes with temporal variation of soil water content defining the theoretical basis can be a problem (Horn & Fleige 2009).Korucu et al (2009) explained the spatial and temporal variation of soil compaction during post harvest period and after stubble burning of wheat by daily penetration resistance measurements.They concluded that moisture loss due to stubble burning resulted in increased penetration resistance, and delayed planting which causes postponing the harvest operations.
Optimizing water applications in drip fertigation systems, Mubarek et al (2009a) conducted some studies to explain the temporal changes of soil properties.Strudley et al (2008) related the dynamics of temporal variability in soil physical properties and processes to tillage management practices; however Mubarek et al (2009a) and Mubarek et al (2009b) stated that temporal variability during a cropping season under drip fertigation remained poorly understood.
The variations in soil characteristics are needed to be monitored to sustain soil quality and enable agricultural production.Thus, the purposes of this study were to characterize spatial and temporal patterns of changes in physical properties of soils with the objective of examining possible differences that might occur when the crop rotation shifted mainly to rice cultivation.

Description of study area
This study was carried out in a 45 ha agricultural field located between the Kelkit River and the Canik Mountains within the middle Black Sea region of Turkey.The physiographic unit of study area is a young stream terrace formed by Kelkit River with 0-2% slope.Coarse particles and gravels found at 50 cm depth of soil profiles in the middle of the study area indicate that a small stream has probably changed its bed several times in the past.Soils were high in clay with an average of 42.9% clay content in surface (0-20 cm) and 43.3% clay in subsurface (20-40 cm).Soils were classified as Typic Ustifluvents (Soil Survey Staff 1999;Özgöz et al 2009).The area has terrestrial climate with a mean annual precipitation of 456.4 mm and majority of precipitation occurs (47.5%) in November, December, April and May; and mean annual temperature is 12.3°C.

Field history
The field probably was not homogenously leveled prior to land leveling studies completed during the 1999 and 2000 growing seasons.A small ridge in the middle of the field was also leveled during this operation.Initial survey on the same land has been conducted on 2005 by Özgöz et al (2009).The rotation used till 2005 was mainly winter wheat (Triticum spp.), sugar beet (Beta vulgaris), maize (second crop) (Zea mays L,), tomato (Solanum lycopersicum) and alfalfa (Medicago sativa).Due to the incentives initiated by the Ministry of Food, Agriculture and Animal Husbandry, the rotation of area started shifting to rice cultivation.At the sampling time, some portions of the study area had been planted with winter wheat, sugar beet and alfalfa, and the rest of the field was prepared for tomato and rice plantation (Figure 1).The rice plantation covers almost 60% of the study area.
The conventional tillage practices were used in crop production.Moldboard plough, cultivator, tooth harrow and drill machine were used for soil preparation and planting.The hoeing operation in sugar beet and tomato was achieved by hand tools, and harvesting was performed by conventional methods.Crops were irrigated by surface flooding.In rice production; conventional tillage is used with a moldboard plough at a depth of 20 cm in October-November, and tooth harrow is used at a depth of 10 cm.Management of rice soils requires ponding water, and usually includes diking, leveling, and puddling.The field is leveled using a scraper to obtain almost level surface enabling even water distribution throughout the field.Rice is cultivated in paddy fields which are flooded to obtain paddy soil.The word puddle means clay that can be worked to a water impermeable stage (De Datta 1981).Premoistened rice seeds are planted via broadcasting.The rice is harvested by combine harvesters.

Sampling design and measurements of soil physical properties
Sampling scheme in this study was the same as in 2005 that the field was divided into 100 x 100 m grid cells with forty four main sampling points.The size of field was 900×500 m, and soil samples were collected from upper nodes of each grid in May, 2011.Ten-fine transects distributed in north-south and east -west directions regularly with 2, 5, 10, 25 and 50 m intervals were sampled to determine the spatial variability of short distance (Figure 1).Soil samples were collected both from topsoils (0-20 cm) and subsoils (20-40 cm).One hundred eighty eight soil samples were collected that was half (94) from topsoil and the rest was from subsoil (Figure 1).Fields have been cultivated to plant sugar beet, tomato and rice, and wheat was at the booting stage during soil sampling.
Bulk density (BD) and gravimetric water content (GWC) were determined as described by Blake & Hartge (1986).Both BD and GWC were measured with duplicate undisturbed samples in 100 cm 3 cylinders.Undisturbed soil core samples were taken both from 10-15 cm and 30-35 cm to represent 0-20 cm and 20-40 cm soil depths, respectively.The wet sieving method of Kemper & Rosenau (1986) was used to determine soil aggregate stability.Penetration resistance (PR) was measured by a hand penetrometer (Eijkelkamp Co.), and the penetrometer was measured up to 5,000 kPa with maximum depth of 80 cm.The measurements were made with a cone with angle of 30° and base area of 1 cm 2 .The penetrometer was pushed into the soil with a speed of 2 cm s -1 , and at every 1 cm (Eijkelkamp 1990).

Statistical analysis
The data collected were grouped in two classes as rice cultivated fields and the rest of the fields cultivated by wheat, sugar beet, tomato and alfalfa.Statistical parameters of mean, standard deviation, maximum, minimum, coefficient of variation, skewness, and kurtosis were calculated for each of the variables.
Data analysis for each of soil properties were conducted in three steps: (i) frequency distribution was examined and the tests for normality were conducted, (ii) exploratory data analysis was carried out calculating minimum, maximum, arithmetic mean, standard deviation, the coefficient of variation (CV), skewness, and kurtosis for each variable, and (iii) correlations between soil properties were calculated.The t-test was used to determine whether parameters differed by years.Classical statistics were performed by SPSS 10.0 statistical software (SPSS 2000).
We checked the skewness of the data, those with skewness greater than  0.5 were subjected to normality test.Normality test was conducted to check the distribution of the variables with the Kolmogorov-Simirnov test, and those variables (gravimetric water content at both depths and penetration resistance at surface) without normal distribution were subjected to log transformation.Spherical and exponential models were fitted to the semivariograms, and selected based on visual best fit and the corresponding coefficient of determination (R 2 ).The data for was checked for anisotropy to calculate the directional semivariograms.The minimum number of data pairs for a lag was limited as 20 for a safe calculation of semivariance, and the maximum and minimum data pairs in each lag varied.Modeling of isotropic experimental semivariogram was performed with GS + (version 7) statistical software Gamma Design Software (2004) and a maximum lag distance of 550 m was applied in experimental semivariogram modeling.The ordinary point kriging method was used with isotropic semivariogram models in preparation of interpolation maps for the soil properties investigated (Isaaks & Srivastava 1989).For ordinary point kriging procedure, at least 12 neighboring points were considered.The maps of soil properties were drawn with Geostatistical extension of ARCGIS 8.1 (Esri 2001) using the semivariogram parameters (nugget, sill, nugget ratio, range) obtained with Software GS + 7.0.
Model parameters (nugget semivariance, range and sill or total semivariance) were also calculated.Nugget semivariance is the variance at zero distance and represents random field and experimental variability that is not detectable at the sampling scale.Sill is the lag distance between measurements at which one value for a variable does not influence neighboring values.Range is the distance at which values of one variable become spatially independent of another (Trangmar et al 1985).Furthermore, the nugget to sill ratio indicates the degree of randomness in the data's spatial variability.This ratio was used to define three classes of spatial dependence for 16 measured soil variables (Cambardella et al 1994), that is, (i) when the ratio was <25%, the measured variable was considered strongly spatially dependent; (ii) when the ratio was between 25% and 75%, the soil variable was considered moderately spatially dependent; and (iii) if the ratio was 75%, or the slope of the semivariogram was 0, the variable was considered random or nonspatially correlated (pure nugget).

Gravimetric water content (GWC) and bulk density (BD)
Soil samples in both years were collected almost at the same time of the year (May), though mean GWC at 0-20 cm (26.7% in 2005 and 43.5% in 2011) and 20-40 cm depths (26.6% in 2005 and 41.4% in 2011) were significantly (P<0.01)higher in 2011 for whole field (Tables 1 & 2).Increase in precipitation and conversion to rice management were the major reasons to observe predominantly high gravimetric water content values in 2011.Based on CV, CV ≤15% is low, 16-30% is medium; and ≥30% is high variability (Wilding et al 1994).Soil heterogeneity affects soil moisture content through variations in soil texture, soil water holding capacity, soil color, soil water retention, and pixel-and pore-scale hydraulic properties (Jacobs et al 2004).Among the stated controlling factors of soil water content, texture particularly clay content of soils affects all the others.Clay content of surface (P<0.01) and subsurface (P<0.05)soils had positive correlations with GWC indicating the significant effects of clay content in GWC (Table 3).The variations of clay content at surface and subsurface soils (21.8% and 23.6%) were higher compared to the CV of GWC (17.5% and 12.8% for surface and subsurface, respectively).The   2004).The homogenization of rotation through rice cultivation probably reduced the variation of GWC in the study area.Paired t-test was used to evaluate the differences in GWC between the years.The results showed that changes in crop rotation reduced the variation in GWC, and yielded significant differences (P<0.01) in GWC between years evaluated (Tables 1 & 2).
Water movement, soil compaction, soil aeration, and plant root development are greatly affected by soil water content and bulk density (Timm et al 2006).Soil bulk density ranged from 1.0 to 1.8 g cm -3 and 1.3 to 1.8 g cm -3 , having a mean of 1.4 and 1.5 g cm -3 in topsoil and subsoil, respectively (Table 1).Similar to the variability's of bulk density values measured in 2005, the bulk densities at surface and subsurface soils had low variability's.The bulk density of surface and subsurface soils had significant negative (P<0.01) effect on aggregate stability indicating that stable aggregates increases the resistance of soil degradation.Although many others have reported  2).In contrast to our results, significant temporal changes were reported by Logsdon & Cambardella (2000) examining the soil bulk densities of surface soils for tillage systems in a sub-humid climate.

Penetration resistance (PR)
The PR values at 0-20 cm in 2005 ranged from 0.4 to 3.9 MPa and for 20-40 cm ranged from 0.7 to 4.8 MPa.After six years of managements, PR at 0-20 cm changed from 1.0 to 2.8 MPa, and for 20-40 cm from 2.3 to 5.6 MPa.The temporal variation in PR values at both depths were significantly different (P<0.01)(Tables 1 & 2 1 & 2).Although PR is known as a function of GWC and bulk density (Ayers & Perumpral 1982;Horn & Lebert 1994;Craig 1998;Dec et al 2011), the variability of PR was higher as compared to the GWC and bulk density.Higher degree of variability in PR is probably related to combined effects of physical, chemical and biological processes that operate at different intensities and scales in the study area (Goovaerts 1998).

Aggregate stability (AgS)
Puddling in rice cultivation naturally destroys soil structure, decreases large pores, and increases small pores in the surface layer (De Datta 1982), and expected to weaken the AgS of soils.Since aggregate stabilities of soils have not been determined in 2005, temporal changes in AgS will not be discussed in this study.The values of AgS in rice cultivated fields (77.9% for 0-20 cm and 74.1% for 20-40 cm) were slightly higher than the rest of the fields (73.8% for 0-20 cm and 72.2% for 20-40 cm) as in organic matter content (Table 1).Organic binding agents promote aggregate stability through binding soil particles to each other more strongly than to adjacent particles to form soil aggregates (Tisdall & Oades 1982).The content of organic matter in paddy soil has probably been increased with the return of organic biomass resources such as rice straw into soil following the harvest of rice (Table 1).Norton et al ( 2006) stated that soil internal (particle size distribution, organic matter, clay mineralogy, CaCO 3 , Fe and Al oxides) as well as external, time dependent properties (climate, tillage, biological activity, wetting-drying) have considerable influence on aggregate stability.Significant correlations (P<0.01)obtained between AgS and clay content (r=0.33),sand content (r=-0.30),GWC (r=0.36),BD (r=-0.45)and PR (r=-0.34)(Table 2).The variations of AgS (CV=11.10 for rice, 12.44% for others at 0-20 cm and 11.54% for rice, 15.91 for others at 20-40 cm) were lower compared to the variations in organic matter and clay content (Table 1) indicating the major controlling factors of AgS are organic matter and clay content in study area.These low CV values could be ascribed to the homogenization effect of rice cultivation and subsequent similar management practices applied within the fields.

Spatial and temporal variability
Measurements of PR under field conditions and GWC and BD in laboratory were used to prepare maps of spatial and temporal variations (Figures 2,3,4 & 5).Although AgS was not determined in 2005, spatial distribution of AgS in 2011 was analyzed and presented in this section.The physical attributes of soils (GWC, BD, PR, and AgS) had spatially dependent structures which described by semivariograms with a defined rank fitting to the spherical and exponential models (Tables 4 & 5, Figures 2, 3, 4 & 5).Geostatistical parameters of physical soil properties for 2005 and 2011 were presented in Tables 5 & 6.Figures 2, 3, & 4 show isotropic experimental semivariograms for gravimetric water content, bulk density, penetration resistance, in topsoil and subsoil for 2011.Although Guedes Filho et al (2010) stated that spherical mathematical is predominantly fitted model in soil science research, the exponential semivariograms were mostly the best fitting models to the experimental semivariances of the GWC, PR, BD and AgS in surface and subsurface soils in 2011.The models used in 2005 were also similar to that of 2011 except for GWC at 20-40 cm (exponential in 2011) and BD at 0-20 cm (exponential in 2011) and 20-40 cm (spherical in 2011) (Tables 4 & 5; Figures 2, 3 & 4).
The measurement error or micro-variability of a property which can't be detected with current scale of sampling causes the nugget variance (Co) (Isaaks & Srivastava 1989).Vieira (2000) stated that the nugget variance/effect (Co) reveals discontinuity of the semivariogram for distances shorter than the distance between samples (Guedes Filho et al 2010).Nugget variance values of soil physical attributes were all higher in 2011 compared the values of 2005 (Tables 4 & 5).
The degree of spatial dependence found for GWC at 0-20 cm and 20-40 cm depths (44.03 and 45.32%), BD (49.82 and 9.12%), PR (43.82 and 36.36%) and AsS (13.04 and 49.96%) were moderate according to the classification scheme of Cambardella et al (1994) (Table 5).Although nugget ratio values of BD at 20-40 cm depth and AgS at 0-20 cm depth indicated strong spatial dependencies, the models fitted were close to the pure nugget effect.The reliability of the estimated data is decreased with the variograms of pure nugget effect and/or very short ranges as occurred with BD for 20-40 cm and AgS for 0-20 cm (Table 5; Figures 3 & 5).

87
The strongly dependant variables may be controlled by intrinsic variations in soil characteristics such as texture and mineralogy (Cambardella et al 1994;Tekin et al 2011;Wang & Shao 2011).Surface AgS and subsurface BD showed strong, and all the other parameters showed moderate spatial dependency (Table 5).The nugget ratio (spatial dependency) values of GWC, BD and PR for surface soils indicated the existence of lower spatial dependency in 2011 as compared to that of 2005 (Tables 4 & 5).Extrinsic variations, such as fertilizer application and tillage, may control the variability of weakly spatially dependent parameters (Ayoubi et al 2007).In addition to this, differences in rotation may also attribute the decline in spatial dependency of a variable.
The range is another important parameter of interpretation of semivariograms and spatial variations in general, indicating the limit distance at which a sample point has influence over another point, i.e., the maximum distance up to which sample points are correlated.All points located within a circle with radius equal to the range can be used to estimate values with smaller spacing (Warrick & Nielsen 1980).Points located at distances beyond the range have no defined spatial dependence but are randomly distributed, behaving independently.Spatial dependence range values at surface varied from 29 m for AgS to 1524 m for BD, at 20-40 cm depth from 40 m for BD to 1182 m for GWC (Table 5) were increased in 2011.The range of BD at surface in 2011 was higher than the range of BD obtained for 20-40 cm that was similar prior to intensive rice cultivation (Tables 4 & 5).
Bulk density and PR are the major indicators of the soil compaction and increased with an increase in GWC (Table 1).Managing large areas as homogeneous, in spite of existing considerable spatial variability inherent to the soil, causes the appearance of distinct soil physical properties and distinct soil fertility zones (Wendroth et al 2003).The maps created are useful in identifying the limiting factors such as fields with high BD and PR within the study area.The BD values in study area ranged from 1.0 to 1.7 g cm -3 at 0-20 cm and 1.2 to 1.8 g cm -3 at 20-40 cm which indicated soil compaction in some of the fields.Similar to the case in BD, PR values particularly at 20-40 cm depth were ranged between 2.7 and 4.3 MPa.Although crop response to soil compaction depends on the interaction among crop, soil type, water content, and compaction degree (Lipiec & Simota 1994), the value of 2 MPa is accepted as restrictive for root growth (Taylor et al 1966).Thus, PR at subsurface and BD at surface and subsurface soils exceeded the threshold values in some of the fields in study area.The spatial distribution of PR and BD maps will give us opportunity to the soil compactions with the agricultural processes (Gunal et al 2012).
Wang & Shao (2011) demonstrated that the spatial distribution and spatial dependence of soil physical properties within a watershed were complex, and influenced by environmental factors as well as land use and topography.The evaluations based on spatial attributes of BD, PR and GWC and visual inspections of distribution maps created in 2005 (Özgöz et al 2009) and 2011 indicated significant spatial and temporal variations (Tables 1, 4 & 5; Figures 2, 3 & 4).Shifting to rice cultivation in the area is probably resulted in temporal and spatial variations of soil attributes evaluated.Dec et al (2011) also reported partial temporal changes in water content and PR, however they have not found any changes in spatial structure of the study area.
The soil's susceptibility to compaction in wet conditions is related to the declining of cohesion between particles.The menisci forces cause soil particles to lose contact between each other when pores are filled water which increases the soil compaction and decreases PR values.The decrease in PR was significantly reflected  (P<0.01) in negative correlation obtained between GWC and PR of surface soils.The correlation between GWC and BD was indicated the increase in compaction with significant positive correlation (P<0.01)(Table 3).In contrast to the wetting, water medisci forces are built up and increase PR with drying conditions (Unger 1996;Beutler et al 2005).At the time of sampling, soils were recently tilled which homogenized the spatial distribution of BD and PR within the study area.Soil management practices for rice cultivation are inherently carried out in wet conditions.Soil tillage practices in rice cultivation are mainly effective at 0-20 cm, thus soil depth below 20 cm is susceptible to severe soil compaction particularly under wet conditions.

Conclusions
Spatial and temporal variations of GWC, BD, PR and AgS in a Typic Ustifluvent soil were determined both classical and geostatistical methods.Reducing the number of crops in rotation homogenized the distribution of soil attributes in study area.The ranges known as spatial correlation distances of soil attributes were greater in 2011 demonstrating the homogenization effect of rice cultivation.Due to the intensive rice cultivation, GWC values were higher in 2011 as compared to 2005.Five years of rice cultivation did not significantly change the bulk densities of surface soils.However, bulk densities of subsurface soils significant (P<0.01)varied in time.Water ponding in rice cultivation probably destroyed the soil structure which increased the penetration resistance and caused a subsequent serious soil compaction.The lower coefficient of variations of AgS compared to the variations in major aggregate binding agents of organic matter

Table 2 -Mean differences and paired t-test results for initial (2005) and final (2011) data
** Significant at P<0.01; ns: non significant CV's of GWC in 2011 were lower than that