Sensitivity analysis of the INRA 2018 feeding system for ruminants by a one-at-a-time approach: Effects of dietary input variables on predictions of multiple responses of dairy cattle

In the INRA 2018 feeding system for ruminants, the prediction of multiple animal responses is based on the integration of the characteristics of the animal and the available feedstuff characteristics, as well as the rationing objectives. In this framework, the characterization of feedstuffs in terms of net energy, digest-ible protein, and fill units requirs information on their chemical composition, digestibility, and degradability. Despite the importance of these feed characteristics, a comprehensive assessment of their impact on the responses predicted by the INRA 2018 feeding system has not been carried out. Thus, our study investigated how variables predicted by the INRA feeding system (i.e., outputs) for dairy cows are affected by variation in feed characterization (i.e., inputs). Five input variables were selected for the sensitivity analysis (SA): CP, OM apparent digestibility (OMd) , GE, effective degradability of nitrogen assuming a passage rate of 6%/h (ED6_N) and true intestinal digestibility ( dr_N ) of nitrogen. A one-at-a-time SA was performed on predicted digestive, productive and environmental output variables for dairy cows with 6 contrasted diets. These 6 diets were formulated to meet 95% of the potential daily milk production (37.5 kg) of a multiparous cow at wk 14 of lactation. Then, the values of the 5 key input variables of each feedstuff were randomly sampled around the INRA 2018 feed table values (reference point). The response of the output variable to the variation of the input variable was quantified and compared using the tangent value at the reference point and the normalized sensitivity coefficient. Among the major final output variables, CP and dr_N had the greatest impact on nitrogen (N) excretion in urine (as a proportion of total fecal and urinary N excretion, UN/TN), OMd and GE had the greatest impact on N utilization efficiency (N in milk as proportion of intake N, NUE), and ED6_N had the greatest impact on milk protein yield (MPY). Additionally, CP, GE, and dr_N had the least effect on methane emission, OMd had the least effect on UN/TN, and ED6_N had the least effect on NUE. The responses of most output variables to ED6_N and dr_N variations were highly dependent on diet, and were related to the ratio between PDI (i.e., metabolizable protein) and UFL (i.e., NE L ) at the reference point of each diet. In conclusion, we were able to analyze the response of output variables to the variations of the input variables, using the tangent and its normalized value at the reference point. The predicted final outputs were more impacted by variations in CP, GE, and OMd. The other 2 input variables, ED6_N and dr_N, had a smaller effect on the final output variables, but the responses varied between the diets according to their PDI/UFL ratio. Among the final output variables affected by ED6_N, MPY was the most impacted, but when quantified this impact was at an acceptable level. Our present study was conducted using 6 representative diets for dairy cattle fed at their potential, but should be completed by the analysis of more diverse conditions.


ABSTRACT
In the INRA 2018 feeding system for ruminants, the prediction of multiple animal responses is based on the integration of the characteristics of the animal and the available feedstuff characteristics, as well as the rationing objectives.In this framework, the characterization of feedstuffs in terms of net energy, digestible protein, and fill units requirs information on their chemical composition, digestibility, and degradability.Despite the importance of these feed characteristics, a comprehensive assessment of their impact on the responses predicted by the INRA 2018 feeding system has not been carried out.Thus, our study investigated how variables predicted by the INRA feeding system (i.e., outputs) for dairy cows are affected by variation in feed characterization (i.e., inputs).Five input variables were selected for the sensitivity analysis (SA): CP, OM apparent digestibility (OMd), GE, effective degradability of nitrogen assuming a passage rate of 6%/h (ED6_N) and true intestinal digestibility (dr_N) of nitrogen.A one-at-a-time SA was performed on predicted digestive, productive and environmental output variables for dairy cows with 6 contrasted diets.These 6 diets were formulated to meet 95% of the potential daily milk production (37.5 kg) of a multiparous cow at wk 14 of lactation.Then, the values of the 5 key input variables of each feedstuff were randomly sampled around the INRA 2018 feed table values (reference point).The response of the output variable to the variation of the input variable was quantified and compared using the tangent value at the reference point and the normalized sensitivity coefficient.Among the major final output variables, CP and dr_N had the greatest impact on nitrogen (N) excretion in urine (as a proportion of total fecal and urinary N excretion, UN/TN), OMd and GE had the greatest impact on N utilization efficiency (N in milk as proportion of intake N, NUE), and ED6_N had the greatest impact on milk protein yield (MPY).Additionally, CP, GE, and dr_N had the least effect on methane emission, OMd had the least effect on UN/TN, and ED6_N had the least effect on NUE.The responses of most output variables to ED6_N and dr_N variations were highly dependent on diet, and were related to the ratio between PDI (i.e., metabolizable protein) and UFL (i.e., NE L ) at the reference point of each diet.In conclusion, we were able to analyze the response of output variables to the variations of the input variables, using the tangent and its normalized value at the reference point.The predicted final outputs were more impacted by variations in CP, GE, and OMd.The other 2 input variables, ED6_N and dr_N, had a smaller effect on the final output variables, but the responses varied between the diets according to their PDI/UFL ratio.Among the final output variables affected by ED6_N, MPY was the most impacted, but when quantified this impact was at an acceptable level.

INTRODUCTION
The French national feeding system for ruminants published in 1978 (INRA, 1978) is widely adopted in France and several other countries and is regularly updated (INRA, 1988;INRA, 2007).The scientific basis of the latest version of the INRA feeding system was published in 2018 (INRA, 2018) and is accompanied by the reference rationing software (INRAtion®V5).The INRA feeding system consists of an energy system (UF Sensitivity analysis of the INRA 2018 feeding system for ruminants by a one-at-a-time approach: Effects of dietary input variables on predictions of multiple responses of dairy cattle Seoyoung Jeon, 1 Sophie Lemosquet, 2 Anne-Cécile Toulemonde, 1 Tristan Senga Kiessé, 3 and Pierre Nozière1* system), a protein system (PDI system), and a fill unit system (UE system).The characteristics of each feedstuff and diet are determined based on these 3 systems.And for that, in addition to their chemical composition, total-tract digestibility, the estimation of ruminal degradability and true intestinal digestibility are needed.Accordingly, it is necessary to evaluate through sensitivity analysis (SA) how output variables of the INRA 2018 feeding system are affected by variation in several input variables, such as chemical composition, degradability, and digestibility values.
Sensitivity analysis of ruminant feeding systems has been already conducted.An SA of the whole nutritional system has been performed for American nutritional models or feeding systems, such as Molly (Baldwin et al., 1995), Cornell Net Carbohydrate and Protein System (Lanzas et al., 2007ab;Higgs et al., 2015), and Nutrient Requirements of Beef (NASEM, 2016) and Dairy cattle (NRC, 2001;NASEM, 2021).Lanzas et al. (2007ab) performed SA to study the effect of carbohydrate and protein fractionation in feed on microbial growth.Bateman et al. (2008) conducted SA on the effect of changes in dietary nutrient composition on milk components predicted by both Molly and NRC.
In NASEM (2016), the effect of the nutrient value of feed on ME requirement, growth, and methane emission was evaluated.NASEM (2021) also recognized that the dairy model was evaluated by an SA, but the method and results were not detailed.Binggeli et al. (2022) conducted SA to evaluate the contributions of DMI, BW, potential milk production, and dietary nutrients (NDF, ADF, CP, and starch) to MP supply and predicted milk protein yield, which are 2 variables included in the INRA 2018 system.However, SA has not been performed on the entire INRA 2018 feeding system applied to dairy cattle, and the variation in key output variables relative to the variation in input variables has not been quantified.
Therefore, we quantified to what extent the main output variables of the INRA 2018 model (including prediction of intake, milk production, nitrogen excretion and methane emission) are affected by variations in the main input variables of feedstuffs.For this purpose, among several SA methods, the one-at-a-time (OAT) method was applied as an easy approach (Groen et al., 2014) which allows evaluation of the individual effect of independent input variables by varying one at a time (Saltelli, 1999;Bar Massada and Carmel, 2008).

MATERIALS AND METHODS
A screening process was first carried out on the feed to select a limited number of key input variables among all input variables of the model, mainly to simplify the SA carried out on the diet in a second step.This screening was conducted using 22 feedstuffs which are widely used for dairy cattle rations; it was performed using OAT analysis and generalized additive models (GAMs), which can consider interactions among the input variables (Hastie and Tibshirani, 1990), unlike the OAT method.Five input variables were selected from this screening process.We varied the values of these 5 input variables of each feedstuff and conducted OAT analysis in terms of diet.Through this, we were able to analyze the relations between feed ingredients and the responses of animals to dietary factors.Abbreviations, meanings, units, and North American terms for the input and output variables included in this study are listed in Table 1.

Feedstuff: screening for selection of key input variables
To use the INRAtion®V5, more than 10 input variables, corresponding to chemical composition, ruminal, intestinal and total-tract digestibilities of feedstuffs, are needed.To avoid a high calculation cost due to an SA performed on all input variables, a previous screening procedure was carried out to select key input variables (Ge et al., 2015).We calculated the effects of all mandatory input variables on the main output variables for 13 concentrates (barley, corn, dehydrated alfalfa, 2 mixed concentrates, faba bean, soybean, beet pulp, wheat bran, corn gluten feed, corn gluten meal, rapeseed meal, and soybean meal) and 9 forages (fresh grass, fresh perennial ryegrass, fresh corn, fresh white clover, alfalfa hay, 1st and 2nd growth grass hay, grass silage, and corn silage) frequently used in France as feedstuffs for dairy cattle (Supplementary Tables 2 and  3).Output variables were the main dietary variables related to energy, protein or fill values and are detailed in each section.All calculations were performed using the PrevAlim application (included in INRAtion®V5 software), which calculates the feed value based on the INRA 2018 feeding system model.
One-at-a-time method.All input variables were considered as being independent.The OAT analysis evaluated output variables by varying each input variable at a time when fixing the other input variables at a nominal value (Iooss and Lemaître, 2015).Output variables were calculated based on random values of each input variable generated through Latin hypercube sampling (LHS, with a sample size n = 50), which used the density function of each input in terms of feedstuff.As a sensitivity index (SI), we calculated the ratio of the CV of each input simulated through LHS and the CV of the output calculated according to the OAT analysis as follows: When the CV ratio was over 50%, an input variable was assumed to be influential on the output variable considered.This cut-off threshold was considered as a compromise that enables selection of a reasonable number of influential input variables.The input variables were gross energy (GE), Ash, CP, NDF, ADF, ether extract (EE), OM apparent digestibility (OMd), effective degradability of nitrogen assuming a passage rate of 6%/h (ED6_N), and true intestinal digestibility (i.e., "digestibilité réelle") of nitrogen (dr_N) for all feedstuffs (forages and concentrates), as well as ADL, starch and the effective degradability of starch assuming a passage rate of 6%/h (ED6_St) for concentrates, and DM for forages.We assumed that all input variables have a normal distribution.Mean values for random sampling were based on the INRA feed table (INRA, 2018), with a range assumed to be ± 10% of mean value, but if the actual inter-laboratory variability was available, the observed value was used (BIPEA, 2022).
The variability of 10% used in the present study is close to the inter-laboratory variability value for a given feedstuff (BIPEA, 2022), and has been considered as reasonable.Additionally, in cases where the theoretical max value cannot exceed 100% (i.e., OMd, ED6_N, ED6_St, dr_N), the range was adjusted so that all random sample values were 99.9% or less.The output variables were GE, ED6_N, dr_N, OMd, protein truly digestible in the intestine (PDI) of alimentary (PDIA) or microbial origin (PDIM), net energy for lactation (UFL, i.e., "Unité Fourragère Lait," 1 UFL = 1.76 Mcal of NE L ), rumen protein balance (RPB), lysine and methionine truly digestible in the intestine (LysDI and MetDI, respectively) for all feedstuffs, as well as the basal fill value (bFVc) and ED6_St for concentrates, and feeding level of reference (FLref ) and fill value for sheep and dairy (UEM and UEL, i.e., "Unité d'Encombrement Moutons" and "Lait," respectively) for forages.When GE, OMd, ED6_N and dr_N are not measured, the INRA 2018 feeding system calculates these values using chemical components, so these 3 items were considered as both input and output variables.If no output variable was impacted by the input variable (i.e., CVratio = 0%), the result was not reported; this concerned only ED6_N and ED6_St for concentrates.

Generalized additive models (GAMs).
Generalized additive models describe an output variable y (response) as a function of some input variables x 0 , x 1 , x 2 , x 3 , ..., x n , (predictors) by using unknown smooth functions.Two kinds of models were fitted: 1) a model that contains all input variables given by y ~ s(x 0 ) + s(x 1 ) + s(x 2 ) + s(x 3 ) + … + s(x n ); [1] 2) a model including all input variables except the target input variable (x 0 in below equation) y ~ s(x 1 ) + s(x 2 ) + s(x 3 ) + … + s(x n ). [2] The output data set y from PrevAlim is restricted to the main output variables (i.e., PDIA, PDIM, and UFL for all feedstuffs, UEM for forages, and bFVc for concentrates); x n are the input variables of the data set obtained using the LHS method (n = 50).The input variables for GAMs included all variables and the same distribution as that used in the OAT analysis in terms of feedstuff (i.e., GE, Ash, CP, NDF, ADF, EE, OMd, ED6_N and dr_N for all feedstuffs, as well as ADL, starch and ED6_St for concentrates, and DM for forages).The GAMs consider interactions among 2 or 3 (or more) input variables by including the terms s(x i , x j ) or s(x i , x j , x k ) (and so on), respectively, in Equation 1, unlike OAT analysis.In addition, GAMs can model complex (e.g., nonlinear) relationships between output and input variables.
The difference between the R 2 of both previous models (Equations 1 and 2) was calculated to evaluate the contribution of the predictor of interest, x 0 , to the variation in the output variable.Thereby, an input variable x 0 that gives a difference higher than 0.4 between the R 2 of both models was considered to be influential on the output variable.Results concerning average differences in R 2 of less than 0.01 were not reported.This concerned Ash, NDF, EE, ADL, starch, and DM as input variables.All statistical analyses were conducted using R software (R Core Team, 2022; version 4.1.3,R Foundation for Statistical Computing, Vienna, Austria).The R package "mgcv" was used for the fitting of GAMs (Wood et al., 2022).A random sample matrix that comprises the values of all input variables was needed.Therefore, 50 samples were randomly sampled using the LHS method.Here, the randomLHS function of R was used (Carnell, 2017), assuming intra-feedstuff variability was independent between input variables, except between dr_N and ED6_N for mixed concentrate and forage (INRA, 2018).Since the intra-feedstuff relationship between dr_N and ED6_N is not linear, the relationship between both values was considered with the Spearman correlation, using the corr.testfunction (Revelle, 2023).

Diet: Quantification of the impact of input variables
In terms of diet, 5 input variables were included, which were selected based on results from the OAT and GAMs approaches in terms of feedstuff.The impact of each input variable on diet was estimated by studying variation in these 5 input variables around the reference point.These variations are described in the next section.
Reference diets for sensitivity analysis.The initial diets selected for SA were called the reference diets.They were conceived according to diets observed in the field, using INRAtion®V5 software, to satisfy intake capacity and reach an objective of 95% of the potential daily milk production (37.5 kg) by a multiparous cow at wk 14 of lactation, at 50 mo of age and 608 kg of BW.Diets were formulated using forages (corn silage, CS; 1st and 2nd growth grass hay, GH1 and GH2, respectively; grass silage, GS; fresh perennial ryegrass, RF), protein concentrates (soybean meal and rapeseed meal), and energy concentrate (barley and corn).The feed table for INRAtion®V5 used for reference diet calculation was constructed with PrevAlim as follows: for each feedstuff, the values for each mandatory input variable (Ash, CP and NDF) and non-predictable variable were the values of the INRA 2018 feed table, and the others (i.e., UEL, PDI, UFL, etc.) were calculated using the PrevAlim model.The diet composition and the predicted responses of the animal (i.e., intake, production, etc.) were calculated by INRAtion®V5, based on the theoretical basis of the INRA 2018 feeding system: the fill value of the diet satisfies the animal's intake capacity (UEL), and that UFL, PDI, MetDI, and LysDI dietary supplies allow to cover the milk yield objective.It is defined as the reference diet or reference point in this study.The composition of the reference diets is shown in Table 2.
Simulation.The feed table for INRAtion®V5 for OAT analysis was constructed with PrevAlim as described for the reference diet, except that the value for each of the 5 input variables (CP, GE, OMd, ED6_N, and dr_N) selected through the above-described (feedstuff) approach varied around the reference values (around ± 10% of mean value; refer to feedstuff).Latin hypercube sampling was performed for each of the 5 input variables, by assuming that they are independent, to construct a data set containing a total of 250 data points (5 variables × 50 replications).For each input variable combination, each feed composing the diet was sorted in the same order based on the value of the target input variable.At this time, the proportions of each feedstuff were fixed to their proportion in the reference diets.Thus, the nutritional content (i.e., GE, CP, ED_N, etc.), affecting dietary UEL, UFL, PDI, MetDI, and LysDI concentration, varied around the 'reference diet' values.The DMI and nutrient supplies (that determines the animal responses) were then calculated using INRAtion®V5 so that the fill value the diet satisfies the animal's intake capacity (UEL).
Sensitivity index.To quantify the influence of each input variable (x i ) on each output variable (y i ), the tangent (T ij ) and normalized tangent (NT ij ) called normalized sensitivity coefficient were used (Hughes et al., 2013), as follows: T ij represents the tangent value around the reference point, NT ij represents the sensitivity of the model output to input variable variation around the reference point; x 2i represents the closest x value larger than the reference value, x 1 represents the closest x value smaller than the reference value, y 2j represents the value of y when x = x 2i , y 1j represents the value of y when x = x 1i .
The input variables x i considered were CP, GE, OMd, ED6_N, or dr_N, and the output variables y j came from INRAtion®V5, including milk protein yield (MPY), energy in methane (ECH4), nitrogen utilization efficiency (N in milk as proportion of intake N, NUE), ratio between urine and total N excretion (UN/TN), global substitution rate (SRg), DMI, PDI, UFL, and PDI/UFL.For a given NT ij , the sign may vary depending on the diet.
To assess both the magnitude and the variability of this indicator, the mean and SE of NT ij has thus been calculated using the absolute value (|NT ij |, Table 6).

Feedstuff: screening for selection of key input variables
Crude protein, GE, OMd, ED6_N and dr_N showed CV ratio higher than 50% for at least one output variable in both forage and concentrate (Tables 2 and 3).In forage, ADF and NDF showed a CV ratio > 50% for UFL and RPB, through an effect on OMd.Thus, CP, GE,   OMd, ED6_N and dr_N, which were influential input variables for both concentrate and forage, appeared as key input variables.As a result of GAMs fitting, the average of R 2 differences in all feedstuffs, when input variables were GE, CP, ADF, dr_N and ED6_St, was not higher than the threshold value of 0.4 (Table 5).Organic matter digestibility was found to have a high influence on protein, energy and fill values: i.e., PDIM for both forages and concentrates, UFL and bFVc for concentrates, UEM for forages.Effective degradability of nitrogen mainly affected protein (PDIA) for concentrate feeds (0.61 ± 0.078), but did not have a significant effect on the others.Therefore, based on an R 2 difference of 0.4 with GAMs fitting, OMd and ED6_N, but not GE, CP and dr_N, appeared as key input variables.Finally, based on results of both the OAT and GAMs fitting approaches, CP, GE, OMd, ED6_N and dr_N were retained as key input variables for SA in terms of diet.
Gross energy.
; i is one input variable, j is one output variable, x 1 and x 2 are the input variable values closest to the reference value, y 1 and y 2 are the values of the output variable when x is x 1 and x 2 , respectively.).The response of UN/TN ratio did not differ across the diets (CV of T dr_N,UN/TN = 14.3%), but the other responses were more affected by diet (CV of T dr_N,MPY = 76.2;CV of T dr_N,ECH4 = −83.6%;CV of T dr_N,NUE = 77.2%; Figure 13), particularly MPY with the RF, GS, and CS diets.The responses of SRg differed depending on diet.With the RF, GS, and CS diets, SRg increased and then decreased; with the other 3 diets, SRg decreased (Figure 14; T dr_N,SRg = −5•10 −4 ± 1.5•10 −4 ).Dry matter intake increased with all diets as dr_N increased, and the response also differed depending on diet (CV of T dr_N,DMI = 77.3%),being greater with the RF, GS, and CS diets.Variation in dr_N affected PDI more than UFL (|NT dr_N,PDI | = 47.5 ± 1.20%; |NT dr_N,PDI/UFL | = 48.6 ± 0.94%; |NT dr_N,UFL | = 1.0 ± 0.36%).The responses of PDI and PDI/UFL were not affected by diet (CV of T dr_N,PDI = 12.5%; CV of T dr_N,PDI/UFL = 17.3%), whereas the low response of UFL differed slightly among the diets (CV of T dr_N,UFL = −86.9%; Figure 15).

Feedstuff: screening for selection of key input variables
The first step of SA is to analyze the contribution of each input variable to the output variable to select key input variables (Ge et al., 2015).The INRA 2018 feeding system calculates the animal's responses according to the characteristics of the diet, which are based on the characteristics of the feedstuff.Consequently, input variables with low impact at the feedstuff level could not have an impact on the diet characteristics and animal responses.Therefore, in this study, the screening of input variables was first performed at the feedstuff level.For that, in the present study, 2 types of analysis methods, OAT and GAMs, were used in terms of feedstuff to screen key input variables for SA in terms of diet.In the OAT analysis conducted using the CV ratio , it is possible to evaluate the individual impact of each input variable on each output variable.However, OAT analysis does not consider the interaction between input variables.To address this limitation, GAMs fitting was used because, unlike OAT analysis, it does not require varying one input variable at a time to evaluate output responses.Although 5 key input variables were selected based on this hypothesis, NT analysis at the diet and animal response levels was performed on all 12 input variables, to verify that no influential input variables are lost through the screening process in terms of feedstuff level.This was confirmed by the result presented in the supplementary table (Supplementary Tables 3  and 4).In the OAT analysis, CP, GE, OMd, ED6_N, and dr_N showed a CV ratio higher than 50%, while in GAMs only OMd and ED6_N showed models with a difference in R 2 higher than 0.4.Although only OMd and ED6_N were key input variables revealed by the GAM approach, it was decided to keep the 5 key input variables revealed by OAT for further SA in terms of   2 CP = crude protein; GE = gross energy; OMd = OM digestibility; ED6_n = effective degradability of nitrogen assuming a passage rate of 6%/h; dr_n = true intestinal digestibility ("digestibilité réelle") of nitrogen.
diet.The significant impact of ADF and NDF on UFL and RPB revealed by the OAT approach for forages was related to their effect on OMd.Since the retained key input variables already include OMd for both forages and concentrates, ADF and NDF were excluded from the list of input variable despite the importance of fiber in ruminant nutrition.
As the operational feed units in the INRA 2018 feeding system are metabolizable proteins (PDI), net energy (UFL), and fill value (UE), the 5 key input variables revealed by SA in terms of feedstuff appeared consistent.Indeed, PDI mainly depends on CP, OMd, ED6_N, and dr_N.And PDI consists of a part originating from diet (PDIA) and microbes (PDIM).It is known that ED6_N and dr_N mainly affect the dietary part (PDIA) and that OMd affects the microbial part (PDIM).Also, sequential losses in energy between GE, DE, ME and UFL depend highly on energy in feces (i.e., on OMd), urine and gas (i.e., on CP and OMd, respectively), and extra heat (i.e., on ME/GE).Lastly, UE is highly dependent on OMd.Sensitivity analysis revealed OMd as a common major input for PDI, UFL and UE feed value, and OMd is also the main criterion to which digestive interactions between feedstuffs are applied (INRA, 2018).c) Nitrogen utilization efficiency (N in milk, g/N intake g), and (d) Urinary nitrogen/Total excreted nitrogen (g/g) when CP (g/kg DM) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
quantify variability when all models are linear (Saltelli et al., 2005).However, as most responses relating input variables to output variables in terms of diet are nonlinear in the INRA 2018 feeding system, a different quantitative approach was needed.Therefore, in the present study, the sensitivity of the output variable to the variation per unit of the input variable was quantified as the tangent value at the reference point.In addition, by calculating the NT value, which is a normalized SI unit, it is possible to compare the respective variations of the output variables with respect to the variation of each input variable (Hughes et al., 2013).
Influence of input variables on PDI and UFL.The 5 main input variables tested in this study can be broadly divided into energy-related (UFL) and proteinrelated (PDI) variables.As a result of the SA, and as expected, OMd and GE had more influence on UFL, while CP, ED6_N, and dr_N had more influence on PDI, and especially on the dietary (PDIA) fraction (Table 6).The contribution of CP and ED6_N to the variation of PDIM was smaller and depended mainly on MCP synthesis, and thus on OMd via fermentable OM (Sauvant and Nozière, 2016).
The most affected output variables are related to nitrogen efficiency.Among the 4 main final output variables retained in the present study (i.e., MPY, UN/TN, ECH4, and NUE), UN/TN and NUE were the most highly influenced by variation in the input vari-  c) Nitrogen utilization efficiency (N in milk, g/N intake g), and (d) Urinary nitrogen/Total excreted nitrogen (g/g) when gross energy (GE, kcal/kg DM) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
ables.These 2 variables are related to N partitioning.Indeed, variations in N excreted through urine mainly refer to the excess in the rumen degradable N and in metabolizable protein (i.e., PDI), and fecal N mainly refers to the apparently undigested N. The variations in the ratio UN/TN were highly and negatively correlated with the variations in NUE (N milk / N intake, i.e., ni-trogen efficiency) when input variables were OMd, CP, or GE (r = −1.00;data not shown).However, when the input variable was dr_N, there was a strong positive linear correlation between both output variations (r = 0.95), and when the input variable was ED6_N, the correlation between both was weak (r = −0.21).This reflects the fact that for a given animal at a given level  of N intake, the digestive utilization of protein (ED6_N and dr_N) highly affects N partitioning between milk, urine, and feces.In other words, the responses of NUE and UN/TN are not always related, depending on the dietary variation factor.Thus, this emphasizes the need to assess indicators of N utilization not only for production (i.e., NUE), but also for its environmental impact (i.e., UN/TN).
The least affected output variable.Among the 4 retained final output variables, ECH4 was the least affected by the fluctuation of input variables, except by OMd.In the model (INRA, 2018), ECH4 is calculated with OM, OMd, proportion of concentrate (PCO), and feeding level (DMI %BW).In our simulation, OM and PCO were fixed values for each diet, and when considering the input variables other than OMd, the variation in dietary OMd was small, so we presume that the variation of ECH4 was small.
Milk protein yield response.In the INRA 2018 feeding system, for a given animal, the MPY is calculated from the PDI and UFL balances (supply -theoretical requirements) and their interactions, adjusted by the LysDI and MetDI contents.The response model, as well as the response model of DMI to the PDI/UFL ratio is nonlinear.In the model, MPY responds at a low level of UFL balance and reaches a plateau at a high level of UFL balance.This explains the variability of some responses to input variables between diets.While  c) Nitrogen utilization efficiency (N in milk, g/N intake g), and (d) Urinary nitrogen/Total excreted nitrogen (g/g) when OM apparent digestibility (OMd, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
MPY increased almost linearly with increasing CP with all diets, it increased with increasing GE with grass hay-based diets, but showed a parabolic curve with the other diets (Figure 4).With grass hay-based diets, it is mainly related to a decrease in the PDI/UFL ratio, without significant change in DMI, in the range where MPY positively responds to change in UFL balance, whereas with the other diets both PDI/UFL and DMI decrease with increasing GE, leading to a less predictable nutritional balance and MPY response, in the range where MPY is less responsive to change in UFL balance.tics (CP, GE, OMd, ED6_N, dr_N, PDI, UFL, PDI/ UFL, and SRg) at reference points were calculated.For ECH4 (with GE, CP, ED6_N, and dr_N as input variables), MPY (with GE, ED6_N, and dr_N), and DMI (with OMd, ED6_N, and dr_N), the correlation was the highest (average R 2 ≥ 0.72) with PDI/UFL than with other reference diet characteristics, such as GE and dr_N (average R 2 ≤ 0.13).In addition, this trend was remarkable when the input variable was dr_N or ED6_N, which may be linked to the fact that these variables have a great influence on PDI and the PDI/ UFL ratio.Moreover, PDI/UFL regulates feed intake  b) Dry matter intake (kg/d) when effective degradability of nitrogen assuming a passage rate of 6%/h (ED6_N, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
Figure 12.Responses of (a) Protein truly digestible in the intestine (g/kg DM), (b) Net energy for lactation (/kg DM), and (c) PDI/UFL (g PDI/UFL) when effective degradability of nitrogen assuming a passage rate of 6%/h (ED6_N, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
in the model (Vérité and Delaby, 2000;Faverdin et al., 2003), and variations in MPY and ECH4 are highly related to DMI.Also, responses of MPY depend on both UFL and PDI balances, through their interactions (Brun-Lafleur et al., 2010;Daniel et al., 2016).Thus, if SI values differ according to diet, this may not be related to the characteristics of each feedstuff, but rather to the dietary PDI/UFL ratio of the reference diet.

CONCLUSIONS
In conclusion, the methods developed through this study allowed us to quantify the responses of the output variables to the variations in the 5 main input variables of the INRA 2018 feeding system, in both absolute and relative values.In addition, it was possible to compare which output variable was more sensitive to variations in input variables through the NT values.The predicted final outputs were more influenced by changes in CP, GE, and OMd.In particular, OMd had a significant effect on MPY, ECH4, and NUE.Conversely, changes in ED6_N and dr_N had little effect on DMI, ECH4, and MPY, on average, but the response varied between diets, mainly according to their PDI/ UFL ratio.However, the ratio UN/TN is expected to be highly affected by variations in dr_N.In this study, the variation in animal response per unit of input was quantified as (normalized) tangent coefficients, and  c) Nitrogen utilization efficiency (N in milk, g/N intake g), and (d) Urinary nitrogen/Total excreted nitrogen (g/g) when true intestinal digestibility (dr_N, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
even if the variability of input variables is not always well known, these results can be used to establish a safety margin when formulating a diet for dairy cattle.However, although the present study included 6 representative diets typically encountered for dairy cows in France, the study of more diverse diets, animal charac-teristics, and rationing objectives is needed to provide more predictive information.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).

ACKNOWLEDGMENTS
Jeon et al.: Local sensitivity analysis of feed unit system

1
MPY = milk protein yield; ECH4 = energy in methane; NUE = nitrogen utilization efficiency (N in milk / N intake); UN/TN = N excretion ratio between urine and N excretion by urine and fecal; SRg = Substitution rate of concentrate; PDI = true digestible protein, UFL = net energy for lactation (1 UFL = 1.76 Mcal of NE L ).

Figure 1 .
Figure1.Responses of (a) Milk protein yield (g/d), (b) Energy in methane (kcal/kg DM), (c) Nitrogen utilization efficiency (N in milk, g/N intake g), and (d) Urinary nitrogen/Total excreted nitrogen (g/g) when CP (g/kg DM) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).

Figure 4 .
Figure 4. Responses of (a) Milk protein yield (g/d), (b) Energy in methane (kcal/kg DM), (c) Nitrogen utilization efficiency (N in milk, g/N intake g), and (d) Urinary nitrogen/Total excreted nitrogen (g/g) when gross energy (GE, kcal/kg DM) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
Figure 5. Responses of (a) Global substitution rate of concentrate (kg DM/kg DM) and (b) Dry matter intake (kg/d) when gross energy (GE, kcal/kg DM) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).

Figure 6 .
Figure6.Responses of (a) Protein truly digestible in the intestine (g/kg DM), (b) Net energy for lactation (/kg DM), and (c) PDI/UFL (g PDI/UFL) when gross energy (GE, kcal/kg DM) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/ input unit).

Figure 7 .
Figure7.Responses of (a) Milk protein yield (g/d), (b) Energy in methane (kcal/kg DM), (c) Nitrogen utilization efficiency (N in milk, g/N intake g), and (d) Urinary nitrogen/Total excreted nitrogen (g/g) when OM apparent digestibility (OMd, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
Figure 8. Responses of (a) Global substitution rate of concentrate (kg DM/kg DM) and (b) Dry matter intake (kg/d) when OM apparent digestibility (OMd, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).

Figure 9 .
Figure9.Responses of (a) Protein truly digestible in the intestine (g/kg DM), (b) Net energy for lactation (/kg DM), and (c) PDI/UFL (g PDI/UFL) when OM apparent digestibility (OMd, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrassbased diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
Figure 11.Responses of (a) Global substitution rate of concentrate (kg DM/kg DM) and (b) Dry matter intake (kg/d) when effective degradability of nitrogen assuming a passage rate of 6%/h (ED6_N, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
Figure13.Responses of (a) Milk protein yield (g/d), (b) Energy in methane (kcal/kg DM), (c) Nitrogen utilization efficiency (N in milk, g/N intake g), and (d) Urinary nitrogen/Total excreted nitrogen (g/g) when true intestinal digestibility (dr_N, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).
Figure15.Responses of (a) Protein truly digestible in the intestine (g/kg DM), (b) Net energy for lactation (/kg DM), and (c) PDI/UFL (g PDI/UFL) when true intestinal digestibility (dr_N, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrassbased diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).

Figure 14 .
Figure 14.Responses of (a) Global substitution rate of concentrate (kg DM/kg DM) and (b) Dry matter intake (kg/d) when true intestinal digestibility (dr_N, %) is varied.The white square with a solid green line is RF (fresh perennial ryegrass-based diet), the black square with a dashed orange line is GH1 (1st growth grass hay-based diet), the white circle with a dotted yellow line is GH2 (2nd growth grass hay-based diet), the black square with a dotted-dashed brown line is GH3 (1st and 2nd growth grass hay-based diet), the white triangle with a long-dashed blue line is GS (grass silage-based diet), and the black triangle with a 2-dashed purple line is CS (corn silage-based diet).Each mark indicates the reference point for each diet, and the accompanying numbers are the tangent values at each point (output unit/input unit).

Local sensitivity analysis of feed unit system Table 1. Glossary of inputs and outputs used in sensitivity analysis
Jeon et al.:

Table 2 .
Diet formulation and nutritive values for multiparous cows 1 12nd parity; body weight 608 kg; 50 mo of age; wk 14 of lactation; 10,000 kg of milk yield per year.

Table 3 .
Jeon et al.: Local sensitivity analysis of feed unit system Average of CV ratio 1 across all concentrate feedstuffs 2 Barley, corn, dehydrated alfalfa, 2 mixed concentrate, faba bean, soybean, beet pulp, wheat bran, corn gluten feed, corn gluten meal, rapeseed meal, and soybean meal.

Table 6 .
Changes in output variables predicted by INRAtion®V5 according to changes in input variables, expressed as a normalized tangent

Table 7 .
Variation of the response of output to input variabilities around the reference situation, expressed as the coefficient of variation (CV; standard deviation / mean value of tangent, in %)