Use of Pedotransfer Functions in the Rosetta Model to Determine Saturated Hydraulic Conductivity (Ks) of Arable Soils: A Case Study

Abstract

A key parameter for the design of soil drainage and irrigation facilities and for the modelling of surface runoff and erosion phenomena in land-formed areas is the saturated hydraulic conductivity (Ks). There are many methods for determining its value. In situ and laboratory measurements are commonly regarded as the most accurate and direct methods; however, they are costly and time-consuming. Alternatives can be found in the increasingly popular models of pedotransfer functions (PTFs), which can be used for rapid determination of soil hydrophysical parameters. This study presents an analysis of the Ks values obtained from in situ measurements conducted using a double-ring infiltrometer (DRI). The measurements were conducted using a laboratory permeability meter (LPM) and were estimated using five PTFs in the Rosetta program, based on easily accessible input data, i.e., the soil type, content of various grain sizes in %, density, and water content at 2.5 and 4.2 pF, respectively. The degrees of matching between the results from the PTF models and the values obtained from the in situ and laboratory measurements were investigated based on the root-mean-square deviation (RMSD), Nash–Sutcliffe efficiency (NSE), and determination coefficient (R2). The statistical relationships between the tested variables tested were confirmed using Spearman’s rank correlation coefficient (rho). Data analysis showed that in situ measurements of Ks were only significantly correlated with the laboratory tests conducted on intact samples; the values obtained in situ were much higher. The high sensitivity of Ks to biotic and abiotic factors, especially in the upper soil horizons, did not allow for a satisfactory match between the values from the in situ measurements and those obtained from the PTFs. In contrast, the laboratory measurements, showed a significant correlation with the Ks values, as estimated by the models PTF-2 to PTF-5; the best match was found for PTF-2.

1. Introduction

Water permeability is a key property of soils, especially with regard to the design of soil irrigation and drainage facilities, modelling of surface runoff and erosion phenomena in land-formed areas, and environmental processes occurring in porous media [1,2,3]. The study of soil physics, and in particular the determination of the filtration coefficient (Ks) based on direct methods, is a very interesting issue; nevertheless, it is both labour-intensive and costly [4,5,6,7,8,9].

Climate change in a region or environment entail changes in the method of soil cultivation and plant production. Changes in the interactions between agriculture (i.e., soil compaction) and natural environment (i.e., weather, soil conditions) are a key feature of the transitions which scientists are trying to explain. Therefore, learning the landscape is an important tool in the proper management of water resources in rural areas. Agricultural development depends on the temperature and rainfall distribution during the growing season and on the fertilities and properties of soils, including their abilities to conduct and retain water. The climatic and soil conditions are often unfavourable for the rational and efficient management of environmental resources. Moreover, global climate change, as characterised by an increased frequency of periods of excesses and shortages of water, may hinder the development of modern agriculture, i.e., by reducing the quantity and quality of crop yields. To counteract these adverse phenomena, water engineering and land reclamation facilities are being used in many regions of the world [10,11]. However, to correctly design the spacing for drains, drainage elements, and irrigation ditches and/or the technical parameters of irrigation devices influencing the intensity of the water supply, it is crucial to identify the water permeability of the corresponding soils, including their infiltration and filtration processes [12,13].

The saturated hydraulic conductivity (i.e., Ks) is a key input parameter for modelling the water flow in a soil [14]. Its value can be determined using in situ and laboratory methods, but this often turns out to be a difficult task, owing to the large spatial and temporal variabilities of soil properties [15], and the need to take special care of samples (e.g., to provide an intact structure [16]. Ks is a basic physical property of the soil that affects all soil-plant-water relationships and processes. However, it is one of the most variable soil properties, as it is related to the soil texture and structure, and is influenced by factors such as the land topography, vegetation, land use, and climate [17]. It is widely accepted that in situ measurements provide the most accurate results for the determination of a soil’s physical and hydraulic properties.

There are different ways to determine the movements of water in unsaturated (infiltration) and saturated (filtration) zones that are of interest to the scientific community and engineers in various industries. The high cost and effort of in situ measurements and the growing demand for such data have led scientists to seek alternative indirect methods [18,19]. In recent years, pedotransfer functions (PTFs) have become popular and have been used to estimate time-consuming and difficult-to-measure soil properties, such as Ks or soil water retention (pF) values [1,20,21,22]. Ahuja et al. [23], Rawls et al. [24], Timlin et al. [25], and Suleiman and Ritchie [26] proposed formulas for estimating Ks values based on the effective porosity. Cosby et al. [27], Wösten et al. [28], Saxton and Rawls [29], and Weynants et al. [30] proposed formulas for determining Ks values based on using the textural composition, density, and other physical properties of soils. Currently, studies demonstrating the possibility of using neural networks to estimate Ks values [5,31] are very popular, along with studies aiming to improve existing formulas with appropriate amendments [32]. Puckett et al. [33] proposed a model to predict Ks based on only clay-sized particles. The authors showed that fine sand, sand, and clay percentages were highly correlated with Ks. Jabro [4] developed a model that used the site-average dry bulk density and grain size as predictive variables of Ks. A simple Campbell’s model, with soil texture data, can be used [34] or Smettem and Bristow’s model from soil clay content using a variety of agricultural topsoil samples [35,36] published a model which were used the relationships between soil texture and soil moisture content at saturation and soil texture and Ks.

Pedotransfer functions, as a five hierarchical models, offered by Rosetta are one of the most applied PTFs, as demonstrated by a number of citations in other research papers more than 1000 times (in Google scholar website). One of the reasons for the Rosetta PTFs popularity is its easy availability and an extensive soil database from North America and Europe containing 2134 soil samples with water retention data and measurements for the saturated hydraulic conductivity (Ks) for 1306 of soils samples and continues to be developed [37]. The majority of work on the performance of PTFs in Rosetta with the aim of determining Ks factor has focused on surface layers of the soil (0–30 cm or 0–50 cm) [13,14,38].

In this work concentrated on the full depth soil pits (150 cm) and including an estimate of Ks for surface layers (0–50 cm) and additionally for a deeper levels of soil (50–150 cm). It was considered that for engineering purposes, such as irrigation or drainage planning, it is advisable to know the Ks parameter.

The aim of this research was to test the possibility of using PTFs in the Rosetta program to determine the Ks values of arable soils based on easily accessible input data, i.e., the soil type, granulometric composition, density, and water content. We hypothesize that the Ks values obtained from PTFs in the Rosetta program for arable soils correspond to the values obtained from (i) direct in situ and (ii) laboratory measurements. To test the hypothesis we examined the arable soils in Central Europe (district of Racibórz, south Poland) in a 150 cm depth soil pits.

2. Materials and Methods

2.1. Description of Study Area

The study area is located in the southwestern part of the Silesian Voivodeship (Poland) in the district of Racibórz (Figure 1). According to Kondracki’s [39] division of Poland into physio-geographical regions, the study area is located in the Central European Plain (31), in the macroregion of the Silesian Lowlands (318.5), and on the border of two mesoregions: the Głubczyce Plateau (318.58) and Racibórz Basin (318.59). The main watercourse flowing through the region is the Oder River, along with its left tributary, the Psina River.

2.2. Meteorological Conditions

In terms of climate, the study area is considered as one of the warmest areas in this region. In general, Poland has a mostly temperate climate, in transition between an oceanic climate dominating in the north and west of the country, and a continental climate in the south and east. In the multiannual period from 1971–2000, the mean annual air temperature was 8.5 °C, and total precipitation was 616 mm (according to the Institute of Meteorology and Water Management—station in Racibórz). Higher mean temperatures and a decrease in rainfall have been observed over the last decade. The characteristics of the meteorological conditions are listed in

Table 1. The characteristics of meteorological conditions of the study area.

Year or Multi-Year	Months of the Growing Season						Period
	Apr	May	Jun	Jul	Aug	Sep	Apr–Sep		Jan–Dec
	Sum of precipitation totals (mm)
2010	69	204	112	95	75	78	633		828
2011	27	69	85	131	62	17	391		508
2012	41	35	75	89	69	58	367		586
2013	21	132	110	14	48	99	424		597
2014	27	137	75	58	92	127	516		645
2015	23	55	28	32	13	17	168		280
2016	45	47	47	116	52	26	333		535
2017	68	36	41	69	44	107	365		562
2018	10	42	64	61	55	39	271		441
2019	31	75	33	30	66	98	333		536
	Mean values
2010–2019	36	83	67	70	68	67	381		552
1971–2000	45	67	79	94	74	56	415		616
	Mean monthly air temperatures (°C)
2010	9.0	12.5	17.2	20.4	18.7	12.7	15.1		7.9
2011	10.6	13.7	17.8	17.4	19.2	15.4	15.7		9.2
2012	9.9	15.3	17.7	19.9	19.1	14.7	16.1		9.2
2013	9.0	13.8	16.9	19.7	19.1	12.6	15.2		9.0
2014	10.8	13.8	16.3	20.4	17.4	15.6	15.7		10.5
2015	8.8	13.1	16.8	20.9	22.3	15.6	16.3		10.4
2016	9.0	14.7	18.4	19.6	18.2	16.4	16.1		9.8
2017	7.8	14.3	18.8	19.1	20.1	13.8	15.7		9.7
2018	14.0	17.0	18.5	20.2	21.5	16.0	17.9		10.6
2019	10.4	11.9	21.9	19.6	20.9	14.7	16.6		10.8
	Mean values
2010–2019	9.9	14.0	18.0	19.7	19.7	14.8	16.0		9.7
1971–2000	8.2	13.5	16.1	17.8	17.7	13.6	14.5		8.5
	Values of HTC (–)
2010	2.6	5.3	2.2	1.5	1.3	2.0		Classification of months:
2011	0.8	1.6	1.6	2.4	1.0	0.4			extremely dry
2012	1.4	0.7	1.4	1.4	1.2	1.3			very dry
2013	0.8	3.1	2.2	0.2	0.8	2.6			dry
2014	0.8	3.2	1.5	0.9	1.7	2.7			relatively dry
2015	0.9	1.4	0.6	0.5	0.2	0.4			optimal
2016	1.7	1.0	0.9	1.9	0.9	0.5			moderately humid
2017	2.9	0.8	0.7	1.2	0.7	2.6			humid
2018	0.2	0.8	1.2	1.0	0.8	0.8			very humid
2019	1.0	2.0	0.5	0.5	1.0	2.2			extremely humid
	Mean values
2010–2019	1.3	2.0	1.3	1.2	1.0	1.6

.

To assess the meteorological conditions of the region where the study was conducted, Selyaninov’s hydrothermal coefficient (HTC) values were calculated, as a measure of the precipitation efficiency in a given month. The calculation was performed as follows [40]:

$$\mathit{HTC}=\frac{10 · P }{\Sigma t} \left(\text{–}\right)$$

(1)

In the above, P is the monthly sum of precipitation (mm), and Σt is the sum of the mean daily air temperature values in a given month (in °C).

Based on the value of the HTC, the months of the growing season (April-September) were classified as follows: HTC ≤ 0.4—extremely dry; 0.4 < HTC ≤ 0.7—very dry; 0.7 < HTC ≤ 1.0—dry; 1.0 < HTC ≤ 1.3—relatively dry; 1.3 < HTC ≤ 1.6—optimal; 1.6 < HTC ≤ 2.0—moderately humid; 2.0 < HTC ≤ 2.5—humid; 2.5 < HTC ≤ 3.0—very humid; HTC > 3.0—extremely humid.

Based on the HTC values calculated for the months of the growing season from the 2010–2019, and as shown in Table 1, it was determined that 58.4% of the months were dry (8.3% extremely dry, 13.3% very dry, 28.5% dry, and 8.3% relatively dry); in 13.3% of the months, the conditions were optimal; and 28.3% of the months were humid (8.3% moderately humid, 6.7% humid, 8.3% very humid, and 5.0% extremely humid). In general, the variable meteorological conditions in the study area, characterised by periods with excesses and shortages of water, indicate the need to use drainage and irrigation facilities for agriculture. To determine their technical parameters, information on the soil permeability is needed. This confirms the advisability of conducting research on the possibility of using mathematical models to determine the soil Ks values, as an alternative to time-consuming in situ or laboratory tests.

2.3. Field Measurement and Soil Sampling

The research was conducted on arable lands during the agricultural season from 2012 to 2015. The field soil tests were conducted at 16 measurement-control points (up to a depth of 150 cm) in Wojnowice (2 points), Bojanów (2 points), Owsiszcze (6 points), Strzybnik (2 points), and Tworków (4 points). Three to five genetic levels were identified for each soil profile.

The soil infiltration was measured in situ, at depths of 10 cm (topsoil) and 35 cm (subsoil), using a double-ring infiltrometer (DRI) method (Figure 2). The DRI comprised an inner ring (9.5 cm diameter) and outer ring (19.5 cm diameter) inserted into the ground at a depth of 10 cm. The DRI was inserted by using a falling weight-type hammer striking on a wooden plank placed uniformly on top of the ring, and without undue disturbance to the soil surface. Each ring of the DRI was filled with a constant head of water level, and the outer ring helped when checking the lateral flow from the inside ring, so as to better estimate the infiltration, reducing losses. The Ks value was estimated when the water flow rate inside the inner ring reached a steady state [1], which in the case of the studied soils, lasted approximately 3–4 h. The infiltration rate was calculated for the respective time intervals ΔT as follows [41]:

$$i=\frac{864 · 4 · V }{3.14 ·{ D}_r^2 · \Delta T} (\mathrm{m}·{\mathrm{day}}^{-1})$$

(2)

where V is the volume of water (cm³) added to the inner ring at time ΔT (s), and D_r is the diameter of the inner ring (cm²).

For steady infiltration (as a constant value for the soil), the Land and Water Development Division [42] developed infiltration classes, as follows: very slow—<0.024 m·day⁻¹; slow—0.024 ÷ 0.12 m·day⁻¹; moderately slow—0.12 ÷ 0.48 m·day⁻¹; moderate—0.48 ÷ 1.56 m·day⁻¹; moderately rapid—1.56 ÷ 4.20 m·day⁻¹; rapid—4.20 ÷ 5.81 m·day⁻¹; and very rapid—>5.81 m·day⁻¹.

Undisturbed soil samples were taken from each genetic horizon, using cylinders with a volume of 100 cm³ (three replicates). In addition, approximately 1.0 kg of disturbed soil from each genetic horizon was used to determine the soil texture, and other laboratory analyses were performed.

2.4. Laboratory Analysis

The properties of the collected soil samples were determined as follows.

• The soil texture was determined using the Bouyoucose-Casagrande areometric method, based on a measurement of the density of the soil suspension during progressive sedimentation, and a sieving method to fractionate the sand. The contents of the particle size classes (sand, 2.0–0.05 mm; silt, 0.05–0.002 mm; and clay, <0.002 mm) were determined according to the Soil Taxonomy system from the United States Department of Agriculture [43].
• The soil bulk density (BD) was determined based on the gravimetric method, based on cylinders (100 cm³) for determining the mass of the dry soil per volume. The weight of each soil core was determined after drying in an oven at 105 °C for approximately 18–24 h. The dry bulk density for each core sample was then calculated as follows [44]:

  $$BD=\frac{M_s}{V_t} (\mathrm{g}·{\mathrm{cm}}^{-3})$$

  (3)

In the above, M_s is the mass of the dry soil weight (g), and V_t is the volume of the total soil sample (cm³).

• The soil water retention was investigated based on determining the soil suction using ceramic plates in a 5/15 bar pressure plate extractor. The pressure plate equipment used in this study was manufactured by the American Soil Moisture Equipment Corporation. In engineering practice, the soil suction is usually calculated in units of pF, as follows:

  $$\mathrm{pF}=\log \mathrm{h} \left(\text{–}\right)$$

  (4)

Here, h is soil suction pressure (in cm H₂O).

In the laboratory, the soil water potentials were measured at 300 cm H₂O ≈ 33 kPa ≈ 2.5 pF (representing the field capacity) and 15,000 cm H₂O ≈ 1500 kPa ≈ 4.2 pF (representing the permanent wilting point) [45]. The soil volumetric water content (θ_V) at pF = 2.5 and 4.2 was determined based on the gravimetric method. In particular, it was calculated as the ratio of the amount of water in the soil sample to the dry weight of the soil, after drying in an oven at 105 °C for approximately 18–24 h, as follows [44]:

$${\theta }_V=\frac{V_{H_{2}O}}{V_t} ({\mathrm{cm}}^3·{\mathrm{cm}}^{-3})$$

(5)

In the above, V_H2O is the volume of water in the soil sample (cm³), and V_t is the volume of the total soil sample (cm³).

• The Ks was measured under laboratory conditions using a laboratory permeameter (Figure 3) and using the Darcy’s law [46] with constant head method (Equation (6)) and falling head method (Equation (7) on undisturbed soil samples for three replications. The constant method can be used with virtually any soil, apart from poorly permeable soils such as clay, whereas the falling-head method is used to measure low-permeability soils, such as f.i. clay or peat samples [47,48]. The soil samples were placed in the laboratory permeameter and then saturated in water for 2–3 days. For the purpose of this study, a permeameter produced by the Eijkelkamp with a closed or open system and 25 holders was used. The Ks was calculated using the Darcy’s [46] equation, as follows:

  $$V=K ·i ·A ·t$$

  (6)

Here, explanations below;

Determining Ks using the constant head method was calculated from the following equation:

$$\mathrm{Ks}=\frac{V ·L}{A ·t ·h}$$

(7)

Here, V is the volume of water flowing through the sample (cm³); K is the permeability coefficient or ‘Ks’ (cm·day⁻¹ or m·day⁻¹); h is the water level difference inside and outside ringholder or sample cylinder (cm); L is the length of the soil sample (cm); i is the permeability rise gradient or h/L (–); A is the cross-sectional area of the sample (m²); and t is the time used for flow through of water volume (day). During measuring the following parameters have been determined: L and A—constants, depending on the type of sample ring used; V—volume measured in the burette (1 mL = 1 cm³); t—length of time lapse; h—calculated with the water levels measured with the water level meter.

Determining Ks using the falling head method was calculated from the following equation:

$${\mathrm{K}}_{\mathrm{s}}=\frac{a ·L}{A ·\left(t_2-t_1\right)} ·ln\frac{h_1}{h_2}+\frac{x ·a ·L}{A · \sqrt{\left(h_1 · h_2\right)}}$$

(8)

Here, Ks is the permeability coefficient (cm·day⁻¹ or m·day⁻¹); h is the water level difference inside and outside ringholder or sample cylinder (cm) h₁ and h₂ water level difference inside and outside the ringholder at respectively t₁(start) and t₂ (end); A is the surface of a cross-section of the sample (cm²); L is the length of the soil sample (cm); t is the time between beginning and end of the measuring t₂ − t₁ (day); a is the cross-section surface of a ringholder or sample cylinder (cm²) for a sample cylinder applies A = a; x is the evaporation factor (literature value): 0.0864 cm·day⁻¹ or 0.000864 m·day⁻¹.

The permeability of the soil is also determined by viscosity of the soil solution. Viscosity depends on the temperature. The laboratory water temperature varies from 18 to 22 °C, whereas the average groundwater temperature is 10 °C. Therefore, for certain applications the permeability will have to be corrected for the viscosity of the soil solution (usually water). The Ks for viscosity was corrected using the following equation [49]:

$$K_{10}=K_T · \frac{h_T}{h_{10}}$$

(9)

Here, K₁₀ is the corrected Ks at 10 °C (cm·day⁻¹ or m·day⁻¹); K_T is the Ks at the applied temperature (cm·day⁻¹ or m·day⁻¹); h₁₀ is the dynamic viscosity of water at 10 °C (Pa·s); h_T is the dynamic viscosity of water at T °C (Pa·s).

The soil classification was established according to the Polish Soil Classification [50], World Reference Base for Soil Resources (FAO and International Union of Soil Sciences (IUSS) Working Group World Reference Base (WRB), [51]), and USDA soil taxonomy (Soil Survey Staff, [43]).

2.5. Rosetta Description

The Ks values were calculated using the published neural network program Rosetta, with hierarchical PTFs based on five levels of input data [5]. The first level (PTF-1) used soil textural classes based on a lookup table providing parameter means for each USDA textural class. The second level (PTF-2) used the sand, silt, and clay percentages as inputs, and in contrast to PTF-1, provided hydraulic parameters that varied continuously with the texture. The third level (PTF-3) included the predictors used in the level PTF-2, along with the soil dry-bulk density. The fourth level (PTF-4) used PTF-3 and the soil volumetric water content (θ) at a water suction of 33 kPa (2.5 pF). The last level (PTF-5) comprised all of the other parameters, plus the θ value at a water suction of 1500 kPa (4.2 pF) [38,52]. The necessary input data for the Rosetta program (e.g., % soil texture group and bulk density) were determined based on the laboratory analyses of soil samples from the different soil horizons (

Table 2. Selected soil physical properties over all depths.

Object Name Profile Number	Depth	Horizon	Soil Texture Group *	Physical Characteristics of Soils
				% Fraction			Bulk Density ρb	θ2.5 pF	θ4.2 pF
	(cm)			Sand	Silt	Clay	(g∙cm−3)	(cm3∙cm−3)
Wojnowice No. 1	0–28	Ap	SiL	30	60	10	1.64	0.273	0.105
	28–62	Etg	SiL	20	69	11	1.62	0.305	0.085
	62–80	Btg1	SiL	15	68	17	1.66	0.337	0.141
	80–120	Btg2	SL	54	35	11	1.91	0.194	0.108
	120–150	Cg	SiL	21	63	16	1.80	0.297	0.137
Wojnowice No. 2	0–30	Ap	SiL	33	56	11	1.61	0.314	0.122
	30–62	Etg	SiL	18	70	12	1.58	0.330	0.102
	62–94	Btg1	SL	74	17	9	1.86	0.260	0.089
	94–116	Btg2	L	39	49	12	1.84	0.387	0.130
	116–150	Cg	SiL	21	63	16	1.72	0.337	0.132
Strzybnik No. 1	0–27	Ap	SiL	26	61	13	1.54	0.210	0.135
	27–100	Bw	SiL	11	70	19	1.69	0.313	0.150
	100–150	C	SiL	12	72	16	1.73	0.335	0.163
Strzybnik No. 2	0–40	Ap	SiL	25	61	14	1.53	0.302	0.132
	40–90	Bw	SiL	12	69	19	1.59	0.325	0.172
	90–150	C	SiL	12	72	16	1.62	0.336	0.124
Owsiszcze No. 1	0–20	Ap	SiL	25	65	10	1.67	0.334	0.129
	20–38	Etg	SiL	21	68	11	1.64	0.320	0.110
	38–65	Btg	SiL	16	64	20	1.68	0.364	0.177
	65–150	Cg	SiL	16	65	19	1.74	0.342	0.153
Owsiszcze No. 2	0–25	Ap	SiL	25	65	10	1.42	0.312	0.129
	25–43	Etg	SiL	21	68	11	1.58	0.343	0.110
	43–70	Btg	SiL	16	64	20	1.60	0.372	0.177
	70–150	Cg	SiL	16	65	19	1.60	0.370	0.153
Owsiszcze No. 3	0–23	Ap	SiL	26	63	11	1.50	0.328	0.029
	23–36	Etg	SiL	19	69	12	1.65	0.320	0.109
	36–68	Btg	SiL	17	63	20	1.59	0.372	0.181
	68–150	Cg	SiL	16	65	19	1.69	0.338	0.171
Owsiszcze No. 4	0–28	Ap	SiL	26	63	11	1.33	0.340	0.029
	28–41	Etg	SiL	19	69	12	1.56	0.352	0.109
	41–73	Btg	SiL	17	63	20	1.58	0.367	0.181
	73–150	Cg	SiL	16	65	19	1.63	0.344	0.171
Owsiszcze No. 5	0–25	Ap	SiL	27	61	12	1.57	0.352	0.103
	25–42	Etg	SiL	20	68	12	1.60	0.331	0.124
	42–65	Btg	SiL	17	63	20	1.60	0.363	0.168
	65–150	Cg	SiL	15	66	19	1.65	0.339	0.135
Owsiszcze No. 6	0–30	Ap	SiL	27	61	12	1.40	0.351	0.103
	30–47	Etg	SiL	20	68	12	1.58	0.358	0.124
	47–70	Btg	SiL	17	63	20	1.56	0.338	0.168
	70–150	Cg	SiL	15	66	19	1.64	0.344	0.135
Bojanów No. 1	0–25	Ap	SiL	34	53	13	1.80	0.303	0.113
	25–55	Etg	SL	54	30	16	1.96	0.246	0.115
	55–117	Btg	SCL	53	26	21	1.95	0.240	0.194
	117–150	Cg	SL	53	28	19	1.95	0.259	0.166
Bojanów No. 2	0–28	Ap	SiL	25	63	12	1.65	0.316	0.152
	28–63	Etg	SiL	16	68	16	1.69	0.299	0.175
	63–106	Btg	SCL	60	20	20	1.92	0.246	0.156
	106–150	Cg	SL	61	20	19	1.97	0.248	0.138
Tworków No. 1	0–30	Ap	SiCL	18	47	35	1.44	0.427	0.269
	30–46	AC	C	31	26	43	1.12	0.526	0.317
	46–63	OCg	C	16	19	65	1.06	0.609	0.331
	63–94	2O	SCL	57	29	14	1.44	0.676	0.309
	94–150	3G	SiL	19	62	19	1.40	0.488	0.192
Tworków No. 2	0–18	Ap	CL	25	35	40	1.32	0.474	0.247
	18–34	AC	C	18	21	61	1.29	0.516	0.289
	34–52	G	L	32	44	24	1.52	0.448	0.256
	52–85	Cg1	SL	55	33	12	1.75	0.354	0.133
	85–150	2G	SL	76	16	8	1.61	0.345	0.103
Tworków No. 3	0–25	Ap	L	43	35	22	1.65	0.377	0.250
	25–47	A/Bw	L	45	33	22	1.81	0.299	0.150
	47–75	Bw	SL	58	24	18	1.66	0.335	0.158
	75–117	Cg	L	40	43	17	1.68	0.348	0.131
	117–150	2Cg	L	41	45	14	1.68	0.342	0.089
Tworków No. 4	0–30	Ap	CL	22	40	38	1.11	0.566	0.294
	30–42	O	SiCL	17	55	28	1.51	0.439	0.287
	42–71	Bw	SiL	27	53	20	1.56	0.423	0.187
	71–100	Cg	SL	54	33	13	1.65	0.269	0.074
	100–150	2Cg	SL	67	25	8	1.62	0.309	0.065

* soil texture group according to USDA (1999): SL—sandy loam; SCL—sandy clay loam; L—loam; CL—clay loam; SiL—silt loam; SiCL—silty clay loam; C—clay.

,

Table 3. Basic descriptive statistics for soil physical properties for topsoils and subsoils (n = 32).

Index Value	Physical Characteristics of Soils
	% Fraction			Bulk Density ρb	θ2.5 pF	θ4.2 pF
	Sand	Silt	Clay	(g∙cm−3)	(cm3∙cm−3)
Minimum	11	21	10	1.11	0.210	0.029
Maximum	54	70	61	1.96	0.566	0.317
Mean	25	57	19	1.55	0.350	0.152
Median	25	62	12	1.58	0.329	0.127
SD	9.21	14.77	12.21	0.18	0.079	0.075
CV (%)	36.9	26.1	66.1	11.7	22.5	49.5
n	32	32	32	32	32	32

and

Table 4. Basic descriptive statistics for soil physical properties over all depths (n = 68).

Index Value	Physical Characteristics of Soils
	% Fraction			Bulk Density ρb	θ2.5 pF	θ4.2 pF
	Sand	Silt	Clay	(g∙cm−3)	(cm3∙cm−3)
Minimum	11	16	8	1.11	0.194	0.029
Maximum	76	72	65	1.97	0.676	0.331
Mean	30	52	18	1.60	0.351	0.155
Median	24	62	16	1.62	0.338	0.138
SD	16.8	17.6	10.5	0.23	0.086	0.065
CV (%)	57.0	33.9	57.2	14.6	24.4	41.8
n	68	68	68	68	68	68

).

In the Rosetta program, the relationship between θ and the water suction (h), i.e., the water retention θ(h), as well as Ks, are described using the well-known Mualem–van Genuchten equations [53] and are given as follows:

$$\theta \left(h\right)={\theta }_r+\frac{{\theta }_s-{\theta }_r }{[1+{\left(\alpha {\left|h\right|}^n\right]}^m}$$

(10)

$$K\left(h\right)={\mathrm{K}}_{\mathrm{s}}\frac{{\left\{1-{\left(\alpha ·h\right)}^{n-1}{\left[1+\left(\alpha ·h{)}^n\right)\right]}^m\right\}}^2}{{\left[1+{\left(\alpha ·h\right)}^n\right]}^{m/2}}$$

(11)

In the above, θ(h) is the soil volumetric water content (cm³∙cm⁻³) at suction h (cm H₂O); θs and θr are the saturated and residual water contents (cm³∙cm⁻³) at h = 0 and 15,000 cm H₂O, respectively; α (>0 in cm⁻¹) is related to the inverse of the air entry suction; n (>1) is a measure of the pore-size distribution, and m = 1 − 1/n; and Ks is the saturated hydraulic conductivity (m∙day⁻¹), as mentioned above.

2.6. Statistical Analysis and Model Performance Evaluation

The data set comprised the analytical results from the soil samples collected from the 16 soil pits. For the statistical analysis, the procedures provided by the program Statistica PL Version 12.5 were used, with a 5% significance level. The minimum and maximum values were determined for each physical parameter of the soil, and the arithmetic mean, median, standard deviation (SD), and coefficient of variation (CV) were computed. Moreover, a Spearman correlation test was conducted for the dataset. The correlation strength and direction of the relationship between two variables were determined based on Spearman’s rank correlation coefficient (rho), which takes values from −1 to 1. A positive sign of the coefficient indicated the existence of a positive correlation, and vice versa. The closer the values were to −1 and 1, the stronger the correlation. When rho = 0, there was no correlation between the examined variables [54]. This statistical method was chosen after finding that most data were not normally distributed (Shapiro-Wilk test). All of data are presented in the tables and graphs, i.e., to produce a visual image that is helpful in interpreting the results.

Based on the available soil data in this study, five widely used PTFs were selected to assess their respective performances in estimating Ks, via comparison with the measured Ks values (via the DRI) for the soil cores collected in the field (using a laboratory permeameter). The calculated values of Ks obtained using Rosetta were compared to the corresponding measured values and evaluated using two statistical parameters: root-mean-square deviation (RMSD) and Nash–Sutcliffe efficiency (NSE). The RMSD gives the mean difference between the measured and calculated values of Ks, and is calculated as follows [38,55]:

$$\mathit{RMSD}=\sqrt{\frac{1}{n} {\displaystyle \sum }_{i=1}^n{\left(m_i-p_i\right)}^2}$$

(12)

where: m_i is the measured value of Ks, and p_i is the corresponding (predicted) value of Ks as obtained using the Rosetta program.

The RMSD values are always non-negative and should be as low as possible; RMSD = 0 indicates a perfect fit of the model to the measurement data.

The NSE compares measured and predicted values (here, for Ks), and is given as follows [56]:

$$NSE=1-\left[\frac{{{\displaystyle \sum }}_{i=1}^n{\left(m_i-p_i\right)}^2}{{{\displaystyle \sum }}_{i=1}^n{\left(m_i-{\overline{m}}_i\right)}^2}\right]$$

(13)

In the above, n is the number of observations, m_i and p_i are the measured and predicted values of Ks, respectively, and ${\overline{m}}_i$ is the mean of the measured Ks values.

The NSE values range from −∞ to 1. A value of NSE = 1 corresponded to a perfect match of calculated values to the measured values of Ks. The ‘efficiency NSE’ = 0 indicated that the calculations of Ks using the Rosetta program were as accurate as the mean of the measured data; in contrast, an efficiency NSE < 0 occurred when the measured mean was a better predictor than the model or, in other words, when the residual variance, as described by the numerator in the equation above, was larger than the data variance described by the denominator.

3. Results and Discussion

According to the PTG [50], FAO and IUSS Working Group WRB [52], and USDA soil taxonomy [43], the examined soils were classified as follows.

• Order 3: Brown forest soils (brown earths, PTG—Polish: Gleby brunatnoziemne; WRB: Cambisols; USDA: Inceptisols—Udepts), Type 3.1. ‘Euthrophic brown soils’ (PTG—Polish: Gleby brunatne eutroficzne; WRB: Haplic Cambisol, Haplic, Stagnic, Endogleyic, or Vertic Cambisol (Eutric); USDA: Typic or Humic or Aquic or Oxyaquic or Vertic Eutrudepts)—occurred in Strzybnik.
• Order 5: Brown forest podzolic soils (Soil lessivé) (PTG—Polish: Gleby płowoziemne; WRB: Luvisols, Albeluvisols; USDA: Alfisols—Aqualfs, Udalfs), Type 5.2. ‘Streak brown forest podzolic soils’ (PTG—Polish: Gleby płowe zaciekowe; WRB: Haplic, Stagnic, Gleyic, Cambic Albeluvisol; USDA: Typic, Arenic, Aquic, Oxyaquic, Haplic, Glossaqiuc, or Vertic Glossudalfs)—occurred in Wojnowice, Bojanów and Owsiszcze.
• Order 7: Chernozemic soils (PTG—Polish: Gleby czarnoziemne; WRB: Chernozems, Phaeozems; USDA: Mollisols—Aquolls, Udolls), Type 7.4. ‘Chernoziemic fluvisols’ (PTG—Polish: Mady czarnoziemne; WRB: Mollic Fluvisol, Endofluvic Phaeozem; USDA: Fluvaquentic Endoaquolls)—occurred in Tworków.

As shown in Table 2 and Figure 4, the examined soils are heterogeneous in terms of their textures. Seven textural classes are observed in the examined soil profiles: silt loam (SiL; n = 42), sandy loam (SL; n = 10), loam (L; n = 6), clay (C; n = 3), sandy clay loam (SCL; n = 3), clay loam (CL; n = 2), and silty clay loam (SiCL; n = 2). The silt loam is the predominant soil texture. The SD values (Table 3 and Table 4) for most soil properties are large, indicating moderate to strong variability; this is not uncommon for soils [12,57,58]. Typically, the CV values are < 45%, indicating the mean variability of the measurement data. The highest CV values are observed for sand (57.0%) and clay (57.2%). The sand content in the soil profiles is between 11% and 76%, the silt content is between 16% and 72%, and the clay content is between 8% and 65% (Table 2 and Table 4). The high swelling-clay content could be a cause of the poor water flow through the samples, as well as a consequence of abnormal Ks results [59].

With regard to the surface soil horizons (n = 32), the mean soil bulk density (ρ_b) is 1.55 g∙cm⁻³, whereas that for overall depths is 1.60 g∙cm⁻³ (Table 3 and Table 4). The soil volumetric water contents at the matric potentials of 33 kPa (2.5 pF) and 1500 kPa (4.2 pF) in both cases have similar mean values (0.350/0.351 cm³∙cm⁻³ and 0.152/0.155 cm³∙cm⁻³, respectively).

Higher values of Ks are found in the upper horizons of the tested soils, owing to the significant macroporosity caused by plant roots, reconsolidation [60], and annual loosening from agrotechnical works and alternating soil freezing and thawing [61]. In addition, zoogenic channels (especially earthworm channels) and plant residues increase soil permeability [62]. For these reasons, the samples taken from arable and sub-arable horizons are less dense than those from lower horizons; the latter are also subject to self-compacting (Table 2, Table 3 and Table 4).

Table 5. Basic descriptive statistics for soil hydraulic conductivity (Ks) for topsoils and subsoils.

Index Value	Ks (m∙day−1)
	Field Measured	Laboratory Measured	PTF-1	PTF-2	PTF-3	PTF-4	PTF-5
Minimum	0.01	0.01	0.08	0.07	0.62	0.24	0.27
Maximum	4.67	1.18	0.38	0.40	1.82	1.90	1.91
Mean	0.85	0.21	0.17	0.25	1.41	0.45	0.45
Median	0.40	0.11	0.18	0.28	1.60	0.33	0.31
SD	1.09	0.29	0.05	0.10	0.37	0.34	0.34
CV (%)	128.8	138.7	29.2	40.5	26.1	75.7	75.2
n	32	32	32	32	32	32	32

presents descriptive statistics concerning the values of Ks for the topsoils and subsoils (n = 32), as obtained from in situ and laboratory measurements, and as estimated with the use of the PTF models in the Rosetta program (variant 1). In contrast,

Table 6. Basic descriptive statistics for soil hydraulic conductivity (Ks) over all depths.

Index Value	Ks (m∙day−1)
	Laboratory Measured	PTF-1	PTF-2	PTF-3	PTF-4	PTF-5
Minimum	0.01	0.08	0.07	0.55	0.24	0.27
Maximum	3.48	0.38	0.66	2.71	4.11	3.93
Mean	0.21	0.20	0.23	1.37	0.64	0.64
Median	0.04	0.18	0.19	1.32	0.32	0.31
SD	0.49	0.08	0.11	0.40	0.76	0.72
CV (%)	230.7	41.8	48.0	28.7	119.2	111.8
n	68	68	68	68	68	68

presents statistical measures considering the Ks values of all of the soil horizons (n = 68), but excludes the in situ tests not conducted on larger and deeper horizons, for practical reasons (variant 2).

In variant 1, the characteristic values of Ks (maximum, mean, and median values) as obtained from the in situ measurements are approximately four times higher than the Ks values obtained from the laboratory measurements. Such large discrepancies between the results obtained from both direct methods result from the manner of conducting the measurements; in particular, small-volume soil samples, which do not allow for the inclusion of all of the factors influencing the soil permeability, were taken for the laboratory tests. The mean Ks value obtained from the in situ measurements, i.e., 0.85 m·day⁻¹, is four times higher than those estimated by PTF-1 and PTF-2, almost half lower than that obtained from PTF-3, and twice as high as the values obtained from PTF-4 and PTF-5. In addition, the median value from the in situ tests (0.40 m·day⁻¹) is from 1.2 to 2.2 times higher than the values obtained from the four PTF models—the exception is model PTF-3, whose median value is four times higher than that from the in situ tests. The maximum and mean values of the Ks values obtained from laboratory tests, compared to those from the in situ method, are usually more similar to those obtained from the PTF models (Table 5). The mean Ks values obtained from the in situ measurements and PTF-3 indicate the mean class, whereas the results obtained from laboratory measurements and the other four PTF models indicate a moderately slow basic infiltration class [42].

Most of the subsurface horizons are characterised by low Ks values, usually not exceeding 0.50 m·day⁻¹. In variant 2, covering the analysis of all soil horizons, the mean value of Ks from the laboratory measurements (0.21 m·day⁻¹) is very close to the results obtained from the PTF-1 and PTF-2 models, almost six times lower than the Ks value from PTF-3, and three times lower than the values estimated by PTF-4 and PTF-5 (Table 6). The mean Ks values obtained from the laboratory measurements and PTF-1, PTF-2, and PTF-3 indicate the moderate class, whereas the results obtained from the other three PTF models indicate a moderately slow basic infiltration class [42].

Regardless of the amount of data taken for analysis (n = 32 or 68), the values of Ks are usually characterised by high or very high random variability (CV = 26.1–230.7%), not only in comparison to the other soil hydrophysical properties, but also compared to the methods used to determine Ks, i.e., field, laboratory, and PTF tests (Table 5 and Table 6, Figure 5). This regularity has been described in various other studies, for example, Rezaei et al. [63]. As reported by Merdun et al. [64] and García-Gutiérrez et al. [65], the high spatial variability between the soil horizons with regard to Ks are the result of the soil heterogeneity, which hampers the prediction of Ks using mathematical models.

The Ks values obtained from the laboratory measurements are statistically significant (p < 0.05), negatively correlated with the clay fraction and density, and positively correlated with the sand fraction (

Table 7. Spearman’s rank correlations coefficients (rho) between measured and predicted saturated hydraulic conductivity (K_s) and other selected soil physical properties.

Variables	% Sand	% Silt	% Clay	Bulk Density	θ2,5 pF	θ4,2 pF	Ks—Field	Ks—Laboratory	Ks—PTF−1	Ks—PTF−2	Ks—PTF−3	Ks—PTF−4	Ks—PTF−5
% Sand	1.00
% Silt	−0.77	1.00
% Clay	−0.31	−0.22	1.00
Bult density	0.23	−0.06	−0.15	1.00
θ2,5 pF	−0.29	−0.04	0.47	−0.61	1.00
θ4,2 pF	−0.30	−0.19	0.87	−0.23	0.53	1.00
Ks—Field	0.29	−0.08	−0.21	−0.24	0.06	−0.04	1.00
Ks—Laboratory	0.29	−0.18	−0.36	−0.30	−0.17	−0.20	0.54	1.00
Ks—PTF−1	0.15	0.07	−0.54	0.35	−0.50	−0.50	−0.13	0.12	1.00
Ks—PTF−2	0.34	0.09	−0.93	0.06	−0.43	−0.78	0.15	0.40	0.53	1.00
Ks—PTF−3	0.37	0.15	−0.96	0.12	−0.47	−0.85	0.21	0.34	0.55	0.95	1.00
Ks—PTF−4	0.81	−0.75	−0.17	0.17	−0.22	−0.20	0.21	0.35	0.03	0.28	0.20	1.00
Ks—PTF−5	0.84	−0.90	0.00	0.04	−0.06	−0.04	0.12	0.25	−0.04	0.13	0.06	0.90	1.00

Red color means, that the determined correlation coefficient is significant (for the significance level α = 0.05).

). Similar observations were noted by Zhao et al. [13], who conducted studies on the Loess Plateau of China. The interdependence of Ks values and soil textures were well documented by Hillel [66], García-Gutiérrez et al. [65], and other researchers [3]. In the case of in situ measurements, the soil hydraulic conductivity is only significantly positively correlated with the sand fraction (rho = 0.29). The PTF-1 model is significantly positively correlated with the soil density, but negatively correlated with the clay fraction and soil volumetric water content (θ) at water suctions of 33 and 1500 kPa (2.5 pF and 4.2 pF). The PTF-2 and PTF-3 models show a significantly positive correlation with the sand fraction and negative correlations with the clay fraction and soil volumetric water content (θ) at 2.5 and 4.2 pF. The PTF-4 and PTF-5 models are significantly positively correlated with the sand fraction, and are negatively correlated with the dust fraction (Table 7).

No significant correlations were detected between the Ks values obtained from the in situ measurements and the PTF models. However, with the exception of the PTF-1 model, the Ks values from the laboratory tests as estimated by the four PTF models are positively correlated (Table 7). For these reasons, only the linear correlations between the values obtained from the laboratory tests and PTF models are presented in Figure 6a–e.

In general, the matching of the Ks results obtained from the direct in situ measurements with the results estimated by the Rosetta PTFs is poorer than that obtained from the laboratory measurements and estimated by the Rosetta PTFs. This is confirmed by the RMSD values, which range from 1.40 to 1.64 in the first case, whereas in the second case, they range from 0.19 to 1.58. In addition, the values of the coefficient of determination (R²), which indicates the extent to which the model explains the gathered measurement data, are usually much higher when comparing the data obtained from the laboratory tests and that from the PTF models (

Table 8. Statistical indices for estimation of saturated hydraulic conductivity (Ks).

Model	Field Measured Ks (m∙day−1)				Laboratory Measured Ks (m∙day−1)
	n	RMSD	NSE	R2	n	RMSD	NSE	R2
Rosetta—PTF-1	32	1.63	−0.42	0.02	68	0.23	0.04	0.01
Rosetta—PTF-2	32	1.54	−0.34	0.02	68	0.19	0.20	0.16
Rosetta—PTF-3	32	1.64	−0.42	0.04	68	1.58	−6.34	0.12
Rosetta—PTF-4	32	1.40	−0.22	0.04	68	0.65	−2.03	0.12
Rosetta—PTF-5	32	1.40	−0.22	0.01	68	0.61	−1.83	0.06
Laboratory measured	32	1.53	−0.33	0.29

).

The RMSD values of the Ks results obtained from the laboratory measurements and estimated with PTFs reveal the best match for the case of the PTF-2 model. The NSE and R² values are the best in this case (Table 8), indicating that the granulometric composition (% of sand, silt, and clay) in the PTF-2 model has the greatest influence on the water permeability of the soils.

4. Conclusions

The search for alternatives to labour-intensive and cost-intensive in situ tests and laboratory methods for determining the value of Ks, at a time of global climate change, seems to be a priority. The increased frequency of periods with excesses and scarcities of water will require the greater use of soil irrigation and drainage facilities, and the modelling of surface runoff and erosion phenomena in land-formed areas. Correspondingly, knowledge regarding soil permeability is required for all technical and non-technical activities.

In general, when using a descriptive comparative analysis, it is difficult to show a significant relationship between Ks values from in situ tests and the diverse hydrophysical properties of arable soils, which are subject to cyclic loosening, reconsolidation, and/or freezing and thawing. The high sensitivity of this physical parameter to biotic and abiotic factors, especially in the soil upper horizons, hampered a satisfactory adjustment to the Ks values obtained from the PTF models. Therefore, in situ tests, although difficult and time-consuming, should continue to be conducted to obtain a reliable measurement database that is free of errors resulting from the soil type variability, and that considers the influences of fauna and flora on the physical properties of soil. This is of particular importance for the modelling of surface runoff and erosion processes and design of surface irrigation systems (e.g., sprinklers), where only 40–50 cm of the upper horizon of the soil is wetted, and the technical parameters of the corresponding drainage facilities are highly dependent on soil permeability.

The in situ measurements of Ks are significantly correlated only with the laboratory tests; nevertheless, the values obtained in situ are much higher, and are influenced by a number of factors, including those of a methodological nature (e.g., the larger diameter of the measuring cylinder). In contrast, the laboratory measurements show a significant correlation with the Ks values estimated using the models PTF-2 to PTF-5. The best match is found (based on the RMSD, NSE, and R² values) when using the PTF-2 model in the Rosetta program. This means that the most reliable Ks values can be obtained from the percentages of the sand, dust, and clay fractions. Therefore, there is a rational alternative, in the form of the PTF, for the costly and labour-intensive determinations of Ks conducted in the laboratory. The use of the PTF-2 model allows for the determination of the hydraulic conductivity of deeper soil horizons, where the variability of the hydrophysical properties is lower and there are large methodological limitations in conducting field work or difficulties in taking intact structural samples for laboratory tests. This type of solution can be used to collect data for designing subsurface drainage and irrigation facilities, some of the technical parameters of which are determined by the Ks values.

In conclusion, the research results confirmed the remarkable difference between the performance of PTFs in the Rosetta and the Ks results obtained from field or laboratory measurements. The study thus provided evidence that the Ks values obtained from local dataset give better results performances in predicting the soil saturated hydraulic conductivity, compared to PTFs in the Rosetta derived from large datasets, which in turn confirm the limitation of applying PTFs developed from one region to other regions.