Representative Elementary Volume as a Function of Land Uses and Soil Processes Based on 3D Pore System Analysis

Abstract

Representative elementary volume (REV) is required for representative measurements of soil physical properties. However, questions may arise whether REV depends on how the soil structure is modified or whether processes in the soil affect REV. Here, we explore REV dependence for contrasting land uses (conventional tillage, no-tillage, and minimum tillage) and applying wetting and drying (W-D) cycles. The effect of different subvolume selection schemes (cube and core) on REV was also investigated. For this study, high-resolution three-dimensional images obtained using the X-ray Computed Tomography (XCT) technique were analyzed. The micromorphological properties measured were porosity (P), fractal dimension (FD), degree of anisotropy (DA), and pore connectivity (C). The results show that REV depends mainly on the land uses for P and C (both selection schemes). The core method showed lower REV due to the larger volume analyzed than that in the cube method. It was not possible to define a REV for DA. The REV obtained using the cube method was more sensitive to changes in the scale of analysis, showing an increasing trend with applied W-D cycles for P and FD. Our results indicate that REV cannot be considered static since land uses and processes influence it.

1. Introduction

Any measurement of soil physical properties requires representative samples. When samples are not representative, they will exhibit fluctuations in the evaluation of their physical parameters due to the microscopic domain region [1]. Therefore, increasing the sample size will cause the analyzed volume movement toward the porous system domain [2,3]. The minimum volume required for the sample to become representative is called the representative elementary volume (REV). This elementary volume indicates that the measured physical property becomes independent of the sample size [4]. Soil is expected to be submitted to different land uses in agricultural areas. These uses can cause significant changes in the soil structure, as occurs in conventional tillage (CT), or minor changes, as in more conservationist systems (e.g., minimum tillage (MT) and no-tillage (NT)) [5]. In CT, the topsoil is revolved by the passage of plows and harrows [6]. In no-tillage, on the other hand, plant residues are kept on the topsoil thus preserving its structure [7,8]. Minimum tillage presents some NT and some CT characteristics since the soil is subjected to disaggregation from time to time [9].

When soil structure is analyzed at the micrometer scale, each land use can cause different changes in its pore architecture [10,11]. Properties related to the pore system, such as porosity, pore size and shape distribution, and pore number, can be significantly modified according to the type of soil management [12,13,14,15]. Other properties of interest related to pore system complexity, such as pore connectivity, tortuosity, degree of anisotropy, and fractal dimension are also sensitive to natural soil structure modification [16,17,18,19]. When these properties are not measured using representative samples, the results will show fluctuations due to the microscopic characteristics of the samples [1,2,20,21]. Thus, for the measured physical properties to be closer to the soil characteristics found in the agricultural area, the use of representative samples is required.

However, as the soil is subjected to different land uses that modify its structure mainly in the topsoil, it is reasonable to question how different types of management will affect REV. The answer to this question is fundamental because, for the same soil, contrasting land uses may require different sample volumes to determine the same physical property. On the other hand, soils usually undergo different processes, such as wetting and drying (W-D) cycles. These cycles can cause modifications in the soil pore system due to swelling and shrinking, aggregate coalescence, changes in the aggregates’ mechanical and structural stability, modifications in the soil effective stress, etc. [22,23,24,25,26]. As the soil is constantly being subjected to dynamic processes, which may affect its structure, it becomes pivotal to ask whether REV can be modified when the soil is submitted to W-D cycles.

The advent of techniques that enable three-dimensional (3D) soil analysis, such as X-ray Computed Tomography (XCT), allowed the study of numerous soil geometrical and morphological properties at the micrometer scale [15,24,27,28]. Computed tomography is considered a non-destructive method applied to the investigation of the internal structure of porous materials such as soils [14]. Non-destructive methods enable the assessment and verification of the soil morphological parameters without affecting its original properties or causing damage to the soil structure [29]. The XCT technique is based on the interaction of X-rays with matter. Material density differences generate images with distinct gray shades after X-ray interaction. The image acquisition involves scanning the sample on a rotating table to obtain a series of two-dimensional (2D) images (radiographs). The combination of the 2D slices allows the acquisition of 3D images of the object of interest. Variations in the size of the sample scanned, detector size and resolution, and position of the sample between the X-ray source and detector also influence the image resolution and accuracy.

It is worth mentioning that the determination of any soil physical property by XCT requires representative samples. Therefore, this technique becomes ideal for defining REV, as it allows analyses inside the samples at different volumes [13,29,30]. Furthermore, computed tomography enables the choice of different subvolumes inside the samples and different subvolume selection schemes without the need for sample destruction [31,32]. Thus, we can analyze possible fluctuations in morphological properties obtained by XCT, such as porosity, pore connectivity, anisotropy, and fractal dimension, as a function of scale changes. This kind of analysis is fundamental since it defines the best sample sizes to be analyzed through XCT. Moreover, appropriate sample sizes enable image resolution optimization allowing the evaluation of different soil functions based on the size and shape of pores.

Thus, this study is based on three hypotheses: (1) REV is affected by different land uses, (2) REV is influenced by processes occurring in the soil, such as W-D cycles, and (3) different selection schemes (cubes or cores) affect REV definition. High-resolution XCT images were employed to verify these hypotheses. The contrasting land uses studied were conventional tillage, no-tillage, and minimum tillage. Soil samples were also submitted to 0 and 6 W-D cycles. However, as described before, our study focuses primarily on REV estimation aiming at measuring geometrical and morphological properties through XCT.

2. Materials and Methods

2.1. Experimental Area and Soil Sampling

This study was carried out with samples collected in the municipality of Ponta Grossa, Paraná, Brazil. The experiment was conducted in a Rhodic Hapludox soil (Soil Survey Staff [33]) on the Research Farm of the IDR (25°06′ S, 50°09′ W; 875 m asl). The rainfall recorded in the region is between 1400 and 1600 mm and the average annual temperatures range from 17 to 18 °C. The climate is classified as humid subtropical (Cfb) according to the Köppen classification system [34]. The soil in the depth studied has clay texture with the following granulometric composition: 58% clay (< 2 μm), 28% silt (2–50 μm), and 14% sand (50–2000 μm). In total, six samples were taken from each of the land uses (3 × 6 = 18 samples). All samples were collected in the topsoil (0–10 cm) after maize harvesting (Zea mays L.) [35,36].

The long-term field experiment (over 35 years old) was initiated in 1981 with the following treatments: (1) conventional tillage (CT), (2) minimum tillage (MT), and (3) no-tillage (NT). The CT plot was prepared by plowing with a 70 cm disk twice a year up to a 25 cm depth followed by two 60 cm harrowing. The NT plot was submitted only to sowing and cleaning operations (no soil disturbance). The MT plot was prepared by one chisel plowing up to 25 cm depth followed by one 60 cm narrow disking. For MT and NT plots, the crop residues were kept on the soil surface. Samples under CT were collected approximately six months after soil preparation.

The undisturbed soil samples were collected following the Kopeck’s ring method using c. 5 cm diameter and height stainless steel core cylinders. First, a trowel was employed to prepare a flat surface of topsoil in the sampling area. Plant litter, living plants, and surface rocks were cleared prior to sampling. Next, the volumetric ring (cylinder) with the beveled edge down was placed onto the soil surface. A block of wood was placed on the top of the ring and using the hand sledge (rubber mallet) the ring was carefully driven down to the desired soil depth. After that, the soil was carefully excavated around the volumetric ring without exerting pressure on it to take out the cylinder from the topsoil. The trowel was placed underneath the volumetric ring to carefully lift it and prevent soil losses. Next, the soil excess outside the steel cylinders was carefully trimmed off with a flat-bladed knife to ensure that the soil volume was equal to the cylinder internal volume [37]. Following sampling, the samples were wrapped in plastic film and placed in a cylindrical plastic sample holder to be transported to the laboratory. To minimize damage to the soil structure due to the force required for the cylinder collection, samples were taken some days after a high-intensity rainfall with the soil near its field capacity. In addition, the samples were collected from the crop rows to avoid compaction due to the traffic of agricultural implements.

2.2. Wetting and Drying (W-D) Cycles

After preparing the samples, a porous cloth was placed on one of their surfaces to carry out the wetting process. The process used to saturate the samples was capillary rise. The saturation process was carried out slowly, placing a water amount of approximately 5 to 10 mm at the base of the samples. After approximately two hours, 5 mm of water was added every hour to approximately 8/10 of the maximum height of the cylinders. The saturation process was completed after a period of 48 h. Following the saturation step, the samples were placed on a tension table (EijKelkamp-model 08.01 SandBox) and partially dried by submitting them to a matrix potential of –6 kPa, equivalent to a c. 50 μm pore diameter (Young–Laplace equation). The time required for thermodynamic hydraulic equilibrium was three days when the water output from the tension table exhaust system ceased. Next, the samples were removed from the tension table and subjected to a new wetting process (W) and subsequent drying (D). This procedure was repeated six times so that the samples under the different land uses were submitted to six wetting and drying cycles (W-D). The cycled samples were then compared to non-cycled samples (0 W-D). The number of cycles was chosen according to the average maximum number of rainfall events in the studied experimental site in the driest month (August).

2.3. X-ray Computed Tomography (XCT)

The soil samples were scanned with X-ray Computed Tomography system Nikon XT V 130C and reconstructed in 32-bit (to avoid greyscale histogram compression) at a voxel size of c. 35 μm using the Inspect-X/CT Pro 3D (image acquisition and reconstruction) (Figure 1). All the soil samples were scanned at voltage, current, and exposure time of 125 kV, 140 µA, and 250 ms, respectively. Aiming to minimize beam hardening artifacts, a 0.25 mm Cu filter was inserted close to the X-ray source. After reconstruction, the images were imported into Volumetric Graphics (VG) StudioMAX^® 2.0 and into ImageJ 1.42 and cropped to a cubic shape with c. 30.1 × 30.1 × 30.1 mm³ (860 × 860 × 860 pixels). Cropping was carried out to avoid possible voids or compacted regions close to the core walls. We employed the region of interest (ROI) tool in ImageJ to crop the imaged samples to the desired volume of interest.

The original grey-level XCT images were processed using ImageJ 1.42 software [38]. First, the raw images were submitted to the following processes before segmentation: (1) a 3D median filter was firstly employed to reduce the noise in the images, with a radius of 2 pixels, and (2) an enhancement contrast (0.5%) tool was utilized to enhance the edges of the pores and improve the image contrast. The segmentation procedure was based on the grayscale histograms of the two-dimensional (2D) images. The images were converted into 8-bit (256 shades of gray) before segmentation. The different materials that make up the soil produce peaks that are related to the attenuation of X-ray photons. For the dry samples analyzed in our study, two peaks occur, with the lower and higher shades of gray in the histograms representing the attenuation by the lower (air) and higher (solid) density materials. Thus, the existence of two peaks in the histogram is related to a material composed of two phases. After that, the segmentation process was based on the nonparametric and unsupervised Otsu method of threshold [39]. Next, we utilized the zoom tool in ImageJ to visually inspect the segmented images to ensure the performance of the segmentation procedure selected. Finally, the segmentation process allowed us to obtain a binary image, in which pores and soil material were represented by black and white pixels (see Supplementary Materials Figure S1).

2.4. Subvolume Selection

In the next step, aiming at the study of the representative elementary volume (REV), the images were cropped again using two different approaches: (1) concentric core—A fixed axis (z = soil depth) was defined and the subvolumes were selected concentric from this axis (Figure 2a), and (2) concentric cubes—A cubic subvolume was selected at the center of the ROI and increasingly larger cubes were selected from the smallest selected to the entire volume of interest (VOI) (Figure 2b). Finally, the crops resulting from all subvolumes were carried out with the REVaux software, developed especially for this study. With the help of this program, we determined the number of sections cropped, defined their sizes in 3D, and used the central section as a reference.

The methods chosen for the selection of the subvolumes were similar to those based on the paper by Costanza-Robinson et al. [31]. Those authors used two approaches, the traditional point-centered cube geometry with the cube growing symmetrically from a central volume, and the core-centered cuboid window geometry. In the latter case, one of the axes has a fixed height equal to the height of the image. Other REV studies for porous media using cube- and core-based selection schemes are those by Koestel et al. [40], Baveye et al. [32], Yio et al. [41], and Wu et al. [42].

Figure 3 shows the subvolume selection obtained for the cubes after the procedure of 3D image cropping. This choice was made, as already mentioned, based on centralized cubes and cores selected along a central axis of symmetry. The subvolume cropping procedure was performed without extrapolating the previously established maximum VOI (27,271 mm³). Nine small subvolumes were acquired by cropping the VOI.

Table 1. Window length-scale (L) and subvolumes (V) selected for the representative elementary volume (REV) analysis for the subvolume selection schemes (centralized cube and centralized core) studied. The subvolume 10 represents the entire volume of interest (VOI).

Subvolumes (V)	Window Length (L)	Core	Cube
	(mm)	(mm3)
1	3.0	271	27
2	6.0	1084	216
3	9.0	2438	729
4	12.0	4334	1728
5	15.1	6863	3442
6	18.1	9861	5930
7	21.0	13,274	9261
8	24.1	17,482	13,998
9	27.1	22,106	19,903
10	30.1	27,271	27,271

shows the window length-scale (L) and the subvolumes selected for both schemes (concentric cube and core).

2.5. Soil Physical Properties

The physical parameters selected and analyzed for the determination of REV were the imaged soil porosity (P), fractal dimension (FD), degree of anisotropy (DA), and pore connectivity (C). These properties were measured using the ImageJ 1.42 software. We define porosity as the ratio between the volume of pores (V_pore) and the total volume of a given sample. We obtained the imaged porosity in this study using the ImageJ voxel counter function. After counting the voxels referring to the pores (V_voxel), the volume occupied by them (V_{total pore}) was divided by the total image volume (V_image) as presented in the following set of equations [35]:

$${\mathrm{V}}_{\mathrm{pore}}={\mathrm{V}}_{\mathrm{voxel}}\times {\left(\mathrm{vr}\right)}^3,$$

(1)

$${\mathrm{V}}_{\mathrm{total} \mathrm{pore}}=\sum {\mathrm{V}}_{\mathrm{pore}},$$

(2)

$$\mathrm{P}\left(\%\right)=\left(\frac{\sum {\mathrm{V}}_{\mathrm{pore}}}{{\mathrm{V}}_{\mathrm{image}}}\right)\times 100,$$

(3)

where vr represents the voxel resolution (mm³) and P(%) the imaged porosity expressed as a percentage.

The fractal dimension provides a systematic method to identify irregular patterns of fractal sizes that contain internal structures repeated at different scales [43]. The measurement of the fractal dimension of any system or porous medium is important because it reflects its heterogeneity. We employed the fractal dimension plugin in ImageJ to determine the FD based on the box-counting method. In this method, the FD is determined by the slope of the straight line of the graph of log(N) (y-axis) against log(r) (x-axis):

$$\mathrm{FD}=\frac{log\left(\mathrm{N}\right)}{log\left(\mathrm{r}\right)},$$

(4)

where N represents the number of boxes that cover the object under analysis while r is the magnification of the inverse of the box size.

The degree of anisotropy provides information on how the given material properties vary depending on the direction in which they are measured. Generally, this property has values ranging from 0 to 1. The first case means that the analyzed property is entirely isotropic; the sample does not have any directionality. In the second case, there is an extreme orientation of specific structures in the image. In our study, we employed the BoneJ plugin to calculate the degree of anisotropy of the soil pore system. The DA calculation is based on a series of vectors of the same length originating from specific positions of the analyzed volume, divided by the number of times they intercept the pores (mean intercept length). Ellipsoids are employed to construct the anisotropy tensors to obtain eigenvalues related to the lengths of the axes of the ellipsoid and eigenvectors giving the orientations of these axes. Thus, the following equation is used to obtain DA [14]:

$$\mathrm{DA}=\left(1-\frac{\mathrm{SE}}{\mathrm{LE}}\right),$$

(5)

where SE and LE stand for the smallest and largest eigenvalues, respectively.

The connectivity is a geometric feature that provides information about pore space distribution. Therefore, it measures the number of independent paths between two points in space, which characterizes the degree of interconnection between pores [44]. The connectivity measurement was performed after the purify filter procedure (connectivity function) in the BoneJ plugin [45]. Next, the BoneJ connectivity tool was used in the C measurements. The pore connectivity quantification was performed based on the measurement of the Euler number (EN) according to the following equation:

$$\mathrm{EN}=\mathrm{N}-\mathrm{C}+\mathrm{H},$$

(6)

$$\mathrm{C}=1-\mathrm{EN},$$

(7)

where N represents the number of isolated objects, C indicates the number of redundant connections, which refers to the connectivity or genus, and H stands for the number of completely enclosed cavities. According to Vogel et al. [21], when the Euler number is analyzed, its value will be positive when the number of isolated unconnected pores exceeds the number of multiple connections between the pores (N > C). For a fully connected network of pores, the Euler number will be negative (C >> N = 1). In this situation, EN counts the number of multiple connections and corresponds to the number of meshes in the pore network. When considering porous systems such as soil, the parameter H can be neglected since an aggregate completely surrounded by pores will hardly be found.

2.6. Representative Elementary Volume (REV) Estimation

The average values obtained for each subvolume, measured properties, and selection schemes were used to estimate REV. Additionally, we chose to use a window length-scale (L) for the two selection schemes, as shown in Figure 2, where L = x. In the cube scheme, L is equal to V^1/3 [32]; however, in the core scheme, L represents an effective length from the center of the subvolume, always keeping a fixed axis (z = depth) [31] (Figure 2a). Next, we present the criteria that were employed in the choice of REV: (1) relative variation (Equation (8)) of the average values of the properties measured (imaged porosity, fractal dimension, degree of anisotropy, and pore connectivity) between adjacent window lengths not over 5%, and (2) at least three or more consecutive window lengths did not present different average values of the properties measured using the variation criterion of the first item. A similar procedure was utilized by Vandenbygaart and Protz [4]. Finally, it is worth mentioning that all average values of the measured properties must have values of Δ less than 5% (existence of a plateau) to estimate REV. This is because the REV value will be estimated from the first length having property average values with Δ < 5%, indicating the existence of a plateau.

$${\Delta }^{\mathrm{i}-1}=\left|\frac{{\mathsf{\omega }}^{\mathrm{i}}-{\mathsf{\omega }}^{\mathrm{i}-1}}{{\mathsf{\omega }}^{\mathrm{i}-1}}\right|\times 100,$$

(8)

where Δ is the relative variation expressed as a percentage, ω is the average value of the property measured, and i (1 to 10) indicates each of the window lengths analyzed (Table 1). The use of Equation (8) to estimate REV was based on the paper by Wu et al. [42]. The value of Δ < 5% was selected considering the spatial variability for the properties studied. In many situations, such as for pore connectivity, the criterion we chose was more restrictive than those employed by Wu et al. [42] based on the coefficient of variation.

Figure 4 illustrates the REV definition considering the regions of the microscopic inhomogeneity domain (left of the dotted line) and the region of the porous medium domain (right of the dotted line). The dotted line indicates the smallest window length (REV) over which a representative measurement can be made. Figure 5 shows the flowchart of the main steps performed during this study.

3. Results

3.1. Soil Porosity (P)

The results of the average imaged porosity showed the subvolume selection scheme (core or cube) influence in the REV determination (Figure 6). We also observed that REV was affected when the soil was submitted to contrasting land use systems and wetting and drying cycles. No clear trend was observed among REVs when the different methods of subvolume selection were employed for P. In the core method, samples submitted to 0 and 6 W-D cycles reached the REV for MT and CT for the same window length (L = 3.0 mm) (Figure 6a,c). For NT (Figure 6e), estimating REV was impossible. In the cube method, we observed that REV could be calculated for all land uses. The application of W-D cycles requires larger window length images in the REV estimation. For MT (Figure 6b), REV was estimated at 12.0 mm (0 W-D) and 18.1 mm (6 W-D). For CT (Figure 6d), REV was reached at 6.0 mm (0 W-D) and 12.0 mm (6 W-D), whereas for NT (Figure 6f), the values found were 18.1 mm (0 W-D) and 21.0 mm (6 W-D). Such a result obtained using the cube method shows that the porosity was sensitive to the chosen land use and the application of W-D cycles.

3.2. Fractal Dimension (FD)

The average fractal dimension, calculated through the box-counting method, was employed to characterize the complexity of the soil pore system. The core method presented lower REV values among the land use systems and W-D cycles when compared to the cube method (Figure 7). Similar to P results, the soil land use showed differences in the REV estimation only in the cube method (except for NT in the porosity case—core method). Using the core method, MT presented the same REV for 0 and 6 W-D (L = 6.0 mm), indicating that for this land use system, the wetting and drying cycles did not influence the REV estimation (Figure 7a). For CT (Figure 7c) and NT (Figure 7e), similar to the results of MT, REV was not affected by the cycles (L = 6.0 mm—0 and 6 W-D). In the cube method, the W-D cycles did not affect the REV definition for MT (L = 18.1 mm) (Figure 7b). However, for CT (Figure 7d) and NT (Figure 7f), the application of the W-D cycles increased the window length needed for the REV estimation. These results follow the same trend as observed for P (Figure 6d,f). For CT, REV increased with the application of W-D cycles from 6.0 mm (0 W-D) to 18.1 mm (6 W-D), while for NT, it moved from 12.0 mm (0 W-D) to 18.1 mm (6 W-D). Interestingly, applying the cycles caused the REV to be achieved for the same length window for all land uses in the core method.

3.3. Degree of Anisotropy (DA)

The degree of anisotropy (Figure 8), a tool used to measure how substructures are oriented in a given analysis volume, showed different results between the subvolume selection methods, land use systems, and W-D cycles. This property identifies whether there are differences in pore distribution throughout contrasting portions of the soil samples. When the core method was used, the average DA values obtained were <0.32, showing more isotropic structures, regardless of the land use adopted. The DA values also presented fluctuations for all the window lengths examined without a clear trend to stabilize among the land use systems studied. For this reason, none of the samples reached the REV for this selection scheme (Figure 8a,c,e). Regarding the cube method, the largest window lengths presented DA reductions but still without a clear trend of stabilization of their values. Similar behavior was observed for all the land uses (Figure 8b,d,f). In this selection scheme, we found that smaller window lengths had more significant anisotropies, showing the dependence of this property on the selection scheme. Thus, we verified the soil land use influence on the DA distribution for the samples analyzed. Finally, as observed in the core method, REV was not achieved for any of the window lengths analyzed.

3.4. Pore Connectivity (C)

Pore connectivity is a physical property often obtained through 3D imaging, providing an idea of pore continuity. This property offers results that directly impact water dynamics in porous systems. Our findings showed the influence of selection schemes as a function of window length-scale on the REV estimation among the land use systems studied (Figure 9). The different land uses affected the connectivity measurement and the REV definition. As observed for P (Figure 6) and FD (Figure 7), the application of W-D cycles impacted the REV estimation for C (except for CT—core method). In the core method, applying the W-D cycles increased the REV of the samples under MT (L = 3.0 mm—0 W-D and L = 6.0 mm—6 W-D) (Figure 9a). On the other hand, opposite behavior was observed for NT (Figure 9e), with REVs of 6.0 mm (0 W-D) and 9.0 mm (6 W-D). The soil under CT (Figure 9c) did not present variations in the REV estimation (L = 3.0 mm) with the application of the W-D cycles. Concerning the cube method, REVs varied when compared to the core method, indicating that this property is susceptible to the subvolume selection procedure. The results show that for MT (Figure 9b), REV was achieved in 15.1 mm (0 W-D) and 21.0 mm (6 W-D). For CT (Figure 9d), REV was reached only for the samples subjected to 6 W-D (L = 6.0 mm), whereas NT (Figure 9f) showed the opposite with REV estimation (L = 21.0 mm) only for the samples not submitted to W-D cycles. These results demonstrate no clear trends in the REV estimation based on the cube method selection scheme considering different land uses.

4. Discussion

This paper focuses on three main assumptions: (1) REV estimation is affected by land uses, (2) REV cannot be considered static, and (3) REV estimation is influenced by the subvolume selection scheme for complex systems such as soils. To analyze these three statements, we selected four soil pore morphological properties, i.e., porosity, fractal dimension, degree of anisotropy, and pore connectivity. We noticed that REV varies as a function of the soil physical property analyzed, as demonstrated by other papers in the scientific literature [2,13,20,21]. The characteristics of each property measured directly impact the REV estimation, mainly those related to the distribution of pores inside the sample volume. This finding reinforces the need for representative samples mainly when XCT studies are carried out.

For the soil porosity assessed using XCT, we observed that except for the sample submitted to the no-tillage system (Figure 6e), the other systems (MT and CT) showed the same REV regardless of the application of the W-D cycles for the core method. This result indicates soil structure similarities for these land uses, which were not significantly affected by the W-D cycles. However, we know that CT causes profound modifications in soil structure and that the soil tends to recover its structure after applying continuous W-D cycles [22]. Conventional tillage is characterized by aggregate breakdown due to the harrowing and plowing operations [6]. In the case of MT, which represents an intermediate land use system presenting characteristics of CT and NT, a more complex structure is expected [9]. However, the smallest subvolume analyzed was already enough for these two land uses to reach REV, which is a surprising result. As Costanza-Robinson et al. [31] indicated, the lower REV values estimated using the core method may be associated with the substantial contribution of the analyzed volume due to the fixed width of one of the analysis axes (z-axis) (Figure 2a). For example, if we look at the first subvolume analyzed, the core method presented a volume almost ten times that of the cube method. This larger volume obtained may explain the stability in the core method’s P values, indicating that the variability of this property with the increased size of the window lengths was compensated by this larger volume (fixed z-axis). Similar results showing minor variability between the measured properties as a function of scale using the core method were also noticed by Wu et al. [42]. Borges et al. [13], using a method similar to that of the centered cube, found more significant variability of P values for the smallest subvolumes when compared to our study. However, the first and second subvolumes analyzed by those authors were approximately 680 and 5.4 times smaller, respectively, than the first subvolume employed in our study. Thus, the first subvolumes investigated using the core method in our study were already large enough to mitigate the effects of P variability.

Regarding the cube method, variations in REV were observed as a function of land uses and W-D cycles. The results show that applying the W-D cycles changed the soil structure, so larger subvolumes are required for representative P measurements [24]. This finding indicates changes in the soil pore system heterogeneity with the application of W-D cycles, as shown by Oliveira et al. [17]. This result may also be associated with the action of organic material and soil clay content in the topsoil, especially for the more conservationist systems (MT and NT). As for MT and NT (Figure 6b,f), which are systems that usually present more complex pore distributions, REV was reached for larger subvolumes than for CT. This result was already expected, especially for the samples collected from the topsoil, because the plant remains on the surface in these land uses [46,47,48]. For CT, the breakdown of aggregates due to the soil preparation process tends to homogenize the soil composition [10]. Similar results to our study were also found by Borges et al. [13] working with CT and NT land uses. Those authors found, in general, higher REVs for the samples submitted to the NT employing the centered cube method. Our results showed that for porosity, even slight variations in this property with depth are enough to influence the REV obtained using the cube method, similar to the findings of Costanza-Robinson et al. [31] and Baveye et al. [32]. Finally, the most significant P variation for the smallest subvolumes was mainly related to the microscopic inhomogeneity domain highlighting the importance of small-scale spatial variability [42,43,49].

Fractal dimension was utilized to identify the complexity of the pore system through the analysis of patterns and shapes [50,51]. Pore system distribution and complexity changes due to soil processes can be quantified through tools such as FD [52]. Our results demonstrated that FD was the unique property with REV for all the land uses, W-D cycles, and selection schemes (Figure 7). The FD values between 2 and 3 are characteristics of 3D structures [53]. The application of W-D cycles did not influence the REV for FD in the core method (Figure 7a,c,e), indicating similarities in the complexity of the soil pore architecture under scale variations [17]. As discussed earlier, the FD homogeneity results verified when using the core method should be associated with the largest analyzed volume, and a higher complexity of the porous system in depth, which was the axis kept fixed in our study for this volume selection scheme. Lopes de Silva et al. [54] recently showed the influence of changes in the region of choice of volumes to estimate REV, especially for properties related to the complexity of the porous system. The same behavior was observed by Borges et al. [13], working with properties such as macroporosity, pore number, connectivity, and tortuosity. We can state that for this selection scheme (core method), it would perhaps be more interesting to use tools that enable the analysis of the fractal characteristics of portions of the sample, such as the multifractal aiming to find subtle differences in the soil pore architecture due to the W-D cycles [50,51].

When using the cube method, as expected, the samples that were not submitted to the W-D cycles presented higher estimated FD REVs for more conservative pore systems (MT and NT) (Figure 7b,f). This result follows the same trend as P and is explained by the greater complexity of the pore networks for these land uses [8,55,56]. The absence (NT) or minimum disturbance (MT) of the soil structure can help to explain this high complexity in comparison to CT [57]. The latter is subjected to severe soil disturbance due to harrowing and plowing, requiring more time for the structure to recover [58]. For MT and NT, the minimum soil disturbance favors the appearance of biopores resulting in branched pore networks [59]. Similar to P results, the application of cycles affected the REV estimate for CT (Figure 7d) and NT, requiring larger subvolumes for these two land uses. This result demonstrates that the W-D cycles increased the irregularity of the shape of the analyzed pores [43]. Pires et al. [60] showed that the application of W-D cycles tends to increase the soil porosity and the pore system complexity, mainly in oxisols such as the one studied in this work, confirming the results found. It is interesting to note that similar results to those obtained in our study regarding the relationship between FD and P were also observed by Wu et al. [61]. Those authors analyzed several possible relations between FD and P, showing that the most common case for coal samples was when the samples reached REV for FD while P declined. They also noticed cases where both FD and P reached REV, similar to our findings.

REV was not reached for the degree of anisotropy (Figure 8). As said before, DA provides information on how the properties of a given material vary depending on the direction in which they are measured [62]. The results of our paper demonstrated the existence of low anisotropic pore structures [14,63]. Based on this, we expected some stability in the DA results, increasing the subvolumes, which was not observed. Only the land uses MT and CT showed some stability trends in the DA for the cubic scheme (Figure 8b,d). Our results indicate that the distribution of the pores is susceptible to the analysis volume [19]. As different volumes comprise distinct portions of the soil pore system, changes in the pore size and shape domains can affect the pore orientation, influencing DA. Thus, the spatial variability of the pore size distribution induces the appearance of clusters of pores in different portions of the sample, affecting DA [18,64]. Based on that result, we concluded that for DA, the choice of higher volumes is necessary for the samples to reach REV for the criterion chosen in our study. Even small changes in the pores configuration (orientation) seem to affect the stability of the DA values when the analysis volume is varied.

The selection scheme affected pore connectivity (Figure 9). REV was reached for all the land uses when employing the core method. This result indicates that this selection method minimizes the pore continuity variation from one volume to the other (more considerable volume of data) when compared to the cube method. The core method maintains the z-axis constant while varying the y- and x-axes as described earlier (Figure 2a). For this selection scheme, samples under MT (Figure 9a) and NT (Figure 9e) land uses required larger volumes to achieve REV when compared to CT (Figure 9c). This result, as already stated, is related to the complexity of the soil pore system under more conservationist uses [65,66]. Generally, soils under MT and NT tend to have large and highly connected pores, as shown by Galdos et al. [8], de Oliveira et al. [10], and Tseng et al. [14] working with Brazilian soils. As regards CT, the plowing and harrowing process breaks the existing connection between the soil pores in the topsoil, homogenizing pore size distribution [67,68]. We observed that REVs were achieved for the smallest volumes analyzed for C, indicating that the maintenance of one of the fixed axes in the selection of the volume in the core method minimized the variability of the pore connectivity values at the scales analyzed similarly to the findings for the other properties [31]. Recently, Borges et al. [13] published a paper evaluating the REV for a type of soil similar to the one studied in our paper. Those authors showed more significant variability in the connectivity values for the different subvolumes analyzed, especially the smaller ones. However, it is worth noting that the first subvolumes investigated by those authors were much smaller than the smaller ones selected in our study, as previously mentioned. Regarding the W-D cycles, they neither showed much influence on the REV estimation for the different land uses, nor a clear trend between them. Again, this result can be explained by highly connected pores inside the samples, mainly in the core method fixed axis direction (sample depth) [12,18,27].

The cube method showed different results from those of the core method for the contrasting land uses studied. Only the soil under MT (Figure 9b) showed REV between samples submitted to 0 and 6 W-D cycles. This finding indicated that the variability of C cannot be disregarded when different subvolumes are analyzed, even though there are dominant and well-connected pores, as shown by other authors [8,10,12,14,18]. Such a result shows that the pore structure and its heterogeneity influence the REV definition, as shown by Xue et al. [69]. When CT was investigated (Figure 9d), the REV was estimated only for the sample subjected to 6 W-D cycles. This result may be associated with the breakdown of aggregates after soil preparation operations under CT. Applying the cycles causes the soil to undergo restructuring favoring the connection between the soil pores [24]. In the case of NT (Figure 9f), the opposite occurs, with REV being estimated for samples not submitted to W-D cycles. For the samples submitted to 6 W-D cycles, we noticed an abrupt variation of C for the largest subvolumes, which may be associated with variations in pore distribution (clusters) within the soil sample [17,62,70,71]. Koestel et al. [40] recently showed that pore connectivity undergoes more significant variations than porosity when different subvolumes are analyzed using the cube method. Such a result is similar to the ones observed in our study for C. Those authors also showed that C is sensitive to the size of the analyzed pores, a situation that was not investigated in our study.

Although our objective in this paper was the REV definition to determine soil properties based on XCT, some results obtained have important implications for many processes in the soil. Water movement, for instance, is greatly influenced by soil permeability, which is closely related to pore connectivity and porosity. Thus, changes in pore connectivity will affect the soil ability to infiltrate and redistribute water due to permeability modifications [72,73]. Adequate soil porosity is also essential for the good development of root systems and to maintain appropriate water levels available to plants. Therefore, reliable determinations of porosity and pore connectivity have several practical implications for agriculture. The water flow direction into the soil also influences crop yield due to its impact on water distribution and retention to plants. Properties such as the degree of anisotropy and fractal dimension give us an idea of the heterogeneity of the pore system. As the water flow direction depends mainly on the pore structure, representative measurements of DA and FD can bring new insights into water distribution inside porous materials. As we can see, there is a close relationship between the morphological properties measured in our study and the soil hydraulic properties [16,73]. Therefore, representative measurements of soil properties through XCT can help to explain many hydraulic processes related to water retention and movement, the destination of soil nutrients and pollutants, and the amount of water needed for crop growth [74,75,76].

Our study demonstrated the importance of selecting representative samples to measure geometrical and morphological properties using XCT. However, depending on the selection scheme, some properties might be representative. This finding indicates that one sample shape type (core or cube) can be sometimes more appropriate than others. Costanza-Robinson et al. [31] indicates that REV estimation using the cube method is useful for literature comparisons, while the core method provides values with practical applicability. Land uses also affect the REV definition, which means it might be necessary to select samples of different sizes for similar soils in the same experimental field [77,78,79,80]. However, our study covered only a few properties measured by XCT. An increase in the number of properties evaluated might also bring exciting results. Although we selected only two schemes of subvolume selection, the choice of other types (e.g., from top to bottom, random, from right to left) might offer essential insights mainly for the properties associated with pore network complexity. Finally, an analysis of REV for different voxel sizes might give us some results about the role of the intra- and inter-aggregate pore structures. Koestel et al. [39] demonstrated that the size of the analyzed pores might influence the REV estimation. Recently Lucas et al. [81] presented results also addressing the issue of scale and image resolutions in representative measurements of morphological and geometrical properties of soil pores. Those authors demonstrated that when changing scales, instabilities can occur in the measured properties as a function of the transition between pore domains and the definition of representative samples for each scale of analysis [82].

5. Conclusions

Based on the three main assumptions investigated, we observed that the selection scheme is generally essential to the representative elementary volume definition, which means that some precaution is necessary before choosing such a scheme. Furthermore, applying the W-D cycles also influences the REV definition among the properties analyzed, indicating that REV cannot be considered a static property. This result indicates that the REV defined for an experimental field in a specific period might change over time. This finding is essential in areas under human intervention, such as those under contrasting land uses. Thus, the land use type is also crucial for the REV definition. Different human activities can change the soil structure differently, so REV might be affected even for the same soil type. We also noticed the REV dependency on the soil pore system complexity. The degree of anisotropy does not reach the REV for any selection scheme (cube or core) studied and land use system, while regarding fractal dimension, the opposite was noticed. Therefore, the results presented in our study demonstrated that for the same soil under contrasting land uses, REV needs to be determined for each management. Additionally, soil processes (e.g., W-D cycles) also affect REV.