Model overview

The e-Dairy model is energy based, dynamic (daily simulation over a 365-day period) and is written in the Visual Basic programming language within Microsoft Excel^®. Dairy farms with either ryegrass- or lucerne-based pastures, with or without summer or winter crops, and with any calving pattern can be simulated. An animal model, e-Cow (Baudracco et al., Reference Baudracco, Lopez-Villalobos, Holmes, Comeron, Macdonald, Barry and Friggens2012), was integrated into e-Dairy to simulate the performance of individual cows at grazing.

The main inputs of e-Dairy are farm area, use of land for either pasture and crops, type of pasture (ryegrass-based or lucerne-based), type of crops (winter or summer crop), monthly pasture growth rate (mean and s.d.; used to predict daily herbage mass (HM) for each paddock), annual crop yield (mean and s.d.), supplements offered, quality of feeds including NDF and metabolisable energy (ME), herd description including herd size, age structure, calving pattern, BCS at calving, probabilities of pregnancy, average genetic merit (potential yields of milk, fat and protein) and economic data for items of income, costs and assets. Stocking rate (cows/ha) is indirectly an input that is the result of the number of cows and the farm area, both set as inputs. The model allows to set management policies such as calving pattern, dry-off policy, target pre- and post-grazing HM, a policy to make hay or silage from pasture and supplementation. The main outputs are herbage and total DM intake, annual pasture utilisation, yields of milk, fat and protein, changes in BCS and LW, economic farm profit and return on assets. A global schematic representation of the model is shown in Figure 1.

Figure 1

Schematic representation of the e-Dairy model. *Can be set to behave stochastically. ^aHerd structure: number of cows per age category (year) and per lactation number. ^bFertility: probability of pregnancy at each service, calving system (seasonal or all-year-round) and number of services. ^cGenetic merit (potential yields of milk, fat and protein, parameters defining milk fat and milk protein content, target BCS parameters and live weight). ^dInitial BCS: mean and s.d. ^ePastures data: number and size of paddocks, pasture type, initial herbage mass and herbage growth rates (mean and s.d.). ^fCrops data: number and size of paddocks, type of crop (summer or winter), amount of DM produced/ha per year, the fate of the crop (grazing, silage or hay) and efficiency of harvesting or grazing. ^gFeed quality: NDF, metabolisable energy (ME) of herbage and supplements, and rate of herbage quality decrease/increase if not used by optimum date. ^hFeed expenses: concentrates price can be set to behave stochastically. ⁱCows randomisation: to allocate a unique genetic merit to each cow, on the basis of average genetic merit and (co)variance matrix. ^jHerd dynamic: number of lactating and dry cows (daily). ^kPaddocks randomisation: to define initial herbage mass (HM) on each paddock. ^lPasture budget sub-routine that matches in advance the amount of pasture available with the amount of pasture required. ^mPaddocks dynamic: daily HM and actual growth rate, grazing dates and grazing time per paddock and actual post-grazing mass. ⁿPasture allowance: actual kg DM offered per dry and lactating cow/day. ^oReserves stock: kg DM reserved as hay and silage. ^pe-Cow model: see schematic representation, inputs and equations in the e-Cow model (Baudracco et al., Reference Baudracco, Lopez-Villalobos, Holmes, Comeron, Macdonald, Barry and Friggens2012) explaining dry matter (DM) and energy intake and energy partitioning within the cow. ^qPaddocks outputs: daily actual pre- and post-grazing herbage mass, daily actual growth rates.

Cows

Each cow has unique values for the parameters defining its genetic merit, and the performance of each cow is individually simulated on a daily basis in e-Dairy, using the animal model e-Cow described in Baudracco et al. (Reference Baudracco, Lopez-Villalobos, Holmes, Comeron, Macdonald, Barry and Friggens2012). The e-Cow model is an energy-based model that predicts intakes of DM and energy, yields of milk, fat and protein and changes in BCS and LW. The main features of the e-Cow model are the combination of physical, metabolic and ingestive constraints in the prediction of herbage DM intake, the homoeostatic and homeorhetic control of body lipid change and its ability to predict performance of cows of different genetic merit.

There is now considerable evidence that genotype affects nutrient partitioning throughout differences in gene expression and enzyme profiles (Friggens et al., Reference Baudracco, Lopez-Villalobos, Holmes, Comeron, Macdonald, Barry and Friggens2012). The genetic merit of the cow, in the e-Cow model, is defined by her potential yields of milk, fat and protein (kg/cow for a 305-day lactation period), LW at calving and genetic targets for BCS at conception and at next calving. The genetic potentials and genetic targets are those achieved by the cow when no feeding restrictions are imposed; however, predictions are later adjusted by the nutritional status of the cow in our model (Baudracco et al., Reference Baudracco, Lopez-Villalobos, Holmes, Comeron, Macdonald, Barry and Friggens2012). This approach of homeorhetic-driven variables and homoeostatic corrections was already used by Martin and Sauvant (Reference Martin and Sauvant2010).

Potential milk yield is calculated in e-Dairy using a mathematical mammary gland model (Vetharaniam et al., Reference Vetharaniam, Davis, Upsdell, Kolver and Pleasants2003), which is based on the interaction of two pools of alveoli (i.e. groups of secretory cells): one active pool and one non-active pool. Further details are given in the description of the e-Cow model (Baudracco et al., Reference Baudracco, Lopez-Villalobos, Holmes, Comeron, Macdonald, Barry and Friggens2012). The default potential milk yields for mature cows (in a 305-day lactation period) are 11 247 and 8011 kg milk per cow, with 836 and 679 kg milk solids (MS; fat plus protein) per cow for North American (NA) and New Zealand (NZ) Holstein–Friesian (HF) strains, respectively. These are default values internally stored in the model. However, cows with different potential milk yield can be simulated by defining other potential yield of milk as inputs. When the input potential milk yield differs from the internally stored values cited above, an iterative procedure is used to find the value of parameters of the mammary gland model that produces a lactation curve with the new potential milk yield (kg milk for 305 days).

The potential yields of fat and protein (kg/cow in 305-day lactation) are the result of the sum of the product between daily potential milk yield and daily milk fat or milk protein concentration. Daily milk fat and milk protein concentrations are calculated using lactation curves produced with the Wilmink (Reference Wilmink1987) function:

$${{y}_t}\,{\rm{ = }}\,a\,{\rm{ + }}\,b{{e}^{{\rm{ - 0}}{\rm{.05}}t}} \,{\rm{ + }}\,ct \eqno(1)$$

where y_t represents the percentages of milk fat or milk protein at day t of lactation and e is the base of natural logarithm (2.718281828), whereas a, b and c are estimated parameters that define the scale and shape of the curve. The values for parameters a, b and c were obtained from a study by Roche et al. (Reference Roche, Berry and Kolver2006), for both NZ and NA HF cows.

The BCS and LW are modelled according to the genetically driven body lipid change model proposed by Friggens et al. (Reference Friggens, Ingvartsen and Emmans2004), with the main concepts of this model further expanded in a study by Friggens and Newbold (Reference Friggens and Newbold2007). This body lipid change model proposes that mobilisation and gain of body reserves are genetically driven to achieve two genetic targets of body fatness: one, at or around conception, and another at the next calving. In the e-Cow model, the genetically driven changes in LW and BCS (i.e. no nutritional constraints) are initially calculated and then adjusted on a daily basis according to the nutritional status of the cow. All equations and parameters used in the e-Cow model are described in detail by Baudracco et al. (Reference Baudracco, Lopez-Villalobos, Holmes, Comeron, Macdonald, Barry and Friggens2012). In the e-Cow model, the genetic merit of the cow is defined by the potential yields of milk, fat and protein, parameters a, b and c of the Wilmink function for milk fat and milk protein concentrations, LW and parameters defining target BCS at conception and at next calving.

Herd

Random generation of cows with correlated variables

Each cow is generated randomly through the following correlated traits: potential yield of milk, LW at calving and parameters of the Wilmink function defining milk fat and milk protein concentration curves (equation (1)). In order to simulate randomly a herd in which each cow has correlated values for all variables, the following matrix operation is performed:

$$\bf H\,{\rm{ = }}\,m'\,{\rm{ + }}\,L\,{\rm{\times }}\,Z$$

where H is the matrix of the simulated herd, with each trait in a column and each cow in a row, m is the vector with the herd mean values for each trait (inputs), L is the lower triangular matrix obtained by Cholesky decomposition of the phenotypic (co)variance matrix (A) between the traits and Z is a matrix with random values derived from a normally distributed function with a mean of 0 and a s.d. of 1. A similar procedure has been used by Gartner (Reference Gartner1981).

Elements of the phenotypic (co)variance matrix A were estimated using phenotypic records from a trial comparing cows of NA and NZ HF strains (Macdonald et al., Reference Macdonald, Verkerk, Thorrold, Pryce, Penno, McNaughton, Burton, Lancaster, Williamson and Holmes2008b). The product between the vector Z that has random values (mean of 0 and a s.d. of 1) and the matrix L (containing information about traits variance and correlation between traits) creates the correlated variance for each cow on each trait. The m vector contains all the traits described as columns of the H matrix and represents the mean herd values for each trait. Thus, summing up m to the product L × Z, a unique value is allocated to each cow on each trait. The average values of the m vector are set as default for each strain (see e-Cow model, Baudracco et al., Reference Baudracco, Lopez-Villalobos, Holmes, Comeron, Macdonald, Barry and Friggens2012), but new values can be set as inputs. Each cow is further allocated an age (see the ‘age structure’ section) and a BCS at calving. The BCS at calving is defined as input (Figure 1), with a mean and s.d. for a normal distribution.

Paddocks

The number and size of paddocks with pasture and crops are inputs (Figure 1). In pasture paddocks, the net herbage accumulation rate depends on a general net herbage accumulation rate curve (input), but is affected by the particular HM of each paddock, which is altered by the events of grazing and its intensity (see the ‘HM, herbage accumulation and grazing dates’ section). This results in a unique growth rate on each paddock.

In paddocks with crops, the type of crop (summer or winter), the amount of DM produced/ha per year (mean and s.d.), the fate of the crop (grazing, silage or hay) and the grazing/harvesting efficiency need to be set as inputs. If the crop is grazed, the grazing efficiency is constant (input) for all crop paddocks and for all the grazing events, in contrast to what happens in paddocks with pasture, where grazing efficiency is predicted for every paddock and every grazing event (see the ‘simulation of grazing’ section). The response to fertiliser is not simulated, but could be indirectly accounted for by changing the input amount of DM produced/ha per year.

HM, herbage accumulation and grazing dates

Daily net herbage growth rates are required as an input (kg DM/day, average of each month; Figure 1). Herbage growth can behave either stochastically or deterministically – it is optional in the model. For deterministic simulations, the growth rates used are those set as inputs. For stochastic simulation, the input herbage growth rates are used as the mean of a normal distribution, and an s.d. is required as input to randomise the herbage growth rates. Thus, the effect of climate on herbage growth rate can be accounted for, by using s.d. obtained from experimental data sets.

A ‘target pre-grazing HM’ (Figure 2), needs to be defined as input, expressed as kg DM for ryegrass-based pastures and as accumulated ‘growing degree days’ for lucerne pastures, which is an indicator of morphological development (Sanderson et al., Reference Sanderson, Karnezos and Matches1994). Daily minimum and maximum temperatures (inputs) are needed to calculate ‘growing degree days’. The date when each paddock reaches the ‘target pre-grazing HM’ is named ‘optimum grazing date’. However, the actual date of grazing, named ‘next grazing date’, may differ from the calculated ‘optimum grazing date’ (Figure 2) according to the use of paddocks. The calculated difference (days) between the ‘optimum grazing date’ and the ‘next grazing date’ is called ‘days away from optimum’, as shown in Figure 2. HM in each paddock will accumulate daily, from ‘previous grazing date’ until ‘next grazing date’, as shown in the timeline of Figure 2. The herbage growth rate set as input is affected by the grazing intensity, using the following equation (Garcia, Reference Garcia2000):

$$ {\rm{Growth}}\,{\rm{rate}} = \left( {{\rm{GR}}\times \frac{{\rm{4}}}{{{\rm{H}}{{{\rm{M}}}_{{\rm{max}}}}}}} \right)\times {\rm{H}}{{{\rm{M}}}_{{\rm{(}}t{\rm{)}}}}\left( {{\rm{1 - }}\frac{{{\rm{H}}{{{\rm{M}}}_{{\rm{(}}t{\rm{)}}}}}}{{{\rm{H}}{{{\rm{M}}}_{{\rm{max}}}}}}} \right) $$

where GR is the input herbage growth rate (kg DM/day, input), HM_{(t )} is the HM of the paddock in the day of simulation, calculated as HM of the previous day (HM_{(t − 1)}) plus GR, HM_max (input) is the HM at which net accumulation rate approaches zero because senescence rate approaches gross growth rate (Bircham and Hodgson, Reference Bircham and Hodgson1983).

Figure 2

Schematic representation of the simulation of grazing dates and herbage mass (HM) for each paddock within the e-Dairy model.

Simulation of grazing

Each paddock is individually grazed. The duration of each grazing event in each paddock ‘grazing time’, is calculated with the following equation:

$$ {\rm{Grazing}}\,{\rm{time}}\,{\rm{(days)}}\, = \,\frac{{{\rm{Pre - grazing}}\,{\rm{HM}}\,\times \,{\rm{paddock}}\,{\rm{area}}}}{{{\rm{Total}}\,{\rm{allowance}}\,{\rm{per}}\,{\rm{day}}}} $$

where ‘total allowance per day’ is calculated as the summed allowance for all cows.

The post-grazing HM is calculated as the difference between ‘pre-grazing HM’ and ‘herbage removed at grazing’, the latter being calculated as the sum of daily herbage DM intake of all cows. Daily herbage DM intake is not assumed but predicted for each cow, as explained in the e-Cow model (Baudracco et al., Reference Baudracco, Lopez-Villalobos, Holmes, Comeron, Macdonald, Barry and Friggens2012). Pasture utilisation (%) is calculated as the kg DM/ha consumed at grazing, divided by the kg DM/ha annual pasture accumulation and multiplied by 100.

Management policies

Pasture budget

The date of grazing and the grazing time (per grazing event) in each paddock are calculated while the model runs; however, it is still possible to define a policy to use paddocks before the simulation starts. The automatic procedure implemented in the model to simulate paddocks dynamics is performed by a sub-routine named ‘pasture budget’ (Figure 3), which runs for a period defined with an input named ‘pasture budget frequency’. The pasture budget in e-Dairy consists of a calculation of future pasture offer (expected pasture offer) and future pasture demand (planned pasture allowance) for the whole farm. The ‘expected pasture offer’ is calculated as the sum of the product between pre-grazing HM in each paddock at the calculated ‘optimum grazing date’ and the area of each paddock (Figure 3). The ‘planned pasture allowance’ is calculated as the target pasture allowance per cow times the number of cows times the number of days in the period to be budgeted.

Figure 3

Schematic representation of the pasture budget performed in the e-Dairy model.

Economics

Two different payment systems can be used, namely: price per litre of milk and multiple component price system (kilograms of fat × A + kilograms of protein × B − litres of milk × C), where A, B and C are the values per kilogram of fat and protein and litre of milk, respectively, as it is used in NZ (Marshall, Reference Marshall1989) and Ireland. The following groups of inputs are required: farm incomes, animal expenses, feed expenses, labour expenses, adjustments, overheads (includes depreciation) and assets. Items included in each section are based on DairyNZ (2009). This allows the calculation of two main outputs: operating profit [farm incomes − (animal expenses + feed expenses + labour expenses + adjustments + overheads)] and return on assets (economic farm profit ÷ assets × 100).

Simulation

Deterministic or stochastic simulations can be carried out with e-Dairy. The stochastic simulation allows multiple runs, which could represent either different farms in 1 year or different conditions for the same farm. A group of variables can be allowed to behave stochastically. Some variables are stochastically triggered at the beginning of the simulation (and stay constant), such as genetic merit of cows, initial BCS, calving dates, cows age, crops yield, milk price and supplement price and pasture growth rates, whereas other variables behave stochastically during the simulation, such as daily herbage allowance and probabilities of pregnancy. Probability distribution functions, and their respective parameters, need to be set as inputs for the stochastic simulation. The normal, flat and gamma functions are available in the model.

Outputs

The main system outputs are herbage DM intake per cow and per hectare, annual pasture utilisation, milk yield per cow and per hectare, annual changes in BCS and LW, economic farm profit and return on assets. However, daily values for each variable, for individual cows or individual paddocks are available after the simulation. Thus, outputs such as daily HM or post-grazing HM per paddock, daily milk yield, BCS and LW per cow are stored after the simulation.

Model validation

Two independent data sets resulting from stocking rate experiments, and not used in the development of the e-Dairy model, were used to validate the model. One data set was obtained from a farmlet trial comparing five levels of stocking rate (from 2.2 to 4.3 cows/ha) with HF cows, winter–spring calving, grazing on ryegrass-based pastures (offered 0.19 t DM supplement per cow/year) during 3 years in NZ (Macdonald et al., Reference Macdonald, Penno, Lancaster and Roche2008a). The second data set was obtained from a farmlet trial comparing three levels of stocking rate (from 1.6 to 2.6 cows/ha) with crossbred HF-Jersey cows, winter–spring calving, grazing on lucerne-based pastures (offered 1.8 t DM supplements per cow/year) during 2 years in Argentina (Baudracco et al., Reference Baudracco, Lopez-Villalobos, Romero, Scandolo, Maciel, Comeron, Holmes and Barry2011).

Measured inputs used to validate the model were stocking rate, monthly pasture growth rates, monthly amounts of supplements used per cow, monthly herbage allowances, monthly ME and NDF of pastures and supplements, weekly calving pattern, monthly pre- and post-grazing HM, farmlet averages for age structure, lactation length, BCS and LW at calving.

Annual outputs for actual and simulated farmlets were compared, with 21 points for validation (five farmlets during 3 years plus three farmlets during 2 years). Outputs compared were: yields of milk and MS (kg/cow per year), BCS and LW at day 365-day of simulation and annual pasture utilisation. The concordance correlation coefficient (CCC; Lin, Reference Lin1989) and the relative prediction error (RPE) (Fuentes-Pila et al., Reference Fuentes-Pila, Ibanez, De Miguel and Beede2003) were used to evaluate the extent of agreement between actual and predicted values. The Landis and Koch (Reference Landis and Koch1977) scale has been used here to describe the degree of concordance, with: 0.21 to 0.40 being ‘Fair’; 0.41 to 0.60 being ‘Moderate’; 0.61 to 0.80 being ‘Substantial’; and 0.81 to 1.00 being ‘Almost perfect’. The CCC reflects both precision, that is, the degree to which the predicted against actual values cluster about the regression line, and accuracy, that is, degree to which the regression line adheres to the 45° line through the origin.

The RPE is defined as the positive square root of the mean square prediction error (equation (5)), the latter expressed as a percentage of the mean of actual values (μA) (Fuentes-Pila et al., Reference Fuentes-Pila, Ibanez, De Miguel and Beede2003):

$$ MSPE\, = \,\frac{{\rm{1}}}{n}\mathop{\sum}\limits_{i\, = \,{\rm{1}}}^n {{{{({{A}_i}{\rm{ - }}{{P}_i})}}^{\rm{2}}} } $$

where P represents the predicted values and A represents the actual observed values for either milk yield, MS yield, pasture utilisation, BCS or LW change. The RPE is used to decide whether the overall level of accuracy could be considered acceptable for practical use. Fuentes-Pila et al. (Reference Fuentes Pila, DeLorenzo, Beede, Staples and Holter1996) suggested that an RPE value <10% is an indication of satisfactory prediction, whereas an RPE between 10% and 20% indicates a relatively acceptable prediction, and an RPE >20% indicates poor prediction.