Long-term monitoring and modelling of pasture regeneration and water quality from a Bothriochloa pertusa site in the Great Barrier Reef catchments

Study site

This study was carried out on a ‘typical’ hillslope within Virginia Park station, a commercial cattle grazing property in the Burdekin catchment, Queensland, Australia (−19.89°, 146.54°) (Fig. 1). The main drainage system on the property is Weany Creek. The area is representative of the highly erodible granodiorite (‘gold-fields’) country between the cities of Townsville and Charters Towers. The site is a 1.2 ha hillslope within a larger 625 ha paddock, and was chosen because it has been studied and characterised in detail from research spanning more than 30 years (Gardener et al. 1990; Leuning et al. 2005; Bartley et al. 2007; Koci et al. 2020, 2021). The key landscape and hillslope attributes are summarised in Table 1.

The dominant land use has been beef cattle grazing for over 100 years. Prior to 1988, pastures at Virginia Park were dominated by tall native perennial grasses, Heteropogon contortus, Dichanthium fecundum, and native legumes, with the upper story consisting primarily of narrow-leafed ironbark (Eucalyptus crebra) and red bloodwood (Corymbia erythrophloia). However, heavy grazing by Bos indicus cattle between 1979 and 1988 resulted in a significant change in the pasture composition, where the dominant species shifted to B. pertusa, annual grasses, non-legume dicots, and legumes (Gardener et al. 1990).

This study covered the period from 2003 to 2021. At the beginning of the project, during a drought, the ground cover and pasture biomass levels were on average 63% and 350 kg of dry matter per hectare (kg DM ha⁻¹), respectively, well below that considered ‘sustainable’ in terms of long term grazing management for this soil type (Ash et al. 2001). Prior to the study (1986–2002), the stocking rate was about 0.25 AE ha⁻¹, whereas in the first 10 years of the study stocking rates averaged about 0.1 AE ha⁻¹, where one animal equivalent (AE) is measured relative to a 450 kg steer at maintenance (McLean and Blakeley 2014). Previous publications documented the data outputs from this site after 6 and 10 years, respectively (Bartley et al. 2010, 2014). After the formal research project finished, stocking rates varied between years, approximated as light stocking rates from 2013 to 2016 (0.1 AE ha⁻¹, except for complete destocking in 2014) and heaving stocking rates in the five remaining years (about 0.5 AE ha⁻¹ in 2012–2013 and 2017–2020).

Vegetation monitoring

Land condition on the flume hillslope was measured using a range of techniques over the 18-year measurement period (Table 2). Measurements were made at the end of the wet season and end of the dry season each year. Between 2003 and 2016, land condition was measured across the hillslope on an 8 × 4 m grid. Between 2020 and 2021, ten 100-m transects were each sampled with ten equidistant 1 × 1 m quadrats. Sampling followed the BOTANAL calibrated visual estimation method of Tothill et al. (1992). Measurements included pasture composition (by species and functional group), pasture standing biomass, percentage pasture cover (associated with standing pasture biomass), percentage ground cover (predominantly pasture cover plus a small (≈5%) area of non-overlapping litter cover) and grass basal area. Between 2016 and 2020 field data were collected from an adjacent hillslope in the same paddock (within 500 m, sharing the same land management and comparable MODIS ground cover). Calibrated MODIS pasture and ground cover estimates (Guerschman 2020) were used when field measurements were not conducted.

Remote sensing of ground cover

MODIS (Moderate Resolution Imaging Spectroradiometer) satellite data were calibrated against field measurements (Table 2) to derive estimates of pasture cover and ground cover (R² = 0.88). High resolution (≤10 m) satellite images (Quickbird, Worldview and Sentinel-2) were analysed and calibrated using the PD54 cover method (Pickup et al. 1993) and used to map the distribution of ground cover for 7 years (2003–2007, 2009–2010).

Rainfall corrected trends in ground cover were derived from Landsat images using the dynamic reference cover method (DRCM) outlined in Bastin et al. (2012). This approach identified minimally disturbed, grazed, non-riparian reference areas in Landsat images to benchmark changes in cover through time. There may be systematic biases in the Landsat estimates of ground cover, particularly in areas dominated by B. pertusa (Bastin et al. 1996; Wilkinson et al. 2014a, 2014b; Bartley et al. 2022).

Rainfall, runoff and sediment yield monitoring

Most rain falls between December and April, but occasional out of season events occur, and therefore rainfall and runoff data extend from July 1 to June 30 (used here to define the ‘water year’, each referenced by its starting year).

To measure water and sediment runoff, a large Parshall flume was installed at the bottom of the hillslope. Details of the data logger, associated instrumentation and data processing can be found in Bartley et al. (2006). Annual percentage runoff was calculated as the total annual runoff divided by total rainfall for that water year. Water sampling was stratified by flow depth, and collected samples were used to represent total suspended solids (TSS) concentrations between 2002 and 2013. In the rare case that samples were not collected for an event, an average TSS concentration from the wet season was applied. In 2014 a turbidity sensor was installed in the flume, and a turbidity-TSS relationship (R² = 0.91) was developed, which allowed continuous sediment loads to be derived when runoff occurred.

Collected water samples were analysed for TSS concentration (considered to represent the silt (0.002–0.063 mm) and clay (<0.002 mm) sediment fractions). Bedload samples (generally between 0.063 and 8 mm) were collected manually from bedload traps but these represent <1% of total load (Bartley et al. 2006), and the data were not used in calculating fine sediment loads. When both concentration and discharge data were available, annual sediment loads were estimated by summing the event loads using the arithmetic mean approach (Letcher et al. 1999). When sediment concentration data were unavailable for an entire event, average values for that wet season were applied.

Modelling approach using GRASP

This study used GRASP (Cedar version 2.1.0), a deterministic, point-based, daily time step, native/naturalised pasture grazing systems model, developed to support research and decision making in semi-arid and tropical grasslands (Rickert and McKeon 1988; McKeon et al. 1990; Owens et al. 2003). The model and parameterisation of its key functions are described by Littleboy and McKeon (1997), with components of the model covering water balance, pasture dynamics, grazing management, runoff and soil loss. GRASP focuses on capturing and representing the gross long-term ‘average’ system behaviours that are of most relevance to property- and catchment-scale management decisions.

Data on the site soil profile (Queensland Government 2022) were used to directly estimate GRASP parameters related to rooting layer depths and water-holding characteristics of each layer (Table 3). Previous pasture measurements were used to directly calculate important parameters related to pasture growth, such as sward nitrogen, maximum pasture growth, and the relationship between pasture cover and biomass. Direct measurements for some other key parameters related to leaf senescence and defoliation were unavailable, and were instead calibrated by selecting pasture parameters that gave a good fit between modelled and observed pasture biomass. The model used daily rainfall measured at the flume site augmented with Charters Towers gap-filled historical weather station data from the SILO online repository (Jeffrey et al. 2001) (Table 3). Model parameters were set to represent the broad shifts in stocking rates described previously in the site description: the intention was not to mimic the full detail of month-to-month changes in stocking rate, but to capture periods of major differences in management. The model was reset to observed pasture states wherever observations were available, which had the effect of limiting the cumulative propagation of small errors that could cause progressive drift of the pasture component of the model from observations over 20 years.

GRASP contains several different sub-model options for representing runoff. Of these, the ‘curve number’ sub-model was used because it provided the best fit in a previous analysis of the first 10 years of Virginia Park data, and because it allowed explicit testing of whether ground cover affected runoff. Curve numbers range between 30 (maximum infiltration) and 100 (an impermeable surface with maximum runoff). The modified implementation of the curve number model in GRASP (Owens et al. 2003) allows specifying different curve numbers for 0 and 100% ground cover (and then interpolating a curve number for intervening levels of cover). Calibration involved running factorials of these two parameters with progressively finer resolution picking the combination that gave the highest Pearson correlation between modelled and observed annual runoff until the deviation of the modelled average was within 2% of the observed value for the study period.

The sub-model for fine sediment yield was parameterised similarly to runoff, using the Rose sub-model option (Rose 1985), which has parameters for slope angle, sensitivity of soil loss to ground cover, and bare soil efficiency of entrainment.

Changes in pasture composition

Although originally dominated by tussock grasses (with no B. pertusa), the site had been transformed by invasion of B. pertusa by 2002 when the flume was installed and a period of reduced grazing pressure implemented (Fig. 2). Such B. pertusa transformations are characteristic of heavy grazing over long periods, and are typical of State 3 and 4 grazing systems within state and transition models for this region (McIvor and Scanlan 1994; Ash et al. 2011). While the B. pertusa invasion was rapid (replacing most tussock grasses within 9 years), the altered state has proven very resistant to improvements in grazing management, with only slow and limited recovery of native tussock grasses over the subsequent 20 years (Fig. 2).

Fig. 2. 

Changes in relative species composition of the pasture layer at the study site measured at the end of the dry season (as percentages of total biomass). Data points for 1979 and 1988 are from previous studies at the site, predating the installation of the flume (Gardener et al. 1990).

Changes in ground cover

Ground cover at the start of the wet season is particularly important for water quality outcomes, as is the spatial arrangement of the vegetation, where cover lower on the hillslope is more important for trapping sediment (Bartley et al. 2010). Both the amount and spatial arrangement of end of dry season ground cover varied during the study (Fig. 3).

Fig. 3. 

Remotely sensed ground cover on the hillslope flume site in 2003, 2007, 2010, 2015 and 2020. The missing area in 2010 was due to cloud interference.

The weather-corrected DCRM data show a slight trend of increasing ground cover over time (not statistically significant, P = 0.134) (Fig. 4). Increases in ground cover occur not only from management practices, but also from changes in pasture composition (as discussed in the next section).

Fig. 4. 

Weather-corrected changes in end of dry season ground cover for the flume site (1990–2020). Data are for Landsat-derived ground cover processed using the dynamic reference cover method.

Modelled vs observed pasture biomass and cover

GRASP modelling reflected the observed fluctuations in pasture biomass at the study site over the duration of the investigation, however, forced resets to observations were required to correct for the gradual accumulation of errors (Fig. 5). Ultimately, the pasture component of the modelling was mainly used to interpolate the vegetation dynamics between observations (Fig. 5) to provide the daily ground cover data required by the runoff and soil loss components of GRASP.

Fig. 5. 

Fluctuations in pasture biomass over the duration of the investigation at the Virginia Park flume site. The dashed line is the unconstrained estimate from the GRASP model with approximated adjustments to stocking rates to represent periods of substantial differences in management. Field observations are marked by diamonds, and the solid line is where the model was used to interpolate between observations.

Pasture cover has a strongly asymptotic relationship to pasture biomass. The structure of pasture vegetation at the site was very different from typical tussock pastures in northern Australia because of the creeping, stoloniferous growth form of the dominant B. pertusa. B. pertusa has much less biomass for the same amount of pasture cover relative to a tussock grass (Fig. 6, where, at 50% pasture cover, the B. pertusa pasture had only a third of the biomass expected for tussock pastures in the region). This structural change in vegetation has important implications for how satellite and field-based measures of pasture are interpreted, as demonstrated in Fig. 7, showing that pasture biomass would be greatly over-estimated if it were assumed to be a typical tussock pasture rather than calibrated for B. pertusa. Although pasture cover is easier to measure, biomass is a better indicator of functional land condition both in terms of water quality (the amount of pasture material that obstructs overland flow of water) (Roth 2004) and productivity (the amount of forage and protein-building nitrogen available for livestock to consume). Inferences from pasture cover alone in B. pertusa dominated landscapes would therefore tend to give overoptimistic estimates about land condition and trend.

Fig. 6. 

Contrasting relationships between pasture cover and biomass for pastures dominated by Bothriochloa pertusa (diamonds: end of dry season observations at Virginia Park; solid line: calibrated curve in the form of the equation used in GRASP) and typical tussock grasses in the region (dashed line: parameterised curve for tussock grasses in the region).

Fig. 7. 

Fluctuations in biomass for Bothriochloa pertusa pastures at the study site for observed data (diamonds) with modelled interpolations (black line) compared to the biomass that would be implied from observed pasture cover alone if pastures were assumed to have the typical structure of tussock grasses (dashed line: calculated using the inverse relationship from Fig. 6 for pasture cover-vs-biomass in tussock pastures).

Hillslope rainfall, runoff, sediment concentrations and sediment yields

The long-term average measured rainfall for the site was 634 mm year⁻¹. The lowest annual rainfall over the study period was 245 mm (2003), the highest 1308 mm (2008) (shown with coincident ground cover data in Fig. 8a). The measurements of runoff (Fig. 8b) and total suspended sediment concentrations (Fig. 8c) over the same period combine to generate fine sediment loads (Fig. 8d).

Fig. 8. 

(a) Annual rainfall and end of dry season ground cover (%) (calibrated MODIS estimate); (b) annual runoff (mm); (c) TSS concentrations (mg L⁻¹); and (d) fine sediment yield (t ha⁻¹) for the flume site from 2003 to 2020. Red bars/plots represent years when cover <70% and blue bars/plots when cover >70%. All years are numbered as the starting year of that ‘water year’ (July–June). Box plots for (c) show upper quartile, median and lower quartile; whisker lines show maximum and minimum; and points show outliers.

The data provide support for a previously documented functional water quality target for end of dry season ground cover of ≈70% (Sanjari et al. 2009; Bartley et al. 2014): summarising the data in Fig. 8 as the averages for years when ground cover was below the 70% threshold vs averages when the 70% threshold was exceeded showed that despite rainfall (506 mm vs 718 mm) and runoff (64 mm vs 115 mm) being substantially higher in years with >70% cover, total sediment yield (0.23 t ha⁻¹ vs 0.13 t ha⁻¹) and TSS concentration (318 mg L⁻¹ vs 119 mg L⁻¹) were reduced. However, in the 2015/2016 wet season, a single event produced 90% of the annual hillslope sediment loss when 101 mm of rain fell in an hour on the 1st of March 2016. The 70% ground cover target for this site is higher than the 40% target recommended in the region for sites with a larger proportion of tussock grasses (Ash et al. 2001).

There was a strong linear correlation between the interannual variations in runoff and rainfall (R² = 0.87) (Fig. 9a), and a moderate relationship between the interannual variations in percent runoff and fine sediment yield (R² = 0.47). Using multiple linear regression, including rainfall and ground cover, the addition of ground cover improved the prediction of annual hillslope runoff marginally (R² = 0.90).

Fig. 9. 

(a) Linear relationship between annual rainfall (mm) and runoff (mm) for the 18-year record at the flume site; and (b) power curve relationship between end of dry season ground cover (calibrated MODIS estimate) and the median annual measured TSS concentration. Outliers in (b) are associated with large rainfall (high intensity) events.

Despite a moderate relationship between median TSS concentration and ground cover (R² = 0.48) (Fig. 9b), the average annual ground cover alone was not particularly useful at predicting year-to-year variation in sediment yield (t ha⁻¹) (R² < 0.2) due to the dominant influence of weather and interactions between the variables within and between years.

Modelled runoff and sediment

The calibrated GRASP curve number sub-model for runoff matched both the long-term observed average annual runoff and provided a reasonable approximation of year-to-year variations (Pearson r = 0.61) (Fig. 10a). The curve numbers that gave the best fit were 39.8 for 100% cover and 82.2 for 0% cover. This demonstrates that, despite its relatively low biomass, there was still a strong relationship between temporal changes in ground cover and runoff for B. pertusa pastures. The calibrated GRASP soil-loss sub-model matched the long-term average annual fine sediment yield but did not capture the fine-scale detailed processes that give rise to annual variations in soil loss (r = −0.09) (Fig. 10b).

Fig. 10. 

Calibrations of the (a) runoff and (b) soil-loss sub-models in GRASP for annual observations vs modelled estimates. These calibrations are for year-to-year variations at a single site. The dotted line is the 1:1 relationship. Averages for annual runoff were 93.5 ± 20.1 (observed) and 95.1 ± 26.7 mm (modelled). Average annual sediment yields were 0.169 ± 0.043 (observed) and 0.168 ± 0.117 t ha⁻¹ (modelled). Annotated biases express the modelled averages and standard errors (s.e.) as percentages of the observed values.

Land condition and water quality trends

Rainfall, and associated runoff, was the dominant, natural perturbation driving year-to-year variation in sediment loads (Fig. 9). A combination of rainfall amount and the number of large (>50 mm h⁻¹) rainfall events accounts for most of the short-term fluctuations in erosion in these systems (Dunkerley 2021), while the modifying influence of grazing management on runoff is more difficult to detect. However, the systems modelling analysis, by accounting for the complex landscape interactions, detected the influence of ground cover on runoff, as evidenced by the large differences in curve number between 0 and 100% cover. While ground cover has a relatively minor influence on the proportion of rainfall that becomes runoff at the annual time scale, it also affects suspended sediment (TSS) concentration, and thus has been shown to be an important parameter for predicting sediment loads at both hillslope (Tiwari et al. 2020) and larger catchment scales (Kuhnert et al. 2012). In fact, Kuhnert et al. (2012) found that TSS concentrations at the Burdekin catchment scale (end of system, Inkerman Bridge) decreased 2.1% for each 1% increase in ground cover, e.g. if ground cover improves by ≈10%, then TSS concentration reduces by ≈20%. Data presented here provide further evidence that ground cover does influence sediment concentrations, particularly when cover falls below 70%. In very wet years, high intensity rainfall events become more influential as the wet season progresses and can generate large amounts of runoff, regardless of the amount of cover, as the soil profile is saturated (Bartley et al. 2006; Jarihani et al. 2017). Cover is likely to exert a stronger influence on sediment loads in lower rainfall years and at the start of the wet season, or when comparing nearby properties that experience similar weather but have contrasting management.

Other plant attributes such as pasture biomass and plant basal area are likely better predictors of infiltration, and thus runoff, particularly in B. pertusa dominated rangelands (Roth 2004). There has been some recent progress in estimating biomass using remote sensing (Chen et al. 2021), however, the methods for hindcasting ground cover data using these methods and linking to historical runoff data are limited, and current approaches cannot adequately differentiate between pasture species (e.g. tussock vs stoloniferous grasses). In this regard, the data show important differences in the pasture cover–biomass relationships for pastures dominated by B. pertusa vs tussock grasses (Fig. 6), where B. pertusa pastures have far less biomass at the same cover than the pastures they replace. This would lead to overoptimistic assessments of land condition and trend (including the inferred sediment yields for a given observed ground cover) if changes in pasture composition and physiognomy are not considered. As B. pertusa expands in many grazed parts of the GBR catchments, traditional cover-based approaches to monitoring the influence of grazing land management on land condition and water quality could be compromised.

Long-term monitoring revealed several other insights into the altered functioning of transformed B. pertusa landscapes, with implications for management. Data show strong asymmetry in the rates of landscape degradation and recovery, where tussock grass regeneration rates have been far slower than their displacement rate by B. pertusa (Fig. 2). At best, there are long time lags between improved land management practices and water quality outcomes for the GBR; at worst, B. pertusa lawns are very resistant to change, and represent a transition to a relatively stable new pasture state (McIvor and Scanlan 1994; Ash et al. 2011) that needs to be managed differently to the native tussock pastures it replaced. The asymmetry in recovery vs degradation also has important implications for managing risks of adjusting stocking rates relative to (uncertain) future pasture growth, since any degradation from setting a stocking rate adjustment too high would require multiple years of intervention (with reduced stocking) for recovery. Conservative stocking rate strategies are therefore particularly favoured in land types with such asymmetric risks associated with irreversible transitions.

Comparison with other rangeland areas in the GBR catchments

The 18-year annual average fine sediment yield was ≈0.16 t ha⁻¹ year⁻¹, which is summarised alongside other hillslope fine sediment yields measured at grazing sites in the GBR catchments in Table 4. Given the variability in slope, rainfall, vegetation cover, grazing intensities and size of the hillslopes, it is not necessarily appropriate to compare these data, as factors such as plot area are known to affect hillslope sediment yield calculations (Ludwig et al. 2007). However, this summary does highlight how variable measured hillslope sediment yields can be between sites and years.

Erosion rates derived using cosmogenic beryllium (¹⁰Be isotope) suggest that the pre-grazing baseline representing the long-term average erosion rate (≈100 to >10 000 years) for this site is likely to be ≈0.27 t ha⁻¹ year⁻¹ (Bartley et al. 2015). Therefore, the contemporary ≈18-year average erosion rate from this flume hillslope site is within the range of tolerable or acceptable erosion rates (Montgomery 2007; Bui et al. 2011). However, data presented in Table 4 demonstrate that many hillslopes in the GBR catchments, including adjacent hillslopes on this property, have average hillslope soil erosion rates several times higher than this baseline value, particularly where there are large bare or scalded areas void of vegetation (Bartley et al. 2014).

Silburn et al. (2011) captured many of the key factors driving hillslope erosion rates in the rangelands draining to the GBR, including the significance of the 50–70% ground cover threshold required to reduce hillslope soil erosion and runoff, the difference between erosion rates predicted from stoloniferous and tussock grasses, the dominance of low-cover bare patches and large rainfall events on erosion rates, and the time lags associated with improved ground cover and reduced erosion rates.

Model calibration and evaluation

Given the large natural variability represented in hillslope erosion measurements in grazed rangelands (Table 4), and the complex interactions and feedbacks among variables, it is difficult to attribute measured changes in hillslope sediment yields to changed grazing management in real time. Robust models are therefore very useful for providing insights and scenario analysis into the potential effects of rangeland management interventions.

The main method for assessing investments in the GBR catchments is the SourceCatchments model (Wilkinson et al. 2014b; McCloskey et al. 2021). The key sub-modelling tool used to estimate sediment loss from grazed hillslopes within the SourceCatchments model is the Revised Universal Soil Loss Equation (RUSLE) (Renard et al. 1994; Rosewell and Loch 2002; Silburn 2011b). Numerous publications outline how the RUSLE is not necessarily suitable for rangeland environments (Weltz et al. 1998; Spaeth et al. 2003); primarily due to the higher levels of heterogeneity in soils and plant properties at the plot scale, relative to cropping landscapes (Nearing et al. 2011). Another challenge is that RUSLE does not explicitly consider runoff. GBR specific studies have also identified that the RUSLE generally over-predicts hillslope erosion (Hughes and Croke 2011; Brooks et al. 2014). Bartley et al. (2006), and other studies, have determined that biomass and the arrangement of ground cover are likely better predictors of runoff than average hillslope cover, although it is unlikely that these attributes will replace average ground cover in predictive models any time soon. Despite these findings, RUSLE continues to be used in the GBR due to its utility and ability to parameterise models that can analyse large spatial areas using sparse data sets. The RUSLE is considered broadly suitable as a prioritisation or ranking tool when used in a relative sense for comparing hillslope erosion rates between regions, acknowledging that there may be high prediction error at any one location (Alewell et al. 2019). However, this means that it may be problematic to apply the SourceCatchments model and associated tools in accounting and financial payment frameworks (e.g. Reef Credits). This is because these models have relatively high uncertainty, compared to on-ground measured data, when estimating runoff or sediment loss at small spatial scales. There is therefore a need to consider incorporating a broader suite of modelling tools for GBR rangelands, such as GRASP (or its component sub-models).

This study is the first to apply the GRASP model to a B. pertusa dominated commercial grazing property. When calibrated using local data on soils and pastures, GRASP performed well at representing long-term average runoff and sediment yield (both within 2% of observed; Fig. 10), but was only able to adequately represent the year-to-year variation in runoff (Pearson r = 0.61), and not fine scale variations in sediment. Erosion models (and complex full-systems models like GRASP) are more suited to representing long-term, broadscale average behaviours. The GRASP calibration in this paper represents a single data point (capturing long-term variation at a single spatial location). A larger number of water quality data sets, from sites with contrasting management and ground cover, would be required to properly evaluate how well GRASP and RUSLE (or other alternatives) can capture the effects of different grazing land management across the GBR catchments (at paddock to catchment scales).

The calibrated model provides valuable insights into what further large-scale encroachment of B. pertusa into existing tussock pastures would mean for GBR catchments. Given the strongly altered physiognomy when tussock pastures are replaced by B. pertusa (as characterised by the contrasting pasture cover–biomass curve relationships), a scaling factor may need to be considered when modelling water quality from un-calibrated B. pertusa dominated landscapes. These changes will have important implications for sites currently being modelled under the Paddock to Reef program (Carroll et al. 2012) given the increasing and pervasive presence of B. pertusa across large parts of the grazing areas draining to the GBR lagoon. The findings presented here will also improve the ability to model how landscape regeneration from changed grazing land management affects runoff and water quality, including being able to evaluate revised guidelines for safe stocking rates (that incorporate the consequences of the asymmetric rates of land degradation and recovery).