Effects of freeze-thaw cycles on soil structure under different tillage and plant cover management practices

m, 300 – 1020 µ m and > 1020 µ m) and analyzed the macroporosity (V t ), mean macropore diameter (d m ) and mean Euclidian dis t ance to the nearest macropore (ED m ). Additionally, we analyzed the effects of tillage and plant cover on several µ CT- derived geometric parameters in Full Range. Overall, NT-B and NT-V resulted in lower macroporosity than in CT-B and CT-V. Similarly, we found fewer, less branched macropores with shorter mean branch length in NT compared to CT for both plant cover treatments. However, we propose that µ CT-derived geometric parameters might be confounded by the overlapping influence of relatively few, complex and voluminous coarse macropores and the more abundant, less complex very fine macropores. Freeze-thaw, in turn, caused crumbling of soil around coarse macropores, reducing V t and d m in Full Range and reducing V t in the > 1020 µ m range. Addi- tionally, FT caused significant increases in V t and reductions in d m and ED m in the < 300 µ m range, indicating creation of new very fine macropores and expansion of previously indiscernible macropores. Overall, the effects of FT were reduced in NT (for equal plant cover treatments) and V (for equal tillage treatments), indicating greater resilience against FT in both cases.


Introduction
Soil structure can be understood as the arrangement and properties of a heterogeneously cemented organo-mineral matrix and a hierarchical system of pores (Vogel et al., 2021). This structure affects the quality of agricultural soil and plays a central role in processes such as water storage and movement, nutrient leaching, gaseous emissions and carbon storage (Ball, 2013;Mueller et al., 2013;Munkholm et al., 2008). Therefore, understanding the processes that drive changes in soil structure is a vital step in understanding how these soil functions can be managed and protected.
Matrix and pores are continuously formed and transformed by a multitude of biotic and abiotic processes in interaction with the spatial properties of the pores and the heterogeneous cohesion of the matrix (Vogel et al., 2021). The cohesion of the soil matrix, in turn, develops as soil particles bound by electrostatic forces are further cemented by deposition of bacterial, fungal, and plant decomposition products as well as microbial and plant-derived binding agents, roots and fungal hyphae (Degens, 1997;Totsche et al., 2018). Conventional tillage fragments the soil and reduces structural cohesion by increasing soil organic carbon (SOC) mineralization (Ye et al., 2020), as well as mixing and breaking of plant and fungal binding structures. This reduction in matric structural cohesion with tillage is directly measurable in intact, functional soils (e. g. Munkholm, Hansen and Olesen, 2008), as well as in clods and aggregates produced in the laboratory (Mondal and Chakraborty, 2022). In contrast, numerous studies (e.g. Munkholm et al., 2003;Celik et al., 2017;Pires et al., 2017;Schlüter et al., 2018) have shown that reduced tillage and no-till positively affect soil structure and macroporosity, resulting in greater soil strength and potentially greater resistance against physical disturbances (Wiermann et al., 2000). Likewise, the use of cover crops promotes short-term stabilization of soil structures via enmeshment of naturally-occurring soil fragments by fine roots, as well as long-term structuring by increasing the carbon input for soil fungi, bacteria and macrofauna in the shape of litter and root exudates (Haynes and Beare, 1997;Rillig et al., 2002). Thus, management choices such as tillage and cover crops are strong determinants of soil structure via their influence over matrix coherence and fragmentation.
One of the less explored abiotic factors affecting soil structure is freeze-thaw. Soils in temperate climates experience freezing air temperatures of different intensity and duration, intermitted with periods of thaw, throughout the winter. It is well-known that a number of physical parameters related to the soil's structure are affected by mechanical stresses caused by the expansion of forming ice lenses and the inhomogeneous wetting and drying that accompany these freeze-thaw cycles (Bullock et al., 1988;Utomo and Dexter, 1981). Furthermore, while most effects take place gradually over tens of freeze-thaw cycles, structural changes have been observed in the laboratory after a few and sometimes a single freeze-thaw event (Henry, 2007). Freeze-thaw affects parameters such as aggregate stability (Dagesse, 2013;Oztas and Fayetorbay, 2003;Skvortsova et al., 2018), bulk density and unsaturated hydraulic conductivity (Leuther and Schlüter, 2021), as well as water retention , saturated hydraulic conductivity and surface runoff (Ma et al., 2019). Thus, due to its potential to alter soil structural properties, freeze-thaw can affect a multitude of biogeochemical processes such as organic matter turnover and stabilization, as well as carbon and nitrogen mineralization (Song et al., 2017).
The effects of freeze-thaw on soil structure are dependent on soil texture, organic matter content, number of freeze-thaw cycles, water content upon freezing and, importantly, the initial structure of the soil and its resilience against mechanical stresses (Oztas and Fayetorbay, 2003). A well-known phenomenon is the fragmentation of large soil clods produced by ploughing during the winter season (Edwards, 2013;Leuther and Schlüter, 2021). Thus, choices in field management have the potential to increase or reduce the effects from freeze-thaw for a given soil, as they render the structure of the soil weaker or stronger to physical strains. However, knowledge about the interaction between freeze-thaw and soil structure as affected by field management, particularly tillage and plant cover, is scarce. For this purpose, X-ray computerized microtomography (µCT) represents an excellent non-invasive method, making it possible to assess the macroscopic properties of soil and its structure at varying spatial scales. Furthermore, technological advances constantly push the limits in image resolution attainable with µCT allowing for the observation of increasingly finer details in the structure of a sample of soil (Vaz et al., 2014).
In this study, we used X-ray µCT to explore the cumulative effect of five freeze-thaw cycles on soil structure (≥58 µm) of intact agricultural topsoil cores from a long-term reduced tillage and crop rotation field experiment. The soil in question was a sandy loam subjected to either conventional inversion tillage (CT) or no-till (NT) and kept either bare or with plant cover post-harvest. We hypothesize that the mechanical effects of water displacement and expansion upon freezing will alter soil structure via the fragmentation of internal soil structures, but that this deterioration will be dependent on the different initial structural stabilities resulting from different field management practices. We reason that soil structure that has developed over a longer time in NT is more stable against the mechanical stresses imposed by freeze-thaw than the less consolidated structure produced by ploughing, which in turn would experience greater fragmentation and crumbling. Likewise, we anticipate that root growth would have a stabilizing effect on soil structure against the effects of freeze-thaw, and thus soil kept under plant cover after harvest would experience smaller structural changes than soil kept in bare fallow. Altogether, we expect the combination of tillage and plant cover to correspond to an increasing gradient in soil structural resilience against freeze-thaw, from bare fallow conventionally-tilled topsoil to no-till with post-harvest plant cover. Conversely, we expect structural changes due to freeze-thaw to be more pronounced in conventionally tilled bare fallow, and to decrease with both plant cover and no-till treatment.

Study site and sampling
Intact soil core samples corresponding to two tillage treatments and two plant cover treatments were extracted from the CENTS long-term crop rotation and tillage experiment at Aarhus University's department of Agroecology in Flakkebjerg (55 • 19 ′ N, 11 • 23 ′ E) (Gómez-Muñoz et al., 2021;Hansen et al., 2015). The soil was classified as a Glossic Phaeozem (Krogh and Greve, 1999) with a sandy loam texture consisting of 14.7 % clay (<2 µm), 13.7 % silt (2-20 µm), 42.6 % fine sand (20-200 µm) and 27 % coarse sand (200-2000 µm), and an average organic C content of 9.8 mg C g − 1 (Gómez-Muñoz et al., 2021). The average bulk density in the top 100 mm was 1.30 Mg m − 3 in CT and 1.42 g Mg m − 3 in NT. The two tillage treatments consisted of conventional inversion tillage (CT) and no-till (NT) for 18 continuous years at the time of sampling, while the two plant cover treatments consisted of bare soil (B) and winter rye volunteers (V) after the current year's harvest. The crop rotation in previous years was the same for all tillage and plant cover treatments (R2 in the field experiment), and consisted of both winter cereals and spring cereals preceded by catch crops. The NT management was carried out by direct sowing with a chisel coulter drill and straw retention after harvest, while the CT treatment was carried out by sowing with a traditional seed drill, with straw retention after harvest and plowing followed by rolling in late autumn (Hansen et al., 2015). The soil at all sampling sites was sprayed with herbicide (1 l ha − 1 Round-up Bio) in late August 2019, killing any established plant cover. Afterwards, the middle portion of all sampling sites remained bare, winter rye spill from the previous harvest established spontaneously towards the edges of the sampling sites, and was allowed to grow freely. Areas of bare soil and areas with established volunteers in the CT treatment were chosen and marked ahead of ploughing in late October 2019. Corresponding areas with bare soil and volunteer cover in the NT Table 1 Experimental design. All cores were scanned twice for microtomography, before and after 5 freeze-thaw cycles. Sampling took place in early December 2019, extracting one intact soil core (Ø=100 mm, h=80 mm) per tillage and plant cover combination at a depth of 20-100 mm. One B and one V sample were taken in an area of bare soil and an area established with rye volunteers, respectively, within each tillage treatment. Aboveground plant biomass was clipped and removed before core sampling in all V tillage treatments and blocks. Sampling was repeated in blocks 1, 2 and 4 of the field experiment and are treated as experimental replicates for the remainder of this study (see Appendix A for an overview of the field experiment and sampling sites). Core extraction was carried out by slowly pressing aluminum cylinders into the topsoil with a handheld soil ring sampler and a small rubber mallet, then carefully excavating around the cores and removing them with support from underneath, trimming the soil at the top and bottom with a sharp blade and capping with plastic lids. Immediately after extraction, all cores were transported to the Department of Agroecology at Aarhus University Foulum, Denmark. Table 1 illustrates the sampling design.

Freeze-thaw experiment
Preparation for both the Control and FT treatments consisted in slowly saturating the cores on ceramic pressure plates with a 0.05 M Ca 2+ simulated soil solution, followed by draining on ceramic pressure plates. To simulate wet winter conditions typical of temperate climates during both FT and Control treatments, the bottom pressure head during draining was adjusted to − 10 hPa using a controlled suction system.
Once drained, all cores underwent the Control treatment, which consisted of storage at 5 • C in a temperature-controlled room for 14 days. Following Control treatment, all cores were drained to a bottom pressure head of − 100 hPa on ceramic pressure plates and taken to the first session of µCT scanning.
After scanning, all cores were prepared for FT treatment by saturating again in simulated soil solution and draining to a bottom pressure head of − 10 hPa. The FT treatment consisted of five cycles, where the samples were placed in a laboratory freezer for 24 h, then allowed to thaw at 5 • C in a temperature-controlled room for 48 h. During FT, the cores were placed in a polystyrene foam container open at the top and filled with vermiculite up to about 10 mm from the upper edge of the cores in order to approximate a realistic boundary between the soil and low-temperature air. The freezer used for FT was kept unpowered and at room temperature during the thaw intervals of the freeze-thaw cycles, and set to minimum power during the freezing portions of the cycles. A temperature probe embedded in a mock soil core during the FT treatment registered the freezing phase transition between 3 and 7 h after the samples were placed in the freezer. A minimum temperature of approximately − 13 • C was measured 11 h after the start of the freezing cycle. The thawing phase transition initiated 4 h after removal of the samples from the freezer, and lasted for an additional 18 h.
After the FT treatment, all samples were drained to a bottom pressure head of − 100 hPa on ceramic tension plates and taken to a second session of µCT scanning.

X-ray µCT and image segmentation
Drained intact cores were scanned using a Nikon XT H 225 X-ray µCT scanner at the Department of Soil System Science, Helmholtz Centre for Environmental Research in Halle, Germany. The scanning was carried out with a beam energy setting of 160 kV and 360 µA, with a 0.7 mm Cu filter to reduce beam hardening. The resulting 2749 X-ray projections were reconstructed using the Nikon CT Pro 3D software into 3D tomograms with a volumetric picture element (voxel) size of 58 µm and a stack size of 1750 * 1750 * 1900 voxels.
3D image processing and analysis was carried out using the FIJI distribution of ImageJ, version 1.53c (Schindelin et al., 2012). The tomograms were cropped into cylindrical volumes of interest (VOIs) 1500 pixels (~90 mm) in diameter and 1166 slices (~70 mm) in height, in order to exclude the top, bottom and edges of the cores, where coring artifacts might have occurred. Visual inspection of recognizable image landmarks was used to ensure that the cropped tomograms of corresponding Control and FT samples encompassed the same volume in order to make quantitative comparisons meaningful.
Classification of macropore features, or segmentation, in the 3D images was carried out by segmentation of the 8-bit grayscale values into binary black and white, where white represents pores, and black represents the soil matrix. Vertical intensity drift (i.e. overall brightness variation between slices of one stack) was corrected by normalizing the mean gray value of each slice to the overall mean gray value of the stack, similarly to the method used by Iassonov and Tuller (2010) and Schlüter et al. (2016). Segmentation was preceded by recalculating the stack histogram so that it occupied the entire range of grayscale values available for an 8-bit image (0− 255), followed by clipping the high and low intensity values to increase contrast. The resulting image was de-noised with a 2D median filter with a radius of 2 pixels, and segmentation was then carried out using the Otsu algorithm (Otsu, 1979).
Post-processing consisted of a 2-cell opening algorithm followed by 1-cell erosion of pore voxels in order to reduce the number of single voxels misclassified as pore space due to the noisy background in the matrix portions of the images. We further controlled for misclassification of matrix voxels as pore voxels by calculating the intersection between the binary and original grayscale images (Fig. 1A) and visually comparing images and histograms against those of the original grayscale images (Fig. 1B), specifically, looking for excessive overlap between pore voxel grayscale values and matrix voxel values.

Pore diameter classification
Using the ImageJ tool 3D local thickness maps, which calculates the diameter of the maximum non-redundant fully-inscribed sphere for every pore voxel in the binary segmented stacks (Hildebrand and Rüegsegger, 1997), pore voxels were further classified according to local pore diameter. We classified the images into three pore diameter ranges ( Fig. 2): between 2 voxels, the minimum volume required by the Local Thickness algorithm, and 5 voxels (<300 µm), between 5 and 17 voxels (300-1020 µm), and larger than 17 voxels (>1020 µm). The pore diameter range < 300 µm represented the very fine macropore volume expected to be filled with water during the freeze-thaw experiments (i.e. at a pressure head of − 10 hPa). It is noted that structural elements smaller than 4-5 voxels are systematically underestimated due to image processing (Vogel et al., 2010). The 300-1020 µm range represented the remainder of the very fine macropores, which were not expected to be filled with water during freeze-thaw. Pores > 1020 µm represented all fine, medium and coarse macropores, which we expected would be most sensitive to consolidation or fragmentation of the matrix. For simplicity, however, we will refer to the > 1020 size range only as coarse macropores. This process caused some edge voxels originally belonging to large diameter pores to be re-classified as < 300 µm, given the irregular shape of many large soil pores. These "nooks and crannies" were not removed from the < 300 µm range, as they potentially remained filled with water during the FT treatment due to capillary adhesion in spite of forming part of a larger soil pore. Finally, we also considered the unclassified images containing pores of all detectable sizes after segmentation, which we refer to as Full Range images.

Structural parameters
Three structural parameters were obtained from the binary images containing all pore sizes (Full Range) as well as the images further classified by pore diameter: total detectable macroporosity (V t , m 3 100 m − 3 ), mean macropore diameter (d m , µm) and the mean Euclidian distance from the matrix to the nearest pore (ED m , µm). The total detectable macroporosity was calculated as the volume of all macropore voxels in a sample's VOI divided by the volume of the VOI in voxels. The mean macropore diameter was calculated as the mean of Local Thickness voxel values weighted by frequency within the VOI of each sample. Finally, ED m was calculated as the mean of the Euclidian Distance Transform voxel values weighted by frequency within the VOI of each sample. Additionally to overall parameters for the full VOIs, profiles of V t , d m and ED m were calculated from Full Range images by restricting the height of the VOI to a single slice of the image and calculating the parameters in the same manner as for the full cores. This allowed us to examine variations in macroporosity and the effects from FT along the vertical axis of each core. This is important as the cores were primarily in contact with the environment at the top during FT, which caused uneven cooling, freezing, warming and thawing during each cycle.

Macropore geometry parameters
To describe overall differences in the shape and complexity of macropores in the soil cores, we collected a series of summarized geometric parameters based on the Full Range segmented images. Individual contiguous macropores in each core were identified and analyzed using the tool Particle Analyzer, which is part of the BoneJ plugin in Fig. 2. Example of macropore diameter classification. Local pore diameters are calculated for every pore-space voxel in the initial binary image (top left). In the resulting Local Thickness image (full range, top right), every voxel has a numeric value equal to the local diameter of the pore (in µm). By classifying local thickness values, images containing only pores with a 3D local diameter of 2-5 voxels (<300 µm), 5-17 voxels (300-1020 µm) and > 17 voxels (>1020 µm) were created to analyze the effect of experimental treatments on different pore diameter ranges.
ImageJ (Doube et al., 2010). From the Particle Analyzer results, the mean volume (mean Vol) and mean surface area (mean SA) of individual pores were calculated for each core. The average macropore density (number of pores per mm 3 ) in each core was calculated by dividing the total number of individual macropores by the volume of the VOI in mm 3 . Additionally, the topographical complexity of the macropores in each core was analyzed using the mean Euler characteristic χ m (dimensionless): where χ i is the Euler characteristic (i.e. one minus the number of loops plus the number of fully-enclosed holes) of each contiguous macropore in an image and n is the total number of macropores, as calculated by Particle Analyzer. In general, greater χ m values reflect a predominance of isolated pores and pore networks with simpler geometries, whereas lower (or negative) values indicate greater topographical complexity with a predominance of branched and interconnected pore networks (Vogel and Roth, 2001). Macropore length and connectivity were analyzed by creating a wire-model of every macropore in each core (tool Skeletonize, Lee, Kashyap and Chu, 1994). From the resulting wire networks, the number of interconnected branch networks, as well as the number and length of their component branches, were extracted using the tool Analyze Skeleton (BoneJ plugin, Doube et al., 2010). The average number of component branches was calculated by dividing the total number of branches by the number of connected pore networks.

Statistical analysis
All data handling, visualization and statistical analysis was carried out using R (R Core Team, 2021) version 4.0.3 (released 10/10/2020). The four field management combinations were arranged as four levels in ascending order of inferred structural strength: CT-B < CT-V < NT-B < NT-V. The effects of FT on V t , d m and ED m for the different levels of field management (and inferred soil strength) were analyzed using linear mixed-effects models (package lme4, Bates et al., 2015) with freeze-thaw and management as main effects. The core identifier was added as random effect given that the FT treatment consisted of repeated µCT scans on the same core, before and after freeze-thaw. Post-hoc multiple pair-wise comparisons using Sidak error corrections (R package multcomp, Hothorn, Bretz and Westfall, 2008) were carried out to aid the visualization of statistical results. Two different contrasts were used for pairwise comparisons, firstly comparing all management levels and freeze-thaw levels within a pore size range, and secondly, restricted to comparing Control and FT levels within each single management level and pore size range. This was meant to aid interpreting statistical differences caused by management, and those associated to freeze-thaw within each level of management, respectively.
The macropore geometry parameters were analyzed collectively via principal component analysis (PCA). The values of all principal component with eigenvalues > 1 were extracted, and the effects of FT and management on principal component values were analyzed using linear mixed models with FT and management as main effects and core identifier as random effect.

µCT-visible structure
Detailed visual examination of all µCT images revealed clear structural differences between field management practices, particularly tillage. The soil in both CT-B and CT-V cores presented many loose fragments a few tens of millimeters in size, especially towards the top of the cores, and the pore space presented striking irregular voids which made up much of the visible coarse macroporosity. In NT-B and NT-V cores, on the other hand, the matrix was less heterogeneous, and the visible macroporosity was more strongly dominated by tube-like pores, likely burrows made by soil macrofauna (earthworms and other invertebrates no larger than a few millimeters in diameter) and plant roots. Fig. 3 (left and center) illustrates some of the more visible structural differences between CT and NT. Structural differences due to plant cover (and the corresponding fine root growth) within the same tillage practice were not easily identifiable upon visual examination. This is in part due to the fact that plant roots tend to follow pre-existing pores in the soil, and that root tissues are usually poorly discernible from the soil matrix in X-ray imaging at this resolution.
Freeze-thaw caused appreciable changes in the soil structure within the cores, which are more striking for CT-B and CT-V than either of the No-Till management practices. Comparing Control and FT images aligned using image landmarks, two effects from freeze-thaw stand out: crumbling of soil along large inter-fragment voids, and formation of visible cracks. Fig. 3 (center and right) presents a few examples of the most common changes observed between Control and FT images. Crumbling was present almost exclusively in CT-B and CT-V cores and most prominently in CT-B. The exception to this occurred in two NT cores, one NT-B and one NT-V, where crumbling within large earthworm burrows was caused, most likely, by movement of the earthworms trapped in the core during freeze-thaw (see Appendix B). Outside the cores affected by earthworm movement during freeze-thaw, FT appears to have caused the least visible changes in the NT-V cores, where both structural deformation and cracking were largely inconspicuous.

Soil structure profiles
Macroporosity profiles (Fig. 4, top) show high variability, with differences between layers of the same core greater than 10 m 3 100 m − 3 , a visible gradient of decreasing V t with depth for CT-B and CT-V and generally lower values in NT-B and NT-V. Mean pore diameter profiles also show a high degree of variability with depth, strongly influenced by the activity of earthworms, as their burrows resulted in large local d m increases in several CT-V, NT-B and NT-V profiles (Fig. 4, middle). Additionally, d m profiles in CT-B are noticeably different from one another, despite belonging to the same management treatment. Finally, there was noticeable variation among ED m profiles of the same management level in CT-B and CT-V, particularly towards the bottom of the cores. In NT-B, this variation among profiles is less noticeable, and greater towards the top of the cores, while in NT-V the variation between profiles is minimal (Fig. 4, bottom).
In spite of visible local differences between FT and Control treatments of corresponding cores, the overall effects of FT are not clearly discernible in the V t , d m or ED m profiles (Fig. 4, different colors). As mentioned above, the movement of earthworms during freeze-thaw can cause disruption on coarse macropore walls. Since earthworms, due to their low density compared to the soil matrix, are classified as voids during segmentation, the movement of an earthworm and subsequent fragmentation of its burrow due to crumbling would have caused a reduction in local thickness for a large number of voxels. This is visible in two cores (NT-B and NT-V), where the fragmentation of large voids due earthworm activity and crumbling resulted in strong local decreases in d m . Importantly, the local changes in question are not apparent in V t or ED m , and seem to have exclusively affected layer-by-layer mean pore diameter values.

Overall soil structure
Total detectable macroporosity in Full Range (diameters 58-10440 µm) was on average around 7 m 3 100 m − 3 for NT and 17 m 3 100 m − 3 for CT, which corresponds to 17 % and 35 %, respectively, of total soil porosity determined from bulk density (Fig. 5). Mean Euclidian distance to the nearest macropore ranged between approximately 500 and 1000 µm in Full Range, with a minimum of 478 µm in the < 300 µm macropore size range. Finally, mean macropore diameter values ranged between 830 µm and 2000 µm in Full Range, with a largest measured macropore diameter of 11530 µm. Total macroporosity, mean macropore diameter and mean Euclidian distance to the nearest macropore data (Fig. 5) all showed a high degree of variation among cores of the same management and freeze-thaw treatment levels. However, certain trends are also visible, with V t decreasing and ED m increasing visibly with increasing expected structural strength (CT-B<CT-V<NT-B<NT-V) associated with field management treatments at all macropore size ranges.
The effects of 5 freeze-thaw cycles on soil structure were small compared to the overall variation among cores. However, differences in V t , d m and ED m between Control and FT repeated measurements were highly consistent among replicates of most management levels. In particular, FT caused consistent decreases in V t in Full Range and the > 1020 µm range and small increases in V t for the size range < 300 µm, as well as a consistent trend of decrease in d m in all size ranges and in ED m in all size ranges except 300-1020 µm.

Total detectable macroporosity
Field management had a significant effect on V t for all pore diameter ranges except > 1020 µm (See Appendix C for detailed linear mixed model statistics). Multiple pairwise comparisons in Full Range showed a Fig. 5. Symbols representing total macroporosity (V t ), mean macropore diameter (d m ) and mean Euclidian distance to the nearest macropore (ED m ) obtained from µCT images of intact soil cores before (cyan) and after 5 freeze-thaw cycles (black). The cores belonged to four management treatments consisting of a combination of either conventional inversion tillage or no-till (CT and NT, respectively), and either bare soil fallow or volunteer plant cover (B and V, respectively). The whiskers represent the 95 % confidence interval of estimates obtained from linear mixed effects models using management and freeze-thaw as main effects and core identifiers as random effects.
trend of decreasing V t with field management in the order CT-B>CT-V>NT-B>NT-V, where V t was 13.3 m 3 100 m − 3 lower in NT-V compared to CT-B in the Control treatment (SE=2.74, p = 0.016). The same was observed in the 300-1020 µm range, where V t was 4.5 m 3 100 m − 3 higher in CT-B compared to NT-V (SE=0.80, p = 0.007), and 3.4 m 3 100 m − 3 in CT-B compared to NT-B (SE=0.80, p = 0.038), all in the Control treatment. This trend is also suggested by the raw data for the < 300 µm or > 1020 µm ranges (Fig. 5, left-most column), however, multiple pairwise comparison did not show any significant differences between management levels.
Pairwise comparisons between Control and FT restricted to corresponding management levels revealed that FT caused a small but significant decrease of 0.4 m 3 100 m − 3 in V t for all management levels in Full Range (SE=0.15, p = 0.026). In the < 300 µm range, however, FT caused significant increases in V t of 0.4 m 3 100 m − 3 for CT-B (SE=0.09, p = 0.002) and 0.2 m 3 100 m − 3 for NT-B (SE=0.09, p = 0.04), with a similar non-significant trend for CT-V. At the range > 1020 µm, FT caused a reduction in V t of 1.1 m 100 m − 3 for CT-B (SE=0.16. 0 <0.001) and 1.0 m 100 m − 3 for CT-V (SE=0.16, p < 0.001) but caused no significant differences for NT-B and NT-V. Finally, FT had no effect on V t in the range 300-1020 µm, regardless of management treatment.

Mean pore diameter
Mean pore diameter, d m , was not significantly affected by field management at any pore size range except for < 300 µm where it followed a trend of decreasing mean diameter in the order CT-B>CT-V>NT-B>NT-V. Pairwise comparisons show d m in the > 300 µm range to be 9.4 µm greater in the Control treatment of CT-B compared to NT-V (SE=2.40, p = 0.050), but no other significant differences within the same freeze-thaw treatment were found.
On the other hand, freeze-thaw resulted in significant decreases in d m at all macropore size ranges, as shown by multiple pairwise comparisons restricted to FT and Control treatments within the same management level. In Full Range, FT resulted in a decrease in d m of 132 µm for all management treatment levels (SE=132, p = 0.005). In the < 300 µm range, FT resulted in decreases of 1.4, 1.5, 1.8 and 1.2 µm for CT-B (SE=0.48, p = 0.020), CT-V (SE=0.48, p = 0.016), NT-B (SE=0.48, p = 0.005) and NT-V (SE=0.48, p = 0.033), respectively. Similarly, FT resulted in a decrease in d m of 4.52 µm for all management levels (SE=0.82. p < 0.001) in the range 300-1020 µm and a decrease of 133 µm for all management levels (SE=54.4, p = 0.033) in the range > 1020 µm.

Euclidian distance to nearest macropore
Mean Euclidian distance from the matrix to the nearest macropore was 695.4 µm for the Control treatment in Full Range and was not significantly affected by field management. Overall, ED m was significantly affected by management only at the 300-1020 µm macropore size range. Multiple pairwise comparisons, however, did not show any significant differences between individual management levels in the 300-1020 µm range.
On the other hand, freeze-thaw had a significant effect on ED m at all macropore size ranges except for 300-1020 µm. Pairwise comparisons show a decrease in ED m of 37.8 µm due to FT for all management levels in Full Range (SE=11.30, p = 0.006), as well as a decrease of 35.3 µm for all management levels in the range < 300 µm (SE=10.5, p = 0.006). At the range > 1020 µm, however, pairwise comparisons showed significant increase in ED m due to FT of 52.5 µm for all management levels (SE=19.1, p = 0.019).

Macropore geometry
Principal component analysis found two linear combinations of macropore geometry parameters (principal components) which described the variance in pore structure among soil cores better than any single measured parameter (i.e. eigenvalues >1). Together, these two principal components (PC1 and PC2) accounted for 95.1 % of the variation among individual cores as described by seven macropore structure parameters (Table 2 and Fig. 6). The values for all macropore structure parameters are shown in Appendix D.
The PC1-PC2 plane highlights the general differences between samples in different treatments, providing an overview of the effects of till, plant cover and freeze/thaw on macropore geometry. Fig. 6 shows that NT cores are clearly differentiated from CT, with higher values in both the PC1 and PC2 axes. This can be interpreted as samples from CT   and NT having substantially different structural properties, primarily in terms of number of pore branches per macropore network, mean macropore surface area and mean Euler characteristic. Additionally, there is a consistent displacement of CT cores in the PC1-PC2 plane due to FT, where FT samples show a tendency to greater macropore density and lower mean branch length, compared to Control cores. This consistent displacement is not as clear for NT samples. Finally, no clustering or consistent displacement can be observed between B and V samples, suggesting that the presence of volunteers did not affect macropore geometry to a large extent. The rotation coefficients of PC1 and PC2 reveal certain correlations between macropore structure parameters among cores (Table 2). PC1, which encompasses 72.5 % of observed variance, shows positive correlations between mean volume and surface area of macropores, mean branch length, and mean number of branches per macropores pore network. PC1 also shows a negative correlation between the mean Euler characteristic of macropores and the mean number of pore branches per pore macropore network, with only a weak relation to macropore density. PC2, describing 22.7 % of the observed variance, was dominated by macropore density, with a weak negative correlation to mean Euler characteristic.
Linear mixed models of principal components (Table 3) show a nearsignificant trend of lower PC1 values in NT cores compared to CT, and significantly lower values in FT compared to Control cores. Additionally, we found significantly lower PC2 values in NT compared to CT, and significantly higher values in FT compared to Control, with a significant interaction between plant cover and freeze-thaw. According to this interaction, the effect of FT on PC2 is significantly lower in V, compared to B. Given the rotation coefficients of PC1 and PC2 (Table 2), lower PC1 values indicate a decrease in the size (volume, surface area, mean length) and complexity (greater Euler characteristic and fewer branches per network) in NT compared to CT, and in FT compared to Control (although with a smaller effect size). Lower PC2 values due to NT, in turn, indicate primarily lower macropore density and overall complexity. Finally, greater PC2 values in FT compared to Control indicate an increase in macropore density and complexity, and this increase is significantly smaller in V cores compared to B cores due to a significant negative interaction between plant cover and freeze-thaw.

No-till vs conventional tillage
It is well investigated that NT affects the structure of agricultural soils, mainly increasing bulk density, mechanical strength, penetration resistance, wet aggregate stability and soil water retention (e.g. Blanco-Canqui and Ruis, 2018;de Oliveira et al., 2021). However, there is less agreement regarding the effects of NT on soil macroporosity, macropore structure and macropore function. It is, for instance, often thought that increased earthworm activity in NT soils leads to greater macroporosity, but Schlüter et al. (2018) found lower µCT-detectable macroporosity, lower saturated hydraulic conductivity and lower pore connectivity after 25 years of NT, in spite of an enhanced earthworm presence. Vogeler et al. (2009), on the other hand, found higher air permeability and saturated hydraulic conductivity in chiseled soil compared to ploughed soil, in spite of very small differences in macroporosity derived from water retention curves. Finally, Abdollahi et al. (2014) compared a series of macropore-related parameters derived from medical-grade tomograms of intact topsoil cores corresponding to direct sowing, harrowing and mouldboard ploughing. They found no significant differences due to field management except when catch crops were used in combination with direct sowing. Here, they found significantly fewer pore network branches, junctions and endpoints, reflecting lower macropore complexity. Our own results show, relative to CT-B, a significant decrease in detectable macroporosity in NT-V in Full Range, and in both NT-B and NT-V in the 300-1020 µm size range, as well as a decrease in mean macropore diameter in < 300 µm in NT-V. This suggests that NT soils contain a greater volume of very fine macropores and a reduced volume of coarse macropores, resulting in lower overall macroporosity compared to CT. However, further explanation of the effect of NT on functional soil properties such as saturated hydraulic conductivity must also take into account the geometry and spatial characteristics of the macropores. Indeed, when studying intact cores of a Brazilian Oxisol Pires et al. (2019) note that in spite of significant differences in total porosity and number of µCT-detectable macropores, pore tortuosity, connectivity and degree of anisotropy were not strongly affected by no-till management, as compared to conventional tillage. Furthermore, the authors didn't find any significant difference due to tillage management in the contributions of pores of different shapes to the total µCT-detectable porosity, suggesting that pore function was similar between no-till and conventional till in their particular case. In contrast, Miranda-Vélez et al. (2022) found important differences in macropore flow and macropore solute transport between CT and NT management in cores similar to those studied here, in turn suggesting that structural (and thus functional) changes caused by tillage are also dependent on soil type.
When examining µCT scans of intact cores from a long-term rotation and tillage experiment, Munkholm et al. (2013) found significantly lower macropore surface area, fewer macropore network junctions and branches in NT compared to CT, with no significant differences in number of pores or pore end-points due to tillage. This resembles our findings from PCA analysis, which indicate that individual macropores in NT, although greater in number, tend to be smaller and less complex than in CT. On the other hand, Garbout et al. (2013) examined medical scanner tomograms of large intact soil cores (200 × 200 mm, voxel size 390 × 390 × 600 µm) and found fewer macropore networks, together with fewer branches and junctions, in NT compared to CT. They also report a higher degree of anisotropy and greater mean macropore branch length in NT compared to CT, indicating longer and more organized macropore segments. These findings support our observations of lower macropore density in NT but seem to contradict our observation of lower mean branch length, as well as the inferences made from V t and d m results, pointing at a greater number of very fine macropores. These contradicting results can be reconciled by considering overlapping effects from NT at different macropore size scales. We propose that the soil structure in NT presents fewer coarse macropores, which are longer and less complex than in CT. Given that surface area and volume generally scale with the square and cube power of a macropore's diameter, respectively, a reduction in the number of coarse macropores results in significantly lower volume and surface area measurements. However, µCT images also contain a much greater number of fine and very fine macropores which are short, simple in shape and mostly randomly oriented. Thus, when examining an array of different µCT-derived measurements as it's commonly done, the effects of tillage on the coarse macropores are likely confounded by different, or nonexistent, effects on the fine and very fine macropores. Remarkably, it is likely that results obtained with medical equipment of lower resolution are more consistent between studies, as they ignore much of the "noise" from the fine and very fine macropores in the sample.
Lack of agreement between macropore geometry results related to different size classes can be addressed in image post-processing by classifying the samples in different size classes as done here for the general macropore parameters. However, such post-processing procedures carry the inherent risk of removing small-diameter features connecting larger macropore networks, and excessively smoothening the surfaces of large irregular features, potentially altering their geometric properties. On the other hand, removing large local diameter voxels to isolate small-diameter macropores, such as done here for the < 300 µm and 300-1020 µm pore size ranges, also isolates the irregularities on the surface or large macropores. This new host of spurious new pores, while a real component of the original macroporosity, does not reflect the geometric properties of actual individual pores (e.g. complexity and length), and thus should be considered as artifacts in measurements of non-extensive properties of the soil macroporosity. Therefore, while potentially helpful, classifying macropores by size is not a without risks when investigating the spatial and geometrical features of the soil structure and should be carried out with caution.
The significance of these and other similar results in the context of macropore function, particularly regarding water movement and solute transport, is not straightforward. Mesocosm laboratory lysimeter experiments on intact cores of the same soils as those analyzed here showed a greatly increased macropore water flow in NT compared to CT, as well as a noticeable bypass effect where native resident nitrate was preserved in NT cores after 200 mm of intense simulated precipitation (Miranda-Vélez et al., 2022). This indicates that there is potential for an increased macropore function in NT, even as the overall macroporosity, complexity and average length are all reduced compared to CT. Here, it is important to consider that relatively few long, well-connected macropores have a strong effect on the hydraulic behavior of a given volume of soil, and at the same time have a small or negligible effect on many summarized parameters commonly obtained from µCT imaging, mainly due to the large number and diversity of elements visible in a microtomogram.

Plant cover
We did not observe any significant main effects from plant cover during PCA analysis, in spite of clear indications that plant cover played a role in the differences observed in V t and ED m due to management. It is possible that the relatively low density of self-sown grain spill, together with a short growing period before sampling, reduced the effect of volunteer plant cover in this experiment. It is also plausible that differences in V t and ED m due to plant cover in the < 300 µm size range were at least in part driven by misclassification of macropores colonized by roots as matrix, a problem that is known in X-ray microtomography (e.g. Phalempin et al., 2021). However, we did observe a significant interaction between plant cover and freeze-thaw in the PCA of macropore geometry parameters, which is discussed later on.

Freeze-thaw
Freeze-thaw had a two-fold effect on the macropore system of the intact soil cores. On the one hand, µCT images taken after 5 freeze-thaw cycles show an overall greater macroporosity and decreased mean pore diameter in the 300 µm size range. Having in mind that the limit of 300 µm is very close to the practical detection limit of the imaging and segmentation process employed here, we suggest that FT caused the expansion of a significant number of very fine macropores at the edge of detection and classification limits, making them visible to µCT imaging and analysis. On the other hand, repeated freeze-thaw cycles caused further fragmentation of large pre-existing fragments that defined the largest macropores in our samples (Fig. 3). This disruption along large pores and voids led to the interruption of coarse macropores matrix fragments, locally decreasing coarse macropore diameter and causing a fraction of the soil macropore volume to fall below the > 1020 µm threshold. This, in turn, caused a significant reduction in detectable macroporosity and an increase in ED m in the > 1020 µm size range.
Together, the two effects from FT resulted in an overall decrease in mean pore diameter while an overall change in macroporosity was more difficult to quantify, given the opposing effects at different pore sizes. It is only when dividing the detectable pore space in different size ranges that the two opposing effects of FT on macroporosity became clearly noticeable. Finally, the expansion of very fine macropores in the < 300 µm range had a diminished effect in NT-B, while the decrease in coarse microporosity (V t in the >1020 µm range) was not significant for either NT-B or NT-V. This supports the hypothesis that soils exhibit decreasing susceptibility to the effects of FT with increasing matric cohesion and structural strength, with the greatest disturbance in CT-B and the least disturbance in NT-V. This interaction is further discussed later on.
Macropore geometry parameters suggest similar effects by freezethaw to those found in macroporosity and mean diameter. Together with an increase in macropore density and total number of pore networks, FT caused a decrease in mean macropore branch length, volume and surface area, which also points at a larger number of very fine macropores and/or collapse and filling of large coarse macropores. These findings are in line with those by Leuther and Schlüter (2021), who found that repeated freeze-thaw cycles caused a significant reduction in both µCT-detectable pore volume and mean pore diameter, and that this effect was dependent on the initial structure of the soil. Indeed, the authors found that the macropore structure of an undisturbed grassland akin to the no-till topsoil examined in our study was less strongly affected by freeze-thaw than that of re-packed clods resembling the surface of a ploughed field. They postulate that fragmentation of re-packed soil clods creates new small pores and cracks at the edge of detectability by µCT techniques, while at the same time reducing the diameter and volume of existing large macropores. When examining smaller regions of undisturbed soil cores at a higher resolution, the authors also note that the opposing nature of these two processes makes it difficult to observe the effects of freeze-thaw reflected in overall porosity and pore structure parameters.
Using synchrotron µCT imagery of considerably higher resolution (3.25 µm per voxel) in extracted soil aggregates 5-7 mm in size, Ma et al. (2021) observed a strong increase in porosity, porosity of equivalent pore size > 100 µm, and pore connectivity over the course of 20 freeze-thaw cycles. Given the increased image resolution and the separation of individual soil aggregates for scanning, these results show that freeze-thaw causes the appearance and expansion of very fine macropores more clearly than images using whole intact cores at lower resolution. Thus, the increase in the volume of macropores > 100 µm within individual millimeter-size aggregates corresponds to the increase in volume of pores close to or below the detection range of µCT proposed here. The authors also report a significant reduction in aggregate stability with repeated freeze-thaw cycles, which is compatible with the mechanism of matrix fragmentation as a cause of reduction in volume and increase in number of large-diameter macropores. This effect, however, is not observable at the very fine spatial scale nor when the soil has been broken into aggregates a few millimeters in size.
We consider important to remark that studies which focus on highresolution µCT imaging of individual aggregates, microaggregates or small re-packed soil samples, are most likely to observe the effects of freeze-thaw cycles in the very fine macropores, and miss effects on the larger-scale structure. Meanwhile, studies focusing on larger intact and undisturbed samples will, due to the accompanying limitations in image resolution, miss changes taking place in the smaller end of the macropore size scale. Thus, it is vital to keep in mind that seemingly contradictory results between different studies might very well correspond to concurrent phenomena observed from different ends of the size spectrum.

Interaction between management and freeze-thaw
The combined effect of management and freeze-thaw on soil structural properties is still poorly understood. However, Soane et al. (2012) conclude in a review that soils gain mechanical strength and aggregate stability in the laboratory with no-till, which would make them less susceptible to the effects of FT. Indeed, Ma et al. (2019) found that well-structured soils are less susceptible to the effects from FT, particularly the fragmentation of the soil matrix. This is supported by our findings, where the effect of FT on macroporosity was reduced by a statistically significant interaction with NT both at the < 300 µm and > 1020 µm ranges. Additionally, V caused a statistically significant reduction to the effect of FT on macropore geometry parameters, as described by PCA. This indicates that both reduced tillage and the presence of volunteers ameliorated the increase in detectable very fine macropore density and the disruption of coarse macropores caused by FT.
The interactions between NT and FT, as well as between V and FT, can be explained by increases in overall matrix cohesion with reduced tillage and plant root growth, due to e.g. increased organic C inputs, decreased SOM mineralization and preservation of root and fungal hyphal networks. Conversely, the matrix in CT and B is likely more susceptible to deformation under the mechanical stresses caused by ice expansion and movement of water and air in narrow water-filled soil pores. Additionally, air and water movement are likely restricted within the irregularly compacted fragments produced by tilling machinery in CT, which would also experience heterogeneous freezing along fragment edges exposed to large air-filled voids. Thus, the reduced matrix cohesion together with pressure build-up from entrapped air and water upon freezing would result in greater expansion of small fissures, crumbling of fragment edges and the filling of the spaces between them. These effects are of particular relevance in the context of reducing nitrogen emissions in agriculture, as both the use of winter cover crops and reduced tillage are considered leading strategies for reducing post-harvest (i.e. autumnwinter) nitrogen leaching (Constantin et al., 2010;Thapa et al., 2018).

Conclusions
Our findings suggest that freeze-thaw has a two-fold effect on the structure of agricultural topsoil. On the one hand, freeze-thaw causes the soil matrix to fragment, partially filling or dividing coarse macropores and reducing their overall volume. On the other hand, new very fine macropores are formed with freeze-thaw, while previously undetectable macropores expand and become increasingly discernible in X-ray µCT imagery. These two effects counteract one another when analyzing the detectable macropore volume as a whole, resulting in small overall decreases in macroporosity, mean macropore diameter and mean Euclidian distance to the nearest macropore in the matrix. It is only when separating the very fine and coarse macropores in µCT imagery, that these opposing effects become apparent. Furthermore, freeze-thaw had a larger effect on more poorly structured soils, i.e. conventionally tilled and/or kept without plant cover after harvest. This suggests that no-till soils, particularly those with plant cover, are less susceptible to disturbance by freeze-thaw.
The effects of tillage and plant cover on soil structure are difficult to interpret via differences in a number of µCT-derived geometric parameters. Here, we found that NT soils generally have lower macroporosity than CT soils and show a significant tendency for smaller diameters in the very fine macropore class. Additionally, we found that macropores in NT are fewer in number, are composed of shorter branches and are less complex than in CT. However, we infer that these results are greatly confounded by the overlapping effects of relatively few, complex and voluminous coarse macropores, and a much greater number of small low-complexity fine and very fine macropores. Therefore, after comparison of our experimental results with other reported findings, we propose that the macropore volume of NT is in general smaller than that of CT, composed of fewer, longer and less complex coarse macropores. Unfortunately, the effect of tillage on the fine and very fine macropores remains more difficult to determine. Further work is required to clarify the contributions from different macropore size classes to summarized geometric µCT-derived parameters and, further, to establish their contributions to functional properties of the soil's structure as a function of tillage.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Appendix B
See Fig. B1 here.   Fig. B1. Earthworms trapped inside two cores (NT-B and NT-V treatments) caused crumbling within their burrows during the freeze-thaw experiment, reducing the local pore diameter of a large number of voxels, as determined by the local thickness method. While these disturbances affected only one very large macropore (the burrow itself), their effects are visible in the mean pore diameter profiles shown in Fig. 4.
ANOVA table of linear mixed model effects for total macroporosity (V t ), mean macropore diameter (d m ) and mean Euclidian distance to the nearest macropore in the matrix (ED m ). Models for each parameter were reduced as much as possible by removal of non-significant interactions, with certain exceptions for near-significant trends.