Applied and Environmental Soil Science
Volume 2017 (2017), Article ID 3597416, 10 pages
https://doi.org/10.1155/2017/3597416

Research Article

Integrated Crop-Livestock Management Effects on Soil Quality Dynamics in a Semiarid Region: A Typology of Soil Change Over Time

J. Ryschawy,¹ M. A. Liebig,² S. L. Kronberg,² D. W. Archer,² and J. R. Hendrickson²

¹Université de Toulouse, INRA, INP-ENSAT, UMR 1248 AGIR, 31324 Castanet-Tolosan, France
²USDA-ARS, Northern Great Plains Research Laboratory, P.O. Box 459, Mandan, ND 58554-0459, USA

Correspondence should be addressed to M. A. Liebig; mark.liebig@ars.usda.gov

Received 13 March 2017; Accepted 22 May 2017; Published 22 June 2017

Academic Editor: Teodoro M. Miano

Copyright © 2017 J. Ryschawy et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Integrated crop-livestock systems can have subtle effects on soil quality over time, particularly in semiarid regions where soil responses to management occur slowly. We tested if analyzing temporal trajectories of soils could detect trends in soil quality data which were not detected using traditional statistical and index approaches. Principal component and cluster analyses were used to assess the evolution in ten soil properties at three sampling times within two production systems (annually cropped, perennial grass). Principal component 1 explained 33% of the total variance of the complete dataset and corresponded to gradients in extractable N, available P, and C : N ratio. Principal component 2 explained 25.4% of the variability and corresponded to gradients of soil pH, soil organic C, and total N. While previous analyses found no differences in Soil Quality Index (SQI) scores between production systems, annually cropped treatments and perennial grasslands were clearly distinguished by cluster analysis. Cluster analysis also identified greater dispersion between plots over time, suggesting an evolution in soil condition in response to management. Accordingly, multivariate statistical techniques serve as a valuable tool for analyzing data where responses to management are subtle or anticipated to occur slowly.

1. Introduction

Integrated crop-livestock systems (ICLS) are recognized globally for their contributions to improve agricultural sustainability [1–4]. Recent emphasis on conservation agriculture, climate-smart agriculture, and sustainable intensification has underscored the potential role of ICLS to create more productive and resilient agricultural systems [5, 6]. An inherent emphasis on multiple enterprises makes ICLS well-suited for future growing conditions, where production synergies between enterprises can serve to enhance adaptability to increasingly variable weather and market conditions, while concurrently minimizing input costs [7, 8].

Integrated crop-livestock systems can be more management and labor intensive than single-enterprise production systems [4]. Accordingly, producers need to know if investments in implementing ICLS translate into improvements in key response metrics. Relevant information would utilize metrics characterizing agroecosystem sustainability and, ideally, do so quickly following adoption of management practices to assess trajectory over time [9]. Such information could serve to justify investments in ICLS or, if necessary, provide quantifiable guidance for adjusting management to more effectively meet production, economic, and/or environmental goals.

Soil quality serves as an important response metric to agroecosystem management given its foundational role in affecting agricultural and environmental outcomes through impacts on ecosystem services [10]. Formally defined, soil quality refers to the capacity of soil to function [11]. Accordingly, soil quality contributes to numerous ecosystem services within agricultural landscapes (e.g., water retention/filtration, climate regulation, and biodiversity conservation) [12].

The central role of soil to deliver key ecosystem services has served as a focal point for the development of numerous assessment methods. Over the past 30 years, soil quality assessments have evolved from the evaluation of a minimum dataset of soil properties [13] to elaborate indices involving scoring functions based on management goals across multiple scales [14–16]. As agricultural practices are intensified to meet anticipated production and energy needs [17], soil quality assessments will be essential to ensure deployed practices are sustainable.

Previous soil quality evaluations of ICLS have focused on treatment effects on single soil physical, chemical, or biological properties or indices of multiple properties [4, 18–21]. Outcomes from these evaluations have helped understand soil responses to ICLS, as results are typically easy to interpret and convey, thereby improving understanding by clientele. Despite this benefit, there are unique attributes associated with ICLS that may make traditional assessment of soil property dynamics, whether by single properties or an index, problematic. Soil properties can change slowly in ICLS, particularly in semiarid regions where soil responses to management occur slowly [22]. Accordingly, long-term research may be necessary to assess ICLS effects on soil quality, yet short-term funding cycles may make it difficult to conduct such research. Secondly, traditional Soil Quality Index methods “merge” outcomes from multiple soil properties to provide an aggregated numerical rating [11, 16]. In doing so, traditional indexing methods may mask nuances inherent to the data and, thus, miss important spatiotemporal trends.

Given this context, we sought to analyze near-surface soil quality indicator data collected over a six-year period from an ICLS in a semiarid region using principal component analysis (PCA). These data were previously analyzed using traditional statistical and index approaches to evaluate soil property responses to residue management, frequency of hoof traffic, season, and production system [18, 23]. We hypothesized a reanalysis of the same data with multivariate statistical techniques would identify trends in soil property dynamics not previously observed.

2. Materials and Methods

2.1. Site and Treatment Description

The ICLS experiment was located at the USDA-ARS Northern Great Plains Research Laboratory near Mandan, ND, USA (46°46′ W, 100°54 N). The experimental site is within the temperate steppe ecoregion of the USA, with a semiarid climate characterized by long, cold winters and short, hot summers. Mean annual precipitation and temperature at the site are 414 mm and 4°C, respectively, and the average frost-free period is 131 days. Gently rolling uplands (0 to 3% slope) characterize the topography of the experimental site, and predominant soil types are Temvik-Wilton silt loams (fine-silty, mixed, superactive, frigid Typic, and Pachic Haplustolls) [24].

A thorough treatment description can be found elsewhere [18]. Briefly, two 6.0 ha crested wheatgrass [Agropyron desertorum (Fisch. ex. Link) Schult.] pastures were sprayed with glyphosate [N-(phosphonomethyl) glycine; 0.7 kg a.i. ha⁻¹] twice in mid-May and converted to an annual cropping sequence of oat/pea (Avena sativa L./Pisum sativum L.), triticale/sweet clover (Triticum aestivum x Secale cereale/Melilotus officinalis L.), and corn (Zea mays L.). Beginning in 2007, the crop sequence was changed to oat/alfalfa (Medicago spp.)/hairy vetch (Vicia villosa Roth)/red clover (Trifolium pratense L.), brown midrib sorghum-sudangrass (Sorghum bicolor L. Moench)/sweet clover/red clover, and corn also using no-till planting techniques. Each phase of the three-year cropping sequence was present in both pastures, which were used as replicates. Management decisions related to seeding, fertilizer, and weed control by herbicides followed recommended practices by area producers.

The oat and triticale/sorghum crop mixtures were harvested for grain from mid-August to early-September with the straw spreader removed from the combine, which created a swath of crop residue for winter grazing. The corn was swathed for forage in mid- to late-September. Each crop strip was split into three residue management treatments that included no residue removal (Control; 0.05 ha), residue removal with a baler (Hayed; 0.05 ha), and residue removal by grazing with livestock (Grazed; 1.69 ha). The Control and Hayed treatments were randomly assigned within each crop strip. For the Grazed treatment, the swathed crop residues from the cropping sequence represented winter forage for ten 4–6-year-old nonlactating bred Hereford cows, due to calve in late March. Plot areas within a crop phase were 27 m × 23 m (0.06 ha) for the Control and Hayed treatments and 54 m × 312 m (1.68 ha) for the Grazed treatment. Stocking rate within a crop phase for the Grazed treatment was 0.2 ha cow⁻¹ (5.9 cows ha⁻¹), corresponding to winter grazing practices used by local producers.

Grazing commenced in mid-November and ended in mid-February, with the oat crop mixture grazed first, triticale/sorghum crop mixture second, and corn last. Access to crop swaths was controlled using electric fences oriented at right angles to the swaths. Fences were moved daily to provide access to fresh forage. A shelter and “frost-free” water fountain were located at the end of each pasture within the Grazed treatment.

Two perennial grass pastures were used for comparison to the annually cropped treatments. The pastures, each 6.0 ha, were composed of a mixture of native and introduced cool-season perennial grasses. Similar to Grazed treatments, swathed grass from the perennial pastures was used as winter forage from 1999 through 2002 for the same number of cows managed in the annually cropped treatments. Within each perennial pasture a nongrazed strip was split into Hayed and Control treatments as outlined above. Drought-induced limitations in forage in 2002 and 2003 resulted in the perennial grass treatments being hayed but not grazed in 2003. From 2004 through 2008, the perennial grass treatment was not swathed but lightly grazed with ten Hereford or Angus cows from mid-October to mid-January.

2.2. Soil Sample Collection, Processing, and Analyses

Soil samples were collected in April of 2002, 2005, and 2008 after the swathed crops had been grazed but not replanted. Within the annually cropped pastures, samples were collected from nine subplots (three per crop phase) in the Control and Hayed treatments, oriented randomly in each treatment but between crop rows. Samples from the Grazed treatment were collected from two transects differing in frequency of hoof traffic, also between crop rows. Nine subplots (three per crop phase) were established in each transect perpendicular to crop swaths approximately 100 m (representing high-traffic; HT) and 200 m (representing low-traffic; LT) from the shelter and water source. Sampling protocol for the western wheatgrass pastures followed that in the annually cropped pastures, with the exception of fewer subplots in each treatment (three in the Control, Hayed, and Grazed hoof traffic transects). Within each subplot, six soil cores were collected from the 0 to 7.5 cm depth using a 35 mm (i.d.) step-down probe and composited. Each sample was saved in a double-lined plastic bag, placed in cold storage at 5°C, and analyzed within 6 wk of collection.

Soil samples were dried at 35°C for 4 d and then ground by hand to pass a number 10 (2.0 mm) sieve. Identifiable plant material (>2.0 mm diameter, >10 mm length) was removed during sieving. Electrical conductivity and pH were estimated from a 1 : 1 soil-water mixture [25, 26]. Soil NO₃-N and NH₄-N were determined from 1 : 10 soil-KCl (2 M) extracts using cadmium reduction followed by a modified Griess-Ilosvay method and indophenol blue reaction [27]. Plant-available soil P was estimated by bicarbonate extraction [28]. Potentially mineralizable N was estimated from the NH₄-N accumulated after a 7 d anaerobic incubation at 40°C [29]. Total soil C and N were determined by dry combustion [30], and as soil pH was <7.2 within the surface 7.5 cm depth, total soil C was considered equivalent to SOC. All data were expressed on an oven-dry basis prior to data analyses.

2.3. Data Organization

Prior to statistical analyses, data were organized to facilitate characterization of all sampled locations as well as sampled locations specific to the annually cropped treatments. Data group 1 (S1) was composed of 10 soil properties from 40 plot locations in both annually cropped and perennial grass treatments (Table 1). Within S1, there were 24 locations in the annually cropped treatments and 16 locations within the perennial grass treatments. Data group 2 (S2) focused on the 24 locations within the annually cropped treatments, with an equivalent number of grazed and ungrazed locations (12 each). For both S1 and S2 data groups, data were included from each sampling time: 2002, 2005, and 2008.

Table 1: Description of soil properties used for data analysis.

Variable
Acronym
Units
Minimum value
Maximum value

Soil bulk density
SBD
Mg m⁻³
0.80
1.61

Electrical conductivity
EC
dS m⁻¹
0.10
0.31

Soil pH
PH
−log[H⁺]
5.48
7.65

Soil nitrate
NO3N
mg kg⁻¹
0
44

Soil ammonium
NH4N
mg kg⁻¹
1
26

Available phosphorus
P
mg kg⁻¹
1
29

Potentially mineralizable N
PMN
mg kg⁻¹
19
166

Soil organic carbon
SOC
g kg⁻¹
13.2
41.9

Total nitrogen
TN
g kg⁻¹
1.2
3.5

Soil carbon to nitrogen ratio
CNRATIO

10.3
12.5

2.4. Statistical Analyses

Individual plot dynamics were quantified by combining a series of multivariate analyses using a method developed by Dolédec and Chessel [31] and adapted by García-Martínez et al. [32] and Ryschawy et al. [33]. The method distinguished impacts of external (the system environment) and internal (system structure and functioning) factors, while considering the evolution of individual plots over time. We applied this statistical methodology to S1 to analyse plot trajectories of near-surface soil quality indicators (Table 1).

To analyse plot trajectories, we used a particular type of PCA (Principal Component Analysis) called Within-Class Analysis (WCA) to consider evolution of individuals (plots) while not considering the general date effect on different plots. We thus considered the deviation of each soil property in a plot from the average across plots, considering time period effect. The first step of the WCA consisted in standardizing the data. For each soil property (column), data were centered by subtracting the mean of the column to the values and scaled by dividing the values of each column by their standard deviation (with a column mean of 0 and a standard deviation of 1). Within-Class Analysis used an R statistical software package dedicated to Analysis of Environmental Data (Exploratory and Euclidean Method; ade4 package) [34]. Within-Class Analysis used orthogonal PCAs where factors were used as covariables to partition rows (time was used as the main factor to partition rows considering each date separately for plots). The analysis characterized differences between plot locations once a general common trend over time was eliminated. We considered the three sampling times (2002, 2005, and 2008) as the factor partitioning rows.

Main factors within data group S1 were summarized by PCA. We then performed a WCA resulting in a data frame W1 containing transformed data for each of the 40 plot locations at the three different dates considered, as the deviation between individual and date-class mean coordinates on within-analysis axes. The results of the WCA therefore represented interplot location trajectories.

To summarize soil trajectories for each individual plot location for all dates, we performed a final PCA on the WCA individuals considered at each date. Per the Kaiser criterion, we thus considered the main factors with eigenvalues >1 as summarizing observed changes. A Hierarchical Cluster Analysis (HCA) (with squared Euclidean distance and Ward’s aggregation method) was carried out on the four main factors of the PCA to establish a typology of plot locations according to their change over time.

To analyse changes in plot locations elucidated by HCA, we calculated means and standard deviation of variables for each measured variable. All statistical analyses were performed using R 3.1.0 software [34]. We repeated this statistical methodology on data group S2 to evaluate soil property dynamics based on grazed (Grazed) and ungrazed (Control, Hayed) areas within the 24 plot locations assigned to annually cropped treatments.

3. Results

3.1. Analysis of All Sampled Locations

The first four principal components accounted for 82.3% of the variability based on WCA of the whole sample (data group S1). Associations between soil properties were evident following the WCA, as three distinct groupings were observed graphically (Figure 1). Soil nitrate (NO3N), ammonium (NH4N), and available phosphorus (P) appeared linked, as did soil carbon to nitrogen ratio (CNRATIO), soil organic C (SOC), total nitrogen (TN), and electrical conductivity (EC). Soil pH (PH) and potentially mineralizable N were also linked. Only soil bulk density (SBD) was not grouped with another soil property.

Figure 1: Spatial representation of measured variables on the first factorial map of the principal component analysis (PCA): SBD, soil bulk density; EC, electrical conductivity; PH, soil pH; NO3N, soil nitrate; NH4N, soil ammonium; P, available P; PMN, potentially mineralizable N; SOC, soil organic carbon; TN, total nitrogen; CNRATIO, soil carbon to nitrogen ratio.

Contributions of measured soil properties to the four principal components are provided in Table 2. Principal component 1 explained 33% of the variability and corresponded to gradients of available N and P (NO3N, P) and CN ratio levels. Principal component 2 explained 25.4% of the variability and corresponded to gradients of PH, SOC, and TN. Principal component 3 explained 13.6% of the variability and corresponded gradients in EC and PMN. Principal component 4 explained 10.3% of the variability and was primarily influenced by EC and NH4N.

Table 2: Principal component scores of measured variables and proportion of variability explained within each axis for data group 1 (S1).

Axis 1
Axis 2
Axis 3
Axis 4

Soil bulk density
0.59
0.42
0.27
0.17

Electrical conductivity
−0.23
−0.32
0.62
0.53

Soil pH
−0.34
0.74
0.42
0.08

Soil nitrate
0.71
−0.53
0.16
0.30

Soil ammonium
0.56
−0.15
0.50
0.54

Available P
0.65
−0.30
0.33
−0.41

Potentially mineralizable N
−0.41
0.48
0.59
0.20

Soil organic carbon
−0.58
−0.65
0.46
−0.04

Total nitrogen
−0.45
−0.70
−0.50
−0.09

C : N ratio
−0.78
0.12
−0.04
0.20

Variability explained
33.0%
25.4.%
13.6%
10.3%

Dispersion between plots showed a trend toward greater dissimilarity in soil properties over time (Figure 2). Based on standard deviation of mean values, plot differences were closer in composition in 2002 and 2005 than in 2008, with the latter showing substantial dispersion at the end of the 6 yr study period (Table 3).

Table 3: Mean values (±1 standard deviation) for soil properties in 2002, 2005, and 2008.

Variable
2002
2005
2008

Soil bulk density (Mg m⁻³)
1.05 ± 0.08
1.09 ± 0.10
1.12 ± 0.11

Electrical conductivity (dS m⁻¹)
0.24 ± 0.04
0.30 ±0.06
0.41 ± 0.10

Soil pH (−log[H⁺])
6.12 ± 0.38
6.23 ± 0.38
6.32 ± 0.55

Soil nitrate (mg kg⁻¹)
4 ± 3
5 ± 4
8 ± 9

Soil ammonium (mg kg⁻¹)
2 ± 1
2 ± 0
5 ± 3

Available phosphorus (mg kg⁻¹)
8 ± 3
5 ± 4
10 ± 6

Potentially mineralizable N (mg kg⁻¹)
74 ± 22
86 ± 24
86 ± 34

Soil organic carbon (g kg⁻¹)
30.1 ± 3.9
26.6 ± 3.6
29.6 ± 5.7

Total nitrogen (g kg⁻¹)
2.7 ± 0.3
2.3 ± 0.3
2.7 ± 0.5

C : N ratio
11.2 ± 0.4
11.4 ± 0.4
11.1 ± 0.4

Figure 2: Spatial representation of three groups of measured variables based on sampling year (2002, 2005, and 2008).

3.2. Analysis of Production Systems

Cluster analysis found distinct groupings based on production system, as annually cropped treatments and perennial grasslands were clearly distinguished graphically (Figure 3). Such trends were evident by comparing soil properties in each treatment between 2002 and 2008 (Table 4). Plots under annual cropping had a greater accumulation of mineral nutrients (NO3N, NH4N, and P) compared to perennial grassland, whereas perennial grasslands had greater accumulation of organic matter (SOC, TN), increased mineralizable N (PMN), and a lower dispersion between plots compared to annually cropped treatments.

Table 4: Mean values of near-surface soil properties within annually cropped and perennial grassland integrated crop-livestock systems in 2002, 2005, and 2008 [18].

Year/change
Annual crop
Perennial grassland

Soil bulk density (Mg m⁻³)

2002

1.02 ± 0.09

2005
1.07 ± 0.12
1.11 ± 0.04

2008
1.16 ± 0.12
1.07 ± 0.06

+0.09
+0.05

Electrical conductivity (dS m⁻¹)

2002
0.24 ± 0.03
0.25 ± 0.04

2005
0.30 ± 0.06
0.29 ± 0.05

2008
0.43 ± 0.12
0.38 ± 0.07

ΔEC
+0.19
+0.13

Soil pH (−log[H⁺])

2002
5.89 ± 0.21
6.46 ± 0.34

2005
6.00 ± 0.33
6.58 ± 0.27

2008
5.95 ± 0.33
6.88 ± 0.24

ΔPH
+0.05
+0.43

Soil nitrate (mg kg⁻¹)

2002
6 ± 3
1 ± 1

2005
7 ± 3
1 ± 1

2008
13 ± 8
1 ± 1

ΔNO3N
+6
0

Soil ammonium (mg kg⁻¹)

2002
2 ± 1
2 ± 1

2005
2 ± 1
3 ± 1

2008
6 ± 4
4 ± 1

ΔNH4N
+4
+2

Available P (mg kg⁻¹)

2002
7 ± 2
9 ± 4

2005
8 ± 4
2 ± 2

2008
14 ± 5
5 ± 1

ΔP
+7
−4

Potentially mineralizable N (mg kg⁻¹)

2002
62 ± 14
92 ± 20

2005
70 ± 14
109 ± 12

2008
66 ± 24
117 ± 21

ΔPMN
+4
+25

Soil organic C (g kg⁻¹)

2002
29.5 ± 3.8
30.9 ± 4.2

2005
27.1 ± 3.5
26.0 ± 3.8

2008
28.1 ± 6.0
31.9 ± 4.4

ΔSOC
−1.4
+1.0

Total N (g kg⁻¹)

2002
2.6 ± 0.3
2.7 ± 0.3

2005
2.4 ± 0.3
2.2 ± 0.3

2008
2.6 ± 0.5
2.8 ± 0.3

ΔTN
0.0
+0.1

C : N ratio

2002
11.2 ± 0.4
11.4 ± 0.5

2005
11.3 ± 0.4
11.7 ± 0.4

2008
10.9 ± 0.3
11.3 ± 0.3

ΔCN
−0.3
−0.1

deviations are given for each variable in 2002, 2005, and 2008. reflects the difference in soil property between 2008 and 2002.

Figure 3: Spatial representation of two groups of soil trajectory types following principal component analysis of the complete dataset.

3.3. Analysis of Grazing Effects under Annual Cropping

In the analysis of grazing effects (data group S2), 83.4% of the variability was explained based on WCA. Principal component 1 explained 29.8% of the variability and corresponded to gradients in total carbon and nitrogen accrual rates (SOC, TN, and CNRATIO), whereas principal component 2 explained 24.7% of the variability and corresponded to gradients of mineral nutrient accrual rates (NO3N, NH4N, and P) (data not shown).

Hierarchical cluster analysis revealed four clusters that were differentiated by grazing and treatment replication (Figure 4). Treatment replications were distinguished graphically along the horizontal plane, with one replication (Field A) left of the vertical plane and the other replication (Field B) right of the vertical plane. Differences in soil properties between treatments were limited to Field B, where soil pH was lower and available P was higher under grazing (Table 5).

Table 5: Mean values of near-surface soil properties within annually cropped fields according to grazing treatment.

Variable
Field A
Field B

Grazed
Ungrazed
value
Grazed
Ungrazed
value

Soil bulk density (Mg m⁻³)

1.04 ± 0.05
0.68
1.18 ± 0.12
1.15 ± 0.06
0.46

Electrical conductivity (dS m⁻¹)
0.22 ± 0.05
0.24 ± 0.04
0.26
0.20 ± 0.04
0.22 ± 0.05
0.40

Soil pH (−log[H⁺])
5.76 ±0.16
5.78 ± 0.15
0.70
5.99 ± 0.13
6.27 ± 0.24
<0.01

Soil nitrate (mg kg⁻¹)
9 ± 4
9 ± 4
0.64
10 ± 10
6 ± 4
0.13

Soil ammonium (mg kg⁻¹)
3 ± 1
3 ± 2
0.67
4 ± 6
3 ± 1
0.31

Available phosphorus (mg kg⁻¹)
9 ± 4
7 ± 4
0.29
12 ± 5
9 ± 4
0.04

Potentially mineralizable N (mg kg⁻¹)
61 ± 16
62 ± 13
0.83
75 ± 20
65 ± 20
0.16

Soil organic carbon (g kg⁻¹)
30.4 ± 5.5
31.2 ± 2.0
0.54
26.5 ± 3.1
24.9 ± 4.0
0.18

Total nitrogen (g kg⁻¹)
2.7 ± 0.5
2.8 ± 0.2
0.34
2.4 ± 0.3
2.3 ± 0.3
0.15

C : N ratio
11.4 ± 0.5
11.2 ± 0.3
0.37
10.9 ± 0.2
10.9 ± 0.3
0.99

standard deviation; value associated with mean comparison within a field between grazed and ungrazed treatments. Field A corresponds to clusters 1 and 2, while field B corresponds to clusters 3 and 4 (Figure 4).

Figure 4: Spatial representation of four groups of soil trajectory types following principal component analysis of annually cropped treatments.

4. Discussion

Traditional statistical and index-based approaches for assessing soil quality may provide a limited characterization of soil property dynamics [35]. Multivariate statistical analyses such as PCA provide a useful method for screening a diverse collection of functionally relevant soil properties to identify potential data trends [36]. These techniques provide an objective approach to extract and weight information in complex datasets [37], thereby providing a potentially valuable tool for researchers seeking to link the status of soil properties, soil function, and agroecosystem management.

Analysis of whole sample dynamics reflected the predominant influence of soil chemical properties on total variance in the dataset. Soil properties most influenced by additions of N and P through fertilization (e.g., NO3N and P) accounted for substantial variance between plots. Annual change in nutrient-related properties would be expected to vary based on weather-driven factors affecting plant growth and nutrient loss [38]. Electrical conductivity, a relative measure of the total quantity of ions in soil solution, also accounted for an appreciable amount of variance in the dataset.

Previous PCA evaluations characterizing soil quality have found organic matter-related properties to account for a large amount of variance in data [35, 36]. A similar outcome was observed in this evaluation, as SOC, TN, and CNRATIO accounted for a moderate amount of variation among soil properties across principal components. The limited measurement timespan (6 yr), coupled with the inherent difficulty in detecting change in properties with a large background signal due to previous management as a perennial grassland [39], may have compromised the usefulness of organic matter-related properties in this evaluation. In general, soil properties associated with organic matter are slow to change in semiarid agroecosystems, sometimes requiring a decade or longer for management effects to be detected [22].

Despite limited change in some soil properties in this evaluation, the aggregate visual assessment of data dispersion proved useful in detecting a general evolution of the dataset. Such dispersion was confirmed upon inspection of standard deviation values for soil properties within each year, as mean values of all soil properties, except CNRATIO, had the greatest standard deviation in the last year of the evaluation. Coupled graphical and tabular data characterizations seem particularly useful for datasets like the one analyzed in this evaluation, as visual assessments infer possible data trends generally, while tabular assessments provide confirmation specifically.

Outcomes from the cluster analysis of the complete dataset showed data groupings based on production system, with perennial grass and annually cropped treatments separated visually. Such separation was evident with cluster analysis despite limited observed differences in mean values of soil properties between production systems (Table 4) [18]. Moreover, multivariate analysis and graphical representation of these previously analyzed data clearly showed limitations associated with the aggregated Soil Quality Index (SQI) approach to assessing soil quality, as differences in the overall SQI between production systems were not observed in any year [18]. From a soil management perspective, outcomes from the cluster analysis underscored the value of perennial grass systems to enhance soil quality through increased accrual of soil carbon and nitrogen, while minimizing levels of nutrients susceptible to loss. Such attributes are important features of sustainable agricultural systems [40] and are central to agroecosystem resilience to climate-induced stressors and improvements in environmental quality from agricultural production [41, 42].

In contrast to whole production system analyses, outcomes were less distinct when data were restricted to annually cropped fields with and without grazing. Cluster analysis revealed an overarching effect of replication on data groupings, of which grazed and ungrazed plots lacked homogeneity. It is possible that the limited number of data points, coupled with a need for greater time-in-treatment to resolve management effects on soil properties, resulted in indeterminate responses to PCA. It is important to note that no differences in soil properties were found between grazed and ungrazed systems in any year when the same data were analyzed using traditional statistical analyses [18]. Future evaluations of grazed and ungrazed treatments in both annually cropped and perennial grassland systems are warranted given the important role of livestock to affect biomass productivity through influences on soil condition [43].

Collectively, statistical methods used in this evaluation contribute to an ongoing evolution of soil quality assessment. As the interpretation of individual soil properties can be complicated by multiple response functions [16], use of aggregate approaches, such as PCA and HCA, potentially offer novel insights into how soil responds to management [14]. Use of statistical approaches capable of resolving changes in soil condition in the near-term can be particularly useful in regions where alterations in individual soil properties require a decade or more to be detected [22]. Similarly, these approaches may be valuable in detecting subtle changes in soil induced by systemic drivers (i.e., climate). The capacity to detect shifts in soil condition quickly under the context of anticipated climate change will be increasingly valuable, particularly for dryland agriculture [44]. Timely management interventions will be essential in the future to mitigate soil degradation and maintain soil function.

5. Conclusions

Use of multivariate statistical analyses revealed trends in data not previously detected using traditional statistical and index approaches for assessment of soil quality. Principal component and hierarchical cluster analyses provided a helpful means to visually discriminate between annually cropped and perennial grass treatments, while concurrently identifying a clear trend toward greater dispersion in the data over time. Based on the typology used in this evaluation, apparent multivariate statistical techniques represent a valuable tool for analysis of data collected at different time periods where changes in response variables are anticipated to occur slowly.

Application of these techniques to datasets in other semiarid regions may contribute toward identifying trajectories in soil responses to applied treatments in the near-term. Such information could justify research expenditures, and therefore continuation of experimental treatments until soil responses can be verified with traditional statistical approaches and/or indices. Conversely, these techniques could be used to inform possible changes in applied treatments if visual discrimination is not detected.

Abbreviations

CNRATIO: | Soil carbon to nitrogen ratio
EC: | Electrical conductivity
ICLS: | Integrated crop-livestock systems
NH4N: | Soil ammonium
NO3N: | Soil nitrate
P: | Available P
PCA: | Principal component analysis
PH: | Soil pH
PMN: | Potentially mineralizable N
SBD: | Soil bulk density
SOC: | Soil organic carbon
TN: | Total nitrogen
WCA: | Within-Class Analysis.

Disclosure

The US Department of Agriculture (USDA) prohibits discrimination in all its programs and activities on the basis of race, color, national origin, age, disability, and, where applicable, sex, marital status, family status, parental status, religion, sexual orientation, genetic information, political beliefs, and reprisal, or because all or part of an individual’s income is derived from any public assistance program (not all prohibited bases apply to all programs.). USDA is an equal opportunity provider and employer. Mention of commercial products and organizations in this manuscript is solely to provide specific information. It does not constitute endorsement by USDA-ARS over other products and organizations not mentioned.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

The authors gratefully acknowledge the contributions of Justin Hartel, Marvin Hatzenbuhler, and Robert Kolberg for planting and swathing crops, Clay Erickson, Gordy Jensen, and Curt Klein for managing cattle, and John Bullinger, Jason Gross, and Holly Johnson for assisting with soil sampling, processing, and analyses.

References

G. Lemaire, A. Franzluebbers, P. C. D. F. Carvalho, and B. Dedieu, “Integrated crop-livestock systems: Strategies to achieve synergy between agricultural production and environmental quality,” Agriculture, Ecosystems and Environment, vol. 190, pp. 4–8, 2014.

J. Ryschawy, A. Joannon, J. P. Choisis, A. Gibon, and P. Y. Le Gal, “Participative assessment of innovative technical scenarios for enhancing sustainability of French mixed crop-livestock farms,” Agricultural Systems, vol. 129, pp. 1–8, 2014.

J. R. Hendrickson, J. D. Hanson, D. L. Tanaka, and G. Sassenrath, “Principles of integrated agricultural systems: Introduction to processes and definition,” Renewable Agriculture and Food Systems, vol. 23, no. 4, pp. 265–271, 2008.

M. P. Russelle, M. H. Entz, and A. J. Franzluebbers, “Reconsidering integrated crop-livestock systems in North America,” Agronomy Journal, vol. 99, no. 2, pp. 325–334, 2007.

K. L. Steenwerth, A. K. Hodson, A. J. Bloom et al., “Climate-smart agriculture global research agenda: scientific basis for action,” Agriculture & Food Security, vol. 3, no. 11, pp. 1–39, 2014.

J. A. Foley, N. Ramankutty, K. A. Brauman et al., “Solutions for a cultivated planet,” Nature, vol. 478, no. 7369, pp. 337–342, 2011.

R. J. Wilkins, “Eco-efficient approaches to land management: A case for increased integration of crop and animal production systems,” Philosophical Transactions of the Royal Society B: Biological Sciences, vol. 363, no. 1491, pp. 517–525, 2008.

J. D. Hanson, M. A. Liebig, S. D. Merrill, D. L. Tanaka, J. M. Krupinsky, and D. E. Stott, “Dynamic cropping systems: Increasing adaptability amid an uncertain future,” Agronomy Journal, vol. 99, no. 4, pp. 939–943, 2007.

C. Bockstaller, L. Guichard, O. Keichinger, P. Girardin, M.-B. Galan, and G. Gaillard, “Comparison of methods to assess the sustainability of agricultural systems. A review,” Agronomy for Sustainable Development, vol. 29, no. 1, pp. 223–235, 2009.

D. A. Robinson, N. Hockley, E. Dominati et al., “Natural capital, ecosystem services, and soil change: why soil science must embrace an ecosystems approach,” Vadose Zone Journal, vol. 11, no. 1, 2012.

D. L. Karlen, M. J. Mausbach, J. W. Doran, R. G. Cline, R. F. Harris, and G. E. Schuman, “Soil quality: a concept, definition, and framework for evaluation,” Soil Science Society of America Journal, vol. 61, no. 1, pp. 4–10, 1997.

K. R. Olson, M. Al-Kaisi, R. Lal, and L. W. Morton, “Soil ecosystem services and intensified cropping systems,” Journal of Soil and Water Conservation, vol. 72, no. 3, pp. 64A–69A, 2017.

W. E. Larson and F. J. Pierce, “Conservation and enhancement of soil quality,” in Evaluation for Sustainable Land Management in the Developing World, vol. 2 of IBSRAM Proc. 12(2), pp. 175–203, Int. Board for Soil Res. and Management, Bangkok, Thailand, 1991.

G. Renzi, L. Canfora, L. Salvati, and A. Benedetti, “Validation of the soil Biological Fertility Index (BFI) using a multidimensional statistical approach: a country-scale exercise,” CATENA, vol. 149, pp. 294–299, 2017.

A. Ferrara, L. Salvati, A. Sabbi, and A. Colantoni, “Soil resources, land cover changes and rural areas: Towards a spatial mismatch?” Science of the Total Environment, vol. 478, pp. 116–122, 2014.

S. S. Andrews, D. L. Karlen, and C. A. Cambardella, “The soil management assessment framework: a quantitative soil quality evaluation method,” Soil Science Society of America Journal, vol. 68, no. 6, pp. 1945–1962, 2004.

J. L. Hatfield and C. L. Walthall, “Meeting global food needs: Realizing the potential via genetics × environment × management interactions,” Agronomy Journal, vol. 107, no. 4, pp. 1215–1226, 2015.

M. A. Liebig, D. L. Tanaka, S. L. Kronberg, E. J. Scholljegerdes, and J. F. Karn, “Integrated crops and livestock in central North Dakota, USA: Agroecosystem management to buffer soil change,” Renewable Agriculture and Food Systems, vol. 27, no. 2, pp. 115–124, 2012.

V. Acosta-Martínez, C. W. Bell, B. E. L. Morris, J. Zak, and V. G. Allen, “Long-term soil microbial community and enzyme activity responses to an integrated cropping-livestock system in a semi-arid region,” Agriculture, Ecosystems and Environment, vol. 137, no. 3-4, pp. 231–240, 2010.

M. W. Maughan, J. P. C. Flores, I. Anghinoni, G. Bollero, F. G. Fernández, and B. F. Tracy, “Soil quality and corn yield under crop-livestock integration in Illinois,” Agronomy Journal, vol. 101, no. 6, pp. 1503–1510, 2009.

A. J. Franzluebbers and J. A. Stuedemann, “Early response of soil organic fractions to tillage and integrated crop-livestock production,” Soil Science Society of America Journal, vol. 72, no. 3, pp. 613–625, 2008.

M. M. Mikha, M. F. Vigil, M. A. Liebig et al., “Cropping system influences on soil chemical properties and soil quality in the Great Plains,” Renewable Agriculture and Food Systems, vol. 21, no. 1, pp. 26–35, 2006.

M. A. Liebig, D. L. Tanaka, S. L. Kronberg, E. J. Scholljegerdes, and J. F. Karn, “Soil hydrological attributes of an integrated crop-livestock agroecosystem: increased adaptation through resistance to soil change,” Applied and Environmental Soil Science, vol. 2011, Article ID 464827, 6 pages, 2011.

USDA, “Natural Resources Conservation Service, United States Department of Agriculture, Web Soil Survey,” 2017, http://websoilsurvey.nrcs.usda.gov.

M. E. Watson and J. R. Brown, “pH and lime requirement,” in Recommended Chemical Soil Test Procedures for the North Central Region, J. R. Brown, Ed., NCR Publ. 221 (revised), pp. 13–16, Missouri Agric. Exp. Stn. Bull. SB1001, Columbia, Mo, USA, 1998.

D. A. Whitney, “Soil salinity,” in Recommended Chemical Soil Test Procedures for the North Central Region, J. R. Brown, Ed., NCR Publ. 221 (revised), pp. 59-60, Missouri Agric. Exp. Stn. Bull. SB1001, Columbia, Mo, USA, 1998.

R. L. Mulvaney, “Nitrogen—inorganic forms,” in Methods of Soil Analysis. Part 3. Chemical Methods, D. L. Sparks, Ed., SSSA Book Series no. 5, pp. 1123–1184, Soil Science Society of America and American Society of Agronomy, Madison, Wis, USA, 1996.

S. Kuo, “Phosphorus,” in Methods of Soil Analysis. Part 3. Chemical Methods, D. L. Sparks, Ed., SSSA Book Series no. 5, pp. 869–919, Soil Science Society of America and American Society of Agronomy, Madison, wis, USA, 1996.

L. G. Bundy and J. J. Meisinger, “Nitrogen availability indices,” in Methods of Soil Analysis. Part 2—Microbiological and Biochemical Properties, R. W. Weaver, S. Angle, and P. Bottomley, Eds., SSSA Book Series no. 5, pp. 951–984, Soil Science Society of America, Madison, Wis, USA, 1994.

D. W. Nelson and L. E. Sommers, “Total carbon, organic carbon, and organic matter,” in Methods of Soil Analysis. Part 3. Chemical Methods, D. L. Sparks, Ed., SSSA Book Series no. 5, pp. 961–1010, Soil Science Society of America and American Society of Agronomy, Madison, Wis, USA, 1996.

S. Dolédec and D. Chessel, “Rythmes saisonniers et composantes stationnelles en milieu aquatique. I- Description d'un plan d'observations complet par projection de variables,” Acta Oecologica, Oecologia Generalis, vol. 8, no. 3, pp. 403–426, 1987.

A. García-Martínez, A. Olaizola, and A. Bernués, “Trajectories of evolution and drivers of change in European mountain cattle farming systems,” Animal, vol. 3, no. 1, pp. 152–165, 2009.

J. Ryschawy, N. Choisis, J. P. Choisis, and A. Gibon, “Paths to last in mixed crop-livestock farming: Lessons from an assessment of farm trajectories of change,” Animal, vol. 7, no. 4, pp. 673–681, 2013.

R Development Core Team, R: A Language and Environment for Statistical Computing, R Foundation for Statistical Computing, Vienna, Austria, 2011, http://www.R-project.org.

M. M. Wander and G. A. Bollero, “Soil quality assessment of tillage impacts in Illinois,” Soil Science Society of America Journal, vol. 63, no. 4, pp. 961–971, 1999.

S. S. Andrews, D. L. Karlen, and J. P. Mitchell, “A comparison of soil quality indexing methods for vegetable production systems in Northern California,” Agriculture, Ecosystems and Environment, vol. 90, no. 1, pp. 25–45, 2002.

S. A. Mulaik, “Objectivity and Multivariate Statistics,” Multivariate Behavioral Research, vol. 28, no. 2, pp. 171–203, 1993.

R. F. Daubenmire, Plants and Environment, John Wiley & Sons, New York, NY, USA, 1974.

M. A. Liebig, J. R. Hendrickson, J. D. Berdahl, and J. F. Karn, “Soil resistance under grazed intermediate wheatgrass,” Canadian Journal of Soil Science, vol. 88, no. 5, pp. 833–836, 2008.

A. J. Franzluebbers, J. Sawchik, and M. A. Taboada, “Agronomic and environmental impacts of pasture-crop rotations in temperate North and South America,” Agriculture, Ecosystems and Environment, vol. 190, pp. 18–26, 2014.

L. W. Bell, R. C. Hayes, K. G. Pembleton, and C. M. Waters, “Opportunities and challenges in Australian grasslands: Pathways to achieve future sustainability and productivity imperatives,” Crop and Pasture Science, vol. 65, no. 6, pp. 489–507, 2014.

R. Lal, J. A. Delgado, P. M. Groffman, N. Millar, C. Dell, and A. Rotz, “Management to mitigate and adapt to climate change,” Journal of Soil and Water Conservation, vol. 66, no. 4, pp. 276–282, 2011.

R. M. Sulc and A. J. Franzluebbers, “Exploring integrated crop-livestock systems in different ecoregions of the United States,” European Journal of Agronomy, vol. 57, pp. 21–30, 2014.

J. Huang, H. Yu, X. Guan, G. Wang, and R. Guo, “Accelerated dryland expansion under climate change,” Nature Climate Change, vol. 6, no. 2, pp. 166–171, 2016.