Influence of Land Use and Land Cover on Hydraulic and Physical Soil Properties at the Cerrado Agricultural Frontier

Western Bahia is one of the most active agricultural frontiers in the world, which raises concern about its natural resources conservation, especially regarding water availability. This study evaluated the influence of five different land uses and land covers on physical and hydraulic soil properties, and developed pedotransfer functions to derive regional hydraulic properties. Significant changes between physical and hydraulic soil properties under agricultural areas and under natural vegetation cover were found, reinforcing that agricultural activity may influence the soil water balance. Cerrado and Forest formation areas have higher infiltration rates (Ksat) compared to managed areas, with average values of 16.29 cm h−1, and 14.47 cm h−1, while irrigated croplands, rainfed croplands and pasture areas have much smaller infiltration rates, with Ksat equal to 3.01 cm h−1, 6.22 cm h−1 and 5.01 cm h−1, respectively. Our results suggest that the agriculture practices do not directly affect the vertical nature of hydrological flowpath, except in the case of intensive irrigated agriculture areas, where Ksat reduction can lead to erosive processes favoring organic matter losses, and decreases in productivity and soil quality. Impacts of land use change on hydraulic and physical soil properties are a reality in the Cerrado agriculture frontier and there is an urgent need to monitor how these changes occur over time to develop effective mitigation strategies of soil and water conservation.


Introduction
In the last decade, the rapid expansion of agribusiness in the Cerrado led to a new Brazilian agricultural frontier known as MATOPIBA (acronym for the states Maranhão, Tocantins, Piauí and Bahia), which have raised concerns about the natural resources conservation [1,2], especially regarding to the water availability [3][4][5].
In MATOPIBA, Western Bahia stands out by the agricultural expansion, representing 49% of the total agricultural area [6] with 1.8 million hectares in 2015, which is equivalent to an increase of 352% since 1985 [7,8].The irrigated croplands alone had an increase of 526% between 1985 and 2002, with a simultaneous decrease of Cerrado, seasonal forest and transition vegetation areas by 881,483 ha, 66,417 ha and 269,592 ha, respectively [9].In fact, the region also stands out for its high productivity records, reaching to 7.4 million tons of soybean, cotton and maize crops in the 2016/2017 harvest [10], confirming the potential of the agribusiness and the systematic dependence of the water availability and natural resources.Although Western Bahia is located above the Urucuia Aquifer and is drained by the Grande, Corrente and Carinhanha basins, there is a major concern about the regional water availability and the impacts of the agricultural activities on Cerrado biome [1,2,5].
The effects of land use change on soil physical proprieties are broadly known, especially when considered conversion of tropical forest to pasture or croplands [11][12][13].Changes in soil bulk density, penetration resistance, porosity, near-surface hydraulic conductivity [14][15][16], infiltrability and saturated hydraulic conductivity [13,17] are described as possible consequences of the land use change.
In Cerrado, these effects have been studied in the last decade [18][19][20], and received a greater focus due to the advancement of the agricultural frontier from Southern Amazonia to MATOPIBA [5,20].Soil compaction [21,22], erosion [23], and decreased permeability and water infiltration rates [22,24] are some effects of land use change that affects the physical and hydraulic soil properties causing damage on soil aeration, and soil water dynamics in Cerrado.However, most of these studies are based on few data and do not allow spatial generalizations.Remote sensing and modeling techniques have been used to broadly study the influence of land use on the water balance [5,[25][26][27], in the dynamic of vegetation [1,7] and climate [28][29][30][31], and show that the replacement of native vegetation by croplands alters in Cerrado the amount of water recycled to the atmosphere at a large-scale [5], affecting the regional climate dynamics [28][29][30][31].The results of these simulations, however, are very sensitive to the soil physical and hydraulic parameterizations in the model.
Thus, the knowledge of the hydraulic and physical soil properties contribute to development of powerful tools as hydrological and dynamic vegetation models, and remote sensing techniques to estimate the water recharge, the soil water availability, the cropland productivity, or the influence of extreme precipitation scenarios in the Cerrado agricultural frontier.Moreover, the hydraulic and physical soil properties are fundamental to make inferences about the soil quality and sustainability, allowing the development of alternatives that can prioritize the water and soil conservancy and preserve the remaining native vegetation while ensuring the increase of agricultural production.
The main goal of this study was to evaluate how land use and land cover affect the hydraulic and physical soil properties at the Cerrado agricultural frontier.In addition, we developed pedotransfer functions to derive hydraulic properties for use in dynamic vegetation models.

Study Area
The study area of this work is the Western Bahia located on the geological formation of the Urucuia Group (Upper Cretaceous), which is one of the main areas of agricultural expansion in the Cerrado biome.This region is drained by three important rivers, Grande, Correntes and Carinhanha, with an area of 131,168.59km 2 (Figure 1).Marked by the arrival of migrant rural producers from the southern region in 1990, Western Bahia is today one of the largest producers of soybean in MATOPIBA and the largest producer of cotton in Brazil.
In this region, the physiognomy vegetation is predominantly Cerrado stricto sensu, with predominance of the tree-shrub stratum [32].The soils are acidic with low fertility, but located at flat or mildly hilly areas.Thus, the use of high level technology and the chemical inputs allows the correction of fertility providing favorable conditions to expansion and intensification of agriculture.The soil granulometry is predominantly sandy and has medium texture, classified as Latossolos (57%), Neossolos (29.6%) and Cambissolos (7%), according to the Brazilian soil classification system [33].
The regional climate is tropical humid according to Köppen [32], and presents two well-defined seasons, dry (October to April) and rainy (May to September).The average annual temperature is 24 • C and the average annual rainfall is 1400 mm in the extreme west, gradually decreasing to 800 mm in the east.
Each class was sampled at 20 sampling points following the criteria: opening year of the agriculture or pasture area between 1990 and 2017, and logistical access to farms, including road conditions and permission to enter the farm by the farmer.In some cases, CDO and the FOR areas were sampled along the road between farms in natural vegetation areas rather than inside the farms.In total, 700 samples were collected considering 5 LULCCs × 20 sample points × 7 depth levels.At each sampling point, undisturbed samples were collected for 0-5, 5-10, 10-15, 15-20, 30, 50 and 70 cm of depth using an Uhland soil sampler, with metal cylinders of 100 cm 3 volume.The 100 samples of 0-5 cm of depth were used to calculate the saturated hydraulic conductivity (K sat ), soil water content at 10 kPa and 1500 kPa (θ f c , θ wp ), texture (coarse sand fraction ( f csa ), fine sand ( f f sa ), and silt ( f s ), total porosity (φ tot ), microporosity (φ m ), macroporosity (φ M ), and soil particle density (ρ s ), while all layers (0-5, 5-10, 10-15, 15-20, 30, 50 and 70 cm) were used to calculate the soil bulk density (ρ a ).

Measuring Methods
The K sat of soils was determined via the constant head permeameter method [34].The constant water flow that flowed through the soil sample was measured and applied into the equation for direct calculation of K sat .For textural analysis and particle density, samples were dried in the open air, and sieved in a 100-mesh sieve.Soil fractions ( f csa , f f sa , and f s ) were separated according to Ruiz [35], using the sieve method for the sand fraction (0.05-2 mm) and the pipette method to determine the silt (0.002-0.05 mm) and clay (<0.002 mm) fraction.The soil water retention curve (SWRC) was measured using a sand tension table at the matric potentials of 1 and −6 kPa, and a Richards chamber with porous plates, at matric potentials of 10, 50, 100, 500 and −1500 kPa [36].The centrifuge method was used to measure the soil water content at matric potential of 30 kPa.The SWRC was measured for 0-5 cm layer for all LULCC, and also for 5-20 cm layer in IRR areas.ρ a and φ tot analyses were carried out using methods described by CLAESSEN et al. [37].The φ m was determined by water content in a volumetric ring under 0.6 m of water column tension, and the soil φ M was estimated by the difference between total soil porosity and soil microporosity.

Adjustment of Soil Water retention Curve
The Campbell and Norman [38] soil model was fitted to soil water content data using Equation (1).
where ψ m is the soil water matric potential in kPa, ψ e is the air entry matric potential in kPa, θ is the volumetric water content cm 3 cm −3 , θ s is the saturated volumetric water content cm 3 cm −3 , and b is the empirical Campbell parameter, related to the particle size distribution.It is strongly dependent on soil texture [39] and is considered an index for soil pore-size distribution [40].This model was chosen due to its minimal set of parameters necessary to describe the soil hydraulic properties, favoring its implementation in regional and global scales, and it has been widely used in modeling studies [41,42].
A linear regression with log-transformed data was used in the equation log y = a + xb, to determine the ψ e and the empirical value of b Campbell parameter, where a is the intercept of the soil water retention in the log-log system and b represents the slope.
The linear least squares method was used to adjust the ψ e and b Campbell constant for each sample.Then, the raw data were also used to develop the pedotransfer functions.

Development of Pedotransfer Functions
The pedotransfer functions (PTFs) were developed using the soil physic proprieties and the SWRC for 100 undisturbed samples collected in the 0-5 cm soil layer.Initially, we tested the normal distribution for f csa , f f sa , f s , φ tot , φ M , φ m , ρ s and ρ a properties using Shapiro-Wilk statistic test considering α = 0.05.The properties that are normally distributed were used to develop the PTF through multiple regression, while the predictors not normally-distributed, such as ρ s , f f sa , f s , φ M and φ m , were excluded from analysis.
To develop and validate the pedotransfer function, the observed data for water retention curve data and physical properties were separated into two groups, where the first one contained 75% of the data (calibration group) and the second 25% of the data (validation group).The samples of calibration and validation groups were randomly chosen.Between these two groups the average of soil physical properties were not significantly different at α = 0.05 according to the Student t test, except for f csa in the forest formations areas FOR.However, for hydraulic proprieties ( ψ e and b), a significant difference was found between b Campbell parameter used in Calibration and Validation Group for CDO and FOR classes.
A multiple linear regression model was adjusted for each hydraulic Campbell parameters (log ψ e , b, K sat , θ f c , and θ wp ) using the soil physical properties measurements as predictors (Equation ( 2)): where y is one of the five hydraulic parameters (log ψ e , b, K sat , θ f c , and θ wp ), ρ a is the soil bulk density in g cm −3 , φ tot is the total soil porosity in cm 3 cm −3 , and f csa represents the coarse sand fraction in percent.The stepwise method with 5% of significance was chosen to select the most important variables for determination of y through the backward and forward mechanism.This stepwise method uses the Akaike criterion to eliminate collinear variables, excluding non-informative variables of the final model [43].

Data Analysis
We tested the homoscedastic and normality of residuals applying the Shapiro-Wilk test [44] for all soil physical properties, and used the Tukey Cramer method to compare the soil physical properties among LULCC classes.The soil physical properties that were normally distributed (i.e., f csa , φ tot , and ρ a , with p-values greater than 0.05 in the Shapiro-Wilk test) were applied to the analysis of variance and were used for the development of pedotransfer function.The t-test was used to compare the average of the 75% of soil properties used for the development and 25% used for the validation of pedotransfer function (0-5 cm layer), and for the comparisons of the differences for each soil layers among LULCC.All data analyses were carried out using R software [43].

Soil Physical and Hydraulic Properties
Most of the soil samples were sandy clay loams (40%), followed by 23.3% of sandy loams and 23.3% of loamy sands.The 13.4% of remaining samples were distributed among sandy clay, sand and silt clay classes (Figure 2).In the Brazilian soil classification system (SiBCS), these samples presented a similar distribution along the soil classes, with predominance of 40% of the samples in MeA (Médio argilosa), and 56% divided equally between MeAr (Médio Arenosa) and ArMe (Arenosa Média) classes, highlighting the predominance of sandy soils in the region (Figure 2).The average sand content ranged between 69.47% and 85.79% within the predominant soil classes, while the average clay content ranged between 11.29% and 26.87%.
The results of average soil bulk density (ρ a ) suggest a compacted soil surface for agriculture land use, ranging between 1.52 and 1.61 g cm −3 (Figure 3a and Table 1).Under natural land cover, FOR and CDO, the ρ a was lower and statistically different from agricultural areas with average ρ a ≤ 1.36 g cm −3 (Figure 3a and Table 1).Along the soil profile, the soil ρ a showed an increase trend from 0-5 cm to 30 cm for all LULCCs (Figure 4).This increase of ρ a in the subsurface layers is higher in areas under managed soil, than in FOR and CDO soils with values above 1.65 g cm −3 (Figure 4).
The particles density (ρ s ) ranged between 2.5 g cm −3 and 2.65 g cm −3 for all LULCCs, showing no statistical differences according to the Tukey Cramer test.The soil total porosity (φ tot ) ranged between 26% and 60%, with average of 43% for all samples collected (Table 1 and Figure 3b).The compaction pattern found in ρ a was also observed in φ tot , with reduction of the φ M in agriculture land use compared to CDO and FOR areas.In natural ecosystems, φ tot ranged between 26% and 57% for CDO and between 38% and 60% for FOR, while in agriculture systems the total soil porosity ranged between 30% and 53% (Figure 3e).Average CDO and FOR φ tot was greater than 45%, while under agriculture land use φ tot were smaller than 41% (Table 1).
The soil compaction decreased φ M in PAST and RAG areas, reducing the soil capacity to drain excess water after a heavy precipitation-lower K sat values.The reduction of φ M contributed to an increase of φ m , which may alter the soil aeration and roots growing conditions.The φ M in managed areas ranged between 9% and 12%, while φ m were above 25%.Although the values of φ M and φ m of the IRR class were not measured, significant differences were found for total soil porosity in relation to CDO and FOR (Figure 3e).Soils under natural vegetation cover had a higher infiltration rates (K sat ) compared to managed areas, with mean values of 16.29 cm h −1 for CDO, and 14.47 cm h −1 for FOR, while among the agriculture land uses, average infiltration rates were much smaller, ranging from 3.01 cm h −1 in irrigated croplands to 6.22 cm h −1 in RAG (Table 1 and Figure 3f).
All LULCCs showed low volumetric water content at the field capacity (−10 kPa), varying from 0.17 cm 3 cm −3 to 0.27 cm 3 cm −3 .Although LULCCs present different values of field capacity (θ f c ) and wilt point (θ wp ), the average difference between θ f c and θ wp for each LULCC is typically around 0.06 cm 3 cm −3 , highlighting the low water retention capacity for these soils (Table 1).
Table 1 shows all hydraulic and physical properties for the Campbell and Norman model.Soils under natural vegetation presented lower air entry potential than soils under agriculture or pasture, with average ψ e values equal to 0.71 kPa for CDO and 0.87 kPa for FOR.The highest values of air entry matric potential were found in soils under pasture (ψ e = 1.80 kPa) and irrigated agriculture (ψ e = 1.46 kPa), which are associated to higher compaction (Table 1 and Figure 3).The Campbell b parameter average values ranged between −4.06 and −5.30, which are similar to literature values for sandy, loamy sand and sandy loam soils.ψ e : soil air potential entry, kPa; b: Campbell parameter; ρ a : soil bulk density, g cm −3 ; K sat : saturated hydraulic conductivity, cm h −1 ; ρ s : soil particle density, g cm −3 ; φ tot : soil total porosity, fraction; θ f c : volumetric moisture at 10 kPa, cm 3 cm −3 ; θ wp : volumetric moisture at 1500 kPa, cm 3 cm −3 ; θ paw : volumetric moisture available to plants, cm 3 cm −3

Pedotransfer Functions
All PTFs used to estimate ψ e , b, K sat , θ f c , and θ paw for soils under PAST and RAG were significant at α = 0.05.The PTFs for IRR, CDO, and FOR showed significance only for a few parameters (Table 2).ψ e was significant in all LULCC, while b was significant only in CDO, PAST and RAG.The PTFs for K sat were not significant for any LULCCs, except PAST and RAG.θ f c PTFs were significant only for CDO, PAST and RAG, while θ wp , PTFs were significant for FOR, PAST and RAG (Table 2).The PTF performance to estimate of the ψ e was generally good, showing the maximum adjustment in FOR areas with r = 0.98, followed by PAST (r = 0.93), RAG (r = 0.84) and IRR (r = 0.84) (Table 2).The lowest correlation was found for CDO areas with r = 0.55.For this LUC, the average log ψ e used for calibration and validation were significantly different, which may have contributed to the lower performance in ψ e estimate.
For the Campbell b parameter, the performance of the PTF estimative showed poor correlation for CDO, IRR, and PAST, with correlations 0.36, −0.69 and 0.24, respectively.However, for RAG and FOR areas, the PTFs showed r > 0.70 agreement between the estimated and observed Campbell b parameter (Table 2).
The soil water content represented by θ f c , and θ wp also had a poor correlation for IRR and RAG agriculture LUs with r < 0.20.In CDO, FOR and PAST, the estimated parameters had a correlation >0.45.The worst correlations were for K sat for all LULCCs with negative r values.

Discussion
The sand fraction is predominant in the granulometric composition of the soils, the clay content varies between 10% and 41% and the silt fraction presents the lower values, typically below 4% (Figure 2).This granulometric composition is the result of the high rate of weathering during the genesis of these soils and the soil source material-in this case, the sandstones of the Urucuia Formation [45][46][47][48].
The general trends showed higher values of ρ a and lower values of φ tot in agriculture areas compared to native vegetation, revealing that soils under agricultural land use are slightly compacted in Western Bahia.This is in line with Fontana et al. [48], who found ρ a in agriculture areas in the region of Luis Eduardo Magalhães equal or above 1.52 g cm −3 .Cunha et al. [46] also found increased ρ a in different crop areas with different periods of cultivation, demonstrating that duration of land use also has an influence on soil structure.Naturally, the Cerrado Western Bahia soils present a cohesive sub-surface horizon [46,[48][49][50], making it more susceptible to compaction by grazing, mechanization and management applied at the soil surface.This transitional or subsurface horizon was observed and characterized by a slight increase in ρ a , generally within 10-30 cm [48].Indeed, our results show the natural cohesive sub-surface in the CDO and FOR areas (Figure 4b), and the intensification of the increase of ρ a in agriculture land uses.The ρ a in RAG was twice as high in the 10-15 cm layer in relation to the 0-5 cm surface layer.In IRR areas, there was also an increase in ρ a in these layers, although with less intensity than observed in RAG areas (Figure 4).
In Western Bahia, where the irrigated and rainfed croplands have a key role in the development of agriculture, this natural sub-surface horizon cohesion might be a concern factor for the maintenance of the rates of infiltration, the permeability and the water availability to the croplands and recharge of the aquifer.In general, a reduction of φ tot , and an increase of ρ a may be occurring due to a combination of agricultural implements used to remove the natural vegetation cover and the applied agriculture management, consequently reducing K sat values.
During the fieldwork, the farmers reported that, after removing the Cerrado, it is mandatory to use the conventional tillage system, with the plowing, complemented by the subsoilers and scarifiers, to break any physical impediments.According to our interviews, after the conversion to agriculture, most farmers adopt a planting cycle system, where for the four years after the removal of Cerrado, the no-till system is applied and in the fifth year the conventional tillage system is used, characterizing a mixed management.However, we do not have enough long-term data to claim that this management is the best practice to adopt in order to avoid cohesive sub-surfaces [51,52], but we evaluated the effects of surface and subsurface compaction in hydraulic properties in IRR areas (Table 3).In these areas, where there are two or more crops planted per year, and the traffic of machinery is more intense, ρ a at the 15-20 cm layer is 1.67 g cm −3 , significantly higher than at the surface (0-5 cm) layer (1.57g cm −3 ) (Table 3), which is a vertical pattern similar to all the other LULCCs (Figure 4).While the difference in ρ a between surface and subsurface layers is significant, it is not sufficient to influence the hydraulic parameters, such as K sat , φ tot , θ wp and θ paw (Table 3), although ψ e and θ f c are significantly different in the vertical.
In this study, the variability of K sat in the natural areas can be explained by the higher heterogeneity of the soils, and, consequently, different accumulation of organic matter, litter, tree density, the soil fauna, source material and root systems acting on the soils structures.Despite the soil compaction, the average K sat for agriculture areas can be considered higher when compared to K sat for latossolos under conventional tillage system in Goiás, K sat = 0.535 cm h −1 [53].For other agriculture areas, in the Cerrado biome, the K sat values presented in the literature range between 5.41 cm h −1 [20], and 15.47 cm h −1 [48].Likely, the use of rotation among maize, soybean, cotton and other croplands, in addition to the mixed management in Western Bahia, contribute to the maintenance of the high rates of K sat in agriculture land uses areas even with the presence of soil compaction.

Conclusions
This study analyzed hydraulic and physical soil properties in Western Bahia at a local scale and in different land uses and land covers.Significant changes were found between some soil properties in agricultural areas and natural vegetation cover, indicating that agricultural activity can influence the soil properties.
The agriculture land use increased soil bulk density at soil surface and subsurface, reducing the K sat by an order of magnitude in relation to Cerrado and Forest areas, and also decreasing the soil porosity.Despite the reduction, K sat still ranges between 30 mm h −1 and 62 mm h −1 , which are still considered high hydrological infiltration rates.Thus, our results reveal that in Western Bahia the agriculture land use areas do not affect directly the vertical nature of hydrological flowpath for visited areas, but, in the case of very intense precipitation events, the K sat reduction may lead to erosive processes favring nutrient and soil losses.In Western Bahia, however, the farmers are very interested in adopting sustainable practices that preserve the soil quality, investing in state-of-the-art technology, increasing the intervals of soil revolving and implementing crop rotation system.
Nonetheless, one must be careful to not extrapolate these results for all Western Bahia region, since we do not know how irrigated agriculture, rainfed agriculture, and pasture influence the hydraulic and physical properties in the long-term, and how the management (conventional tillage or no-tillage) can influence these properties in loamy/sand soils.
In the literature, several studies monitoring physical soil properties show that reduced tillage practices have promising results for soil moisture conservation and for crops growth.However, the number of scientific studies remains low for the Cerrado agricultural frontier.We emphasize that the impacts of land use change on hydraulic and physic soil properties are a reality in MATOPIBA and there is an urgent need to monitor how these changes occur over time in order to develop effective mitigation strategies of soil and water conservation.

Figure 1 .
Figure 1.Study area and location of the soil samples collected considering different land use and land cover in Western Bahia.

Figure 2 .
Figure 2. Average soil texture fractions for 0-5 cm layer in Western Bahia according to the USDA soil classification and to the Brazilian soil classification system-SiBCS.

Figure 3 .
Figure 3. Soil hydraulic and physical properties considering 0-5 cm depth for different LULCCs in Western Bahia.In the box plots, the lower limit of the box indicates the 25th percentile, the black line within the box marks the median, the red point within the box marks the mean, and the upper limit of the box indicates the 75th percentile.Bars above and below the box indicate the confidence interval.The samples are distributed in Cerrado areas (CDO), n = 22; Forest formations (FOR), n = 19; Irrigated agriculture (IRR), n= 20; Pasture (PAST), n = 21; and Rainfed agriculture (RAG), n = 20.Different letters means that averages are statistically different according to Tukey Cramer test at α = 0.05.

Table 1 .
Average soil physical parameters for 0-5 cm layer under different land use cover in the Western of Bahia used in this study.

Table 3 .
Hydraulic and soil physical properties for irrigated land use in different depths.−4.26 a 1.574 a 3.00 a 2.612 a 0.397 a 0.169 a 0.0858 a 0.0840 a 15-20 0.03 b −4.16 a 1.665 b 3.07 a 2.662 a 0.374 a 0.172 b 0.0902 a 0.0818 a ab Values significantly different according to the t Student test at α = 0.05 are followed by different letters.