Simulating grazing beef and sheep systems

• Animal-plants-environment interactions in complexed livestock systems need process-based models to be investigated. • Functions in an existing model were extend and validated based on a datarich monitoring farm platform. • The extended model accurately simulated liveweight changes of beef and lamb within the natural variations expected. • Simulated emissions (N2O, NH3, CH4 and CO2) and components in the nitrogen budget are compatible with the published values. • The model has the potential to investigate the consequences of agronomic practices and grazing strategies. A R T I C L E I N F O Editor: Mark van Wijk


H I G H L I G H T S G R A P H I C A L A B S T R A C T
• Animal-plants-environment interactions in complexed livestock systems need process-based models to be investigated. CONTEXT: Ruminant livestock make an important contribution to global food security by converting feed that is unsuitable for human consumption into high value food protein, demand for which is currently increasing at an unprecedented rate because of increasing global population and income levels. Factors affecting production efficiency, product quality, and consumer acceptability, such as animal fertility, health and welfare, will ultimately define the sustainability of ruminant production systems. These more complex systems can be developed and analysed by using models that can predict system responses to environment and management. OBJECTIVE: We present a framework that dynamically models, using a process-based and mechanistic approach, animal and grass growth, nutrient cycling and water redistribution in a soil profile taking into account the effects of animal genotype, climate, feed quality and quantity on livestock production, greenhouse gas emissions, water use and quality, and nutrient cycling in a grazing system. METHODS: A component to estimate ruminant animal growth was developed and integrated with the existing components of the SPACSYS model. Intake of herbage and/or concentrates and partitioning of the energy and protein contained in consumed herbage and/or concentrates were simulated in the component. Simulated animal growth was validated using liveweight data from over 200 finishing beef cattle and 900 lambs collected from the North Wyke Farm Platform (NWFP) in southwest England, UK, between 2011 and 2018. Annual nitrous oxide (N 2 O), ammonia, methane and carbon dioxide emissions from individual fields were simulated based on previous validated parameters.

Introduction
We are at a critical juncture for global livestock production as competing requirements for maximal productivity and minimal pollution have driven the requirement for sustainable intensification (Springmann et al., 2018). Ruminants make an important contribution to global food security by converting feed that is unsuitable for human consumption into high value food protein, demand for which is currently increasing at an unprecedented rate because of increasing global population and income levels (Tilman and Clark, 2014). Reduction in red meat and dairy intake is increasingly seen as a pathway to improving human and environmental health (e.g. Westhoek et al., 2014), but globally, ruminant livestock will be important for the foreseeable future and demonstrating methods for sustainable production will become increasingly important. Sustainable intensification of ruminant livestock may be applied to pastoral grazing, mixed-cropping, feedlot, and housed production systems. All these systems have associated environmental impacts such as water and air pollution where greenhouse gas (GHG) emissions, soil degradation and erosion are all of particular concern. In addition, factors affecting production efficiency, product quality, and consumer acceptability, such as reduced animal fertility, health and welfare, also impact on the sustainability of ruminant production systems. These challenges necessitate multidisciplinary solutions that can only be properly researched, implemented and tested in real-world production systems (Eisler et al., 2014). As a consequence, there is a call to 'redesign' livestock systems, including the integration of both crops and livestock (Dumont et al., 2014). These more complex systems can be developed and analysed by using models that can predict system responses to environment and management.
Several reviews of grassland-based ruminant production models have been published (Bateki et al., 2019;Bryant and Snow, 2008;Snow et al., 2014). In order to simulate ruminant livestock systems, the components of animal genetics (breed), nutrition (forage), management practices and their subsequent impact on the surrounding environment (emissions to air and water) must be considered as a whole in computational models. Several mechanistic process-based simulation models have attempted to simulate the whole system, e.g. the Hurley Pasture Model (Thornley, 1998) and its subsequent revisions -PaSim (Graux et al., 2011), WFM (Neil et al., 1999), GRAZPLAN (Donnelly et al., 2002), GrazeIn (Faverdin et al., 2011), SEDIVER (Martin et al., 2011), e-Dairy (Baudracco et al., 2013) and LiGAPS-Beef (van der Linden et al., 2019). Challenges remain in modelling ruminant systems, due to the symbiotic relationship between rumen microbial anaerobic fermentation and subsequent mammalian metabolism of a combination of derived rumen microbial biomass (microbial protein), fermentation byproducts (volatile fatty acids and ammonia) and dietary components which by-pass rumen fermentation. As well as associated microbial activity which influences lipid profiles (biohydrogenation), atmospheric pollutants (methanogenesis) and which ultimately drives the partitioning and retention (milk, live weight, faeces, urine) of dietary nutrients. A systems approach to investigate ruminant production through modelling and simulation is therefore recommended (Dougherty et al., 2019;Hirooka, 2010).
The SPACSYS model (Wu et al., 2007) is a weather-driven dynamic simulation model at a field scale with up to a daily step. Since it was first published in 2007, it has been developed to provide added functionality, e.g. the impact of vernalisation on overwinter crops (Bingham and Wu, 2011), biological nitrogen (N) fixation by legumes (Liu et al., 2013), microbial-based N 2 O emissions (Wu et al., 2015) and soil phosphorus (P) cycling (Wu et al., 2019). The model can simulate the interactions of soil carbon (C), N and P, plant growth and development, water redistribution and heat transformation in agricultural fields. The model has been applied to grassland systems in the assessment of GHG emissions (Abalos et al., 2016), responses to environmental change (Ehrhardt et al., 2018;Wu et al., 2016) under various climatic and soil conditions, nutrient cycling (Carswell et al., 2019b) and C fluxes (Sándor et al., 2020). However, as there is no component to describe animal growth, simulations involving animals have required pre-processing and direct input of data on grass intake rate and nutrient returns in animals, rather than deriving directly from animal performance, constraining model application.
This study presents a framework that dynamically models animal and grass growth, nutrient cycling and water redistribution in a soil profile taking into account the effects of animal genotype, climate, feed quality and quantity on livestock production, GHG emissions, water use and quality, and nutrient cycling in a grazing system, using a processbased and mechanistic modelling approach. Simulated animal growth was validated using liveweight data collected from over 200 finishing beef cattle and 900 lambs collected from the North Wyke Farm Platform (NWFP) in southwest England, UK, between 2011 and 2018. The framework could potentially integrate economic, environmental and social factors to provide decision makers with the ability to forecast, interpret and respond to potential threats to UK livestock production systems.

SPACSYS model
In this study, a component to estimate ruminant animal growth, AnimalCom, was developed, implemented and integrated with the existing components of SPACSYS (Fig. 1). Existing model components are published in detail elsewhere (Wu et al., 2007(Wu et al., , 2015(Wu et al., , 2019, while the new AnimalCom is described here.

AnimalCom
The AnimalCom component consists of two parts: intake of herbage and/or concentrates and partitioning of the energy and protein contained in consumed herbage and/or concentrates. Herbage intake by grazing ruminant livestock is assumed to be regulated by one of three factors (Loewer et al., 1983): a) the physiological limit on intake (or thermodynamic limit), b) the feed availability and c) the physical ability of the animal to consume feed (Fig. 2). The first factor is partially determined by the energy requirement of the animal. There are several systems developed for nutritional evaluation (Tedeschi et al., 2005), where the description below is mainly based on the UK Agricultural and Food Research Council metabolizable energy (ME) and protein (MP) system (Agricultural and Food Research Council, 1993) in which the dynamics of the rumen microbial population plays a vital role.
Metabolisable energy intake (MEI) through grazing and concentrate feeds is partitioned among that required for maintenance, pregnancy (for cow and dry ewe only), growth and fattening, and milk production (for cow and lactating ewe only). In general, the requirement for animal maintenance is given the highest priority, then pregnancy, and the lowest for milk production and liveweight gain.

Energy requirements
Physiological ME requirement (MJ head − 1 day − 1 ), as shown in Fig. 2 as "animal requirement" is defined by a generic term: where E req-main is the ME requirement for maintenance (MJ head − 1 day − 1 ), E req-growth is the ME requirement for growth and fattening (MJ head − 1 day − 1 ), E reg-preg is the ME requirement for pregnancy (MJ head − 1 day − 1 ), and E req-milk is the ME requirement for milk production (MJ head − 1 day − 1 ).
2.2.1.1. Energy requirement for maintenance. Following Agricultural and Food Research Council (1993), E req-main including fasting metabolism and activity allowance for the animal is estimated by: where LWT is the live weight of the animal (kg), a, b and c are empirical parameters and e e-main is the efficiency of utilisation of ME for maintenance. For sheep, it has been documented for some time that the AFRC (1993) model may underestimate maintenance energy requirement (e.g. Yang et al., 2019). The equation adopted here is therefore that used in the UK inventory model for agricultural GHG emissions, developed by Steven Anthony (pers comm. ADAS, 2021). Consequently, the requirement is estimated as: where k 2 , k 3 and k 4 are constants and setting 0.26, 0.03, and 1.0 (female and castrate male) or 1.15 (entire male), respectively; D w is the weaning length of sheep (days); A sheep and A lamb are the age (days) of sheep and the lamb, respectively; min and max are math functions for a minimum and maximum value of two values, respectively. This adaptation to AFRC (1993) added a further 9% of MEI to the requirement.
2.2.1.2. Energy requirement for pregnancy. E reg-preg is estimated as (Agricultural and Food Research Council, 1993): where a e , b e , c e and d e are parameters, D preg is the pregnancy period (days) and e e-preg is the efficiency of utilisation of ME for the conceptus.
2.2.1.3. Energy requirement for liveweight change. E req-growth is based on the potential live weight growth rate that is expressed as the Gompertz function (Lewis et al., 2002;Taylor, 1968): where g f is a Gompertz constant that tends to be smaller as the mature size becomes larger (Emmans and Kyriazakis, 2000), LWT m is the average weight of the animal at maturity (kg), e g is the ME requirement per unit live weight increase of the animal (MJ kg − 1 ), and e e-growth is the efficiency of utilisation of ME for liveweight change.
For finishing beef cattle, however, the potential energy requirement for growth and fattening is determined by the potential gain in protein (ΔP) and fat content (ΔF) of the empty body weight (Topp, 1999): where P e (MJ kg − 1 ) and F e (MJ kg − 1 ) are the energy values of protein and fat, respectively.

Energy requirement for milk production.
Energy requirement for milk production from a lactating animal is estimated by: where e m is the energy requirement per unit milk produced (MJ kg − 1 ), e emilk is the use efficiency of ME for milk production and Y milk is potential milk yield (kg head − 1 d − 1 ), that is controlled by a lactation curve described by Wood (1980) and then corrected by the period of milk production and the weeks of calving (and lambing) from the beginning of a year (Mainland, 1985).
where f w and f c are parameters to reflect seasonal and calving (lambing) date effects on milk production, which are the tabulated functions of weeks from the beginning of a year. Y init is the initial yield (kg head − 1 d − 1 ) of milk and affected by lactation number. a w and b w are parameters. The efficiencies of ME utilisation are determined by a linear function of the metabolisability (M e ) of gross energy at maintenance. For grazed grasses, it is estimated by: where ME g and GE g are the ME and the gross energy content of the forage (MJ kg − 1 ), respectively.
where C d-M is the conversion coefficient from digestible to metabolisable energy, with a default value of 16 MJ kg − 1 (Agricultural and Food Research Council, 1993), and d g is digestibility of forage, i.e. D-value (g kg − 1 DM), that is estimated in the model. GE g was calculated as (Murray, 1991): where CP is the crude protein content (g kg − 1 DM) of the grass and estimated by the N content of the grass multiplied by 6.25.

Protein requirements
Following Hulme et al. (1986), the protein requirement for maintenance (P reg-main ) is estimated as: where e p-main is the conversion efficiency of metabolizable protein to net protein.
The protein requirement for pregnancy (P reg-preg ) is estimated as (Agricultural and Food Research Council, 1993): where a p , b p , c p and d p are parameters and e p-preg is the efficiency of utilisation of protein for the conceptus. Protein requirement for growth is estimated as: where f p-growth is the fraction of protein in faeces. The protein requirement for milk production (P reg-milk ) is estimated by: where P per is the protein percentage in milk, f true is the fraction of true protein in milk and e p-milk is the efficiency of utilisation of protein for milk production.

Herbage intake
Mechanisms for the long-term regulation on feed intake are still unclear and will differ between grazing and stall feeding. It was assumed that actual daily intake for the animal (DMI, kg DM head − 1 d − 1 ) is determined by the most limiting factor among feed availability, physical ability and physiological requirement for intake: where DMI PHYSICAL is the physical ability on herbage intake (kg DM head − 1 d − 1 ), DMI PH is the herbage intake (kg DM head − 1 d − 1 ) based on energy requirements, and DMI G is the intake rate (kg DM head − 1 d − 1 ) based on herbage availability in the field.

Physical ability.
Feed intake by the animal is controlled by the rate of passage through and the amount of undigested material in the digestive tract. For cattle, this is used (Kahn and Spedding, 1984): and for sheep (Blaxter et al., 1961): where d max and dg DM are the average faecal DM output rate per unit liveweight (kg DM day − 1 ) and digestibility of feed, i.e. D-value (g kg − 1 DM), respectively.

Physiological requirement.
DMI PH is regulated by the daily ME requirement of the animal (Topp, 2001), and given by: where E conc is the daily ME intake rate of concentrates if supplied (MJ head − 1 day − 1 ), and M Fod is the ME (MJ kg − 1 DM) of ingested herbage, defined by Moran (2005):

Feed availability.
When the quantity of herbage available for consumption is less than that required for 95% of maximum daily intake, the daily allowance of green herbage regulates intake. The green herbage allowance is taken to be the green herbage mass above the minimum herbage mass required for grazing. DMI G was estimated as (Zemmelink, 1980): where p shape is a constant, H is the daily allowance of green herbage for the animal (kg DM head − 1 d − 1 ) and I max is the maximum daily intake of herbage (kg DM head − 1 d − 1 ) and is described by: where F max is the maximum DM intake rate per kg of metabolic weight (kg DM (liveweight) -0.75 head − 1 d − 1 ).

ME intake partitioning
There are four possible scenarios to partition ME intake depending on ME supply and animal requirements (Tess and Kolstad, 2000;Topp, 1999), meeting: 1) the physiological requirements of the animal (MEI ≥ E PH_C ); 2) the maintenance and pregnancy requirements but not the potential energy requirements for milk production and growth and fattening (E PH_C > MEI ≥ E req-main + E req-preg ); 3) the maintenance requirements but not pregnancy and the potential energy requirements for milk production and growth and fattening are not fulfilled (E req-main + E req-preg > MEI ≥ E req-main ); and 4) no requirement (MEI < E req-main ).

Scenario 1.
In this case, all the requirements can be met, and potential milk production (Eq. (9)) and growth (Eqs. (6) or (7)) will be achieved.

Scenario 2.
The energy requirements of the animal for maintenance and pregnancy are met. The energy deficit (ME d , MJ head − 1 d − 1 ) for milk and liveweight change is: It was assumed that the energy deficit is partitioned in equal amounts to reductions in potential energy requirements for milk and growth, i.e.
If E a_growth ≥ 0, milk production and growth are calculated based on E a_growth and E a_milk .
If E a_growth < 0, then maternal body tissue will be catabolised for milk production (ΔE m ) by: with the rate of change in body weight as: and milk production estimated as: where N l is the net energy produced per unit of catabolised liveweight (MJ kg − 1 ) and k bm is the efficiency of utilisation of maternal body tissue for milk production.

Scenario 3.
The ME requirement for pregnancy (ΔE p ) and milk production (ΔE m ) are met from maternal tissue catabolism: and where k bc is utilisation efficiency of maternal body tissue for pregnancy. Actual milk production and liveweight change rate are: 2.2.4.4. Scenario 4. The ME required from maternal body tissue to meet the maintenance (ΔE ma ), pregnancy (if needed) and milk production that is estimated by Eq. (31).

Protein degradation
Protein degradation in the rumen and efficiency values for the degradation followed the metabolisable protein system proposed by AFRC (1993). During the degradation of intake crude protein, fermentable metabolisable energy from MEI is incorporated to estimate microbial crude protein supply (MCP), shown in Fig. 2  emission rate (g C head − 1 d − 1 ) from an adult dairy or beef cow was estimated by: However, for a lamb or ewe, the rate was estimated (CIGR, 2002;Haque et al., 2014) by: where P a is atmospheric pressure (Pa), T a is air temperature ( • C) and HP is heat production (Watt): where GE, ME and DE (MJ kg − 1 DM) are the gross energy, ME and digestible energy in the dry matter intake, including forage and concentrate, respectively, and GEI is the gross energy intake (MJ head − 1 d − 1 ). Values for GE, ME and DE were estimated as the weighted average across forage and concentrate. Following Stergiadis et al. (2015), digestible energy from grass was estimated as: For sheep and lamb, the equation proposed by Blaxter and Clapperton (1965) was adopted: where 0.05565 (MJ g − 1 ) is the energy generated by CH 4 (de Haas et al., 2011).

Case study grazing system
Simulated animal performance was validated with data collected from the NWFP from 2011 to 2018 (50 • 46 ′ 10 ′′ N, 3 • 54 ′ 05 ′′ W and 120-180 m a.s.l.). North Wyke has a temperate climate with average annual precipitation of 1030 mm and mean daily minimum and maximum temperatures of 7.0 and 13.6 • C, respectively, from 1989 to 2018. The site overlays clay shales and the predominant soil type is a Stagni-vertic Cambisol under the FAO classification (Harrod and Hogan, 2008), which comprises a slightly stony clay-loam topsoil, overlying a mottled stony clay derived from the carboniferous culm measure. The platform is a 63 ha systems-based experimental facility divided into three 21 ha farmlets (small farms) with five hydrologically isolated subcatchments in each. Over the simulation period, the farmlet treatments (pasture-type) were one of permanent pasture (PP) predominantly perennial ryegrass (Lolium perenne L.), monoculture reseed with high sugar perennial ryegrass (L. perenne cv Aber Magic) and a reseed mixture of high sugar perennial ryegrass and white clover (Trifolium pratense L.). From April 2011 to March 2013, the baseline period, all three farmlets were as one (PP) with no separate treatments in operation. From April 2013 to September 2015, the two reseed farmlets transitioned into a post-baseline phase with the third continuing as PP. Thus, some subcatchments entered a post-baseline phase much earlier (say in 2013) than others (say in 2015). Given this and to furnish a long time series of consistent / coherent data for a robust calibration, validation and interpretation of the SPACSYS simulations, only outputs from the PP farmlet ( Fig. 3) were reported in this study. The size of each the five subcatchments and the seven fields / paddocks for the PP farmlet together with management activities are shown in Table 1.
For the study period from 2011 to 2018, each farmlet was grazed by 30 finishing beef cattle and typically 75 ewes with their lambs (typically 135 assuming a lambing rate of 1.8). Cattle were introduced to the farm platform after weaning, at 6 months of age, and were initially housed over the winter period (typically from October to March) and fed silage harvested from their respective farmlet, and then grazed on their respective farmlet at turnout until removed for slaughter on achieving a target weight of 555 kg (heifers) / 620 kg (steers) and fat class (4 L). Ewes typically grazed into the winter season (late November to early January) and were then housed and fed off the platform prior to lambing, which normally occurs between mid-March and early April; they were subsequently returned to the platform (typically March/April) with their lambs, which were finished at a target weight of 43 kg and fat class (3 L). All animal movements were recorded using unique identifier tags. Prior to 2017, a Hereford x Friesian herd provided predominantly Continental x calves, with heifers first calved to a Hereford bull. The breeding herd was subsequently transitioned to Stabilisers™. In total, seven breeds dominate: Charolais cross (CHX), Hereford (HEX), Limousin (LIMX), Stabiliser cross (STX), Stabilisers (ST), Simmental cross (SMX) and Belgian Blue cross (BRBX). Lamb were predominantly progeny of Suffolk x Mule (SUFMU) ewes crossed with Charolais (CHA) or Lleyn (LLE) rams.
In total, data for over 200 finishing beef cattle and 900 lambs were used in this study for the period 2011 to 2018. This resulted in 1383 periodic liveweight beef cattle measurements (across the above seven cattle breed combinations) and 3997 periodic liveweight lamb measurements (across the above two sheep breed combinations) for use in model performance assessment.

Simulation configurations
The SPACSYS model has previously been validated using the NWFP data in terms of water fluxes, N 2 O emissions, grass biomass accumulation and soil C and N budgets (Carswell et al., 2019b;Li et al., 2017;Liu et al., 2018;Wu et al., 2016). Hence, all the initialised states and parameters for soil water redistribution, heat transformation, grass growth and soil C and N cycling were adopted from previous studies. No further validation on these variables was made for this study. Information on agronomic management, animal movement and liveweight was freely accessed and downloaded from the NWFP data portal (http://resources. rothamsted.ac.uk/farmplatform). In most study years, liveweight was measured once every two weeks, while in the latter years this was reduced to once every four weeks.
The growth rate of each animal was simulated from its first grazing day in a field of the PP farmlet to its last day in the PP farmlet. The record of the birth date and liveweight at the beginning of grazing for each animal was used as model input for the initial weight and age of the animal. For each animal, dates of movements between fields were used to determine the grazing period within a given field. It was assumed that feed supply during the off-paddock periods was sufficient to meet the animal requirements. Weaning date for lambs each year is set at the end of June. For ewes, if there was no initial weight recorded, a default value of 70 kg was applied.
Annual NH 3 , CH 4 and CO 2 from both animals and soils and N 2 O from soils were simulated and N cycling in each field was analysed. A hydrological year from April to March was used to calculate annual values.

Statistical analysis
To assess the performance of the finishing beef cattle and lamb liveweight simulations, six accuracy indices were found (the mean error (MErr), the percentage bias (PBIAS), the mean absolute error (MAE), the normalized root mean square error (NRMSE), the Nash-Sutcliffe efficiency (NSE), and the Kling-Gupta efficiency (KGE)), which can be respectively defined as: Where N is the total paired number, ẑ i are simulated liveweight values, z i are measured liveweight values,ẑ i is the mean of the simulated values, z i is the mean of the measured values, z max and z min are the maximum and minimum values among the measured data, respectively. Further, r is the Pearson product-moment correlation coefficient (between simulated and measured) and σẑ and σ z are the standard deviations for the simulated and measured data, respectively. The ideal value of the four error-based indices (MErr, PBIAS, MAE and NRMSE) is zero such that the closer to zero, the more accurate the model simulation. Negative MErr and PBIAS values indicate a tendency to under-prediction, while positive values indicate a tendency to overprediction of liveweight by the model. NSE takes values from − ∞ to 1, where unity corresponds to an exact match between simulated and measured data, zero indicates that the simulations are as accurate as the mean of the measured values and a negative value indicates that the simple arithmetic mean of the measured is a better predictor than the model. KGE incorporates the correlation coefficient r, the ratio between the means of the simulated and of the measured data and the variability ratio. As with NSE, KGE takes values from − ∞ to 1. MErr and PBIAS provide complementary error indices for over-and under-prediction, MAE and NRMSE provide complementary error indices for prediction accuracy, while NSE and KGE provide complementary indices for levels of agreement between simulated and measured data. Performances indices are calculated using the 'hydroGOF' R package. A seventh model performance index is also reported with the usual R 2 value (the coefficient of determination) for a regression fit to the simulated and measured data. Performance indices are found across different animal ages, breeds and grazing years.
Using simulated outputs only, one-way ANOVAs were used to test differences in liveweight, growth rate, CH 4 and CO 2 emissions between cattle and sheep breeds. Note for CH 4 and CO 2 emissions, no model validation data exist.

Model performance assessment with measured data
Model performance indices per individual liveweight measurement are presented in Tables 2 and 3. Performance indices (MErr, PBIAS, MAE, NRMSE, NSE, KGE and R 2 ) were found conditional to age, breed and grazing year. Graphical depictions of model performance according to breed are presented in Fig. 4, where animal age was plotted against average liveweight for simulated and measured data. In each case, the plots for simulated and measured data were fitted with polynomial functions so that simulated (on average) growth curves could be visually assessed against measured (on average) growth curves. The first model assessment (Tables 2 and 3) is more detailed as it is conducted on each individual liveweight measurement, while the second assessment (Fig. 4) is broader as it is conducted on average liveweights.
For beef cattle, all accuracy indices (Table 2) suggest model performance moves from a high to a low level of accuracy as animals age. There was no consistent over-or under-prediction given MErr and PBIAS could be both positive and negative. Model accuracy was poor for animals that were aged around 600 days and over (e.g., with NSE dropping to 0.19 and R 2 dropping to 0.29 for the 600-to 660-day class). Although this threshold coincided with a sharp decrease in observations as animals reached their target weight ready for slaughter. For cattle breed, all accuracy indices suggest little difference in model performance, where prediction accuracy was commonly strong, and where under-prediction was more likely than over-prediction (as MErr and PBIAS tended to be negative). Strongest levels of model accuracy were found for the SMX breed (with NSE, KGE and R 2 values all close to unity), while weakest levels of accuracy were found for BRBX cattle. Both SMX and BRBX breeds were relatively small in number, where the predominant breeds (LIMX, HEX and CHX) were all predicted with strong levels of accuracy. Similar to cattle, all indices for lambs (Table 3) suggest model performance moves from high to low levels of accuracy as animals age. However, for lambs, prediction accuracy became weak at a relatively   early age, with NSE dropping to 0.09 and R 2 dropping to 0.47 for animals in the 140-to 160-day class (and this performance became weaker still for all remaining ages). There was also a consistent over-prediction of liveweight for each age class given that MErr and PBIAS were always positive. For lamb breed, the accuracy indices were more diverse and harder to interpret, where moderate to strong levels of accuracy were found for the dominant SUFMU and CHA breeds (with R 2 and NSE values of 0.80 and 0.64, and 0.70 and 0.47, respectively), while moderate to weak levels of accuracy were found for the LLE breed (i.e. an R 2 = 0.67 coupled with a poor NSE = − 1.45). Again, there was a consistent over-prediction of liveweight given that MErr and PBIAS were always positive. For lambs by grazing year, most accuracy indices suggest moderate to strong levels of model performance across all eight years with, R 2 values >0.60 and KGE values >0.74. However, the NSE index only suggested moderate to strong levels of model performance for some years (say 2011, 2012, 2014, 2017 with NSE > 0.51) and not others (say 2013, 2015, 2016, 2018 with NSE < 0.51). In summary, overall model accuracy for simulating a lambs liveweight (regardless of age, breed or grazing year) was moderate to strong with NSE = 0.50, KGE = 0.83 and R 2 = 0.72 (i.e. one out of three indices were close to unity).
For the graphical descriptions of model performance, where average liveweights are assessed against age (Fig. 4), liveweight for all seven cattle breeds was simulated with strong levels of accuracy (as the fitted polynomials were highly similar). However, such levels of accuracy could weaken as cattle get older, confirming that found above. For the BRBX breed, the model tended to under-predict liveweight (as the fitted polynomial to the simulated data mostly lies below that for the fitted polynomial to the measured data), while conversely the model tended to over-predict liveweight for the LIMX breed. Such clear under-or overprediction was not present for the remaining five breeds. Liveweights  for the lambs tended to be over-predicted for all three breeds, where this over-prediction was stronger as the lambs aged for both CHA and LLE breeds. Overall, model simulation accuracy for a lamb's liveweight was weaker than that found for a cow's liveweight.

Simulated performance for different breed combinations
Simulated liveweight gain and gaseous emissions rates during the grazing period from individual beef and sheep (lamb) breed combinations are shown in Table 4. For cattle, STX emitted the least CH 4 per head compared with other cattle breed combinations, while SMX had the highest emission, although when expressed on a per LWG basis there were no significant differences. There was no significant difference in CO 2 respiration among the cattle breed combinations. There was a significant difference in the growth rate between sheep breed combinations, with LLE at the greatest rate and CHA at the least. Across the sheep breed combinations, CHA showed the lowest emissions of both CH 4 and CO 2 .

Simulated gaseous emissions
Averaged annual emissions of GHGs and ammonia from different sources over the simulation period are shown in Table 5. There was less variation in N 2 O and CO 2 emissions from plants and soils between fields than for NH 3 , CH 4 and animal-derived CO 2 emissions, which relate to animal type (cattle or sheep), stocking density and duration in each field. For example, annual stocking density is 340 head⋅d ha − 1 for cattle and 412 head⋅d ha − 1 for sheep in Burrows but in Golden Rove annual stocking density is 224 head⋅d ha − 1 for cattle and 1500 head⋅d ha − 1 for sheep.

Nitrogen cycling
Averaged annual N inputs to and outputs from the individual fields over the simulation period are shown in Table 6. Total N input ranged from 190 to 260 kg ha − 1 . Between 37 and 60% of the N added to the fields was removed through harvested biomass (silage) or animal intake, and 15-26% of N was lost through surface runoff or lateral drainage. Annual averaged gaseous losses of N 2 O and NH 3 over the simulated period were 2.77 ± 0.24 and 8.95 ± 1.54 kg N ha − 1 , respectively.

Model performance of beef finishing cattle and sheep growth
For both beef cattle and sheep, individual animals differ in their growth rate and their health status naturally within any livestock enterprise. Growth rates will similarly vary between breeds and the change in meteorological conditions during each grazing season. Given this, when objectively assessing model performance on simulating liveweight of cattle and sheep, for different age ranges, breed combinations and grazing seasons (Tables 2 and 3 and Fig. 4), the extended SPACSYS model could accurately simulate the dynamics of animal liveweight within the natural variations expected. Relatively, liveweight simulations for cattle were shown to be more accurate than those for sheep, where in both instances, simulation accuracy weakened as animals aged, which possibly is the result of the estimation of the potential growth rate. Further, levels of accuracy differed more across sheep breeds than it did across cattle breeds. Grazing year could also influence simulation accuracy, although reasons for this are not entirely clear.
The extended SPACSYS model is capable of simulating not only animal growth but also other elements of livestock (either beef finishing cattle or sheep) production at a systems level. Therefore, the model has the potential to investigate the responses of the system on and consequences of a range of agronomic management and grazing strategiesi. e., not only those as analysed across the farmlet (small farm) of this research with its specific (single) management and (single) grazing approach.

Gaseous emissions from cattle and sheep
The simulated averaged CH 4 emission rate was between 242 and 289 g head − 1 d − 1 for beef cattle and between 22 and 31 g head − 1 d − 1 for sheep aged between three and seven months (Table 4). There are few measurement datasets available for UK grazing systems, but the simulated data are within the expected range according to those datasets that have been published and, more broadly, with the default values provided by the IPCC Guidelines for national GHG inventories (IPCC, 2019). Meo-Filho et al. (2021) reported average emission rates of 183-213 g head − 1 d − 1 for growing beef cattle grazing the same PP farmlet of the NWFP during late summer of 2019, measured using the SF 6 tracer gas technique (Berndt et al., 2014). Fraser et al. (2014), also using the SF 6 tracer gas technique, measured emissions from upland and lowland grazing beef cattle and reported emissions in the range 173-217 g head − 1 d − 1 . For sheep, using an emission chamber methodology and a cut and carry system for feeding fresh herbage, Moorby et al. (2015) measured emissions from mature ewes fed permanent pasture of 11-15 g head − 1 d − 1 ,  reported emission rates in the range 12-17 g head − 1 d − 1 for growing lambs. More generally, default emission rates provided by the IPCC (2019) equate to 142 and 25 g head − 1 d − 1 for finishing beef cattle and productive sheep in Western Europe, respectively. While there were significant differences between breed Table 4 Simulated average daily liveweight gain (kg d − 1 ) and methane and carbon dioxide emissions (g head − 1 d − 1 ) during the grazing period for beef and sheep (lamb) breed combinations (different letters in a column either for cattle or sheep indicate a significant difference among breed combination, p < 0.05).   combinations in the simulated emissions per head for both beef cattle and sheep (Table 4), the literature evidence is that breed is a far less important variable (generally non-significant) influencing CH 4 emission than other factors such as diet characteristics and feed DMI (Duthie et al., 2017;Fraser et al., 2014;Fraser et al., 2015;Moorby et al., 2015). Any differences in emissions per head between breeds are generally accounted for through differences in body size, productivity or feed intake and, therefore, on an emission intensity basis (CH 4 kg − 1 LWG) breed is considered relatively unimportant. There was less variation in respiration rate between different beef and sheep breed combinations (Table 4) suggesting that breed plays only a minor role and that body size is the major determinant of the respiration rate (data not shown). Although a direct comparison with measurement data is lacking, relative errors of less than 10% between the simulated and reported values for animals of the same size (Chaves et al., 2006;Gunter and Beck, 2018) support the model output. There have been few measurements, reported to date on CO 2 emissions from sheep. In an early study, Whitelaw et al. (1972) reported that an average of 232 g CO 2 -C head − 1 d − 1 was produced by sheep weighing 56-78 kg at 12 • C ambient temperature, which is slightly lower than we estimated for an animal with an average body size of 40 kg.

Nitrogen cycling
Averaged annual N input to the individual farmlet fields ranged from 190 to 260 kg ha − 1 , which mainly reflected variations in stocking density and duration across fields ( Table 6). The estimated output components in N balance are within the range of the reported values. For example, an annual average loss rate through surface runoff or lateral drainage of 42 (±5) kg N ha − 1 over the farmlet, which was close to the estimate from the NWFP in a previous study (Shepherd et al., 2019).
An annual average of 2.77 (±0.24) kg N ha − 1 as N 2 O over the simulated period was emitted to the atmosphere. Although as a proportion of the total input this is small and agronomically of little consequence, it is of environmental significance because of the high global warming potential of N 2 O. Sources for this emission include the atmospheric N deposition, the applied fertiliser N and farm-yard manure (FYM) N as well as the in-field recycled N being deposited as dung and urine by the animals (making the N content of the grazed herbage available to the soil microbial processes of nitrification and denitrification) and N from senescent above-and below-ground plant material. The simulated N 2 O emission was equivalent to 1.14 ± 0.05% of the N input for these sources. This estimate is a composite of the various N 2 O sources and therefore difficult to compare with emission factors reported elsewhere for individual N sources. It is in the range of 0.1-1.8% given as the default emission factor (EF 1 ) by IPCC (2019) for fertiliser and FYM N additions to the soil, and the range for the default IPCC emission factor (EF 3 ) of 0-1.4% for cattle excreta returns during grazing (IPCC, 2019). It is also of a similar order of magnitude to empirical data from recent UK studies. Cowan et al. (2020) reported an average value of 1.33% for synthetic N fertiliser (ammonium nitrate) based on 202 observations for grassland soils in the UK and Ireland. Thorman et al. (2020) reported an average emission for FYM applied to grassland of 0.37%, based on three experimental sites, with a value of 0.13% specific to the North Wyke site. Chadwick et al. (2018) analysed available UK data for N 2 O emissions from cattle dung and urine returns to soil, developing average emission factors of 0.69 and 0.19% for urine and dung, respectively, based on five sites and applications at three times of the year across the grazing season. There are large uncertainties in these estimated emission factors for agricultural soils because of many influencing environmental and management factors (Cowan et al., 2020).
Agriculture is the major source of NH 3 emissions to the atmosphere, primarily deriving from livestock excreta, including manure and urea / NH 3 -based fertiliser applications (Behera et al., 2013). In SPACSYS, NH 3 volatilisation from chemical fertilisers is not yet considered. Ammonium nitrate was applied in this study, which is associated with much lower NH 3 emissions than other fertiliser types, e.g. urea (Forrestal et al., 2016), typically of less than 5% of the applied fertiliser N (e.g. Misselbrook et al., 2004). We simulated NH 3 volatilisation from applied FYM and excretal grazing returns at an average annual value of 8.95 ± 1.54 kg N ha − 1 , equivalent to 12.9% of the FYM and excreta N. Ammonia emissions from applied FYM can be low, as the ammonium-N content of the FYM is typically low (Chambers, 2003), particularly for FYM that has previously been stored for some months, because of volatilisation losses and immobilisation processes during storage. Nicholson et al. (2017) quoted a mean emission for livestock FYM based on UK experiments of 4.5% of the total N applied while Misselbrook et al. (2005) reported a loss of 69% of the available N at spreading, which equated to approximately 8% of the total N applied. Emissions from excretal returns at grazing derive primarily from the urine (Laubach et al., 2013) and previous experiments in the UK and Netherlands give emissions typically in the range 5-10% of urine N deposited (Bussink, 1994;Jarvis et al., 1989;Jarvis et al., 1991;Lockyer, 1984), although Laubach et al. (2013) reported somewhat higher values (c. 25%) from trials in New Zealand. As with N 2 O emissions, NH 3 emissions can vary considerably according to application techniques, N forms, soil texture, soil wetness and weather conditions at the times of application to the field. However, the rate might be underestimated in the model and should be further investigated, including an implementation of the NH 3 volatilisation process from chemical fertilisers.
On average, the study farmlet annually received 208 kg N ha − 1 and took 150 kg N ha − 1 from the system (Table 6), which resulted in a surplus of 58 kg N ha − 1 . The imbalanced N budget suggested that the N application rate could be reduced to a certain extent or the livestock density might be increased to graze more forage during the grazing season. However, average simulated annual N uptake is 264 kg N ha − 1 (data not shown). Considering the contribution from soil N mineralisation, the N budget could be balanced. Although volatilisation from FYM application or animal excreta has been considered, the loss from chemical N fertiliser through the process has not been included. It was reported that NH 3 emissions represented 7 and 21% of the total applied N for ammonium-nitrate and urea, respectively, on grassland in the UK (Carswell et al., 2019a). Not including this loss in the model adds uncertainty to the N cycle.

Future development
As shown in this study, the extended SPACSYS model can dynamically simulate animal and grass growth, nutrient cycling and water redistribution in a soil profile considering the effects of animal genotype, climate, feed quality and quantity on livestock production, GHG emissions, water use and quality, and nutrient budgets at a field scale. It is novel to link animal, plant, soil and atmosphere together into a whole system model to quantitatively investigate the dynamics of animal and grass production and nutrient fate, and their interactions under varied environmental conditions. Through this study, the configuration for a permanent pasture grazing system has been validated. All PP farmlet fields were reasonably homogenous and dominated (>60%) by perennial ryegrass, with a smaller biomass of creeping bent (Agrostis stolonifera), Yorkshire fog (Holcus lanatus) and marsh foxtail (Alopecurus geniculatus) also contributing to the sward; legumes, on the other hand, comprised <1% of the overall composition (Takahashi et al., 2018). As more diverse, multi-species swards with higher proportions of legumes and forbs in intensive grasslands are becoming more common in practice, the modelling of these more diverse botanical composition swards needs to be validated as a subject of future work. Such modelling could also include the dynamics of individual species in the swards and their impact on animal growth and nutrient flows. Furthermore, animal intake preference and deliberate selection for specific species within the sward should be considered as it affects both the intake parameters and the subsequent development and composition of the sward, although this is challenging to implement in a process-based model. To date, no components have been implemented to simulate the impacts of extreme events such as temperature and rainfall, systematic animal-mediated nutrient transfers, pests, weeds and plant and animal genetic characteristics -environment interactions (GxE) on an agricultural ecosystem, which is highly desirable (Bryant et al., 2011). Evidence suggests that current guidance (Agricultural and Food Research Council, 1993) on nutritional requirements needs to be updated, where the ongoing research project (https://www.cielivestock.co.uk/improve-beef-feed-gu idelines/) may lead to revisions to the energy requirements of beef cattle.

Conclusions
In this study, the extended SPACSYS model was shown to accurately and dynamically model finishing beef cattle, lamb and grass growth, nutrient cycling and water redistribution in a soil profile considering the effects of genotype, climate, feed quality and quantity on livestock production, GHG emissions, water use and quality, and nutrient cycling in a permanent pasture grazing system consisting of seven fields. Averaged annual N input to the individual fields ranged from 190 to 260 kg ha − 1 , of which 37-60% removed from the fields in terms of biomass cut or animal intake, and 15-26% through surface runoff or lateral drainage and 1.14% emitted to the atmosphere as N 2 O. About 13% of the FYM and excreta N in the farmlet volatilised from the soil. There are significant differences in animal growth rate, CO 2 and CH 4 emissions between different sheep breeds. However, there are less differences between the cattle breeds. Although the extended model was validated with data specific to Southwest England and for a permanent pasture grazing system, the model has clear potential to explore more innovative practices to maintain / increase livestock production whilst reducing adverse environment impacts across different livestock breeds, climates and soil types.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.