the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 12 Jun 2025
Soil physico-chemical indicators for ecosystem services: a focus on water regulation
Abstract. This study investigates the intricate relationship between soil properties and water-related processes, with a focus on their collective impact on ecosystem service provision, particularly water regulation. Conducted in three diverse regions Marchfeld (Austria), Bologna (North Italy) and Rmel (Tunisia), the research aims to identify key soil properties that influence water infiltration (INF), groundwater recharge (GWR), and crop water stress indexes (CWSI). Key soil characteristics such as saturated hydraulic conductivity (KS ), available water content (AWC), bulk density (BD), saturated water content (θs ), organic matter (OM), clay content and soil depth were analyzed for their role in regulating water movement and the overall hydrological balance. Pairwise correlation and multiple linear regression analyses were used to assess the interactions among soil water balance processes and soil properties. The results reveal significant variations between regions in terms of the factors that control infiltration, groundwater recharge, and CWSI. For example, in Marchfeld infiltration showed a strong positive correlation with BD (r = 0.74, p < 0.001), while CWSI had the most significant negative correlation with soil depth (r = -0.35, p < 0.001). Futhermore, multiple linear regression models were developed to assess the relevance of the different soil properties and of their interactions on the components of the soil water balance. As an example, in Marchfeld, the model for infiltration (r = 0.79, p < 0.001) was highly predictive, incorporating Clay, OM and soil depth. These results emphasize the critical role of key soil properties KS , AWC, BD, OM, clay content, θs and soil depth in controlling soil water processes. The study highlights the value of using these properties in predictive models to inform water management practices to optimize crop performance and soil conservation in different agricultural settings.
- Preprint
(1424 KB) - Metadata XML
- BibTeX
- EndNote
Status: open (until 03 Oct 2025)
-
RC1: 'Comment on egusphere-2025-1927', Anonymous Referee #1, 28 Jul 2025
reply
General comments
The authors studied the link between soil physical and chemical parameters and indicators of hydrological processes and ecosystem services. They performed an extensive simulation study with the hydrological FLOWS model for 3 contrasting regions regarding soil variability and climate. They showed that the most relevant soil properties for soil functions and ecosystem services depend on the environmental conditions in a region, which is an important innovation.
The paper reads fluently, the methodological setup is sound, and the results are well presented in Tables. However, I still have questions and remarks, as described below.
Particular comments
Figure 2: For each site, a different soil legend is used. This makes comparison hard. Can you create a figure with a single uniform soil legend?
Figure 3: In the Figure of FLOWS, you show the groundwater. Does FLOWS consider capillary rise?
Also related to this Figure: “lateral drainage” is drawn in the unsaturated zone. As the majority of lateral drainage occurs in the saturated zone, better to depict “lateral drainage” there.
Line 130: The FLOWS model plays an important role in your analysis. Can you mention studies that showed the reliability of FLOWS?
Line 147: At which soil depth did you apply the free drainage option?
Line 159: You define here the potential soil evaporation, which in general is much larger than the actual soil evaporation, an important part of your water balance. How did you calculate the actual soil evaporation?
Line 165: Is potential root uptake also assumed to be proportional to root density?
Lines 174-183: You select for each site a median year based on the amount of rainfall. Did you also check whether the rainfall distribution over the year was representative of the site's climate?
Precipitation input: At which interval did you provide input of rainfall fluxes? In order to simulate surface runoff, daily rainfall intensities are generally not suitable.
Lines 184 and 185: The data of soil texture and organic matter: to which soil layer do they belong? Did you assume these data to be uniform over the entire root zone depth (max 80 cm)?
Line 187: The crop coefficient is defined depending on the stage of crop development. How did you derive these stages?
Table 2: The water balances of Bologna and Rmel do more or less close. However, the fields in Marchfeld, on average gained 94 mm of water. This is a substantial amount compared to other water balance components and a bit strange for a median year. Please clarify this in the text.
Line 274: Minimum and maximum values are not depicted in Table 2.
Lines 417 – 433: Important that you mention the ambiguity on the role of the saturated conductivity. You might add here that the simulated amount of runoff not only depends on the saturated conductivity but also on the time intervals with which rainfall is provided. In general, daily intervals underestimate the runoff amounts.
4 Discussion: Mention important limitations of your study:
- Your hydrological results are based on median weather years. Commonly, in these years, hydrology functions well. The problems might arise in extreme dry or wet weather years. A thorough analysis of soil functions and ecosystem services should include more extreme weather years.
- You focus on the durum wheat crop. In the case of other crops or multiple crops per year, the results will be different.
- Soil layering, which is explicitly mentioned in the introduction, is not included in the analysis.
Detailed comments
Lines 110 – 114: Current text is wrongly phrased and repetitive. Define the climate of site Rmel more concisely and clearly.
Line 134: “Genucthen” should be “Genuchten”.
Line 138: The symbol “h” was already defined in line 134. Define each variable once, upon first usage in the text.
Lines 189 and 190: Include references.
Line 191: Include units.
Line 134: “groupf” should be “group”.
Table 2: In the caption, clarify CWSI.
Table 3: In the caption, clarify SD.
Citation: https://doi.org/10.5194/egusphere-2025-1927-RC1 -
AC1: 'Reply on RC1', Marialaura Bancheri, 29 Aug 2025
reply
Referee #1, 28 Jul 2025
General comments
The authors studied the link between soil physical and chemical parameters and indicators of hydrological processes and ecosystem services. They performed an extensive simulation study with the hydrological FLOWS model for 3 contrasting regions regarding soil variability and climate. They showed that the most relevant soil properties for soil functions and ecosystem services depend on the environmental conditions in a region, which is an important innovation.
The paper reads fluently, the methodological setup is sound, and the results are well presented in Tables. However, I still have questions and remarks, as described below.
A. We would like to sincerely thank the Reviewer for their careful reading of our manuscript and for their constructive comments and suggestions. Their feedback has been very helpful in improving the clarity, quality, and overall rigor of the work.
Particular comments
Q. Figure 2: For each site, a different soil legend is used. This makes comparison hard. Can you create a figure with a single uniform soil legend?
A. Although we agree with the suggestion of the reviewer, the adoption of a single legend for the three sites would result in a loss of detail in the spatial distribution of the different soil types presented in the three areas. This is due to the fact the three soil maps have a different reference scale: 1:25,000 for Marchfeld (detailed soil map, where different anthropic soils units are identified and mapped), 1:50,000 for Bologna (semi-detailed soil map, where mapping units are pedolandscapes with distinct soil series) and 1:100,000 for Rmel (reconnaissance soil map, where mapping units are in most cases associations or complexes of contrasting soil groups). This results in different levels of taxonomic details for the mapping units of the three sites, which in practice hamper the adoption of a single legend without reducing the detail of the information and the soil variability depicted in the maps.
Q. Figure 3: In the Figure of FLOWS, you show the groundwater. Does FLOWS consider capillary rise?
A. Yes, FLOWS consider capillary rise of both water and dissolved solutes. FLOWS is based on the numerical solution of the Richards equation and water fluxes are calculated according to the hydraulic gradient. In the case of a shallow water table, thus close the bottom of the simulated soil domain and interacting with the soil profile, upward water fluxes (capillary rise) will develop anytime the hydraulic head at the water table will be higher than that at the simulation node just above the water table.
Q. Also related to this Figure: “lateral drainage” is drawn in the unsaturated zone. As the majority of lateral drainage occurs in the saturated zone, better to depict “lateral drainage” there.
A. There appears to have been a misunderstanding regarding the use of the term drainage. In the FLOWS model, “drainage” does not refer to lateral drainage processes naturally occurring anytime there is a lateral hydraulic gradient. Actually, FLOWS is a 1D model and would not allow for calculating these lateral fluxes. Instead, in FLOWS this term has been used to identify fluxes to artificial drains, calculated according to the Hooghoudt’s theory, which are included in the Richards equation as a sink term.
For clarity, in the figure the term drainage was changed to artificial drainage, to distinguish it from the vertical drainage and to avoid any confusion with lateral Darcian flow.
Q. Line 130: The FLOWS model plays an important role in your analysis. Can you mention studies that showed the reliability of FLOWS?
A. The model has been used for several years in several papers related to water and solute balance, even in complex soil profiles, both at the theoretical (see for example, Coppola et al., 2009; 2012; 2015a) and applicative levels (see for example Coppola et al., 2015b; 2019a, 2019b; Hassan et al., 2024; Wang et al., 2013; Porru et al., 2024).
Coppola, A., M. Abdallah, G. Dragonetti, P. Zdruli, and N. Lamaddalena. 2019a. “Monitoring and Modelling the Hydrological Behaviour of a Reclaimed Wadi Basin in Egypt.” Ecohydrology 12, no. 4: e2084. https://doi. org/ 10. 1002/ eco. 2084.
Coppola, A., G. Dragonetti, A. Sengouga, et al. 2019b. “Identifying Optimal Irrigation Water Needs at District Scale by Using a Physically Based Agro- Hydrological Model.” Water 11, no. 4: 841. https:// doi. org/10. 3390/ w1104 0841.
Hassan, S.B.M., Dragonetti, G., Comegna, A., Lamaddalena N., Coppola A, 2024. Analyzing the role of soil and vegetation spatial variability in modelling hydrological processes for irrigation optimization at large scale. Irrig Sci (2024). https://doi.org/10.1007/s00271-023-00882-7
Coppola, A., Basile A., Comegna A., Lamaddalena N., 2009. Monte Carlo analysis of field water flow comparing uni- and bimodal effective hydraulic parameters for structured soil, Journal of Contaminant Hydrology (2008), doi:10.1016/j.jconhyd.2008.09.007.
Coppola, H. H. Gerke, A. Comegna, A. Basile, V. Comegna, 2012. Dual-permeability model for flow in shrinking soil with dominant horizontal deformation. Water Resources Research, Vol. 48, W08527, doi:10.1029/2011WR011376.
Coppola A., A. Comegna, G. Dragonetti, H. H. Gerke, A. Basile, 2015a. Simulated Preferential Water Flow and Solute Transport in Shrinking Soils. Vadose Zone J. doi:10.2136/vzj2015.02.0021
Wang X., G. Quan, Y. Pan, R. H.,Y. Zhang, A. Tedeschi, A. Basile, A. Comegna, A. Coppola, R. de Mascellis, 2013. Comparison of hydraulic behaviour of unvegetated and vegetation-stabilized sand dunes in arid desert ecosystems. ECOHYDROLOGY (ISSN:1936-0584). DOI: 10.1002/eco.1265;
Coppola A., N. Chaali, G. Dragonetti, N. Lamaddalena and A. Comegna, 2015b. Root uptake under non-uniform root-zone salinity. Ecohydrol. (2014). Wiley. DOI: 10.1002/eco.1594
Porru M.C., Hassan S. B. M., Abdelmaqsoud M., Vacca A., Da Pelo S., Coppola A., 2024. Using index and physically-based models to evaluate the intrinsic groundwater vulnerability to non-point source pollutants in an agricultural area in Sardinia (Italy). Front. Water. 6:1399170. doi: 10.3389/frwa.2024.1399170
In addition, two manuscripts specifically addressing an inter-code comparison with the HYDRUS model are in preparation and will be submitted shortly to Computational Geosciences. Their provisional titles are:
FLOWS: A model to simulate water and solutes fluxes in agricultural and environmental systems. 1. Theory and applications with inter-code benchmarking
FLOWS: A model to simulate water and solutes fluxes in agricultural and environmental systems. 2. Modules for specific applications in agro-hydrological modeling
Q. Line 147: At which soil depth did you apply the free drainage option?
A. The free drainage boundary condition was applied at the bottom of each soil profile, with the depth corresponding to the specific profile considered in the simulations. We will clarify this point in the revised manuscript.
Q. Line 159: You define here the potential soil evaporation, which in general is much larger than the actual soil evaporation, an important part of your water balance. How did you calculate the actual soil evaporation?
Q. Line 165: Is potential root uptake also assumed to be proportional to root density?
A. The two questions are somehow related. FLOWS uses the reference evapotranspiration and a crop coefficient to calculate the potential evapotranspiration, ETp. This, in turn, is separated in the potential evaporation, Ep, and potential transpiration, Tp, components by using the Beer’s law, which includes an attenuation factor for radiation penetration through the vegetation.
As for actual evaporation, Ea, FLOWS calculates the maximum upward flux at the soil surface, Emax, that the soil conditions may allow, based on the hydraulic gradient between the top boundary and the first node in the soil profile, as well as on the hydraulic conductivity in that numerical compartment. Accordingly, the code uses the following condition for calculating Ea: if Ep<=Emax, Ea=Ep; if Ep>Emax, Ea=Emax.
As for actual transpiration, Ta, FLOWS uses a macroscopic root uptake (see for example, Feddes, Kowalik, and Zaradny 1978; Feddes and Raats 2004). This approach distributes the potential transpiration, Tp, along the different simulation (numerical) nodes within the rooted depth (which is an input for the model) according to a root density distribution function, g(z). The user may opt for different root distributions, such as uniform, triangular, logistic, etc. The fraction of Tp in a node multiplied by the node compartment depth, delta_z, provides the potential sink, Sp, at that node. In other words, the integral of Sp over the rooted depth corresponds to Tp.
To obtain Ta, FLOWS firstly calculates an actual sink term, Sa, by multiplying the Sp by an uptake reduction function, α(h,hos). The sink reduction factor accounts for the possible water stress (related to soil water pressure head, h) and osmotic stress (related to the osmotic potential, hos) experienced by roots at different simulation nodes within the root zone. As the reduction uptake factor is related to the local (the specific numerical node) conditions within the root zone, there may be nodes where roots are drawing water at the potential conditions (α=1 and Sa=Sp) and stressed nodes (α<1 and Sa<Sp). The integral of Sa over the rooted depth corresponds to Ta, which may thus be <= to Tp, depending on the stress conditions along the root zone.
Accordingly, details on how the model calculates actual evaporation and transpiration will be provided in the revised text.
Q. Lines 174-183: You select for each site a median year based on the amount of rainfall. Did you also check whether the rainfall distribution over the year was representative of the site's climate?
A. In addition to selecting the median rainfall year for each site, we also examined the intra-annual rainfall distribution to ensure that the selected year was representative of the typical seasonal rainfall pattern of each site.
Q. Precipitation input: At which interval did you provide input of rainfall fluxes? In order to simulate surface runoff, daily rainfall intensities are generally not suitable.
A. The model simulations were carried out at a daily time step, with precipitation provided as daily amounts. We acknowledge that this temporal resolution does not allow for the explicit simulation of surface runoff processes driven by short-duration, high-intensity rainfall events. However, the purpose of our work was to simulate the soil–water balance at the daily scale, focusing on infiltration, evapotranspiration, and soil moisture dynamics, rather than a single event-based runoff generation. Moreover, precipitation data at shorter time steps are rarely available at the large spatial scale of our three case studies, which further justifies the use of daily inputs.
Q. Lines 184 and 185: The data of soil texture and organic matter: to which soil layer do they belong? Did you assume these data to be uniform over the entire root zone depth (max 80 cm)?
A. The data on texture, bulk density, and organic matter content were available for each horizon of the investigated soil profiles. For each horizon, the hydraulic parameters were then derived by applying the PTFs. Finally, for the purpose of the statistical analysis, we computed (i) the weighted mean of these properties over the horizon depths, and (ii) the arithmetic mean across all soil profiles within each site, in order to obtain overall spatial statistics.
Q. Line 187: The crop coefficient is defined depending on the stage of crop development. How did you derive these stages?
A. The crop development stages were determined by integrating the FAO56 guidelines with expert knowledge, literature, and site-specific field observations. This approach ensured a robust and site-appropriate parameterization of crop development for each case study. For example, based on the typical crop calendar and the growth stages described in FAO56, flowering occurs earlier at RMEL than at Marchfeld and Bologna; while sowing takes place earlier at Marchfeld and Bologna compared to the Rmel site. To further strengthen transparency, we will provide the representative development stages for each site in the supplementary material and revise the main text accordingly.
Q. Table 2: The water balances of Bologna and Rmel do more or less close. However, the fields in Marchfeld, on average gained 94 mm of water. This is a substantial amount compared to other water balance components and a bit strange for a median year. Please clarify this in the text.
A. From Table 2, the higher average water storage (91.8 mm) observed in Marchfeld soils compared to RMEL (14.4 mm) and Bologna (-3.3 mm) can be attributed to their greater soil depth (Table 1) and the resulting right-skewed distribution of soil profile depths, which allows more infiltrated water to remain within the profile. Moreover, the use of a single climatic year prevents proper closure of the water balance as would be expected when extending the simulations over longer periods.
Q. Line 274: Minimum and maximum values are not depicted in Table 2.
A. We have now updated the table to include the minimum and maximum values for each variable, providing a more complete picture of the variability across sites.
Q. Lines 417 – 433: Important that you mention the ambiguity on the role of the saturated conductivity. You might add here that the simulated amount of runoff not only depends on the saturated conductivity but also on the time intervals with which rainfall is provided. In general, daily intervals underestimate the runoff amounts.
A. Thank you for your comment. We will add the following text to clarify this aspect:
“This analysis was conducted at a daily scale; however, for certain terms of the water balance, such as potential runoff production, a finer temporal resolution would be more appropriate. This is particularly relevant in the context of increasing high-intensity, short-duration rainfall events associated with climate change (Prein et al., 2017).”
Prein, A. F., Rasmussen, R. M., Ikeda, K., Liu, C., Clark, M. P., & Holland, G. J. (2017). The future intensification of hourly precipitation extremes. Nature climate change, 7(1), 48-52.
4 Discussion: Mention important limitations of your study:
- Q. Your hydrological results are based on median weather years. Commonly, in these years, hydrology functions well. The problems might arise in extreme dry or wet weather years. A thorough analysis of soil functions and ecosystem services should include more extreme weather years.
A. We fully agree that assessing soil functions and ecosystem services under extreme weather conditions is essential. While the current analysis is based on median weather years to provide a baseline characterization, we are already working on extending the assessment to include both dry and wet extremes. This will allow us to better capture the variability and potential vulnerabilities of soil functions under climate variability and change.
Accordingly, the following text will be added to the manuscript in the discussion section at line 455:
“It should be noted that, in assessing the complex interplay between meteorological variables and soil properties, this study was conducted using a median year, which provides baseline information on climate–soil interactions. Clearly, extending such analyses over longer time spans, under more extreme climatic conditions, and at finer (e.g., hourly) temporal resolution would provide additional and complementary insights.”
- Q: You focus on the durum wheat crop. In the case of other crops or multiple crops per year, the results will be different.
A: Certainly, the results would vary with different crops or multiple cropping systems. We focused on durum wheat because it is a widely cultivated rainfed crop in the three study areas, and its growth cycle is well-documented. This allowed us to rely on robust phenological information and expert knowledge to accurately define its development stages. Moreover, the absence of irrigation simplifies the analysis by reducing variability due to water management, making it a suitable reference case for our assessment framework.
- Q: Soil layering, which is explicitly mentioned in the introduction, is not included in the analysis.
A: Soil layering was in fact considered in the analysis, but we acknowledge that this was not sufficiently explained in the manuscript. As also noted in the previous comment, we will revise the text to clearly indicate how soil layering was incorporated into the model setup and analysis.
Detailed comments
Q: Lines 110 – 114: Current text is wrongly phrased and repetitive. Define the climate of site Rmel more concisely and clearly.
A: Corrected
Q: Line 134: “Genucthen” should be “Genuchten”.
A: Corrected
Q: Line 138: The symbol “h” was already defined in line 134. Define each variable once, upon first usage in the text.
A: You’re generally right. We are going to revise and delete the repetition of the variables when it is possible, according to the specific contexts.
Q: Lines 189 and 190: Include references.
A: Corrected
Q: Line 191: Include units.
A: Corrected
Q: Line 134: “groupf” should be “group”.
A: Corrected
Q: Table 2: In the caption, clarify CWSI.
A: Corrected
Q: Table 3: In the caption, clarify SD.
A: Corrected
Citation: https://doi.org/10.5194/egusphere-2025-1927-AC1
-
RC2: 'Comment on egusphere-2025-1927', Anonymous Referee #2, 27 Aug 2025
reply
General remarks:
The authors present a study on the use of the FLOWS model for identifying the correlation between soil indicators and ecosystem services. The manuscript is properly written, describing how the study is performed adequately. However, there are a number of issues that before this study can be published.
- At lines 39 – 43 you state that the current literature is in a simplified manner showing integration of soil indicators and their overall effect on soil functions, soil processes and ecosystem services. This is a rather simplified summary of the literature review, since there are many studies actually showing an integrative, holistic approach of soil indicators (e.g. Drobnik et al., 2018, Schoenholtz et al., 2000, Dominati et al. ,2010, Dominati et al., 2016, Constantini et al., 2016, to mention a few). In fact, you are not extending your analysis to a holistic approach. You use a limited set of physical soil indicators (soil hydrological parameters, partly even taken from PTFs form soil texture data) to explain their effects on a limited set of ecosystem services (groundwater recharge, infiltration, and the calculated rest products of these). What is your contribution to your statement that current studies show a simplified manner in integrating soil indicators and ecosystem services? Clarify this in your introduction and results/discussion accordingly.
- On the methodology part: you have used the FLOWS model to calculate your identified ecosystem services (RUN, INF, GWR, CWSI). Basically, this study is a sensitivity analysis of the model, without a proper description of the model itself. You indicate that you use a ‘representative median year’ for running the model (line 170), but not what the calculation timesteps of your runs are. Coppola et al.(2019) used rain events for calibrating/validating the model, with the, among others, Ks as calibration factor. You state at line 353 that Ks has a moderate to weak correlation for infiltration (and even more pronounced at line 418). But in what way is this due to the configuration of your model-design? How does the timestep affect this, or the number of profile layers in the model setup? Is there a top-layer and a more compact bottom layer in your profile? In addition, the Ks you use is not the actual Ks, since macropore flow is not taken account for by the MvG equations. There is a limited number of publications on the use of the FLOWS model, so more information is essential to understand its concept. Moreover, you do not calibrate or validate your model runs, so it remains unclear if the magnitude of the results is realistic. Explain your model configuration, validation and perform a sensitivity analysis.
Other remarks:
Line 68: you are not evaluating chemical soil parameters, only physical. OM is often considered both chemical and physical, but still leaves you with only one (partly) chemical.
Line 75: refer to the original paper of the model (Coppola et al., 2019), although in that paper the model is also not properly described. This makes understanding the model difficult.
Line 111 – 112: Both rainfall and precipitation are mentioned, use the same phrasing.
Line 134: Typo in Genuchten
Line 199 (and elsewhere): you state at places that the runoff is calculated. But the model does not calculate runoff. Actually, at line 200 you mention the correct terminology: ‘potential runoff’. Basically, this is the outcome of the equation:
Potential runoff = Precipitation – Interception – infiltration.
The remainder is ponding, potentially available for runoff (but also for evaporation or infiltration in the next timestep of calculations). Change this term accordingly.
Line 207: Typo in ‘group’
In Results: Section 3.1 should go to ‘methods’, it describes the study area and input data.
Table 1: You present the MvG parameters, as a result of the PTFs given by Wôsten et al., 2001. I assume you are aware of the very low R2 of most of those values, making your study even more difficult without having the option for validating your model results. Describe this in your discussion, how does this fit affect your results.
Line 289: The discrepancy is not due to runoff, since runoff (better: ponding) is the final component of the calculations. So, your statement that runoff removes substantial portion of the precipitation is not true, since infiltration is calculated first.
Citation: https://doi.org/10.5194/egusphere-2025-1927-RC2
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 313 | 51 | 18 | 382 | 12 | 22 |
- HTML: 313
- PDF: 51
- XML: 18
- Total: 382
- BibTeX: 12
- EndNote: 22
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1