Effect of Nutritional Variation and LCA Methodology on the Carbon Footprint of Milk Production From Holstein Friesian Dairy Cows

The UK livestock industry urgently needs to reduce greenhouse gas (GHG) emissions to contribute to ambitious climate change policy commitments. Achieving this requires an improved understanding of emission sources across a range of production systems to lower the burden associated with livestock products. Life cycle assessment (LCA) methods are used in this study to model milk production from two genetic merits of Holstein Friesian cows managed in two novel and two conventional UK dairy systems. Select merit cows sired by bulls with high predicted transmission for fat plus protein yield are compared with Control merit animals sired from UK average merit bulls. Cows were managed in conventional housed and grazed dairy systems with novel Byproduct and Homegrown feeding regimes. A LCA was used to quantify the effect of allocation and management of feed components on the carbon footprint of milk production. Natural variation in nutritional quality of dairy system rations was investigated to quantify uncertainty in the carbon footprint results. Novel production system data are used to assess the effect of introducing home grown legumes and co-product feeds. Control merit footprints across each of the management regimes were significantly higher (p<0.001) in comparison with a high production Select merit, on average by 15%. Livestock emissions (enteric, manure management and deposition) and embedded emissions (purchased feeds, fertiliser, and pesticides) were also significantly higher from control merit cows (p<0.01). Mass and economic allocation methods, and land use functional units, resulted in differences in performance ranking of the dairy systems, with larger footprints resulting from mass allocation. Pairwise comparisons showed GHG's from the systems to be significantly different in total and source category emissions, with significant differences in mean embedded emissions found between most management systems (p<0.05). Monte Carlo simulated system footprints considering the effect of variation in feed digestibility and crude protein also differed significantly from system footprints using standard methods (p < 0.001). Dairy system carbon footprint results should be expressed using multiple units and where possible calculations should incorporate variation in diet digestibility and crude protein content.

The UK livestock industry urgently needs to reduce greenhouse gas (GHG) emissions to contribute to ambitious climate change policy commitments. Achieving this requires an improved understanding of emission sources across a range of production systems to lower the burden associated with livestock products. Life cycle assessment (LCA) methods are used in this study to model milk production from two genetic merits of Holstein Friesian cows managed in two novel and two conventional UK dairy systems. Select merit cows sired by bulls with high predicted transmission for fat plus protein yield are compared with Control merit animals sired from UK average merit bulls. Cows were managed in conventional housed and grazed dairy systems with novel Byproduct and Homegrown feeding regimes. A LCA was used to quantify the effect of allocation and management of feed components on the carbon footprint of milk production. Natural variation in nutritional quality of dairy system rations was investigated to quantify uncertainty in the carbon footprint results. Novel production system data are used to assess the effect of introducing home grown legumes and co-product feeds. Control merit footprints across each of the management regimes were significantly higher (p<0.001) in comparison with a high production Select merit, on average by 15%. Livestock emissions (enteric, manure management and deposition) and embedded emissions (purchased feeds, fertiliser, and pesticides) were also significantly higher from control merit cows (p<0.01). Mass and economic allocation methods, and land use functional units, resulted in differences in performance ranking of the dairy systems, with larger footprints resulting from mass allocation. Pairwise comparisons showed GHG's from the systems to be significantly different in total and source category emissions, with significant differences in mean embedded emissions found between most management systems (p<0.05). Monte Carlo simulated system footprints considering the effect of variation in feed digestibility and crude protein also differed significantly from system footprints using standard methods (p < 0.001). Dairy system carbon footprint results should be expressed using multiple units and where possible calculations should incorporate variation in diet digestibility and crude protein content.

INTRODUCTION
The UK government is committed to reducing GHG emissions to net zero by 2050 with an even more ambitious target date of 2045 in Scotland (Committe on Climate Change, 2019). Agriculture is estimated to be responsible for 10-12% of global anthropogenic greenhouse gas (GHG) emissions (Smith et al., 2014). Emissions stemming from milk production in developed dairy regions such as the UK are estimated at between 1.2 and 1.4 kg CO 2 e/kg, respectively, which is lower than the global average of 2.5 kg CO 2 e/kg of fat and protein corrected milk (FPCM) (FAO, 2019). Dairy products have been processed and consumed in Britain since the Neolithic era and grassland, including rough grazing, covers ∼80% of the land area in Scotland (Charlton et al., 2019). However, agriculture is now the second largest source of GHG emissions in Scotland and there is an urgent need for this sector to contribute to national GHG emission reductions (SG, 2018). GHG emissions from livestock need to be reduced at a time when global demand for these commodities is increasing (Opio et al., 2013;FAO, 2019). UK GHG emissions from agriculture have declined by ∼14% since 1990, and reductions have largely arisen from a change to the Common Agriculture Policy (CAP), which ended a link between subsidy amounts and animal numbers (DBEIS, 2019, AHDB, 2014. Fewer total livestock numbers have led to lower stocking densities, less manure, and thus lower emissions (Rotz, 2004;del Prado et al., 2010;DBEIS, 2015). Formulating policies to enable further emission reductions on dairy farms will require an understanding of mitigation measures appropriate for specific production systems. Models used to quantify GHG's are important tools to aid the understanding of mitigation pathways that lie within the intricate footprints of livestock systems (Opio et al., 2013). Estimates of GHG emissions from livestock systems contain uncertainties from model boundaries and allocation, variation in input values, or epistemic uncertainty arising from modelled biological processes, all of which present challenges for researchers and decision makers (IPCC, 2006;Flysjö et al., 2011;Opio et al., 2013;Röös and Nylinder, 2013). Epistemic uncertainty analyses of modelled dairy and beef livestock systems have shown that, overall, nitrous oxide (N 2 O) and methane (CH 4 ) emissions from manure, fertiliser N input, and enteric CH 4 contribute most to variability and can inflate GHG emissions in livestock footprints (Ross et al., 2014;Zehetmeier et al., 2014;Sykes et al., 2019). Variation of N 2 O emissions in dairy system footprints was found to mainly stem from the IPCC emission factor for volatilisation and atmospheric deposition of N (Ross et al., 2014). Sensitivity analysis provides a deeper technical understanding of complex systems and is recommended to clarify potential impacts, however Baldini et al. (2017) show that only 20% of 44 reviewed LCA studies of milk production published between 2009 and 2015 carried out sensitivity analysis.
Gradual improvements in methodology have allowed more precise estimates of emissions arising from agricultural systems. However, simultaneously annual global GHG emissions have continued to increase and the interval available to implement reduction strategies narrows (Boden et al., 2017). Numerous examples can be found in the literature comparing carbon footprints arising from various dairy production methods in Scotland (Ross et al., 2014) and across the world (Cederberg and Mattsson, 2000;Basset-mens et al., 2005;Flysjö et al., 2011;O'Brien et al., 2014) including suggestions for establishing a system emitting less CO 2 per unit product or management type (O'Brien et al., 2012). However, differences in LCA methodology, allocation methods (for milk and meat) or functional unit hinder comparability (Baldini et al., 2017) and a meta-analysis of 30 published LCA's with 87 carbon footprints from pasture, mixed, and confinement dairy systems found no average footprint differences per kg of FPCM (Lorenz et al., 2019). Choice of functional unit is also important when expressing results from LCA studies assessing environmental impacts arising from differing dairy systems and when considering effects of intensification (Ross et al., 2017, Salou et al., 2017. As far as the authors are aware, allocation methods assessing feed components and their effect on dairy system carbon footprints are not available in literature. Studies quantifying uncertainty and assessing sensitivity of milk production LCA's have investigated management changes, C sequestration (O'Brien et al., 2012), manure storage (O'Brien et al., 2012;Battini et al., 2014) and changes in energy consumption (Roer et al., 2013). Methodologies to model CH 4 and N 2 O emissions from manure, and enteric CH 4, require measurements of diet digestibility and crude protein (CP) (Dong et al., 2006). The consequences of dietary and other variabilities should be considered, and uncertainties communicated when quantifying dairy farm carbon footprints (Zehetmeier et al., 2014;Milne et al., 2015). Of the nine LCA studies reviewed by Baldini et al. (2017) that incorporated sensitivity analysis an assessment of variability of crude protein and digestibility was not reported, this information is not available in literature for multiple dairy systems as far as the authors are aware.
The digestibility of a dairy cow diet relates to the chemical composition of feed components and also the ration as a whole. Digestibility can be described as the fraction of a food that is absorbed, and this is affected by fibre content of feeds, of which the forage components tend to exhibit wider variation (McDonald et al., 2011). Predictions of enteric CH 4 emissions are lower from diets with high digestibility (Röös and Nylinder, 2013) and diet digestibility has been shown to influence uncertainty in beef production footprints (Sykes et al., 2019). Rations containing optimum digestibility and balanced CP can lead to lower GHG emissions because a cow would require less feed to meet nutritional requirements. Conversely, too much CP leads to higher N excretions, which can cause nutrient surpluses that contribute to air and water pollution, as well as climate change. Edouard et al. (2019) showed that IPCC GHG estimates were not as accurate for high levels of dietary N because of increased NH 3 emissions. When compared with soybean meal, the addition of legumes such as faba beans and peas in a ruminant ration were shown to have higher digestibility, CP and energy, which can be beneficial for dairy cow nutrition (Volpelli et al., 2012), as long as rations are balanced to ensure higher levels of CP do not reduce nitrogen use efficiency (NUE).
Legumes are found in a wide range of ecosystems and the majority are genetically distinct form other plant species due to a symbiotic relationship with rhizobia. These are soil bacteria located within root nodules with the ability to fix nitrogen from the atmosphere (Kenicer, 2005). Within crop rotations, legumes can displace the need for imported nitrogen fertilisers, as well as nurture and condition the soil, which can have positive environmental and resource security consequences along with disease suppressing qualities (Lüscher et al., 2014;Stagnari et al., 2017). Home-grown protein feeds for animal production are increasingly being encouraged in the EU to reduce the protein deficit that relies upon soya imports which can fluctuate in price and can be associated with rainforest loss (European Parliament, 2018;Taherzadeh and Caro, 2019). Introducing legumes such as spring beans into crop rotations has the potential to reduce emissions through displacement of fertilisers, which in Scotland would translate to 100 to 180 kg/ha less N per year for spring and winter cereals, respectively (Iannetta et al., 2019). Legumes are estimated to generate <20% of the emissions associated with synthetic fertilisers, however N 2 O emissions can still occur from leguminous crop residues (Senbayram et al., 2015;Stagnari et al., 2017). Increasing the ratio of corn to alfalfa silage on large and small dairy farms in northern US has been shown to raise farm gate GHG emissions per kg of FPCM (Kim et al., 2019). An increase in the use of forage legumes within dairy systems should therefore be considered to improve outcomes for livestock and the wider environment.
This study adds to literature on Scottish dairy farm carbon footprinting carried out by Ross et al. (2014) by presenting novel Homegrown and By-product feeding systems in comparison with more conventional management techniques. This paper seeks to clarify the impact of dairy systems on the environment using a modelling approach to address specific questions; (i) what effect does the method of allocation of emissions from animal feeds in dairy systems have on the carbon footprint of milk produced, (ii) what effect does alternative system inputs, such as legumes and co-products have on the composition of carbon footprints and (iii) what is the effect of variation in feed digestibility and CP content on the global warming potential of milk produced in dairy systems under a range of management scenarios. Carbon footprints from Holstein Friesian cows managed in novel and conventional UK dairy feeding regimes are presented using life cycle assessment (LCA) methods. Monte Carlo simulations were applied to describe uncertainty brought about by variation in nutritional quality of the diets. Carbon footprint results using mass and economic allocation of feed components, land use functional units and considering variation in CP and digestibility were compared by ranking performance. Knowledge surrounding effects of UK farm grown forage legumes on dairy carbon footprints and insight into the influence of variation of feed digestibility and CP content in LCA uncertainty are novel aspects of this study. Additional impact categories are not presented in this manuscript, in order to focus on GHG emissions and climate change.

Dairy Systems Description
Data in this study originates from the Langhill herd of Holstein Friesian cows which form one of the worlds' longest running genetic line × feeding systems experiments (Pollott and Coffey, 2008). Data were used from all cows belonging to the herds based at the SRUC's Crichton Royal Farm, Dumfries, Scotland between 2006 to 2010 and 2012 to 2015. A Select (S) group of cows sired by bulls with high predicted transmitting abilities (PTA) for fat plus protein yield are compared with a Control (C) group sired from UK average merit bulls (Pryce et al., 2001). System experiments were managed according to the same rules and each regime was designed to allow animals to express their potential for milk production, within the limitations based on the feed rations offered. All experimental cows were milked three times per day and if not grazing were housed in the same building, with cubicles and concrete passageways that were cleaned with automatic scrapers. A complete diet was offered as a total mixed ration (TMR), irrespective of milk yield and stage of lactation. The four dietary treatments compared in this analysis were, i. a high forage (HF) composite system which can be defined as a conventional regime; grazing cows when availability of grass is adequate and housing during inclement winter months when animals are fed conserved forage and concentrate through a total mixed ration (TMR), ii. a novel home grown (HG) partially housed system, defined here as a regime where all feed is grown on farm using legume-based protein sources with no purchased feeds except minerals, and where animals are housed for one period each day and fed a conserved forage TMR, iii. a low forage (LF) conventional housed system with animals confined all year round and being fed a diet of conserved forage and concentrate through a TMR, iv. a novel by-product (BP) fully housed system that required no on-farm land by feeding mainly non-human edible coproducts from the food industry with no forages except straw.  March et al. (2017). Youngstock were managed as one group, with rations attributed by age, for heifer calves 0-12 months, and from 12 to 24 months.

Data Collection
Each sub-group was treated as a whole farm in the life cycle inventory. Four dietary treatments and two genetic merits of cow allowed a comparison of eight diverse dairy production systems summarised annually. Milking and dry cow populations were evaluated for each sub-group and replacement animal populations were categorised by age for calves and heifers as 0-12, 12-24, >24 months. Milk yield was measured for individual animals at each milking with a sample taken once per week and analysed for fat and protein content by infrared spectroscopy. Samples were analysed using a Milkoscan minor spectrophotometer (Foss Ltd., Denmark) and calibrated by methods of AOAC (2000). Outputs of milk were summed for the systems annually and weekly fat and protein constituents were averaged. Liveweights were recorded three times per day after milking and weekly liveweights were recorded for dry cows and replacements.
In each dietary regime Select and Control merits were housed together on one side of the building in two management groups which rotated every 3 days either being fed using automatic gates or behind a strap. A formulated TMR was offered daily, and average rations for each diet are presented in Table 1. Individual feed intakes for milking cows were measured using HOKO gates (Insentec BV, Marknesse, The Netherlands). LF and HF cows were fed 0.75 kg/day fresh weight of a standard 16% CP complementary parlour concentrate whilst milking. Dairy system inputs and outputs were determined annually using data extracted directly from SRUC's Langhill herd database and averages of production indicators are presented in Table 2 for Select and Table 3 for Control systems. Intakes of grass were not measured, however, periods of time spent grazing were recorded and samples of fresh grass were taken and analysed. Dry cows consumed a straw based diet and were then fed a transition diet 3 weeks before calving which consisted of 33% of the average daily dry matter intake of the appropriate milking cow ration plus 5 kg straw. Housed youngstock were managed as a group and offered a ration that included straw, distillers grains, and molasses.
Forage components of the LF, HF and HG system diets were grown on the farm. Maize, wheat, and spring beans were sown annually, lucerne every 2 years and grassland for pasture and silage every 5 years. Up to three cuts of ryegrass silage were harvested each year and any instances of double cropping of fields were noted with the lengths of time attributed to each crop allocated accordingly. For example, a field to be sown for maize may have been cut for silage before ploughing and there were instances where a grass silage cut was taken from a field sown for red clover bales. Harvested crops were stored on farm in covered clamps. Land required annually for those systems consuming crops grown on the farm was determined from amounts and DM's of each crop component fed to the herds and the crop DM and yield at harvest. Dry matter losses occurred at harvest, during ensiling or baling with estimated losses from grass silage, wheat alkalage, red clover bales and maize silage applied when considering land requirements for sub-groups because crops were not grown or ensiled separately for each of the dairy systems. An estimate of surface and effluent DM losses in the ensiling clamp were added to an estimate of wilting, leaching, respiration and mechanical field losses of crop to determine total land requirement by crop (Bastiman and Altman, 1985;Xiccato et al., 1994). Total land for each system year was calculated by  adding on-farm land required to an estimate of off-farm land. Off-farm land was approximated using economic allocation of feed components within each of the diets, national data for crops SAC Consulting (2016) and Feedprint (Vellinga et al., 2013) for processed feeds. Table 4 presents estimated feed component land use requirements and GHG emissions associated with mass and economic allocation methods for purchased feed inputs. TMR's were sampled monthly and wet chemistry analysis provided measurements of metabolisable energy (ME) content, dry matter (DM), digestibility and crude protein (CP) content of the ration (AOAC, 2000). Within the BP system only non-human edible purchased concentrates and straw were consumed ( Table 1). Leguminous by-products, soya bean meal and soya hulls were included in the LF and BP housed system TMR's at proportions of 14 and 9%, respectively. Legumes grown on the farm for the HG system represented 25% of the winter TMR and there were no legumes or leguminous by-product components fed within the HF regime ( Table 1). Applications of nitrogen (N), phosphorus (P), and potassium (K) fertilisers and organic fertiliser to crops grown on the farm were determined by the farm manager using a long-term nutrient management plan. Fertiliser use data recorded for each crop area with application rate and fertiliser type was compiled annually and apportioned to each system by crop land requirement for that system. Organic fertiliser was applied as solid manure or as liquid slurry using a splash plate, trailing shoe, or by shallow injection. Manure was managed as solid or liquid storage or deposited at pasture. Liquid manure was pumped from the building and stored in uncovered slurry tanks and solid manure was stored uncovered outdoors. All youngstock and dry cow manure was managed as solid storage. All manure from housed milking cows was stored in liquid storage. Types, amounts and application rates of insecticide, fungicide and herbicide sprays applied to each crop were obtained from the Langhill herd database and directly from the supplier (pers. comm. Richard Bray). Forestry and other land managed within the farm boundary such as broadleaf woodland and biodiversity strips were apportioned using annual data prepared by the farm manager for farm subsidy applications (Pers. Comm. H. McClymont, SRUC) depending on the age and type of woodland. Use of petrol and diesel, including the fuel needs of contractors, for each required activity was recorded in the Langhill database and then attributed to a feeding system by task, such as fertiliser application, and then by genetic group. Activities on the farm that required fuel related to crop management included fertiliser applications and herd management, such as feeding. Electricity use (kWh) was estimated from average milk yield per cow in each system using the method of Sheane et al. (2010) which applies 0.051 kWh/kg milk for yields of 6,500-8,500 litres per cow and 0.045 kWh/kg milk for yield >8,500 litres. Electricity was estimated from milk yield because power consumed was not recorded separately for each of the systems. Average annual energy use data is provided in the Supplementary Material. Annual diet digestibility and CP for each of the systems were determined from proportional intakes from monthly TMR and feed component sample analysis.

Goal and Scope
Life Cycle Assessment (LCA) is a key approach to determine environmental or other impacts along a product chain and is carried out according to ISO 10440 and ISO 10444 standards. This research applies an attributional "cradle to farm gate" LCA method as defined by the British Standards Institute (BSI) PAS2050 (BSI, 2011) for the assessment of the life cycle GHG emissions of goods and services. Boundaries applied in this study were "cradle to farm gate" which included all stages from production of farm inputs and raw materials until the milk or animals left the farm. This boundary included supply chain input resources and emissions that are generated off farm, such as those associated with growing and manufacturing purchased feeds, transport, fertiliser manufacture as well as fuel and electricity production. On farm system components included applications of fertilisers, sprays, fuel, crops and field activities, animal feed, livestock of all ages and the management of their manure. Minor emission sources excluded on the basis of materiality PAS2050 (BSI, 2011) in this study were indirect emissions such as staff travel, maintenance of farm buildings, disposal of dead animals and ancillary purchases such as medicine and disinfectants used to clean infrastructure.
A standard functional unit (FU) related to dairy LCA's of fat and protein corrected milk was applied using the following equation ( FU's applied in this study are one kilogramme of FPCM milk leaving the farm gate, total hectares of land required and FPCM production per hectare of total land required. When considering intensification land use is important, because globally land availability is a limiting factor for agriculture (Salou et al., 2017). Allocation describes how GHG's are attributed to the products, and possible co-products that leave the farm and the methods applied affect the estimation of emissions. In this case, as no crops were sold, co-products included animals culled and manure exported from the system. Methods of allocation available in LCA studies include biological causality, system expansion, economic allocation, mass allocation and no allocation (Audsley et al., 1997). In this study a whole farm approach is taken and emissions are allocated between milk and meat using no allocation to meat and 100% to milk.
Emissions from co-product feeds were allocated proportionally by component as described by Vellinga et al. (2013). Greenhouse gas emissions attributed to feeds purchased for the LF, HF and BP systems followed an economic allocation method by feed component in the first instance and a mass allocation method as a comparison, shown in Table 4. This comparison is reported for Select merit herds in the BP, LF, HF and HG feed systems. Additional impact categories are not presented in this study to focus on types and sources of GHG emissions from a range of dairy systems.

Inventory Analysis
A life cycle inventory of annual data from five system years 2006-2010 for HF and LF diets and four system years 2012-2015 for HG and BP treatments was compiled for both genetic merits. Emissions from livestock were calculated using monthly herd dynamic data that was prepared for each of the systems for all years to determine livestock within each of the age categories, those culled, died, or sold, as well as dry and transition cows. Manure management emissions for each of the systems were allocated by the length of time cattle spent at either liquid storage, solid storage or depositing at pasture. Data showing proportions of time that cattle spent in each manure management category are provided in the Supplementary Material. Liquid slurry stems from the housed milking cows and the proportions of time spent grazing were determined and used to allocate deposition directly at pasture. Dry, transition cows and young stock generated solid manure. Manure generated by the dairy systems that was not applied to the crops was exported.
Livestock emissions from dairy cattle included those stemming from manure and enteric CH 4 , direct and indirect N 2 O from manure management, and additionally leaching of N from the soil, and ammonia (NH 3 ) volatilisation arising from the application and deposition of manure. Amounts of N excreted were estimated from N consumed minus N utilised for production, growth and maintenance, which were derived from dry matter intake and CP content of the diets (Dong et al., 2006). Tier II emission factors (EF's) were applied for livestock and manure management and Tier I for fertiliser and crop residue N 2 O (de Klein et al., 2006;Dong et al., 2006). GHG emissions from production, processing, and distribution are embedded in purchased feeds brought onto the farm. Embedded emissions from fertiliser included those associated with manufacture and distribution, which, in the case of the Haber process for N, can be energy intensive. The global warming potentials associated with each feed component within the TMR's of the four diets are presented in Table 4. The LF and HF diets included a proportion of distillery products and the BP system comprised of purchased co-products from the bakery, distillery, brewing and confectionary industries. Estimated GHG emissions per kg for crops in the HG diet are presented in the Supplementary Material and the GWP of minerals in the LF, HF and HG diets was 261 kg CO 2 e/kg using mass and economic allocation (Vellinga et al., 2013). The CT (2010) database was used to source emission factors applied to the production of fertilisers and pesticides ( Table 5). Land category emissions arising from fertiliser application include N 2 O from soil, volatilisation, leaching and run-off as well as N 2 O emissions from crop residues. Carbon sequestration of farm woodland (by age and type) is modelled using Tier 1 IPCC (2006) methodology and reported separately as a reduction of net emissions. Selected emission factors and equations applied within this study are shown in Table 5.

LCA -Impact Assessment
The impact category of interest in this study is climate change, which was assessed by estimating total GHG emissions expressed in kg CO 2 e stemming from the annual inventories of the dairy systems. The inventory prepared for each of the eight systems provided annual farm inputs and outputs from 36 distinct annual inventories which refer to nine calendar years in total for subsequent analysis using SRUC's Agricultural Resource Efficiency Calculator (Agrecalc) v1.4 (SRUC, 2014). Agrecalc (SRUC, 2014) is a carbon foot-printing and resource efficiency tool designed to model emissions at farm level using IPCC methodology (Dong et al., 2006). A PAS2050 (BSI, 2011) accredited model is available online and the tool is used by consultants, farmers, and livestock researchers (Toma et al., 2013;Sykes et al., 2017). Factors applied in Agrecalc (SRUC, 2014) to convert GHG emission flows to CO 2 e were 25 and 298, for emissions of CH 4 and N 2 O, respectively.

Statistical Analysis
Statistical analyses were carried out in R version 3.5.2 (R Core Team, 2013) using lme4, car, and lattice packages (Bates et al., 2015), to determine the effect of dairy production system upon product GHG emissions. Footprints applied in the statistical analysis were estimated using economic allocation of feeds and a FPCM FU for comparability with other studies. A linear mixed model was fitted and included fixed effects of feeding regime, genetic merit, and a random effect of year. An ANOVA, and a Tukey pairwise comparison test was carried out to determine significance of the production systems using the following model, Where, y ijk GHG emissions using economic allocation and expressed per kg FPCM µ = grand mean F i = feed type (i = 1 to 4) fixed effect M j = genetic merit j = 1 to 2 fixed effect T k = year k = 1 to 9 random effect ε ijk = residual error

Sensitivity and Uncertainty Analysis
Stochastic simulations were carried out using ModelRisk5 (Vose Software) to assess the effect of annual variation in neutral cellulase gammanase digestibility (NCGD) and CP content of the rations on dairy system GHG emissions. Agrecalc (SRUC, 2014) was used to estimate baseline carbon footprints from all feed systems and both genetic merits using economic allocation of feeds and average annual values for NCGD and CP content. A FPCM FU is used for consistency and ease of comparability. Variation applying mass allocation of purchased feed emissions is not assessed in this study. Monte Carlo simulation using repeated random sampling was used to generate distributions of footprints for the dairy systems that accounted for uncertainty stemming from variability in NCGD and CP content for each treatment group. Descriptive statistics for the NCGD and CP distributions are shown in Table 6. Exponential and Log Laplace distributions were fitted to NCGD and CP analysis results, respectively, and Monte Carlo simulations with 10,000 iterations (seed = 2,605) were carried out. To determine if there was a significant difference in sensitivity the footprint outputs from the linear mixed model detailed in the previous section were refitted and modelled and an Anova was carried out to test for a significant difference between the models.

Statistical Analysis
An ANOVA was conducted to compare GHG's and test for significant differences between the four feeding regimes and two genetic merits using an economic allocation of feed emissions and results ( Table 7). Normality checks and Levene's test were carried out and the assumptions were met. The effect of the interaction was significant (p < 0.01). There was a significant difference in mean GHG's per kg FPCM between the feed groups [F (3, 28) = 15.6, p < 0.001] and the genetic merits [F (1, 28) = 46.5, p < 0.001]. Post hoc comparisons using the Tukey test showed mean Select and Control merit GHG totals to be significantly different (p < 0.05) in LF, BP and HF feed types but no significant difference was found between the HG and HF diet. Tukey test results showed that the LF diet was significantly different from BP (p < 0.05), the HF and the HG (p < 0.001) diets for GHG totals. Significant differences in mean embedded emissions were found between all management systems (p < 0.05) and livestock emissions were all significantly different (p < 0.05) apart from HFS and HGS regimes ( Table 7).

Effect of Allocation Method
For Select merit herds the average annual carbon footprint for milk produced in each of the dairy systems was calculated using both economic and mass allocation of feed components which resulted in large differences in the ranking of the different systems (Figure 1). In the housed systems, with economic allocation, the BP diet is associated with higher emissions per kg product at 1.07 kg CO 2 e/kg FPCM compared with the LF system, which averaged 0.95 kg CO 2 e ( Table 7). Economic allocation of feed components in the HF and HG grazed system footprints led to similar product emissions per kg of FPCM at 1.15 and 1.16 kg CO 2 e, respectively, as a result of lower embedded emissions in the HG system ( Table 7) that were outweighed by higher emissions from energy and fuels. The HG system is connected with higher fuel use associated with crop production on farm. Economic allocation of emissions for the HF and HG TMR's were  206 and 252 kg of CO 2 e, respectively. Mass allocation of feed components led to increases in system footprints per kg FPCM because emissions were higher for all diet TMR's, except HG. TMR emissions were 473, 2,072, 757 and 252 kg CO 2 e/tonne, for the LF, BP, HF and HG systems, respectively. On average, in the housed systems, BP diet emissions increased per kg product from 1.07 to 3.79 kg CO 2 e/kg FPCM using economic and mass allocation, respectively. LF system product emissions averaged 0.95 kg CO 2 e/kg FPCM using economic allocation and 1.3 kg CO 2 e/kg FPCM using mass allocation (Figure 1). Control merit results followed the same rank order and are reported in the Supplementary Material.

Sensitivity and Uncertainty Analysis
Simulated footprints were generated using economic allocation of feeds to obtain distributions of dairy management system results if variation in NCGD and CP levels were considered. Mean dairy system footprint simulations considering the effect of NCGD and CP variation differed significantly from mean system footprints estimated using average annual NCGD and CP values (p < 0.001). Mean milk footprints were increased in the BP and HF systems and decreased in the LF and HG systems, in comparison with applying an average annual figure for digestibility and CP. Accounting for nutritional variation of the rations throughout the year had widened footprint ranges and altered comparative dairy system performance ranking. Select merit cows in the BP regime incurred greater average emissions per kg FPCM, at 1.21 kg CO 2 e, when compared to the LFS, HGS and HFS systems which averaged 0.92, 1.15, and 1.17 kg CO 2 e/kg FPCM, respectively (Figure 2). Average Control merit footprints followed the same rank and were higher than the Select merit in the same system at 1.09, 1.26, 1.35, and 1.42 kg CO 2 e/kg FPCM for LFC, HGC, HFC and BPC systems, respectively (Figure 2). Higher average diet digestibility combined with a lower average CP content in the LF system ( Table 6) led to lower mean emissions, in comparison with BP, HF and HG diets in both genetic merits (Table 7) and the standard method applying an average annual figure for digestibility and CP. When compared to the LF, BP and HF regimes Select merit cows managed in the HG system attracted a wider range of potential carbon footprints (Figure 2). Sources of GHG's within the carbon footprints varied by dairy management regime, therefore farm mitigation strategies may prove more effective if applied by system type. Land and crop GHG emissions stem from crop residues, manure and fertiliser application and these ranged from zero in the BPS system to 0.11 kg CO 2 e in the HGS system (Table 8). Embedded emissions are generated by energy consumed in the manufacture of feeds, fertilisers and pesticides and also in the use of bedding. Embedded emissions were greatest in the BPS housed system, at 0.46 kg CO 2 e /kg FPCM, because all feed and bedding were imported. The HGS grazed system attracted higher embedded emissions than the HFS system, as a larger area of on farm crop land replaced purchased feeds (Table 8).
Livestock emissions that arise from enteric fermentation and manure management were greater in the HF and BP systems, at 0.79 and 0.68 kg CO 2 e/kg FPCM, compared with 0.57 and 0.66 kg CO 2 e/kg FPCM in the LF and HG systems, respectively (Figure 3). Higher emissions arose from greater amounts of manure stored in the BP system and from depositions while at pasture in the HF system. Emissions related to energy use were greater in the HG system, as this stemmed from the fuel used for crop related activities. Sequestered carbon estimated to occur within the woodland in the LF, HF and HG systems, lowered  Select merit footprints by 0.01, 0.03, and 0.05 kg CO 2 e/kg FPCM, respectively ( Table 8). Control merit footprints were reduced by 0.02, 0.03, and 0.06 kg CO 2 e/kg FPCM, for the LF, HF and HG, systems, respectively ( Table 8).

Effect of Increased Legume Forages
Economic allocation of feed components generated similar average product emissions for HF and HG systems at 1.15 and 1.16 kg, respectively ( Table 7). Mass allocation increased the HF milk footprint to 1.67 kg CO 2 e/kg due to the proportion of distillers' grains in the ration. Accounting for nutritional variation slightly reduced the HG average to 1.15 kg CO 2 e per kg FPCM and increased the HF to 1.17 kg CO 2 e/kg FPCM (Table 8). If C sequestration was not included, the footprints were, on average, equivalent at 1.19 kg CO 2 e/kg FPCM (Table 8). Trade-offs between livestock manure emissions and energy use to grow crops has led to similar milk total emissions being returned from the HG and HF systems ( Table 7).
Total on-farm land use per milking cow increased from an average of 0.86 ha to 1.23 ha when comparing the HF and HG systems. The HG system attracted greater embedded emissions than the HF system, these stemmed from the use of fungicide and herbicide applications to the wheat and spring beans.

Effect of Genetic Merit
Control merit total product footprints across each of the management regimes were significantly higher (p < 0.001) in comparison with high production Select merit cows, on average by 15% (Table 7). Livestock and embedded emissions were also significantly higher from control merit cows (p < 0.01). On average, comparing raw milk quality across each of the management regimes, Control merit cows yielded less milk volume and produced lower percentage fat and protein content ( Table 3) than Select merit animals consuming the same diet ( Table 2). When encompassing nutritional variation lowest to highest mean system ranking for Control merit was LF, HG, HF, BP and this rank order was preserved for Select mean footprints (Figure 2). Control merit carbon footprints were higher than the Select merit animals apart from in the LF system, where the Control merit resulted in slightly lower emissions than Select merits in the BP, HG and HF management. The housed LF regime incurred less GHG's / kg FPCM than the BP system irrespective of merit or allocation method mainly because of emissions embedded in the production of feeds. Raw milk at the farm gate produced by Control merit cows, within the BP system, were associated with greater product emissions at 1.42 kg CO 2 e/kg FPCM, than other systems and merits.

Effect of Land Use as a Functional Unit
On average, considering both on and off-farm land requirements for Select merit cows (Table 9), the BP system required the least amount of land in total, due to the high proportion of human inedible crop products and industry co-products. Land as a functional unit showed the HG system as least GHG intensive. When output of product is included with total land, it was the BP dairy system that emitted fewer GHG per hectare (Table 10).

DISCUSSION
Using an LCA approach, this study demonstrates the importance of allocation method used to attribute GHG emissions of animal feeds and, in addition, the effect of nutritional variation on the carbon footprint of novel and conventional UK dairy systems. The ranked performance of dairy management types depends on the approach used to calculate impact and whether or not uncertainty is included (Table 11). Economic allocation resulted in mean dairy system emissions that ranged from 0.95 to 1.16 kg CO 2 e/kg FPCM and 1.13 to 1.28 kg CO 2 e/kg FPCM for S and C merit, respectively ( Table 7). Mean system emissions were lower than the UK average of 1.25 kg CO 2 e (AHDB, 2014), except for C merit cows in HG and HF regimes ( Table 7). Results expressing emissions totals per kg FPCM and per hectare produced in this study are in line with similar research such as Salou et al. (2017) and O'Brien et al. (2014). Carbon footprint results presented for the housed and composite systems are also in line with confinement and mixed systems reported in a review of 30 LCA studies from 15 different countries (Lorenz et al., 2019).

Effect of Allocation Method and Uncertainty Analysis
Carbon footprints were, on average, higher using mass allocation EF's, and accounting for uncertainty, stemming from changes in diet CP and digestibility, altered the dairy system ranking. Mass allocation of feed component emissions raised product emissions on average by 41%, for the LF and HF rations, which comprised of a mixture of grown crops and purchased concentrates. The BP ration was formulated from mainly nonhuman edible food and drink industry co-products and the TMR produced only marginally less GHG's than the LF diet. With economic allocation, ration EF's per tonne in the LF diet was 256 kg CO 2 e and 249 kg CO 2 e in the BP diet. Milk quality in the BP diet was, on average, lower in fat and protein content in comparison with the LF system. In the grazed HF system, mass allocation of feed components increased product footprints because distillers' grains and rapeseed meal elevated emissions. The effect of the elevation in TMR emissions using mass allocation in the HF system is not as pronounced as in the BP system due to the time spent grazing by HF cows. The HG grazed system footprints were least effected by allocation method as purchased feed was limited to minerals and diet component EF's were all equivalent apart from wheat grain. Ration EF's increased using mass allocation because feed components such as industry co-products tend to be associated with higher emissions when additional processing into animal feed is required. Nutritional quality of animal feed varies, and in this study feed analyses showed a higher mean CP and lower mean digestibility in the BP ration when compared to the LF system. The BP system had the lowest ranges in digestibility and CP content, possibly because there was no effect of local climate on farm grown crops in this ration. Reducing the CP intake of the dairy cow diet would help in reducing GHG emissions (particularly N 2 O) and UK research has shown that loss of production can be lower than expected (Reynolds et al., 2016). Other environmental and financial strategies to improve nitrogen use efficiency, such as home-grown legumes, should have positive consequences for GHG emissions through increased protein supply and the reduced need for N application from imported inorganic fertilisers.

Effect of Increased Legume Forages
The HF system ration included feeds requiring crop rotations of grass silage, maize and wheat, which required N fertilisers and were ensiled on farm. Crop products were combined with purchased concentrates and on average the HF diet consisted of ∼75% forage on a DM basis and 1.3 tonnes of concentrate per cow (March et al., 2017). In comparison with the HF diet, the HG ration required less maize crop, as the purchased distillers' grains and rapeseed meal were replaced with farm grown proteins, such as, spring beans and lucerne. The HG herds were grazed for an average of 26% of the year, whereas the HF cows were grazed for an average of 30% annually and attracted greater emissions from dung and urine deposition at pasture. For Select merit cows, however, feed intakes on a DM basis were similar in both the HF and HG systems ( Table 2). Average milk yield reduced slightly, by 98 kg per cow from 7,575 kg in the HFS system to 7,477 kg in the HGS system, although, milk quality was similar in both the systems ( Table 2).
The HG system is a comparatively high emitter using economic allocation however, Table 11 shows this regime outranked all the other systems using mass allocation because no additional emissions are generated by imported products. In this case mass allocation methods and sensitivity analysis of nutritional variability highlight the benefits of a self-sufficient agricultural system, which may contain positive consequences when incorporating mitigation measures or when moving to more circular economic methods of farming. The HG system also had the lowest area-based emissions, a consequence of the replacement of synthetic N fertilisers by N inputs through biological fixation in comparison with the HF system. Legumes altered the composition of the footprints however the longterm effects of soil conditioning or crop disease prevention were not quantified by carbon footprinting and carbon sequestration modelling needs to be improved (Sykes et al., 2017) to reflect these and other desirable consequences. Mitigation of emissions related to inputs could be achieved in the HG system by reducing pesticide use and using renewable energies on farm.

Effect of Genetic Merit
On average across all diets, product emissions from Control merit cows using an economic allocation were 15% higher, indicating that improving genetic merit offers an immediate emissions reduction strategy, mainly through increased milk yields. In addition, GHG emissions can be reduced by selecting for feed efficiency (Bell et al., 2011). This could be accelerated using techniques such as genomics in the herd to enhance overall feed efficiency and in vitro fertilisation (Gifford and Gifford, 2013;Pryce and Bell, 2017;Hailu, 2018). Financial analysis of the LF and HF regimes found a Control merit cow in a housed regime to be least profitable because milk yields were not sufficiently high to justify the feed costs (March et al., 2017). Considering the diets, the total emissions differed significantly (p < 0.05), apart from the HF and HG rations, however, livestock, energy and embedded emission types did differ significantly between these systems. This highlights that dairy systems mitigation potentials, and measures implemented, should be quantified and designed by considering the production method and the emission source.

Effect of Land Use as a Functional Unit
Novel rations such as those used in the BP system required less land, however, incorporating high percentages of co-product based animal feeds can lead to greater GHG emissions as a consequence of upstream processing, such as drying or milling, which can be energy intensive (Vellinga et al., 2013). Ruminant diets for high yielding cows can be formulated to achieve lower emissions, and to make more efficient use of human inedible co-products (Wilkinson and Garnsworthy, 2017), however, not all co-product feeds are low carbon, and feeding TMR's all year round usually requires cows to be housed in adequate modern animal housing facilities with slurry storage systems. In Scotland, industry co-products have traditionally been used as animal feeds, however, feeds such as distillers' grains contain added water, which stems from the mashing stage of the whisky making process. This can hugely inflate mass balance emissions, and for some products drying grains requires a substantial input of energy, as the water content has to be reduced from ∼75% to under 10% (Bell et al., 2012). Product quality in the BP system was also reduced (Tables 2, 3), this was reflected through lower milk fat and protein which have financial consequences for farm income.
Comparisons of low input grass based, mixed, and fully housed intensive dairy systems are valuable to explore uncertainty and mitigation pathways, rather than to justify efficacy of a particular method of farming. Between and within countries agricultural practises vary and livestock farming is to some extent governed by history, culture, and tradition (Boogaard et al., 2011). Overall focus should be turned to mitigation of emissions, adaptation to changing climates, improving comparability of LCA's, and communicating uncertainty. GHG emissions from dairy farming can be mitigated through multiple pathways such as increasing the longevity of cows within a herd, improving fertility, lowering initial calving age and by reducing enteric methane and improving digestibility of cow rations (Garnsworthy, 2004;Wilkinson and Garnsworthy, 2017). Reduced enteric CH 4 and increased digestibility could be achieved through the reformulation of the diet or through feeding additives and supplements (Knapp et al., 2011(Knapp et al., , 2014. In less intensive dairy systems, enteric CH 4 emissions can be reduced by increasing yields (Yan et al., 2010). Renewable energies and technologies such as anaerobic digestion can be effective in reducing emissions from manure storage, with one study reporting reductions of up to 36% (Weiske et al., 2006;Battini et al., 2014).

CONCLUSIONS
Mass and economic allocation methods, and land use functional units, are shown to generate alternatively ranked footprint results. Monte Carlo simulated system footprints considering the effect of variation in feed digestibility and crude protein differed significantly from system footprints using standard methods. Implications are that dairy farm footprint calculations should incorporate the variation in diet digestibility and crude protein content where possible. Using an economic allocation, a localised home-grown feeding regime had the highest C footprint, however, this more self-sufficient system was associated with the lowest footprint using mass allocation and attracted the lowest area-based emissions, when not considering milk output. This result suggests the need for dairy system footprint results to be expressed in multiple units and to be mindful that methods used to allocate inputs can affect carbon footprint results. It is expected that in developing economy-wide reductions in greenhouse gas emissions, mass and area-based assessments of mitigation are most likely to guide the delivery of policy objectives.

DATA AVAILABILITY STATEMENT
The datasets presented in this article are not readily available because raw data is not publicly available. Requests to access the datasets should be directed to Margaret D. March, maggie.march@sruc.ac.uk.

ETHICS STATEMENT
The animal study was reviewed and approved by SRUC Animal Experiments Committee (AEC).

AUTHOR CONTRIBUTIONS
MM, RR, and AS contributed to the concept and methodology. MM and AS carried out the analysis. All authors contributed to the article and approved the submitted version.

FUNDING
This research was funded by the TRUE project, under the EU Horizon2020 Research and Innovation Programme (Grant