Spatial Variability of Preferential Flow and Infiltration Redistribution along a Rocky-Mountain Hillslope, Northern China

: Rock fragments in soil strongly increase the complexity of hydrological processes. Spatial variability of preferential flow and infiltration characteristics, especially along a rocky-mountain hillslope are poorly understood. In this study, five rainfall – dye tracer experiments were performed in the rocky Taihang Mountains, northern China, to investigate the spatial variability of preferential flow and infiltration redistribution on different hillslope positions. Tracers were used to distinguish macropore flow and actual water flow patterns, and preferential flow indices and spatial non – uniformity of the infiltration redistribution were calculated using image analysis. Results showed increasing trends in the dye coverage, maximum infiltration depth, and steady infiltration rate with increased hillslope position, with a preferential flow fraction of 0.10, 0.11, 0.15, 0.29, and 0.26 for the bottom – , down – , mid – , upper – , and top – slope positions, respectively. With increased hillslope position, the spatial non – uniformity of the infiltration redistribution gradually increased in orthogonal and parallel directions to the stained section, and was supported by the fractal dimensions. Positive (gravel mass ratio, saturated water content, altitude, hydraulic conductivity and roots) and negative (bulk density and clay content) impacts on preferential flow and infiltration redistribution were quantitatively emphasized. The characteristic and mechanism of infiltration process were further identified along a rocky-mountain hillslope.


Introduction
Infiltration is an important process within the watershed water cycle, in which the soil absorbs, regulates, and transfers water to the soil layer under the combined actions of gravity and soil water forces [1,2]. Preferential flow plays an important role in the infiltration process. The non-equilibrium flow phenomenon, which penetrates through and bypasses most of the soil matrix and pores along some specific paths, widely exists in heterogeneous, anisotropic permeable media, especially in mountain areas with a complex geological structures [3,4]. In comparison to matrix connective flow, preferential flow is mainly driven by gravity and greatly decreases the contact time and area between the seepage liquid and solid media [5]. It is vital to have an effective buffer for the peak and duration of surface runoff [6,7] for the immediate replenishment of groundwater and required chemical composition [8,9], and to reduce the risks of soil erosion and landslides [10].
Many factors contribute to the formation of preferential flow, thus leading to various patterns of manifestation. The main and current researches are related to macropore flow and unstable finger flow, with other preferential flow patterns also given attention due to their accompanying hydrological and environmental geological problems [11,12]. The corresponding studies on preferential flow in the vadose zone can be traced back to 1864, while Darcy's law, Richard's equation, etc. have dominated subsequent theories of seepage and infiltration. It was not until the 1870s that the movement of heterogeneous water in macropores challenged the traditional homogeneous theory [13]. Beven et al. [4] discussed the relationship between soil pores and preferential flow, and initiated the subsequent preferential flow research boom. Much research has been undertaken regarding the basic theory of flow development, morphological characteristics, influencing factors, and model simulation of priority flow [13][14][15][16]. Peters [7] studied the hillslope hydrological process of subsurface flow in Shield forest basin using a hydrological test and geochemical tracer method, and found that most water fluxes in the forest soil were priority flow that did not comply with Darcy's law. Some scholars believe that most of the priority channels would not be activated until the soil reaches saturation [9,17,18]. Alaoui [19] found that large pores could transport approximately 74-100% of soil water with 0.23-2% of the total soil volume, and that different combinations of preferential types led to different water transfer rates. Studies on the influencing factors of preferential flow have mainly focused on vegetation type and plant roots [5,6,[20][21][22], soil structure [23,24], soil physical properties, and insect activities [25,26].
Rock fragments in stony soil also increase the complexity of soil hydrological processes [27]. In the past 20 years, scholars have studied how rock fragments in soil influence soil moisture and solute transport from the aspects of gravel type, particle size, shape, content, and their distribution in the soil matrix [28][29][30][31]. Rock fragments can increase the curvature of water flow and restrict water infiltration, although they can also promote water infiltration by increasing the number of macropores depending on their coverage and content in soil, physical and chemical properties, and spatial relationship with the soil matrix [32]. Dadkhah [22] found that the infiltrated water quantity was positively correlated with the coverage of rock fragments without compaction, and negatively correlated with that compacted in soil. However, results by others have not been consistent due to difficult sampling conditions, a variable spatial distribution, and limited field measurements [33]. Moreover, relatively few experimental studies have focused on the identification and quantification of preferential flow in soil with rocky fragments, especially on a hillslope or mountain scale [17,[34][35][36][37][38][39].
Rocky mountains are widely distributed in northern China and cover a total area of ~754,000 km 2 , which includes the Haihe River basin, Huaihe River basin, and part of the Yellow River basin. In the present study, the rocky Taihang Mountains-where the Jing-Jin-Ji Economic District (one of the most developed metropolitan regions in China) is located-was chosen as a representative study region. Rapid socioeconomic development has greatly aggravated the shortage of resources and the deterioration of the ecological environment in this region. As a result, accurate simulation of the hydrologic cycle is vital for local water resources management and ecological health restoration. However, currently, the contribution of preferential flow to the actual water flow infiltration and its spatial variation were generally ignored by the mainstream distributed hillside hydrological model [40,41]. Observations of δ 18 O and soil water potential dynamics have demonstrated the contribution of preferential flow for infiltration and groundwater recharge in the Taihang Mountains [42]. Similarly, heterogeneous infiltration and preferential flow have been indirectly proved by studying the process of hillslope precipitation-runoff and the relationship between surface water and groundwater for the Taihang Mountains [39,43,44]. It is therefore necessary to study the spatial variation characteristics and contribution of preferential flow as well as its nonuniform redistribution on a rocky-soil hillslope.
The objective of the present study is to investigate the spatial variation of degree of preferential flow and heterogeneous redistribution of infiltration in soil media with rock fragments on different hillslope positions. For this purpose, a series of dye tracer experiments were performed at the bottom, lower-, mid-, upper-, and top-slope positions of a typical rocky slope in the Taihang Mountains.
Two kinds of dye tracers, Brilliant Blue FCF (hereafter BBFCF) and KI with starch (hereafter KIstarch) were used to quantitatively determine (i) the degree of preferential flow, (ii) the spatial redistribution of infiltration, and (iii) factors impacting preferential flow and infiltration redistribution. Our findings will help researchers to further understand the infiltration process of preferential flow in rocky mountain areas, and also contribute to the improvement of simulation accuracy of montane hydrological models.

Study Area
The study area is the Chongling catchment area, which is located in the eastern part of the Yunmeng Mountains-A branch of the Taihang Mountains approximately 105 km southwest of Beijing, China ( Figure 1). The Chongling watershed (115°21′-115°23′ E, 39°22′-39°25′ N) belongs to the secondary tributary of the northern Yishui River, and covers a total area of ~6 km 2 . There are five main river gullies in the watershed and the average river bed gradient is 23%. The landscape is characterized by low rolling hills and valleys with soil-stone mixtures. The altitude decreases from northwest to southeast and range from 85 m to 300 m. Geologically, there are limestone and marble in the northwest, purple-red conglomerate in the southeast, granitic gneiss weathering with a thickness of 1-3 m in the central area, and quaternary sediments composed of sand, clay, loess, and gravel in the valley area where the groundwater table is shallow at a depth of 2.5-4 m [43,44]. The Chongling watershed has a temperate continental monsoon climate. The annual precipitation in the study area ranges from 217.0 mm to 1004.3 mm, with an average of 641.2 mm. Affected by the East Asian monsoon, rainfall during the rainy season (usually between June and September) accounts for >70% of the total annual rainfall. The region is dry and cold during the winter with an average temperature of −7.2 °C, and hot and humid during the summer with an average temperature of 32.2 °C. The annual evapotranspiration (measured by evaporating dishes) is 1905.9 mm/year, >45% of which occurs between April and June due to the aridity and high temperatures [43,44].
The vegetation in the watershed is relatively well protected with high coverage and abundant vegetation types because the region was chosen as the Qing Dynasty Royal Tombs. Due to measures taken to enclose the hills for natural grassing or afforestation in the upstream hillside area, the vegetation coverage is generally 70-80% and mainly includes communities of Eriophorum brachyantherum and Angelica keiskei with shrubs such as Vitex negundo L. and Ziziphus jujuba Mill. The downstream hillside areas are mainly planted forests and nursery fields that include Platycladus orientalis, Pinus tabulaeformisCarr., Robinia pseudoacacia Linn., and Populus tomentosa Carr. The area of cultivated land is lesser and is mainly distributed in the valley area.

Experimental Procedure
The field campaign was conducted in May 2018. A representative hillslope in the southeast of the Chongling watershed near its outlet was selected as a typical hillslope to conduct dye experiments with the aim of evaluating the degree of preferential flow and infiltration redistribution at different hillslope positions.  Table 1) and dye experiments were performed using BBFCF and KI-starch in each plot. The experimental plots were located at least 1.5 m from trees to reduce the interference from roots. There had been little rainfall in the two weeks before and during the experiments. Before the dye experiments, the hillslope surface was flattened and vegetation (e.g., grass, shrubs and fallen leaves) on the surface were carefully removed. On each prepared plot, a rectangle frame made of polymethyl methacrylate (80 cm length × 60 cm width × 30 cm height) was embedded 10 cm into the soil, as illustrated in Figure 2. The soil on both sides of the frame was tamped by a wooden hammer in case of preferential infiltration along the inner wall of the frame. BBFCF and KI-starch (colorized by spraying a starch solution) were selected as the staining tracers due to their high visibility in soils, ease of dissolution, and relatively low toxicity to the environment [15,45]. BBFCF is only accessible to infiltration paths with a larger diameter because the particles have a relatively large diameter and are easily blocked or absorbed by the soil matrix. However, negative monovalent iodide ions of the KI dye are easily soluble in water and are hardly absorbed by the soil matrix nor cause an obstruction. Hence, in this study, the soil sections stained by BBFCF demonstrated a pattern of preferential flow that was channeled via large pores, whereas those stained by KI-starch demonstrated the actual pattern of water infiltration. Thirty liters of mixed solution containing BBFCF (4 g/L) [15] and KI (10 g/L) [46] was sprayed evenly through a rainfall sprinkler for each plot. According to pre-tests, the simulated rainfall/dye intensity was set at 18 mm/h to avoid ponded or running water on the ground. After the simulated rainfall, the quadrat area was covered with waterproofing and insulation materials to prevent moisture from evaporating. After 24 h, the cover materials and polymethyl methacrylate frame were carefully removed.
Immediately after the remove of the cover and frame, seven vertical soil sections were excavated at each plot to a depth of 80 cm in 10 cm intervals starting from the lower edge along the hillslope ( Figure 2); hence, a total of 35 vertical soil sections were obtained ( Figure 3). After the excavation of each section, a rectangle scale frame made of aluminum alloy and a cubic white balance calibrator were placed on the section for the correction of geometric distortion and color analysis of images. Firstly, the soil section stained by BBFCF was photographed by a Canon digital camera under a shading plate ( Figure 2). Then, the color developing agent containing ferric nitrate (25 g/L) and water-soluble starch (50 g/L) was evenly sprayed on the vertical soil section [46]. After the complete development of colors (~20 min), the soil section stained by KI-starch was also photographed from the same position as for the BBFCF staining. Earthworms and other soil-dependent animals were not observed during the excavation. Bulk density (BD), saturated hydraulic conductivity (Ks), field capacity (FC), saturated soil content (θsat) and soil particle size were measured with the soil samples collected by cutting ring with a volume of 100 m 3 in different layers of soil. Meanwhile, about 10 kg mixed soil samples (including lager gravels) from different depths were collected to measure the gravel mass ratio (GMR), which is equal to the ratio of the mass of the dried gravels (>2 mm) to that of the total dried soil sample, to reduce the underestimation due to the limited volume of cutting ring ( Table 2). The stable infiltrability was also estimated by means of curving fitting to Philip's equation using double ring infiltrometers near the experimental plot during the act of excavation [47]. Note. Data in the table represent means and standard deviations (SD). Lowercase letters ( a,b,c,d ) indicate significant differences (p < 0.05) among hillslope positions whereby the letter indicates the relevant plot. ST: soil thickness; BD: bulk density; GMR: gravel (>2 mm) mass ratio, which is defined as the ratio of the mass of the dried gravels (>2 mm) to that of the total dried soil sample; Ks: soil saturated hydraulic conductivity; FC: field capacity; θint: initial water content; θsat: saturated water content.

Image Analysis Techniques
Geometric distortion occurred when the stained sections were photographed with a digital camera. Therefore, the geometrical correction of the images was performed firstly by establishing a relationship between the location of the pixels in the image and the corresponding coordinates on the soil section [15], and the size of the image was set as 8 cm × 8 cm to embody the vertical soil section with a practical size of 80 cm × 80 cm. The image resolution was set as 100 pixels/cm and the total pixels per image was 800 × 800. A series of image analysis procedure for an image of a vertical soil section ( Figure 4) was conducted to extract available data for calculation and evaluation. The image recordings of the soil sections that were stained by either BBFCF or KI-starch were treated in the same way, and the processing methods of the stained section images obtained under the two different dying treatments were also the same. The abundance of roots were measured with the method of image processing. The root distributing in the soil vertical section were recognized, color-replaced from the digital image and the pixels of roots were then extracted from the binary image. The abundance of roots (Φroot) were calculated as the proportion of root pixels in the whole section. The maximum depth of the roots (Droot) were also measured ( Table 1).

Figure 4.
Flow chart of the image analysis procedure for an image of a vertical soil section, taking the 2 nd section in Plot C for example. S1 (step 1) was to reduce the noise of the original image from the field section. S2 (step 2) was to maximize the color saturation so as to subsequently recognize target color as the stained area. S3 (step 3) was to binarize the image by replacing the staining color with white and replacing the unstained area with dark, then to extract data to calculate the statistics of the stained pixels in the binary bitmap as well as the preferential flow indices and evaluate the infiltration redistribution. S4 (step 4) was to merge the binary image of Brilliant Blue FCF (BBFCF) with that of KI-starch, replacing the stained area of BBFCF with blue and that of KI-starch with purple for more directly visual comparison.

Preferential Flow Indices
Dye coverage, DC (%), is defined by Flury et al. [45] as the percent coverage of the stained area. The DC at some depth, h, and the total DC (TotDC) of the soil section were calculated using Equations (1) and (2) respectively: where Dh is the number of dye-stained pixels at the depth of h, and UDh is the unstained pixels on the depth of h. The maximum infiltration depth, MaxID (cm), is the maximum penetration depth of dying tracer solution on the soil section. The uniform infiltration depth, UniID (cm), was defined by van Schaik [48] as the depth where the stained area strongly decreases, thus presenting the depth of the uniform infiltration front above which the infiltration process is dominated by uniform flow [26]. The depth of 80% of the DC is usually selected as the critical depth to divide the uniform flow area and the preferential flow area [21,26,48].
The preferential flow fraction, PF-fr, is the fraction of the total infiltration that flows through preferential flow paths [48]. The calculation method for the PF-fr was widely used in previous studies that only relied on the UniID as the boundary between the matrix flow area and the priority flow area [20,21,26,48,49]. However, the traditional calculation can lead to the inability to truly and completely distinguish the priority flow from the matrix flow. The first limitation of previous approaches is that they all used BBFCF as a single tracer for the analysis and calculation. As mentioned, during the process of infiltration with water flow, BBFCF particles are easily adsorbed and blocked by the soil matrix, and the redistribution process results in obvious hysteresis when compared with the actual water flow. The use of KI-starch could greatly reduce hysteresis and more accurately simulate the actual distribution of water flow. Another limitation is the imprecise division of the matrix flow region and priority flow region into upper and lower parts. A considerable portion of water flow is transversely diffused through the matrix during rapid infiltration via macropore channels. Although BBFCF dye can mark the macropore channels, it can hardly trace the part of water flowing horizontally through the macropore channel. The use of KI-starch tracer can, however, compensate for the above limitations. Hence, the present study modified the method of predecessors by simultaneous tracing with BBFCF and KI-starch. The matrix flow region is defined here as including the piston matrix flow region above the UniID (DC = 80%) and the transverse diffusion region, which represents the portion stained by the KI-starch while unstained by the BBFCF dye below the UniID in the KI-stained image.
When MaxIDBBFCF ≥ UniIDKI-starch, the modified calculation of the PF-fr is calculated using Equation (3); when MaxIDBBFCF < UniIDKI-starch, the modified calculation of the PF-fr is calculated using Equation (4).
where TotStARKI-starch is the total area stained by the KI-starch (cm 2 ).

Non-Uniformity Analysis of Infiltration Depth
In the tracer experiment, the vertical soil sections were not completely stained above the maximum staining depth due to flow inhomogeneity. Therefore, the actual staining depth can replace the maximum staining depth as recommended by others [17]. However, in the vertical section, the difference between the above two depths on the vertical section is mainly due to lateral infiltration, which is one of the important manifestations of preferential flow [32]. The replacement may fail to truly and completely describe the non-uniform characteristics of the preferential flow. Therefore, the non-uniformity coefficient (Cµ) of the preferential infiltration depth defined in this study is defined by Equation (5): where ℎ max ( , ) is the maximum infiltration depth at ( , ), ℎ is the average of the maximum infiltration depth on a section, and N is total number of pixels parallel to (or orthogonal to) the section.

Fractal Dimension of the Preferential Flow Wetting Front
Fractal theory was first proposed by Mandelbrot [50], and was subsequently introduced by Tylor and Wheatcraft [51] for the study of soil science. The preferential flow wetting front trace recorded by the digital camera in the dye tracer experiment had distinct fractal features with obvious irregular, split, and disordered complex structures. In this study, the box-counting fractal dimension is used to describe the non-uniformity of the wetting front distribution in the two staining patterns. The method of calculating the box dimension involves first covering the wetting front by square pixel grids with a side length of δ. The number of boxes covered by the wetting front is N(δ), and then the fractal dimension, Df, is defined by Equation (6) Df = 1 indicates that the wetting front is a straight line, thereby meaning that there is no preferential flow. Df > 1 indicates an irregular curve with the existence of preferential flow, and the value of Df is positively correlated with the non-uniformity of preferential flow. In this study, MATLAB was used to automatically calculate the Df of the wetting front, which is verified by the Koch curve and the Sierpinski triangle with high accuracy (relative error < 5%).

Correlation and Error Analysis
Pearson's correlation coefficient analysis was applied to investigate the correlations between the indices of preferential flow and spatial redistribution and environmental variables by using SPSS 18.0 software. The differences of the results were compared using one-way analysis of variance with SPSS 18.0 software.

Degree of Preferential Flow
The two different dye tracers were applied to mark the pattern of soil water infiltration and redistribution on the experimental rocky hillslope. Thirty stained images of vertical soil sections on the five hillslope positions were investigated. Five representative stained section examples and their corresponding change in dye coverage with depth are presented in Figure 5.
The dye was distributed consistently throughout the area in a longitudinal direction at the bottom-and the lower-slope positions (Figure 5a,b). However, the stained area was distributed as the fingerlike, fragmented, and heterogeneous pattern on the mid-, upper-, and top-slope positions as evidenced by the existence of fingerlike extensions or the unstained fraction surrounded by a stained area (Figure 5c-e). In general, the purple KI-starch stained areas in Figure 5 are larger than the blue BBFCF stained areas, which implies that the infiltration flow distribution could not be recorded with a single dye tracer to estimate the degree of priority flow. The purple area is the extension and spread of the blue area, which indicates that soil water mainly infiltrated through some macropore channels and gradually diffused to the surrounding soil matrix. As mentioned, the negative monovalent iodide ions of the KI-starch dye are easily soluble in water and are hardly obstructed or absorbed by soil particles, which could be approximated as a representation of the actual flow [53]. Our results were consistent with those of Germán-Heins [54] and Jiang [55]. In addition, with increased elevation of the hillslope position, the blue area in Figure 5c-e occupies an increasing (previously purple in Figure 5a,b) area, which indicates that the contribution of priority flow to the actual water infiltration gradually increased with the elevation, and that preferential flow was the main contributor to soil water infiltration on the higher hillslope positions. The variables of preferential flow treated with the two dye tracers on the different hillslope position are shown in Figure 6. The degree of preferential flow that was determined by the two dyes in addition to the soil steady infiltrability both increased with hillslope position (except for UniID with the BBFCF dye). All of the TotDC treated with KI-starch were significantly higher than those treated with the BBFCF dye, irrespective of the hillslope position. However, differences in the MaxID between the dye treatments among each of the higher hillslope plots were not statistically significant, while those between the bottom-and lower-slope plots were significant (p < 0.01). This indicates that a large fraction of the infiltrating dye tracer was channeled through preferential pathways, which is supported by the generally increasing trend of the modified PF-fr with increased hillslope position (bottom: 0.10; lower-slope: 0.11; mid-slope: 0.15; upper-slope: 0.29; top: 0.26). The mean value of the UniID, above which is generally defined as the region dominated by matrix piston flow, was 7.79, 8.36, 5.18, 2.68, and 3.19 cm for the BBFCF treatment and 18.02, 18.14, 16.53, 26.63, and 28.73 cm for the KI-starch treatment, from the bottom to the top of the hillslope, respectively. The mean values of stable infiltrability from the mid-slope to the top of the hillslope were 4.06, 5.01, and 6.13 mm min −1 and were significantly higher than that of the bottom (2.51 mm min −1 ) and lower-slope (2.65 mm min −1 ). The variation range of each variable that is related to preferential flow also exhibited a gradually increasing trend with increased elevation, whereas the range of steady infiltrability was opposite. This indicates that infiltration channels, especially preferential channels, were more irregular in the upper-and top-slope positions and more uniform in the bottom-, lower-, and mid-slope plots.
The degree of preferential flow directly affects the infiltration capacity of water, and presents the more complex and chaotic redistribution pattern of water in soil [14,17,19]. Topographic factors can affect the infiltration process of soil water by influencing the physical and chemical properties of soil, soil thickness, and gravel content amongst others [23,24,26,34]. Even small topographic changes can cause drastic changes in soil conditions over a short distance. Hillslope position is an important terrain factor. With the increased elevation of the terrain, the soil layer became shallower, and the content of gravel increased with slower soil development and enhanced weathering ( Table 2). Gravels with a diameter >30 mm were found to account for a larger proportion of the soil in the upper-and top-slope plots. The aforementioned increase in the degree of preferential flow with the rise of hillslope position agrees with the findings of Mei et al. [21], who suggested that the role of rock fragments and plant roots play an emphasized role in leading macropore flow. Macropores were determined by Alaoui et al. [6] to be more efficient at transporting water in comparison to piston flow.

Spatial Redistribution of Infiltration
In order to reflect the spatial characteristics of soil water redistribution after infiltration, the ground surface of each experimental plot was discretized into a square grid of 0.1 cm × 0.1 cm, thereby obtaining 4.8 × 10 6 square grids for each experiment. MATLAB was applied to draw a threedimensional space contour of the maximum infiltration depth of each experimental plot (Figure 7). Great differences were observed in the soil water redistribution after infiltration between experimental plots/hillslope positions. The Cμ of the infiltration depths were analyzed for each three-dimensional spatial distribution map along the parallel and orthogonal directions of the stained experimental soil cuboid (Figure 8). The solid line with a triangle in Figure 8 represents the variation of Cμ with the distance to the origin in the orthogonal direction (InOrD hereafter) to reflect the degree of heterogeneity in the infiltration depth between the dye-stained sections. The dashed lines with solid circles represent the variation of Cμ in the parallel direction (InPaD hereafter) to reflect the degree of internal heterogeneity in the dye-stained sections. The Cμ-distance curve in Figure 8a for the bottom of the hillslope was similar for the two dye treatments, whereby the orthogonal variation ranged from 0.02 to 0.18 and the parallel variation ranged from 0.11 to 0.28. In comparison to the Cμ InOrD, the Cμ InPaD yielded a higher value and less fluctuation. In Figure 8b for the lower-slope, the variation of the Cμ in the different directions and with the two dyes was consistent with that of Figure 8a, where the overall variation range was between 0.05 and 0.34. The increase of the Cμ fluctuation InOrD was observed for the BBFCF treatment. The Cμ-distance curve of the two dying treatments in Figure 8c show an obvious separation. The curve of the BBFCF dye fluctuates more (ranging from 0.14 to 0.52) than the curve of KI-starch, which maintained a relatively stable trend in both directions, and was closer to that of Figure 8a,b. In Figure 8d, the value and range of the Cμ for the BBFCF dye increases further, and both increase greatly InPaD and InOrD with increased distance to the origin. Meanwhile, the fluctuation of the curve for the KI-starch treatment also increased significantly, but was still more stable than that of the BBFCF dye. The maximum values of Cμ for the BBFCF and KI-starch treatments were found in Figure 8e to be 0.86 and 0.62, respectively. The fluctuation of the Cμ-distance curve also peaked while that for the different dying treatments remain separated.
To accurately compare the Cμ of the water flow infiltration depth at different hillslope positions and in different directions under different dying treatments, the statistical characteristic values of the Cμ curve in Figure 8 are summarized in Table 3.  With increased hillslope position, the statistical characteristic values of the Cμ (i.e., the average, Ave. Cμ, the range, Ran. Cμ, and the standard deviation, Std. Cμ) mainly showed a gradual increasing trend. Although some minimum values were evident for plot B (lower-slope), they were close to the corresponding values determined for plot A (bottom of the hillslope), for example, the Ave. Cμ of the BBFCF dye InOrD and the Ave. Cμ for the KI-starch treatment InPaD. This indicates that the nonuniformity of the spatial redistribution of the soil water infiltration gradually increased with the rise in the hillslope position. The stable curve of the Cμ for the bottom-(plot A) and lower-slopes (plot B) might have been due to a more homogeneous seepage medium. In this case, soil water would have been stored in the soil pores under molecular and capillary forces, thus making any transfer to the surrounding environment difficult [30]. The spatial redistribution of infiltration in the mid-and upper-slopes revealed more potential preferential channels that were caused by the heterogeneous soil structure. These trends were supported by the findings of Liu et al. [56], van Schaik et al. [48], and Öhrström et al. [57].
Comparison of the Cμ under different treatments in the same direction showed that all of the statistical characteristic values of the Cμ for the BBFCF dye were significantly higher than those of the KI-starch dye, which indicates that the spatial heterogeneity of the infiltration pattern recorded by the BBFCF dye was significantly larger than that of the KI-starch. Sheng et al. also found the infiltration of an ionic solution tracer revealed the spatial redistribution more accurately than that of a macromolecule solution [47], while the latter was better for detecting macropores channels. Moreover, the Cμ in different directions under the same treatment were compared, and showed that the Ave. Cμ InOrD was larger than that of the InPaD, while the Ran. Cμ and Std. Cμ presented the opposite trend, which might have represented the non-uniformity lateral flow pathway [32].
The fractal dimension (Df) values of the wet front curve of each plot for the two different dying treatments are shown in Table 4. The Df increased with increasing elevation regardless of the dye type. The larger the Df value, the higher the irregularity of the wet front curve, which further supported the influence of the hillslope position on the redistribution of infiltration [17]. In addition, the Df values for the BBFCF dye were generally larger than those for the KI-starch dye, which emphasized the role of water diffusion in the infiltration process with priority flow [1]. It is necessary to take the convection and diffusion of soil water movement into consideration in the study of matrix and priority flow.

Factors Impacting Preferential Flow and Infiltration Redistribution
The main environmental factors affecting the preferential flow and infiltration redistribution were assessed using the correlation coefficients between the environmental variables and indices of preferential flow and infiltration redistribution (Table 5).
In the BBFCF treatment, altitude, slope gradient, gravel mass ratio (GMR), soil saturated hydraulic conductivity (Ks), initial water content (θint), saturated water content (θsat), sand and the abundance of the roots in stained area (Φroot-stained) generally had a negative correlation with UniID and a positive correlation with other preferential and redistribution indices. Soil thickness (ST), bulk density (BD), clay, silt, the total abundance of roots (Φroot) and the maximum depth of roots (Droot) exhibited a positive correlation with UniID and a negative correlation with other indices. The relationships between altitude, BD, GMR, θsat, soil particle size and Φroot-stained and MaxID, UniID, Cμ, and Df showed varying levels of significantly correlation. Ks had a significant impact on the Cμ between sections. In the KI-starch treatment, ST, BD exhibited a negative correlation with preferential and redistribution indices, whereas some other environmental variables (altitude, gradient, GMR, Ks, θsat, Φroot-stained) exhibited positive correlations with those indices. The relationships between altitude, BD, GMR, Ks, soil particle size and Φroot-stained and TotDC, MaxID, and Df were significantly correlated. Furthermore, the correlation between BD, θsat, soil particle size and PF-fr was also significant ( Table  5). In addition, the results showed that there is little correlation between θint, Φroot and Droot and the preferential flow and infiltration redistribution.
Preferential flow and infiltration redistribution were statistically affected by the GMR and BD. Gravels with a diameter of >30 mm were found in the soil layer on the middle-and upper-slope plots. Katra et al. [30] found that gravel with a diameter of >4.76 mm had a significantly greater impact on the soil infiltration rate than gravel with a diameter of 2.38-4.76 mm. In the uncompacted soil layer, an increased gravel content has been found to significantly increase the number of soil macropores, while the decrease of the overwater section and the increase in the water flow curvature have been shown to have lesser effects on the soil water inhibition [22,32]. In addition, the increased gravel content may also provide potential channels for the lateral soil midstream [32], thus manifesting the increase of the non-uniformity InOrD.
The GMR has also been found to indirectly affect the degree of priority flow and infiltration redistribution by changing the BD, Ks, and θsat amongst others [30]. The BD in the present study exhibited strong negative correlations with the degree of preferential flow and the spatial infiltration redistribution, which were consistent with the findings of Mei et al. [21], who demonstrated that the impact of BD was related to the total porosity. In addition, we found that Ks and θsat were strongly related to the infiltration process, especially with respect to spatial redistribution, which was supported by the results of Laine-Kaulio [58].
The correlation between the roots in stained area (Φroot-stained) and the indices of preferential flow and infiltration redistribution demonstrated that the impact of roots can not afford to be neglected. The plant roots were generally proved to change the soil properties by reducing bulk density and increasing soil properties, which is critically responsible for the depth and heterogeneity of infiltration [6,21,26]. However, there is little correlation between the total abundance of roots (Φroot) and the infiltration indices. The large differences between Φroot and Φroot-stained in Plots A and B were found, contrasting with those in Plots C, D, and E. As mentioned in Table 1, the main difference in surface vegetation types in Plots A and B from those in other higher plots is the Robinia pseudoacacia Linn. The root systems of Robinia pseudoacacia Linn were more deeply spread than those of shrubs [21], and there were coarse roots found in the soil vertical section (like the example of Plot A in Figure  3). However, stained spots were not found around the most of the roots in Plot A, especially around the coarse roots. In contrast, the vertical soli section exhibited more fine roots, particularly in the upper zone, and most of them were surrounded by dye tracer. These results agreed with those obtained by Mei [21], who indicated that root systems in the upper soil have more strong influence on preferential flow. However, Yunusa [59] indicated that roots pathways are deeper and wider under woody perennials compared to annual crops and pasture, contributing to enhance the preferential flow, which conflicted with our results. Limited amounts of simulated precipitation and denser soil texture might led to this confliction.   Note. ST: soil thickness; BD: bulk density; GMR: gravel (>2 mm) mass ratio, defined as the ratio of the mass of the dried gravels (>2 mm) to that of the total dried soil sample; Ks: soil saturated hydraulic conductivity; FC: field capacity; θint: initial water content; θsat: saturated water content; Φroot: the abundance of roots; Φroot-stained: the abundance of the roots distributed in the stained area; Droot: the maximum depth of roots; TotDC: total dye coverage; MaxID: maximum infiltration depth; UniID: uniform infiltration depth; Cμ: non-uniformity coefficient of of the preferential infiltration depth; Ave.Cμ in PaD: average value of the Cμ in parallel section direction; Ave.Cμ in OrD: average value of the Cμ in orthogonal section direction; Std.Cμ in PaD: standard deviation of the Cμ in parallel section direction; Std.Cμ in OrD: standard deviation of the Cμ in orthogonal section direction; Df: fractal dimension PF-fr: preferential flow fraction.* Indicates a significant correlation (p < 0.05). ** Represents significant correlation (p < 0.01).

Conclusions
To evaluate spatial variation of the degree of preferential flow and infiltration redistribution on different hillslope positions, dye tracer experiments were performed at the bottom-, lower-, mid-, upper-, and top-slope positions on a representative hillslope in the rocky Taihang Mountains. The dye tracer solution was mixed with BBFCF and KI-starch to distinguish the flow patterns of macropore flow and actual water flow. The three main conclusions of this study are as follows: (1). With increasing hillslope position, the total dye coverage, maximum infiltration depth, and saturation steady infiltration rate all showed obvious increasing trends; the contribution of preferential flow to the actual water infiltration gradually increased with plot elevation, whereby the value of the mean preferential flow fraction (PF-fr) was 0. Additionally, the degree of non-uniformity of the spatial redistribution for macropore flow revealed by the BBFCF dying pattern was larger than that of the actual water flow revealed by the KI-starch dying pattern. (3). The GMR, θsat, altitude, Ks, Φroot-stained were found to be the main positive factors impacting the degree of preferential flow and spatial non-uniformity of the infiltration redistribution. A lower BD and lower clay content could effectively promote preferential flow and the spatial distribution of heterogeneity.
The preferential flow and infiltration redistribution exhibited significant spatial variability along a hillslope in the rocky mountain area. It is necessary to take the contribution of preferential flow and heterogeneous infiltration into consideration, particularly in the application and modification of numerical modeling which is worthy further studying. Our findings therefore provide information for other researchers investigating the infiltration process of preferential flow in rocky mountain areas and reference for parameter application and selection, contributing to the improvement of simulation accuracy of montane hydrological models.