The study area

The study was conducted across a collection of farms to the north–west of the Ashley Irrigation Area (centered approximately on latitude 29°20′45″S, longitude 149°49′7″E), which is located north of Moree in the Gwydir Valley of northern NSW, Australia (Fig. 1). The total study area encompassed ∼81 877 hectares (ha) and included grazing areas, and dryland and irrigated cropping fields. Irrigated cotton (Gossypium hirsutum) and grains including wheat (Triticum aestivum) are the most common crops grown in the study region. The climate of the region is semi-arid with long, hot summers, cool to mild winters and an annual mean temperature range between 26.9 and 12.6°C (Bureau of Meteorology (BOM) 2021). Rainfall is summer dominant, averaging 561.3 mm annually with a mean monthly maximum precipitation of 79.1 mm in January and minimum of 22.7 mm in April (Bureau of Meteorology (BOM) 2021). Grey and Brown Vertosols are the dominant soil types of the region according to the Australian Soil Classification (Isbell 2016).

Fig. 1.

Location of the Ashley Irrigation Area relative to New South Wales, and the location of legacy soil surveys across the Ashley Irrigation Area.

Soil datasets

The soil data used in this study were obtained from two surveys conducted across the Ashley Irrigation Area (Fig. 1). Soil samples for Survey 1 were collected from a single field to the north–west of the Ashley Irrigation Area. A detailed description for the soil sampling design for Survey 1 can be found in Triantafilis et al. (2001). Briefly, soil samples were collected at 81 sites and subsampled at depth increments of 0–0.3, 0.3–0.6, 0.6–0.9, 0.9–1.2 m, resulting in a total of 405 individual samples available for Survey 1. Samples for Survey 2 were collected across the Ashley Irrigation Area from a mixture of soil cores subsampled at 30 cm depth increments down to 1.2 m and via an electromagnetic (EM) induction survey. A total of 145 sites and 1156 individual samples were available for Survey 2. A detailed description of data collection, soil sampling and laboratory analysis for Survey 2 can be found in Zhao et al. (2019). In each of the surveys, measurements were taken for , EC_1:5, ESP, and particle size fractions (Triantafilis et al. 2001; Zhao et al. 2019).

The number of sampling sites, and the number of and depth of sampling intervals varied between the two legacy soil surveys. Only soil data from the top 120 cm of the soil profile was used in this study. The total soil dataset included 1561 samples from 226 individual sites (Fig. 1).

Spatial covariates for modelling and mapping

The covariates used for modelling and mapping consisted of satellite imagery, gamma radiometrics, a silica index map, elevation, and derived terrain attributes. These covariates were at different spatial and temporal resolutions.

Landsat 7 satellite imagery was accessed through the Google Earth Engine (GEE) (Gorelick et al. 2017). All tier 1 surface reflectance imagery from 1 January 2000 to 31 December 2017 was accessed at a 30 m spatial resolution. Landsat 7 satellite imagery was chosen as it provides data over a longer time period compared to other products (e.g. Sentinel 2), and can be used to represent long-term trends using consistent data with a ∼16-day revisit time. The banding error generated by the failure of the Scale Line Corrector (SLC) in the Landsat 7 imagery was negated as statistics were obtained from many images over the available time period (Filippi et al. 2020b). A cloud-masking filter was applied to all images to remove pixels affected by cloud cover, and the 5th, 50th and 95th percentile statistics were calculated for each pixel for the entire data time-series to reflect the full range of the data. The Landsat 7 imagery was used to derive Red Band and Normalised Difference Vegetation Index (NDVI) imagery. The 5th, 50th and 95th percentiles of the NDVI imagery were selected to reflect crop growth and vigour through time and to indicate patterns of production potential. The Red Band imagery (5th, 50th and 95th percentiles) was used to represent topsoil colour. The timing of the remote sensing imagery and soil samples was not considered to be an issue as soil constraints are expected to remain relatively stable over time (Dang et al. 2011). The imagery was used to represent long-term trends to filter out more temporary, seasonal effects, rather than the surface reflectance at the time of sampling.

A digital elevation model (DEM) was obtained from the ELVIS (ELeVation Information System) platform at a 30 m resolution (Department of Finance, Services and Innovation 2020). Terrain attributes, derived from the Shuttle Radar Topography Mission (SRTM), were obtained from CSIRO’s Data Access Portal at a 30 m spatial resolution, including slope percent, degree, relief, topographic wetness index (TWI), multi-resolution ridge top flatness (MrRTF), multi-resolution valley bottom flatness (MrVBF), aspect, plan and profile curvature, and topographic position index (TPI).

Airborne gamma radiometric dose rate, uranium (Rad U), thorium (Rad Th) and potassium (Rad K) data, and their respective ratios (Rad U:Th, Rad U:K, Rad U2:Th, Rad Th:K), were obtained through the Geophysical Archive Data Delivery System (GADDS), Geoscience Australia, at a ∼100 m spatial resolution. All radiometric products were processed with low-pass filtering (Minty et al. 2009).

A silica index map at a resolution of 1:250 000, as described by Gray et al. (2016), was also used as a covariate for modelling. This index map represents the silica content of a soils’ parent material, relating to soil texture and influencing other important soil properties.

Spatial modelling of soil properties

Across the study area, soil pH, EC, and ESP were modelled. These were considered the most agronomically-important subsoil constraints to irrigated and dryland cropping in the region. Clay, sand, and soil organic carbon (SOC) content was also mapped across the area as inputs for the pedotransfer function (PTF; Padarian Campusano 2014) to estimate AWC. All data analyses were performed in the open-source software R, ver. 4.0.4 (R Core Team 2021).

Modelling soil properties

At each of the 226 sites, the associated spatial covariates described above were extracted for the study area. The 1 cm splined soil data for all soil properties were then stacked and combined with these covariates to create a single dataset. The mid-depth between the upper and lower depth of each splined interval (i.e. 0.5 cm for the 0–1 cm depth interval) was also calculated and added as a predictor variable.

Random Forest models (Breiman 2001), which are a suite of decision trees, were used to build a relationship between each soil property and the predictor variables (Table 1). Random Forests were chosen due to their efficiency in terms of identifying relationships between the response and predictor variables that can be non-linear and involve interactions (Breiman 2001). Spatial covariates were selected based on their variable importance measure, discussed below, and the combination of covariates that resulted in the best model quality. All depths were modelled together instead of modelling each depth increment individually.

Table 1.

Spatial covariates used in the spatial models for the different soil properties.

For the Random Forest models, the ‘ranger’ package, ver. 0.12.1 (Wright and Ziegler 2017), in the software R was used, which enables the fast implementation of Random Forest models that is well suited to high-dimensional data (Wright and Ziegler 2017). For each Random Forest model, the number of trees was set to 500 and the number of variables to possibly split at in each node (mtry) was set to the default option, which was the rounded-down square root of the number of variables within each model (Wright and Ziegler 2017).

Modelling available water capacity

The AWC of the top 120 cm of the soil profile was calculated using PTFs from modelled clay and sand content, and bulk density derived from SOC (Padarian Campusano 2014). The PTFs were applied at each 1 cm depth increment down to 120 cm to calculate field capacity (FC) and permanent wilting point (PWP), and AWC was calculated as the difference between the two for each layer. In the PTF, the FC was defined as the volumetric water content of an initially saturated soil following 2–3 days of drainage, also known as the drained upper limit (DUL; Veihmeyer and Hendrickson 1949). The PWP was defined as the volumetric water content remaining in the soil after a healthy crop has reached maturity in water-limited conditions with uninterrupted root development, also referred to as the crop lower limit (CLL; Hochman et al. 2001).

Interactions between depth-to-constraint and available water capacity

The interactions between the depth-to a constraint (pH, EC, and/or ESP) and AWC was assessed. For constraints reached within the top 30 cm of the soil profile, the depth-to-constraint was fixed at 30 cm as constraints were assumed to be non-limiting to plant access to water above this depth. ‘Unconstrained’ AWC was calculated as the AWC without considering the presence of constraints. The AWC for each 1 cm soil layer was summed down to the depth any constraint was reached so that the ‘constrained’ AWC map excluded soil water beyond the depth-to a constraint. The total unavailable water due to the presence of constraints was calculated as the difference between unconstrained and constrained AWC. The unconstrained and constrained AWC, and unavailable water due to the presence of any constraint, was mapped across the study area. The impact of individual pH, EC, and ESP constraints on AWC was also calculated to estimate potential yield loss (see Results - Potential yield loss due to the presence of constraints) and most-limiting constraint to yield (see Results - Most-limiting constraint to yield) but are not presented.

Potential yield loss due to the presence of constraints

The map of unavailable water due to the presence of all constraints combined, and of each individual constraint (described above), was used to determine the potential lost yield for wheat and cotton crops. Water use efficiency (WUE) equations developed by French and Schultz (1984) and Roth et al. (2013) were used to estimate potential lost yield for wheat and cotton, respectively. The French and Schultz (1984) equation estimated WUE for wheat grain yield to be 20 kg/ha.mm water. The WUE equation developed by Roth et al. (2013) estimated cotton lint yield to be >3 kg/ha.mm water, so in this study 3 kg/ha.mm water was used as a conservative estimate of yield loss. For cotton, potential yield loss was presented in bales/ha as the standard industry measurement of lint yield is 1 bale = 227 kg (Roth et al. 2013). A land-use map of NSW, obtained from the Office of Environment and Heritage (OEH), NSW Government (Office of Environment and Heritage – OEH 2017), was simplified into two categories – cropping (arable) and other. The maps of potential yield loss were masked to only include arable land across the study area.

Most-limiting constraint

The potential yield loss maps for the depth-to pH, EC, and ESP constraints were stacked and the greatest potential yield loss was calculated at each grid point. This value was assigned to the corresponding constraint (pH, EC, or ESP) and was determined to be the most-limiting constraint to yield. This was then mapped across the study area.

Soil analytical results

The range, mean, and median values were calculated for soil pH, EC, ESP, clay, sand, and SOC contents (Table 2). Soil pH was generally constant with depth (Fig. 2) with a range of 6.60–10.44 and a mean of 8.58 (Table 2). On average, soil EC_1:5 was 0.60 dS/m and increased with depth up to a value of 3.30 dS/m (Fig. 2, Table 2). Soil ESP averaged 10.99% across the study area and steadily increased with depth (Fig. 2, Table 2). Most samples were classified as medium-heavy clays, with a mean clay content of 53.97% that remained relatively constant down the profile (Fig. 2, Table 2). Sand content was steady with depth, averaging 19.56% (Fig. 2, Table 2). SOC content ranged between 0.18% and 1.06%, averaging 0.47% overall and decreased down the soil profile (Fig. 2, Table 2).

Table 2.

Soil analytical results for soil properties of splined data used for modelling and mapping.

Fig. 2.

Boxplots of soil (a), EC_1:5 (b), ESP (c), clay (d), sand (e), and SOC (f) attribute values with depth. Median values represented by solid black lines inside boxes, and the dashed vertical red lines indicate the threshold values deemed to be constraining to crop growth for pH, EC, and ESP ( 9, EC_1:5 1 dS/m, ESP 15%). No data for SOC were available for a depth > 90 cm.

Model quality

After performing the LOSOCV by whole soil profiles, the LCCC and RMSE were used to assess model predictions. The quality statistics showed that clay content was predicted with greatest accuracy (LCCC = 0.64) (Table 3). Soil EC, sand, and SOC were predicted with comparable accuracy (LCCC = 0.60), followed by ESP (LCCC = 0.50). The model used to predict soil pH performed most poorly (LCCC = 0.38) (Table 3).

Table 3.

Statistics of the quality of soil attribute models for the LOSOCV.

Variable importance

The top five most important covariates for each soil property and model are presented (Table 4). Generally, a mixture of covariates representing terrain attributes, radiometrics, and satellite imagery were important within each of the six models (Table 4). All available covariates were required to build the clay model (Table 1). Conversely, the pH model was built using the least number of covariates, using predominately satellite imagery and terrain attributes (Tables 1 and 4). Overall, radiometrics and satellite imagery, including both NDVI and Red Band imagery, were the most important groups of covariates used within all models (Table 4). Mid-depth, MrRTF, Rad U:Th, and RED 95 imagery were used consistently across each of the six models (Tables 1 and 4).

Table 4.

Ranked importance of the top five variables for each Random Forest model.

Depth-to-constraint maps

Overall, there was little similarity in the spatial distribution of the depth-to a constraint when comparing each of the soil constraint thresholds (i.e. pH, EC, and ESP) (Fig. 3). Approximately 54% of the study area was affected by at least one constraint in the top 120 cm of the soil profile. Less than 1% of the total study area was constrained by a of 9 (Fig. 3, Table 5). Over 25% of soils were found to be constrained by an EC_1:5 constraint of 1 dS/m somewhere in the top 120 cm of the soil profile, mostly between 30 and 90 cm (Fig. 3, Table 5). These were found throughout the east of the study area, as well as a smaller central region and areas to the north–west (Fig. 3). The map of the depth-to an ESP constraint of 15% showed a high degree of spatial variability and revealed that over 35% of soils were constrained by sodicity somewhere in the top 120 cm of the soil profile (Fig. 3, Table 5). Most soils constrained by an ESP of 15% were within the top 30 cm of the soil profile and were widely distributed throughout the study area (Fig. 3, Table 5).

Fig. 3.

The depth-to three constraint thresholds:  > 9, EC_1:5 > 1 dS/m, ESP > 15%.

Table 5.

The depth (at 30 cm increments from 0 to 120 cm) to which predicted pH, EC, and ESP constraints at defined threshold values were first reached in the soil profile, and the proportion of the study area this represents, including total area constrained and unconstrained.

Interactions between depth-to-constraint and available water capacity

The interactions between the depth-to shallowest constraint (i.e. pH, EC, or ESP) and AWC was mapped across the study area. Without the presence of constraints, there was little spatial variability in AWC and much of the study area was predicted to have a moderate-to-high AWC that averaged 161.89 mm (Fig. 4). As expected from the methodology, total constrained AWC averaged 111.15 mm across the study area and a lower AWC due to the presence of constraints was apparent across the north–east of the study area (Fig. 4). Soils with a higher constrained AWC were located predominately to the south of the study area and within narrow bands running throughout central regions (Fig. 4). Soil AWC with or without constraints ranged between 14.13 and 282.15 mm. The total unavailable water due to the presence of constraints was highly spatially variable and averaged 50.73 mm across the study area, ranging between 0 and 206.44 mm. More water was unavailable due to constraints in clearly defined fields to the far west and north–east of the study area, as well as a small central region (Fig. 4).

Fig. 4.

Unconstrained and constrained AWC, and total unavailable water due to the presence of all constraints within the top 120 cm of the soil profile.

Potential yield loss due to the presence of constraints

Based on calculations from WUE equations (French and Schultz 1984; Roth et al. 2013) using AWC estimates (Fig. 4), the distribution of potential yield losses for wheat and cotton due to the presence of all constraints (i.e. pH, EC, and ESP) across arable regions of the study area (Fig. 5) were spatially similar to the distribution of unavailable AWC due to the presence of constraints (Fig. 4). On average, potential yield loss due to all constraints for wheat was 1.1 tonnes (t)/ha and ranged from 0 to 3.8 t/ha (Table 6). For cotton, potential losses in lint yield averaged 0.8 bales/ha with a range between 0 and 2.5 bales/ha (Table 6).

Fig. 5.

Potential yield loss due to all constraints, and due to the presence of individual constraints pH > 9, EC_1:5 > 1 dS/m, and ESP > 15% for cotton and wheat on arable land.

Table 6.

Potential yield loss due to the presence of pH, EC, and/or ESP constraints for wheat and cotton calculated using WUE equations developed by French and Schultz (1984) and Roth et al. (2013), respectively.

As for the depth-to-constraint maps (Fig. 3), there was little similarity in the spatial distribution of yield losses due to the presence of each individual pH, EC, and ESP constraint (Fig. 5). Therefore, the greatest potential yield losses were predicted to be attributed to an ESP constraint of 15%, averaging 0.4 t/ha for wheat (range 0–3.8 t/ha) and 0.8 bales/ha for cotton (range 0–2.5 bales/ha) (Table 5). Maximum potential yield losses due to an EC_1:5 constraint of 1 dS/m were 3.5 t/ha for wheat and 2.3 bales/ha for cotton, averaging 0.4 t/ha and 0.3 t/ha, respectively (Table 5). Potential yield losses due to a constraint of 9 ranged between 0 and 3.5 t/ha for wheat and 0–2.3 bales/ha for cotton, however these impacts were negligible as <1% of the total study area was constrained (Table 5, Fig. 5).

Most-limiting constraint to yield

The most-limiting constraints to yield were identified from the potential yield loss maps, and the predictions were mapped across arable regions of the study area (Fig. 6). Over 36% of arable soils were most-limited by a sodicity constraint (ESP > 15%), followed by salinity (EC_1:5 1 dS/m) which affected approximately 23% of arable soils (Fig. 6). Less than 1% of arable soils were most-limited by an alkalinity constraint ( > 9) (Fig. 6). Yields across the remaining 38% of arable land were not limited by any constraint within the top 120 cm of the soil profile (Fig. 6).

Fig. 6.

Identification of the most-limiting constraints to yield for arable land across the study area.

Model quality and modelling approach

Random Forest models were developed to predict soil properties in 3D. As whole soil profiles were removed from the model during LOSOCV, potential interactions with other samples at different depths from the same site were eliminated. A comparison of LOSOCV techniques performed by Filippi et al. (2020a) found that validating models with splined data, as was applied in this study, is a robust approach that is comparable with using un-splined data, therefore giving us confidence in our predictions and model quality assessment used here.

Various soil properties have been modelled on comparable soils across northern NSW. Soil ESP was predicted with moderate accuracy and model performance in the current study was similar to the Random Forest-based estimates made by Filippi et al. (2020a) at the catchment scale. The quality of predictions for pH in this study were better than those obtained by Ma et al. (2021) across the Edgeroi District but poorer than those from Filippi et al. (2019a). The substantially lower sampling density here and uneven distribution of samples across the study area may explain the poorer model performance, with over 35% of all legacy soil survey samples obtained from a single field in the north–west of the study area (i.e. Survey 1, Fig. 1) potentially contributing to a more accurate prediction of soil properties in this region compared to the south of the study area. This is a common problem when using legacy soil data, as new data are expensive to obtain.

Good predictions were obtained for clay and were comparable to estimates across the Cox’s River catchment made by Bishop et al. (2015) and better than those obtained by Ma et al. (2021) across the Edgeroi District, both located just south of the Ashley Irrigation Area. Similar model quality statistics were obtained by Filippi et al. (2018b) when estimating soil EC across similar cotton-growing soils in southern NSW. Predictions of SOC content made across the Edgeroi District by Malone et al. (2009) and Ma et al. (2021) were both poorer than those obtained in this study, as were predictions obtained for the sand particle size fraction by Buchanan et al. (2012) further west across the Bourke Irrigation District.

The most important covariate for soil pH, EC, ESP, and SOC models was mid-depth, which is unsurprising considering these properties were most variable with depth (Fig. 2). The relatively high importance of radiometric covariates within almost all models makes sense, as these factors provide information about soil parent material and relate to soil forming processes (Bishop et al. 2015).

Depth-to-constraint

While subsoil alkalinity is considered an inherent property of the alluvial cotton-growing valleys of northern NSW, less than 1% of the study area was constrained by alkalinity within the top 120 cm of the soil profile. This is substantially less than predictions across the Namoi Valley, located south of the Gwydir Valley (Filippi et al. 2019a), but is not unexpected as soil pH values for the Namoi Valley area are higher than those for the lower Gwydir Valley (Singh et al. 2003).

Both soil EC and ESP values increased with depth and EC threshold limits were mostly exceeded in the subsoil between 30 and 60 cm (Fig. 3), which is typical of soils across the plains of northern NSW (Dang et al. 2006). The distribution of salinity and sodicity constraints across the study area were largely contrasting, despite soil salinity and sodicity constraints being commonly found together across the shrink-swell clay soils of northern NSW (Dang et al. 2006). The majority of soils across Australia’s cotton growing regions are considered sodic (Northcote and Srene 1972; Naidu et al. 1995), and similar proportions of an ESP constraint of 15% have been found across the Namoi Catchment (Filippi et al. 2020a).

This study mapped the depth-to fixed constraint threshold values for pH, EC, and ESP based on known agronomic impacts to yield potential (Hazelton and Murphy 2016; Isbell 2016). While fixed constraint thresholds are useful in theory, in practice they do not reflect the complex interactions between different soil properties or the different conditions under which potential yield may be constrained. This study assumed that the presence of a constraint at a fixed threshold value was an immediate barrier to root exploration down the profile (van Gool et al. 2005). However, a plants’ physical response to a soil constraint is not defined by a set value, and the magnitude of impact on root growth may vary for a range of constraint values. The use of ‘fuzzy’ boundaries that capture a range of potentially constraining threshold values and have relative impacts on root growth (and therefore yield), as opposed to rigid limits, may be more appropriate in future studies to capture the variation in response by different crops. Likewise, not all crops may incur a growth or yield penalty at the same threshold value, and future work should consider different constraint threshold limits for a range of crop types and growth conditions.

Interactions between depth-to-constraint and AWC

Uncertainty about the availability of water is considered one of the most significant limitations to irrigated and dryland crop production in Australia (Roth et al. 2013). Estimates of AWC made in this study were comparable with others across northern NSW soils (Cull et al. 1981; Malone et al. 2009). However, these estimates were either limited to single points (Cull et al. 1981), or did not explicitly consider the influence of constraints on AWC (Cull et al. 1981; Malone et al. 2009). Sadras et al. (2003) and Dang et al. (2010) developed conceptual models to estimate plant available water and plant available water capacity, respectively, considering the influence of soil constraints. However, these relationships were derived from point-based observations and were not considered spatially.

Soil AWC has been mapped across large spatial extents and at the sub-field scale. Based on in-field observations of DUL and CLL, Padarian et al. (2014) used DSM techniques to map AWC across the Australian wheat belt. However, direct measurements of FC/DUL and PWP/CLL are time-consuming and expensive (Rab et al. 2009). The application of PTFs enables the estimation of soil properties using more readily available or more easily measured soil information that can provide more accurate predictions of AWC compared to DSM techniques (Austin et al. 2019). Both Rab et al. (2009) and Malone et al. (2009) applied PTFs to estimate AWC at the district level, however neither study considered, quantified, or mapped the impact of constraints on AWC.

The PTFs applied in this study to estimate DUL and CLL for AWC were generated and calibrated in an Australian context on similar (predominately Vertosol) soils (Padarian Campusano 2014). Thus, the AWC estimates produced here (Fig. 4) address precautions around model transportability (Padarian Campusano 2014), giving us confidence in our estimations.

Lost yield potential

Many agricultural soils across Australia feature physical or chemical characteristics that constrain crop production (Bot et al. 2000). Although interactions between subsoil constraints impact yields in different ways, as discussed, the deeper a constraint is reached in the soil, the higher crop yields will be (Filippi et al. 2019a, 2020a). Dalal et al. (2002) showed that a decrease in rooting depth from 105 to 45 cm resulted in a yield reduction for wheat from 3.5 to 1.5 t/ha across southern Queensland cropping soils. However, we cannot accurately state the degree to which subsoil constraints reduce AWC across the region without appropriate validation with observed yield data. Future work should validate potential crop yield losses with real yield maps to assess the validity of WUE equation-based estimates of constrained AWC and potentially identify other drivers of yield variability.

The interactions between subsoil constraints and yield have been extensively studied. Filippi et al. (2019a) and Filippi et al. (2020a) mapped the depth-to pH and ESP constraints, respectively, and assessed their impact on crop yields at the sub-field scale. While both studies found clear relationships between the depth-to a constraint and crop yield, the impact of constraints on AWC was not quantified and much of the available literature does not consider these interactions with yield. Although the conceptual model developed by Dang et al. (2010) considered realistic yield potential as a function of the interactions between constraints and AWC, the realistic yield potential was not explicitly calculated, and this relationship was also not considered spatially.

The map of unavailable water due to the presence of constraints produced in this study was integrated with WUE equations developed from French and Schultz (1984) and Roth et al. (2013) to estimate the potential yield loss due to individual and all constraints. Potential yield loss for cotton was lower than for wheat, although this is to be expected as cotton has a lower WUE of the harvestable component compared to wheat. Although similar constraint thresholds were considered for both cotton and wheat in this study, it should be acknowledged that these theoretical threshold values may not capture the full extent to which constraints impact AWC for different crop species.

In their simplest form, the French and Schultz (1984) and Roth et al. (2013) WUE equations do not account for the influence of soil factors on available water, nor do they consider the timing of rainfall, length of growing season, or variation in soil evaporation. This study attempted to account for interactions between soil constraints and AWC, while still maintaining the original simplicity of the WUE frameworks. By calculating potential crop yield with respect to constraint-limited AWC, the impact of subsoil constraints on yield and the soil environment has been accounted for and provides a more simplistic estimation of the extent of these interactions. Future work should use soil water balance models, for example that described by Wimalathunge and Bishop 2019 to assess the true impact of reduced AWC on yield rather than assuming lost AWC equals lost water over the entire season.

Practical applications

The effective amelioration and management of subsoil constraints relies on useful decision-making tools. Mapping the spatial distribution of subsoil constraints in 3D is not a novel approach, however, the integration of this information with PTFs to estimate AWC has not been explored in detail. This study explored a potential new workflow for identifying the most-limiting constraint to yield by considering the interactions between depth-to subsoil constraints, AWC, and potential crop yield. This single map presents these combined interactive effects and simplifies multiple sources of information into a single decision-making tool. This can then be used by growers and farm managers to guide the targeted amelioration or management of soil constraints.

The identification of limiting constraints to yield has long been a goal in research. Shatar and McBratney (2004) applied the boundary-line analysis (BLA) method to compare actual and potential yields and identify the most-limiting constraints to yield at the field scale. While the map produced by Shatar and McBratney (2004) is visually similar to the map identifying the most-limiting constraint produced in this study (Fig. 6), a BLA requires high-quality and accurate yield and soil information that is often not attainable by growers on-farm. Further, many BLA studies fail to consider the depth at which constraints were reached in the soil profile, nor do they consider the impact of constraints on AWC.

From this study, it is evident that ameliorating or managing subsoil sodicity across much of the study area should be a priority, while particular fields to the north should be managed to address salinity constraints (Fig. 6). While chemical and mechanical amelioration strategies, such as the application of gypsum in combination with cultivation, have been shown to be effective in treating soil sodicity (Page et al. 2018), understanding the depth-to and nature of a constraint is important for determining the feasibility of amelioration strategies and minimising unnecessary expenditure on resources.