2.1 Farm description

The study was conducted based on a research farm in the University of Helsinki at Viikki. The farm is located along the coast of southern Finland about 7 km from the center of Helsinki. The area is classified as a cold hemi-boreal humid continental climate. The research farm has a dairy barn with 61 dairy cows, whose average milk production is around 9 000 kg per cow per year. The farm administrates 155 hectares of arable land, and the dominant soil texture is clay loam. Feed is mainly produced on the farm, and it consists of grass silage and feed grain. Field plots closest to the cowhouse are mainly used as pastures, and a 3-4-year silage grass crop rotation is applied to the rest of the fields. A grassland is renewed by an autumn ploughing followed by application of dry manure, produced at the farm, and cultivation of spring cereal, which is undersown by perennial grasses. The grasses are thereafter harvested and fertilized 2–3 times during each growing season.

The Station for Measuring Ecosystem - Atmospheric Relations Agriculture (SMEAR-Agri) was established in 2021 and is located at the research Farm (60°13’11"N 25°01’13"E). The monitored field is 6.8 ha under a rotation of perennial grass with one year of cereal cultivation in between 3 and 4 years of grassland cultivation. It is representative of the typical grass-crop rotation and management practices used across the entire farm. The soil at the site is characterized by a clay texture and classified as artificially drained Stagnosols and Planosols. Given that historical land use and management play a key role in shaping the long-term dynamics of SOC change, historical data are essential for simulating the long-term trends and initializing soil models. The field’s grass-crop rotation history from the Viikki Research Farm’s founding in 1996 until 2023 was presented in Figure A1. Field management during the measurement years 2022–2023 included grass cultivation, harvest and fertilization, grassland renewal by ploughing, manure application and cereal crop (barley) cultivation, and a final newly established grassland (Table A2). In detail, the field was plowed down to 20 cm for the termination of perennial grass in November 2022. The field remained bare without vegetation cover until April 2023 when it was harrowed, fertilized with dry manure and sown with barley undersown with a mixture of perennial grasses. Mineral N fertilizer, ammonium nitrate, was applied on 10th May 2022 (108 kg N ha− 1) and 25th May 2023 (67.5 kg N ha− 1), and barnyard manure in 1701 kg ha− 1 dry mass containing 12.9 kg N ha− 1 and 5.4 kg P ha− 1 was applied on 27th April 2023. Grass was harvested twice in July and August 2022. Barley seeds were sown on 25th May 2023 and harvested on 17th August 2023.

Fig. 1 The alternative text for this image may have been generated using AI. Full size image The SMEAR-Agri Station of Viikki Research Farm in Helsinki on the southern coast of Finland. The target field with the flux measurements is marked in light green

2.2 Life cycle assessment of CFP

2.2.1 Scope, system boundary and functional unit

The CFP of milk production in Viikki Research Farm in 2022 and 2023 was assessed. The Solagro Carbon Calculator (Tuomisto et al. 2015) was used to calculate the total GHG emissions over the two years except for CO 2 emissions from SOC stock changes. The Solagro Carbon Calculator is a farm-level tool designed to quantify GHG emissions for the most common farm types in the EU-27, based on international standards for LCA and carbon footprinting (GHG Protocol 2011; ISO 14040 2006; ISO 14044 2006; ISO/TS 14067 2013). Its methodology for quantifying emissions aligns with the 2006 IPCC Guidelines for national greenhouse gas inventories. Characterization factors for background processes, such as emission factors for agricultural inputs (including mineral fertilizers, feedstuff, pesticides, and seeds), are sourced from Weiss and Leip (2012), Wood and Cowie (2004), ADEME (2012), GESTIM (2010) and Brentrup and Palliere (2008). Energy use emission factors are country-specific for electricity and account for both upstream and combustion emissions for fuels, using data from ELCD (2001).

The assessment used a “cradle to farm-gate” system boundary (Fig. 2), which included the raw milk life cycle from the production of inputs until milk leaves the farm. Materials and processes on-farm were considered, including forage and crop production, fuel and electricity use, manure and livestock management and associated emissions. The direct farm emissions included cows’ emissions from enteric fermentation, emissions from manure management for storage and treatment of the manure, N 2 O emissions from both organic and inorganic fertilizer application to the soil, indirect N 2 O emissions via NH 3 volatilization and N leaching and runoff, and CO 2 emissions from energy use. The quantification methods for those direct farm emissions were based on IPCC 2006; while emission factors for electricity and fuel type were from data source ELCD (2001). In addition to the direct farm emissions, the assessment considers the emissions occurring during the production and transportation of resources used on the farm, including fertilizers, seeds and herbicides, purchased feedstuffs, rearing of purchased animals, farm buildings and materials (Source of emission factors: ADEME 2012), and fuels manufacturing and transportation (Source of emission factors: ELCD 2001).

GHG emissions are expressed as tons of CO 2 equivalent (t CO 2 eq). The CFP of milk production is expressed as total GHG emissions in ton CO 2 equivalent per 1 ton of fat and protein corrected milk (FPCM). FPCM is calculated based on fat and true protein content, following IDF (2022) shown in Eq. (1):

$${{\rm{M}}_{{\rm{milk}}}}{\rm{ = P \times }}\left({{\rm{0}}{\rm{.1226 \times }}{{\rm{M}}_{{\rm{fat}}}}{\rm{ + 0}}{\rm{.0776 \times }}{{\rm{M}}_{{\rm{prot}}}}{\rm{ + 0}}{\rm{.2534}}} \right)$$ (1)

Where M milk is FPCM (kg year− 1), P is production of raw milk (kg year− 1), M fat is fat content in raw milk (%) and M prot is true protein content in raw milk (%).

Fig. 2 The alternative text for this image may have been generated using AI. Full size image System boundary of the studied milk production system. GHG: greenhouse gas

2.2.2 Co-product allocation

Co-product allocation for milk and meat was performed independently of the Carbon Calculator, employing a biological causality-based approach. This method, described by Thoma et al. (2013), partitioned emissions according to the energy invested by the animal in milk and meat production, following the Eq. (2):

$${\rm{A}}{{\rm{F}}_{{\rm{milk}}}}{\rm{ = 1 - 6}}{\rm{.04 \times BMR}}$$ (2)

Where the allocation factor for milk (AF milk ) is determined by beef: milk ratio (BMR), expressed in kilograms of meat per kilograms of FPCM. The calculated AF milk value used in this study was 84.2% and 87.3% in 2022 and 2023, respectively.

2.2.3 Life cycle inventory data

Farm characteristics for life cycle inventory are listed in Table 1. These inputs were used to produce the annual milk production. The GHG emission outputs in each year was subsequently normalized to the functional unit in t CO 2 eq t FPCM− 1.

Table 1 Farm characteristics in 2022 and 2023 Full size table

2.2.4 Sensitivity analysis

The sensitivity of the CFP assessment of milk production using the Solagro Calculator was tested to evaluate how variables of the input parameters affect the results. The critical parameters regarding their contribution to the uncertainty of results include the live body weight of animals which affects the calculation of methane (CH 4 ) emitted from enteric fermentation and the N fertilizer application which is a key input factor for calculating N 2 O emissions from soil management. To test the sensitivity of the assessment, we adjusted the critical input parameters of average cow live weight and N chemical fertilizer application based on assumptions of their uncertainty level. The uncertainty range of average cow live weight and N fertilizer used for the sensitivity analysis was ± 2% and ± 10%, respectively, according to Monni et al. (2007).

2.3 Integrating SOC stock changes into CFP of milk production

Emissions from SOC stock changes for on-farm feed production were excluded from the CFP calculation within the Carbon Calculator. Instead, the emissions were quantified independently, using three approaches: the IPCC Tier 1 method, process-based DNDC modelling, and direct CO 2 flux measurements. For purchased feedstuffs, emissions for processing and transportation were regarded as indirect emissions and were calculated using emission factors from the Carbon Calculator. The SOC stock changes to produce purchased feedstuffs were not included in the study. All three methods for SOC stock changes used agricultural management and soil data from the field monitored by SMEAR-Agri to ensure comparability. The SMEAR-Agri was established in 2021, and flux measurement data were available across 2022 and 2023. The monitored field was in its third year of perennial grass in 2022 and was subsequently cultivated with spring barley in 2023. It represented the typical mode of grass-crop rotation and management practices on the whole-farm level for feed production. Therefore, the CO 2 emission from SOC stock changes of the monitored field in 2022 was associated with grass silage production and that in 2023 was associated with crop grain production. To integrate SOC stock changes into CFP of milk production, results on monitored field in 2022 were used to calculate the CO 2 emissions from SOC stock changes for on-farm grass silage production in both 2022 and 2023. Similarly, SOC stock changes on monitored field in 2023 were used to calculate the CO 2 emissions from SOC stock changes for on-farm crop production in both 2022 and 2023.

The change in SOC stock (ΔSOC) was expressed in kg C ha− 1 a− 1 with negative values referred to loss of SOC stock and positive values referred to SOC storage (Tables 2 and 5). To integrate it into the CFP of milk production, the ΔSOC was first converted to CO 2 equivalent by applying the conversion factor 3.664 (the ratio of molecular weight of CO 2 to the atomic weight of C). To be consistent with LCA results, the loss of SOC stock was converted to positive CO 2 emissions and the SOC storage was converted to negative CO 2 emissions. The resulting CO 2 emissions from ΔSOC per ha per year (expressed in kg CO 2 ha− 1 a− 1) were then scaled up to the farm level by using the actual area (ha) utilized for either crop grain or grass silage production in the respective year. The calculation yielded the total emissions from SOC stock changes (in kg CO 2 a− 1). Since these CO 2 emissions attributed to the production of on-farm feeds (grain and silage), they were included in the farm’s total GHG emissions for animal production and subsequently allocated between milk and meat products, then normalized to t CO 2 eq t FPCM− 1.

2.3.1 Calculation of SOC stock changes via IPCC Tier 1 method

The calculation of SOC stock changes followed Tier 1 method in Chap. 2 of volume 4 (IPCC, 2006). The method is based on the assumptions that: (i) over time, SOC stock reaches a stable value specific to the soil, climate, land-use and management practices; and (ii) SOC stock change during the transition to a new equilibrium occurs in a linear fashion over a period of 20 years. According to the time-based steady SOC stock assumptions, Tier 1 gives a default reference for the initial SOC stocks (SOC REF ) in top 30 cm, in terms of the soil type and climate region. The new equilibrium SOC stock after 20 years is calculated by scaling the SOC REF with a land-use factor, a management factor, and a factor for C input levels, according to each of their default categorizations.

In this study, the default method for mineral soils was applied, and the change is computed using Eq. (3) and Eq. (4).

$${\rm{\Delta SOC = }}\left({{\rm{SO}}{{\rm{C}}_{{\rm{Equilibrium}}}}{\rm{- SO}}{{\rm{C}}_{{\rm{Initial}}}}} \right){\rm{/20 \times T}}$$ (3)

Where ΔSOC is the change in SOC stock (t C ha− 1) for the transition of land use and management practices, SOC Equilibrium is the assumed new equilibrium SOC stock 20 years after transition, and SOC Initial is the SOC stock before transition. T indicates the actual operating years of the certain land use and management practices after the transition. Nine transitions occurred due to the rotation between perennial grass and cereal crops during 1996–2023 on the monitored field. The SOC Equilibrium and SOC Initial were thus calculated for each of the transition. The SOC Initial before the first transition in 1996 was obtained based on assumptions that the land was initialized as grassland, natural land and cropland, respectively, according to the IPCC default values.

$${\rm{SO}}{{\rm{C}}_{{\rm{Equilibrium}}}}{\rm{ = SO}}{{\rm{C}}_{{\rm{REF}}}}{\rm{ \times }}{{\rm{F}}_{{\rm{LU \times }}}}{{\rm{F}}_{{\rm{MG}}}}{\rm{ \times }}{{\rm{F}}_{\rm{I}}}$$ (4)

Where SOC Equilibrium indicates the new equilibrium SOC stock (t C ha− 1) 20 years after transition, SOC REF means the reference SOC stock in the topsoil from 0 to 30 cm (t C ha− 1), F LU is the stock change factor for land-use systems or sub-systems for a particular land-use, F MG is stock change factor for management regime and F I is the factor for the input of organic matter. The SOC REF value used in this study was 95 t C ha− 1 for cold temperate moist climate region and soils with high-activity clay minerals. The F LU , F MG and F I used for grassland were 1, 1.14 and 1, respectively, according to practices during the grass year on the target field (Chap. 6 of volume 4 of NGGI-IPPC-2019 “Grassland”). The F LU , F MG and F I factors used for cropland were 0.7, 1.04 and 1, respectively, according to practices during crop year on the target field (Chap. 5 of volume 4 of NGGI-IPPC-2019 “Cropland”).

2.3.2 Flux measurements and calculation of C balance for SOC stock change

To quantify changes in SOC stocks, we measured the NEE of CO₂ between the ecosystem and the atmosphere using the EC technique. This micrometeorological approach quantifies the turbulent exchange of CO₂ by measuring high-frequency fluctuations in vertical wind and gas concentrations, providing direct, continuous, and area-averaged estimates of ecosystem-scale CO₂ fluxes. Continuous CO 2 and H 2 O flux were measured using EC technique at the SMEAR-Agri Station from 1st January 2022 by employing an enclosed infrared gas analyzer (LI-7200, LI-COR Biosciences, USA) and a 3-D sonic anemometer (uSonic-3 Sci HS, Metek GmbH, Germany) at 10 Hz. The sensors were positioned above the soil surface at 2.4 m and a horizontal separation of 0.2 m. The sample air was transported to the gas analyzer through a 91 cm long heated tube with a 12 L min− 1 flow rate and sampled at 10 Hz. Auxiliary meteorological and continuous soil measurements at the field included: Relative humidity (HMP110, Vaisala Oyj, Finland), photosynthetically active radiation (LI-190R quantum sensor, LI-COR Environmental GmbH, Germany), and soil temperature and moisture at 10 cm depth (Teros-12, METER Group, USA). Air temperature was acquired from the official Helsinki Malmi airport weather station operated by the Finnish Meteorological Institute, located 2 km north-east of the measurement field. The 10 Hz raw data were processed to produce 30-min turbulent CO 2 fluxes, which were calculated with Eddypro software (v. 7.0.9, LI-COR Biosciences, USA) applying standard quality control and corrections protocols. The effects of temperature and water vapor fluctuations were accounted for with a long enough sampling line and internally within the LI-7200 analyzer, respectively. The 30-min fluxes were screened for relative stationarity and integral turbulence characteristics and discarded if the combined flag of these tests indicated the lowest quality on a three-level system (Foken and Wichura 1996). Also, the periods of weak turbulence were removed by applying a friction velocity threshold (0.19 m s− 1), which was estimated with the moving-point-transition method (Reichstein et al. 2005). Finally, a source area analysis was carried out to remove periods when the data did not originate from the target area. This was performed by calculating the cross-wind-integrated footprint distribution from the measurement tower to the edge of the field by using the model developed by Kljun et al. (2015). The 30-min period was discarded if less than 70% of the measured flux originated from the target area. Missing data were filled using gradient boosted decision trees (XGBoost; Chen and Guestrin 2016) following Vekuri et al. (2023).

Based on the CO 2 flux measurements, the SOC stock change (ΔSOC) on an annual basis was calculated by Eq. (5):

$${\rm{\Delta SOC = NEE - }}{{\rm{C}}_{{\rm{Harvest}}}}{\rm{ + }}{{\rm{C}}_{{\rm{Manure}}}}{\rm{ - \Delta }}{{\rm{C}}_{{\rm{Residues}}}}$$ (5)

Where NEE is the net ecosystem exchange from flux measurements during the reference year and positive value indicates carbon capture from the atmosphere (kg C ha− 1), C Harvest indicates the aboveground biomass C removed by harvest (kg C ha− 1), C Manure is the C input to soil from manure application, and ΔC Residues means the change of C stored in residues pool (kg C ha− 1) in the year. Positive values in ΔSOC and ΔC Residues refer to increase in the corresponding C pool. Detailed values applied for the calculation were shown in Table 2.

Table 2 The net ecosystem exchange (NEE), carbon removed by harvest (C Harvest ), carbon input from manure (C Manure ) and change of carbon in the residue pool (ΔC Residue ), as well as estimated soil organic carbon stock changes (ΔSOC) according to Eq. (5) in 2022 and 2023 Full size table

The C Harvest was calculated by multiplying the aboveground biomass in dry matter basis with the C content, which was assumed as 45%.

The C Manure was calculated based on the dry mass of applied manure (1701 kg ha− 1 applied in 2023), organic fraction of dry mass in the manure (79.9%) and C in organic fraction (50.4%), the latter two parameters were estimated according to Kätterer et al. (2011).

The C stored in undecomposed root biomass was recognized as the C in the residue pool. In this study, the grassland was ploughed in November 2022 and April 2023 and was converted to cropland for spring barley in 2023. The residue C pool in 2022 was assumed to be unchanged due to the consistent land use for perennial grass, as the belowground residues were kept stable until the end of this year. Therefore, ΔC Residues in 2022 was assumed as zero. Whereas the C in root biomass of grass was assumed to decay in 2023 due to the land use conversion after ploughing in the spring. Therefore, a negative value (-1014 kg C ha− 1) for ΔC Residues in 2023 was obtained. The calculation was based on the C content (45%) and root biomass of grass, which was estimated from aboveground biomass (3000 kg ha− 1, from cumulative harvest) and root: shoot ratio (1.33). The later was obtained from measured average grass shoot and root dry mass in Southern Finland according to Palosuo et al. (2015)

2.3.3 DNDC modelling for SOC stock changes

The DNDC model simulates C and N dynamics in agro-ecosystem considering the interaction between ecological drivers including climate, soil, vegetation and management (Gilhespy et al. 2014; Brilli et al. 2017). To simulate SOC changes in agricultural land, DNDC requires climate data, soil properties and cropping practices as input. The historical daily air temperature and precipitation used in the model were obtained from the Helsinki Vantaa weather station (1996–2006) and Helsinki Kumpula weather station (2007–2023) operated by the Finnish Metrological Institute. Site specific soil properties i.e. soil bulk density, texture, clay content and pH were measured based on samples collected in autumn 2021 (Ahrends et al. 2023). Initial SOC content in the topsoil layer was measured from samples collected in 1996, data were acquired from Mokma et al. (2000). Cropping practices regarding crop rotation, tillage, fertilization and manure application, cultivation and harvesting were used in the model. Historical land-use data from 1996 to 2020 (Table 1) were applied for model initialization. Recognizing the significant influence of crop growth on soil C, N and water dynamics, and consequently on biogeochemical processes, calibrating the crop growth parameters was prioritized following the user’s guide for the DNDC model. The DNDC model was therefore parameterized using crop yield data in 2021, 2022 and 2023 to ensure reliable representation of SOC dynamics in 2022 and 2023.

2.4 Sources of uncertainty