# A new method to estimate residual feed intake in dairy cattle using time series data

Pauline Martin, Vincent Ducrocq, D.G.M. Gordo, Nicolas N.C. Friggens

► To cite this version:

Pauline Martin, Vincent Ducrocq, D.G.M. Gordo, Nicolas N.C. Friggens. A new method to estimate residual feed intake in dairy cattle using time series data. Animal, 2021, pp.100101. 10.1016/j.animal.2020.100101. hal-03064668v2

HAL Id: hal-03064668

https://hal.inrae.fr/hal-03064668v2

Submitted on 6 Jun 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A new method to estimate residual feed intake in dairy cattle using time series data

<figure><img src="image_2.png" /></figure>

Check for
updates

P. Martin<sup>a,*</sup>, V. Ducrocq<sup>a</sup>, D.G.M. Gordo<sup>b</sup>, N.C. Friggens<sup>c</sup>

## ARTICLE INFO

**Article history:**

Received 1 January 2020

Received in revised form 23 September 2020

Accepted 24 September 2020

Available online 10 December 2020

**Keywords:**

Correlations

Dynamic approach

Feed efficiency

Lactation parameters

Multitrait random regression

## ABSTRACT

In dairy, the usual way to measure feed efficiency is through the residual feed intake (RFI) method. However, this method is, in its classical form, a linear regression, which, by construction, does not take into account the evolution of the RFI components across time, inducing approximations in the results. We present here a new approach that incorporates the dynamic dimension of the data. Using a multitrait random regression model, the correlations between milk, live weight, DM intake (DMI) and body condition score (BCS) were investigated across the lactation. In addition, at each time point, by a matrix regression on the variance-covariance matrix and on the animal effects from the three predictor traits, a predicted animal effect for intake was estimated, which, by difference with the actual animal effect for intake, gave a RFI estimation. This model was tested on historical data from the Aarhus University experimental farm (1 469 lactations out of 740 cows). Correlations between animal effects were positive and high for milk and DMI and for weight and DMI, with a maximum mid-lactation, stable across time at around 0.4 for weight and BCS, and slowly decreasing along the lactation for milk and weight, DMI and BCS, and milk and BCS. At the Legendre polynomial coefficient scale, the correlations were estimated with a high accuracy (averaged SE of 0.04, min = 0.02, max = 0.05). The predicted animal effect for intake was always extremely highly correlated with the milk production and highly correlated with BW for the most part of the lactation, but only slightly correlated with BCS, with the correlation becoming negative in the second half of the lactation. The estimated RFI possessed all the characteristics of a classical RFI, with a mean at zero at each time point and a phenotypic independence from its predictors. The correlation between the averaged RFI over the lactation and RFI at each time point was always positive and above 0.5, and maximum mid-lactation (>0.9). The model performed reasonably well in the presence of missing data. This approach allows a dynamic estimation of the traits, free from all time-related issues inherent to the traditional RFI methodology, and can easily be adapted and used in a genetic or genomic selection context.

© 2020 The Authors. Published by Elsevier Inc. on behalf of The Animal Consortium. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/).

## Implications

This paper proposes a new methodology to estimate feed efficiency in a continuous manner across the lactation. The efficiency is evaluated from a residual feed intake derived from a multitrait random regression model, which allows the coefficients to vary over time. This approach allows a dynamic estimation of the traits and their correlations, free from all time-related issues inherent to the traditional methodology, and can easily be adapted and used in a genetic or genomic selection context.

## Introduction

With feed costs representing above 50% of the total costs of dairy production (European Commission, 2018), the issue of feed efficiency

has become a priority for the sector. The notion of feed efficiency refers to improving the balance between output (production) and input (feed intake). The most common way to determine feed efficiency in dairy cattle is through residual feed intake (RFI). First proposed by Koch et al. (1963), RFI is the difference between the actual feed intake of an animal and its predicted feed intake based on its performance, i.e. the intake necessary to cover the demands of the different energy sinks, estimated by regression. By construction, this method benefits from the RFI being phenotypically independent (or genetically in case of genetic regression) from its predictors, which theoretically allows RFI to reflect digestive and metabolic variabilities (Archer et al., 2002; Berry and Crowley, 2013).

Although this definition of RFI is widely used in dairy cattle (Connor, 2015), some issues remain with this approach. First of all, the time influence on RFI is a key question. RFI is usually measured between two given days. If the trial duration is too short, the number of measures will be low and the results of the prediction will not be accurate. On

<footnote><sup>a</sup> UMR GABI, INRAE, AgroParisTech, Université Paris-Saclay, 78350 Jouy-en-Josas, France</footnote>

<footnote><sup>b</sup> Center for Quantitative Genetics and Genomics, Department of Molecular Biology and Genetics, Aarhus University, DK-8830 Tjele, Denmark</footnote>

<footnote><sup>c</sup> UMR MoSAR, INRAE, AgroParisTech, Université Paris-Saclay, 75005 Paris, France</footnote>

<footnote>* Corresponding author.</footnote>

<footnote>E-mail address: pauline.martin@inra.fr (P. Martin).</footnote>
the other hand, different biological processes are involved at the different lactation stages, in particular body reserve mobilization in early lactation and reserves accretion associated with pregnancy. Therefore, regression coefficients associated with the different predictors are likely to vary when taken at different lactation stages, leading to a biased assessment of RFI when estimated from point measures over the whole lactation (Li et al., 2017). Moreover, with the advent of precision farming technologies, high-frequency time series measures are becoming available and thereby increasingly offer the opportunity to track efficiency. Key questions in this context are how to deal with changes over time, when to measure RFI and for how long. A second issue comes from the nature of RFI itself. Being a residual, in addition to the actual animal efficiency, it contains all the modeling and measurement errors. Within this context, Fischer et al. (2018) used random regression on the individual level to isolate the cow-specific part of RFI. However, they faced some difficulties due to strong correlations between their predictors (correlations that changed over time during the lactation), which, combined with a limited data set size, restricted their modeling and conclusions.

With the aim of solving these two issues, we investigate in this paper the relationship between intake and its predictors across the lactation using a multitrait random regression model on experimental data. With such a model, based on the variance-covariance functions, it is possible to predict intake from other traits and to compare this prediction with actual intake. This leads to a consistent definition of RFI, with regression coefficients free to vary over time.

## Material and methods

### Population resources and feeding management

Data were collected between 2002 and 2016 at the Danish Cattle Research Centre (Foulum, Denmark) and shared as part of the GenTORE project (https://www.gentore.eu). Animals were Holstein cows with lactation rank ranging from first to third lactation. Cows were housed in a freestall barn with cubicles and slatted floor and milked in an automatic milking system (AMS; DeLaval, Tumba, Sweden) allowing free cow traffic. Data were collected during various trial periods, and cows were not nested within trial. They had ad libitum access to a partial mixed ration (PMR) varying in nutritional content in accordance with the particular trial in which the cows were involved. The PMR diets used during this period contained the following ingredients (typical g/kg DM): rapeseed meal (106), barley (135), dried sugarbeet pulp (51), grass silage (264), maize silage (422), urea (7) and minerals (18). The corresponding chemical composition was (typical g/kg DM): CP (158), NDF (340), starch (195), sugar (46), crude fat (39), ash (74), with a Net Energy of Lactation (NEL) of 6.9 MJ/kg DM. Across the different feeding trials carried out during this period, diet differences did not exceed +8% of the average value for CP and likewise +15% for NEL. The nutritional values of the PMR were all within the range of typically recommended requirements, formulated to support the milk yield level of the herd and allocated in amounts allowing approximately 10% orts to ensure ad libitum intake. In addition, cows were supplemented daily with a maximum of 3 kg of concentrate to ensure voluntary access to the AMS. These feeding data were already described in previous studies (e.g. Li et al., 2016; Byskov et al., 2017).

### Phenotypes and data editing

Data from 1 469 lactations of 740 cows were collected. A weekly measurement of average daily milk yield per cow was obtained from the average of daily milk yield records per cow in each week. Milk samples were taken weekly for analyses of fat and protein. On the same basis, the DM contents in PMR and concentrates were analyzed regularly and the compositions were aligned and merged with feed intake

records to obtain weekly DM intake (DMI) values for individual cows. Animals were also automatically weighed at each milking so that BW records were averaged to obtain a weekly record of BW per cow in each week. Body condition score (BCS) was evaluated every 2 weeks and scored on a scale from 1 to 5.

A corrected milk (cmilk) trait was created following the FAO formula which defines as a standard milk with 4.0% fat and 3.3% protein (FAO [Food and Agriculture Organization of the United Nations], 2010):

$$ \text{Corrected milk (kg)} = \text{raw milk (kg)} * (0.337 + 0.116 * \text{Fat content (\%)}) + 0.06 * \text{Protein content (\%)} $$

To avoid non-sensical performances, a filter was used to discard records that differed too much from the previous record registered in the same lactation for that animal. Therefore, records implying differences between two consecutive records higher than 12 l for cmilk, 50 kg for the liveweight, 7 kg for the daily DM intake and one unit of BCS were discarded. These threshold values were defined after studying a previous data set with similar performance level. Less than 2% of the data was discarded for cmilk and BCS and about 4% for weight and DMI. A second step of filtering was made on the duration of data collection for each lactation independently: only animals that had been recorded for a minimum duration of 200 days during their lactation were kept. This data set, that was used to conduct the main analyses, contains 40 619 records of cmilk, 40 662 records of weight, 42 177 records of DMI and 19 661 records of BCS, which represents about 73% of the overall records. Additionally, the overall records were used to assess the robustness of our model in case of missing data, and the details of this analysis are presented in Supplementary Material S1.

### Statistical analysis

#### Multitrait random regression model

To analyze the relationship of the four traits along the lactation, a multitrait random regression analysis was performed using the Wombat software (Meyer, 2007). The model used was as follows:

$$ y_{ilmr} = c_i + \sum_{n=0}^{3} \beta_{ln} \varphi_{nr}(t) + \sum_{n=0}^{2} \alpha_{mn} \varphi_{nr}(t) + e_{ilmr} $$

where $y_{ilmr}$ is the observation of the rth trait, $c_i$ is the fixed effect of the ith month-year combination corresponding to the record date, $\beta_{ln}$ is the nth fixed regression coefficient specific to parity class l, $\alpha_{mn}$ is the nth random regression coefficient of the animal m, $\varphi_{nr}(t)$ is the nth coefficient of Legendre polynomial of degree ($d=2$ for the animal effect and $d=3$ for the parity class), evaluated at day in milk (DIM) t, and $e_{ilmr}$ is the random residual effect. Residual effects were assumed to have a homogenous residual variance. Successive lactations of the same animal were considered separately, as if they were from different animals, implicitly ignoring any permanent environmental effect across lactations.

From the random part of the equation, we get as outputs the three coefficients of the Legendre polynomial at each time point, as well as the animal solutions. For each trait, we can then obtain at each time and for each cow an estimated animal effect. These animal effects will be denoted $a_{cmilk}(t)$, $a_{weight}(t)$, $a_{DMI}(t)$ and $a_{BCS}(t)$ in the article.

Variance-covariance components of the random regression were obtained and gathered into a 12 by 12 matrix, corresponding to the four traits with three Legendre coefficients each. Variances ($V_{trait1}(t)$) and covariances ($C_{trait1, trait2}(t)$) for each time point were estimated by pre- and post-multiplying the variance-covariance matrix by the corresponding time coefficients of the Legendre polynomials using the following formulas:
$$V_{trait1}(t) = [l_0(t) \quad l_1(t) \quad l_2(t)] \begin{bmatrix} \sigma_{a_0 trait1}^2 & \sigma_{a_0 a_1 trait1} & \sigma_{a_0 a_2 trait1} \\ \sigma_{a_0 a_1 trait1} & \sigma_{a_1 trait1}^2 & \sigma_{a_1 a_2 trait1} \\ \sigma_{a_0 a_2 trait1} & \sigma_{a_1 a_2 trait1} & \sigma_{a_2 trait1}^2 \end{bmatrix} \begin{bmatrix} l_0(t) \\ l_1(t) \\ l_2(t) \end{bmatrix}$$

$$C_{trait1,trait2}(t) = [l_0(t) \quad l_1(t) \quad l_2(t)]$$

$$\begin{bmatrix} \sigma_{a_0 a_0 trait1, trait2} & \sigma_{a_0 a_1 trait1, trait2} & \sigma_{a_0 a_2 trait1, trait2} \\ \sigma_{a_0 a_1 trait1, trait2} & \sigma_{a_1 trait1, trait2} & \sigma_{a_1 a_2 trait1, trait2} \\ \sigma_{a_0 a_2 trait1, trait2} & \sigma_{a_1 a_2 trait1, trait2} & \sigma_{a_2 a_2 trait1, trait2} \end{bmatrix} \begin{bmatrix} l_0(t) \\ l_1(t) \\ l_2(t) \end{bmatrix}$$

These variance-covariance components were then used to calculate correlations between traits across the lactation.

### Derivation of residual feed intake

Using elements from the 12 by 12 variance-covariance matrix and the animal effects for cmilk, weight and BCS, it is possible to estimate a predicted intake from the other variables by a matrix regression. If we use the following notation for the variance-covariance matrix at a given time:

$$\begin{bmatrix} V_{\text{DMI}} & C_{\text{DMI,cmilk}} & C_{\text{DMI,weight}} & C_{\text{DMI,BCS}} \\ C_{\text{DMI,cmilk}} & V_{\text{cmilk}} & C_{\text{cmilk,weight}} & C_{\text{cmilk,BCS}} \\ C_{\text{DMI,weight}} & C_{\text{cmilk,weight}} & V_{\text{weight}} & C_{\text{weight,BCS}} \\ C_{\text{DMI,BCS}} & C_{\text{cmilk,BCS}} & C_{\text{weight,BCS}} & V_{\text{BCS}} \end{bmatrix} = \begin{bmatrix} \mathbf{B}_{11} & \mathbf{B}_{21} \\ \mathbf{B}_{12} & \mathbf{B}_{22} \end{bmatrix}$$

then a predicted animal effect for DMI can be calculated from a regression similarly to what was proposed by Kennedy et al. (1993), but in its matrix form with:

$$\hat{a}_{e\text{DMI}}^{*}(t) = \mathbf{B}_{12}(t) * \mathbf{B}_{22}^{-1}(t) * \begin{bmatrix} a_{\text{cmilk}}(t) \\ a_{\text{weight}}(t) \\ a_{\text{BCS}}(t) \end{bmatrix}$$

where $\hat{a}_{e\text{DMI}}^*$ is the predicted animal effect for DMI at time t. Finally, we can construct a RFI estimate which is the difference between the actual animal effect for DMI and the one predicted from the three other variables:

$$\text{RFI}(t) = a_{\text{DMI}}(t) - \hat{a}_{e\text{DMI}}^*(t)$$

As this RFI is dependent on time, we also defined RFI<sub>tot</sub> as the averaged RFI of each animal over the whole lactation.

### Use of the outputs to study changes between times

One could also consider that changes of BCS are more important than BCS itself in the estimation of RFI. To examine this point, we also estimated RFI<sub>ΔBCS</sub>(t) where all the occurrences of BCS-related variables are replaced by ΔBCS-related variables. ΔBCS(t) is defined as the difference of BCS between two consecutive time points. In the above formulas, it translates as a difference of Legendre polynomial coefficients as follows for the animal effects, the variances and covariances, respectively:

$$a_{\Delta\text{BCS}}(t) = [l_0(t) - l_0(t_{-1})] * \hat{a}_{c_0\_BCS} + [l_1(t) - l_1(t_{-1})] * \hat{a}_{c_1\_BCS} + [l_2(t) - l_2(t_{-1})] * \hat{a}_{c_2\_BCS}$$

$$V_{\Delta\text{BCS}}(t) = [l_0(t) - l_0(t_{-1}) \quad l_1(t) - l_1(t_{-1}) \quad l_2(t) - l_2(t_{-1})]$$

$$\begin{bmatrix} \sigma_{a_0 BCS}^2 & \sigma_{a_0 a_1 BCS} & \sigma_{a_0 a_2 BCS} \\ \sigma_{a_0 a_1 BCS} & \sigma_{a_1 BCS}^2 & \sigma_{a_1 a_2 BCS} \\ \sigma_{a_0 a_2 BCS} & \sigma_{a_1 a_2 BCS} & \sigma_{a_2 BCS}^2 \end{bmatrix} \begin{bmatrix} l_0(t) - l_0(t_{-1}) \\ l_1(t) - l_1(t_{-1}) \\ l_2(t) - l_2(t_{-1}) \end{bmatrix}$$

$$C_{\Delta\text{BCS},trait2}(t) = [l_0(t) \quad l_1(t) \quad l_2(t)]$$

$$\begin{bmatrix} \sigma_{a_0 a_0 \Delta\text{BCS}, trait2} & \sigma_{a_0 a_1 \Delta\text{BCS}, trait2} & \sigma_{a_0 a_2 \Delta\text{BCS}, trait2} \\ \sigma_{a_1 a_0 \Delta\text{BCS}, trait2} & \sigma_{a_1 a_1 \Delta\text{BCS}, trait2} & \sigma_{a_1 a_2 \Delta\text{BCS}, trait2} \\ \sigma_{a_2 a_0 \Delta\text{BCS}, trait2} & \sigma_{a_2 a_1 \Delta\text{BCS}, trait2} & \sigma_{a_2 a_2 \Delta\text{BCS}, trait2} \end{bmatrix} \begin{bmatrix} [l_0(t) - l_0(t_{-1})] \\ [l_1(t) - l_1(t_{-1})] \\ [l_2(t) - l_2(t_{-1})] \end{bmatrix}$$

And therefore, the predicted animal effect for DMI is calculated from:

$$\hat{a}_{e\text{DMI}}^{*}(t) = \mathbf{B}_{12}(t) * \mathbf{B}_{22}^{-1}(t) * \begin{bmatrix} a_{\text{cmilk}}(t) \\ a_{\text{weight}}(t) \\ a_{\Delta\text{BCS}}(t) \end{bmatrix}$$

where the **B** matrix contains the ΔBCS-related variances and covariances.

The same transformation was also performed on changes of BW.

## Results

### Descriptive statistics

Descriptive statistics, averaged over the lactation, are reported in Table 1 for the raw data of the four estimated traits, the four animal effects and the RFI. The four traits are variable, with a large range of values. Some values are especially low (e.g. 3.06 kg/d of cmilk or 6.02 kg/d of DM), but they are unlikely measurement errors because, first, data were averaged on a weekly basis and, second, they are relatively close to the previous value for the same animal (due to the filtering step). These extreme values are possibly due to health issues.

The second half of the table presents the animal effects for the four traits estimated from the model and RFI. The animal effect corresponds to the deviation of this animal from the curve fitted to describe the time course of the performance of an average animal for a given trait × parity group, after correction for the other fixed effects. The fixed effects curves are presented in Supplementary Figure S1. For the four traits, the averages of the animal effects obtained from the model are exactly zero, at each time and overall. The statistics presented in Table 1 were obtained only from animal effects for which corresponding raw data existed. By construction, the model can extrapolate and estimate animal effects for each animal at time points beyond the range of data for that animal, even though there was no corresponding raw data. As these extrapolations

<caption>**Table 1**</caption>

Means, SD, minima (Min.), maxima (Max.) and 5% and 95% centiles for the raw data of the four considered traits of dairy cattle, animal effects ($a_{\text{trait}}$) and RFI.

<table><thead><tr><th>Type of data</th><th>Trait</th><th>Mean</th><th>SD</th><th>Min.</th><th>Max.</th><th>5% centile</th><th>95% centile</th></tr></thead><tbody><tr><td rowspan="4">Raw data</td><td>Milk (kg of corrected milk)</td><td>33.4</td><td>7.8</td><td>3.1</td><td>68.9</td><td>21.1</td><td>46.9</td></tr><tr><td>Weight (10 kg)<sup>1</sup></td><td>64.3</td><td>7.5</td><td>39.9</td><td>98.3</td><td>52.7</td><td>77.2</td></tr><tr><td>DMI (kg)</td><td>21.7</td><td>3.4</td><td>6.0</td><td>38.7</td><td>16.3</td><td>27.4</td></tr><tr><td>BCS (scale of 10 to 50)<sup>1</sup></td><td>31.3</td><td>3.40</td><td>15.0</td><td>47.5</td><td>25.0</td><td>37.5</td></tr><tr><td rowspan="5">Computed data</td><td>a<sub>cmilk</sub></td><td>0.20</td><td>5.89</td><td>-28.61</td><td>24.05</td><td>-9.12</td><td>10.27</td></tr><tr><td>a<sub>weight</sub></td><td>-0.11</td><td>6.75</td><td>-22.02</td><td>29.91</td><td>-10.27</td><td>11.65</td></tr><tr><td>a<sub>eDMI</sub></td><td>0.08</td><td>2.74</td><td>-14.72</td><td>12.61</td><td>-4.16</td><td>4.69</td></tr><tr><td>a<sub>BCS</sub></td><td>-0.03</td><td>2.43</td><td>-10.88</td><td>13.48</td><td>-3.78</td><td>3.88</td></tr><tr><td>RFI</td><td>0.04</td><td>1.37</td><td>-18.26</td><td>7.16</td><td>-1.99</td><td>2.21</td></tr></tbody></table>

RFI: Residual Feed Intake.

DMI: DM Intake.

BCS: Body Condition Score.

<footnote n="1">¹ The original scale was divided (for Weight) or multiplied (for BCS) by 10 in order to have raw data in the multiple trait analysis with similar orders of magnitude for the 4 traits.</footnote>
were discarded to compute these statistics, it is expected that the means in Table 1 are not exactly zero. The exact same phenomenon occurs for the RFI estimates, with a mean of 0.04.

## Correlations between animal effects across time

The correlations between the animal effects for the four traits across the lactation are presented in Fig. 1. The curves of correlation between DMI and cmilk and between DMI and weight showed an increase in early lactation and a decrease in late lactation. The correlations between DMI and cmilk were always higher by 0.1 or 0.2. However, the correlation between cmilk and weight steadily decreased along the lactation up to a point where it became slightly negative at the extreme end of lactation. Animals that produce more milk than average in early lactation are heavier than average, benefiting from a higher intake capacity, but these higher producing animals are less persistent. The correlation between cmilk and BCS was negative throughout the lactation, meaning that those who produce more milk always have a lower BCS than average, and more so in late lactation (they may lose more in early lactation and then have difficulty to rebuild body reserves at the end). The correlation between weight and BCS remained stable over time at around 0.4, reflecting the weight equivalent of a unit BCS (animals having more reserves are heavier at every time point). Finally, the correlation between DMI and BCS started around zero and slowly decreased down to -0.2 in late lactation. At the Legendre polynomial coefficient scale, the correlations are estimated with a high accuracy (averaged SE of 0.04, min = 0.02, max = 0.05), which correspond to variance estimation errors (averaged on the lactation) of below 5% for milk, 1% for BW, 7% for intake and 4% for BCS.

Correlations between animal effects at different time points within each trait are presented in Supplementary Figure S2. The pattern is the same for the four traits with only positive correlations and with higher correlations between animal effects that are close in time. Correlations were always high for weight (> 0.8), which is a cumulative trait, while they were a little more variable for the three other traits (correlations > 0.4). The largest changes were mostly at the beginning or end of the lactation, and the middle part was very stable (correlations > 0.8). These correlations are very accurate with sampling errors being below 4% for intake, 2% for milk and BCS and 1% for BW.

## Correlations between the predicted animal effect of DM Intake and the other traits

Correlations between the predicted animal effect for intake ($\hat{a}_{eDMI}^*$) and the four original animal effects were calculated at each time and are presented in Fig. 2. Among the predictors, $\hat{a}_{eDMI}$ was very highly correlated

<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Time trends of the correlations between the animal effects (i.e. the animal differences from the average at any given time point, across the lactation of dairy cattle) of the corrected milk (cmilk), the body weight, the DM intake (DMI) and the body condition score (BCS).</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Correlations between the predicted animal effect for intake ($\hat{a}_{eDMI}^*$) and the four original animal effects (corrected milk (cmilk), weight, DM intake (DMI) and body condition score (BCS)) in the considered dairy herd.</figcaption></figure>

with the animal effect for cmilk, the correlation being always above 0.8 and almost 1 at the beginning of lactation. The second largest correlations were with the animal effect for weight, which starts at 0.7 and slowly decreases down to 0.4 along the lactation. Finally, the correlation with the animal effect for BCS was always very low, starting positive but below 0.2 in early lactation and decreasing down to -0.4 at the end of lactation. The correlation between $\hat{a}_{eDMI}^*$ and $a_{DMI}$ was between 0.6 and 0.8 in early lactation and then remained above 0.8 after 50 days of lactation.

Correlations of the animal effects for all the traits with RFI were also calculated. As expected, by construction, correlations between RFI and the three predictors, as well as with the predicted intake, were zero all along the lactation. The correlation between RFI and $a_{DMI}$ started at almost 0.8 and decreased down to 0.5 during the first to third of lactation and then remained stable until 250 DIM where it slightly increased to 0.6 at the end of lactation.

## Evolution of residual feed intake across the lactation

Various individual profiles of RFI across lactation were observed, with some animals being efficient/non-efficient during the whole lactation, while others were efficient in early lactation and non-efficient at the end or the opposite. In order to evaluate the possibility to measure only a part of the lactation and predict the overall RFI, we calculated correlations between RFI at each time and the RFI averaged over the entire lactation (RFI<sub>tot</sub>). Results are presented in Fig. 3. The correlation starts at 0.45, increases up to 0.9 at 92 DIM where it remains stable until 197 DIM, decreasing down to 0.6 afterwards. These results suggests that RFI in the middle of lactation is a good predictor of the mean RFI over the whole lactation.

<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Correlation between the averaged residual feed intake (RFI) over the whole lactation and the RFI calculated at any given time point in the considered dairy herd.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Fig. 4.** Correlation between the residual feed intake (RFI) calculated using the trait as a predictor and the RFI calculated using changes of this trait as a predictor, for both weight and body condition score (BCS) in the considered dairy herd.</figcaption></figure>

### Using changes rather than the trait itself

Correlations over time points between RFI and $RFI_{ABCS}$ and between RFI and $RFI_{ABW}$ are presented in Fig. 4. Between RFI and $RFI_{ABCS}$, the correlations were high ($>0.87$), reaching 0.99 both at the beginning and the end of the lactation, indicating that RFI and $RFI_{ABCS}$ are almost the same trait. Between RFI and $RFI_{ABW}$, the correlations are also always high ($>0.7$), especially in early lactation. This suggests that changes are already taken into account by the use of the three Legendre polynomial coefficients (the changes being the derivatives of the trait). The dramatically high correlation at the beginning (and end of the lactation for BCS) could also indicate that, during these lactation stages, changes are more important than the trait itself, as we do not lose information when considering only changes as predictor. However, this is not true in the middle of lactation where changes are less important and where the traits bring additional values per se.

## Discussion

### Contribution of the proposed model to the current methodology

In this study, we proposed a modeling strategy which allows the integration of time series data to compute RFI, not only as repeated differences between the end and the beginning of arbitrary periods, but also by considering all traits dynamically. The traditional way of defining RFI – by linear regression – is not appropriate to jointly describe a dairy cow's non-linear trajectories of relevant traits during her entire lactation (Li et al., 2017). Indeed, the allocation of energy varies among the different energy sinks across the lactation, leading to changes in regression coefficients over time. As a result, feed efficiency in dairy cows is usually computed over short time periods (Prendiville et al., 2011; Fischer et al., 2018). This implies that results are based on a small number of data points and therefore are highly sensitive to measurement errors or one-time events (e.g., mastitis). In addition, even on short time periods, RFI estimates may still be subject to approximations due to fixed regression coefficients not being able to reflect biological changes, for instance a switch between loss and gain of body reserves. With the approach presented here, the number of measurements included in the system limits its sensitivity to one-time errors or disturbance. Moreover, it allows a permanent re-adjustment of the relevant coefficients over time, leading to estimations of RFI that avoid the biases inherent to the traditional RFI estimation methodology. As the model accounts for changes in variance for each trait over time, it also enables the correlations between traits to change over time, giving us a better overview of the relationships between the four traits across lactation.

Various studies already tried to explore further the methodology of RFI or the relationships between its components. For example, a

multitrait random regression model was previously developed by Manzanilla Pech et al. (2014) to investigate the relationship between DMI, milk and live weight over time but RFI was not estimated. Lu et al. (2015) used a modified Cholesky decomposition from a multitrait linear model that allowed greater accuracy in genetic merit prediction in case of partially missing data. Strathe et al. (2014) proposed a RFI derived from a bivariate random regression model for BW and cumulative feed intake in pigs. More recently, Islam et al. (2020) used a Bayesian multivariate random regression to analyze DMI, energy-corrected milk, BW and BCS and derived a genetic RFI from it. The approach presented here is an additional step on the way to improving the modeling of feed efficiency.

### Evolution of the correlations between traits across the lactation

The approach used here gives us an overview of the evolution of correlations between traits over the lactation. Such results are still rarely available in the literature where most studies calculated correlations on specific time points, on a fixed period of the lactation, and not in a continuous, dynamic way. The works of Veerkamp and Thompson (1999), Spurlock et al. (2012), Liinamo et al. (2012) and Manzanilla Pech et al. (2014) are some of the rare exceptions. In addition, most of the studies describing correlations include pedigree information and therefore are able to split the animal effect between a genetic component and a permanent environmental effect, which was not our case. Even though the animal effects studied here are not exactly the same as additive genetic effects, some similar patterns with genetic correlations reported in the literature can be observed.

Similarly to what was found here, Manzanilla Pech et al. (2014) reported that genetic correlations within traits across the lactation were generally positive and maximum during mid-lactation and that correlations taken between times further apart were smaller. The very high positive correlations were also in accordance with previous studies (Koenen and Veerkamp, 1998; Veerkamp and Thompson, 1999; and Liinamo et al., 2012 for live weight; Tetens et al., 2014 for DMI; and Veerkamp and Thompson, 1999 and Hüttmann et al., 2009 for milk production). However, Manzanilla Pech et al. (2014) reported slightly negative genetic correlations (−0.2) between milk in early lactation and milk during the rest of the lactation and similarly for DMI, while our correlations between animal effects were always positive.

Correlations observed here for animal effects of the four traits at the same time point were also in general accordance with the literature. High positive genetic correlations were reported by Veerkamp and Thompson (1999) and Spurlock et al. (2012) between DMI and weight. Hüttmann et al. (2009) found that the genetic correlation between these two traits was changing over time with an almost null correlation between 31 and 60 DIM and a correlation of 0.4 between 121 and 180 DIM, while Manzanilla Pech et al. (2014) reported the correlation to be maximum at 34 DIM (0.56) and minimum at 153 DIM (0.29). Our results are showing a pattern similar to what Hüttmann et al. (2009) found, but our correlations are always higher and closer to the values proposed by the other studies. In the literature, general positive genetic correlations (rg) were also reported between milk and DMI (van Elzakker and van Arendonk, 1993; rg = 0.46; Veerkamp and Brotherstone, 1997; rg = 0.34; Vallimont et al., 2010; rg = 0.52). However, Manzanilla Pech et al. (2014) found a negative correlation between milk and DMI in early lactation that become positive and high ($>0.7$) only after 38 DIM and reached a maximum at 195 DIM (rg = 0.91). Our results are essentially different at the onset of the lactation and similar after. It may be due to differences in animal management around calving or to genetic differences between the data sets. Our evolution of the correlation between milk and weight over time is also very different to what was reported by Manzanilla Pech et al. (2014). These authors found a slightly negative correlation (−0.1) both at the beginning and the end of lactation and a positive correlation of about 0.3 in
mid-lactation. In contrast, our results were similar to what those re-
ported by Karacaören et al. (2006) with a correlation decreasing with
time and becoming negative at the end of lactation.

Body condition score is usually not among the traits considered in these kind of studies, and therefore, estimated correlations are scarce. It is expected that cows mobilize body reserves in early lactation, while being in a situation of negative energy balance (Tamminga et al., 1997; Grummer, 2007). In the present study, the negative correlation between animal effects on milk and BCS is in accordance with this as higher than average milk production is associated with lower than average BCS. The stable correlation between weight and BCS supports the suggestion made in previous studies that live weight change can be a good indicator of body reserve mobilization (Thorup et al., 2013; Manzanilla Pech et al., 2014).

Use of the residual feed intake

With the methodology proposed in this paper, we obtain a RFI
changing through time for every single animal. This affords new possi-
bilities. First, we identified that the averaged RFI over the lactation is
highly correlated with RFI measured in mid-lactation, in accordance
with Prendiville et al. (2011) and Connor et al. (2012). If the objective
is to identify animals which are on average the most efficient, then the
costly measure of intake could be done only during a few weeks during
mid-lactation (ideally between 115 and 175 dim based on our results)
with only a small loss of information, rather than during the entire
lactation.

However, the question of what type of animal is desirable to se-
lect for is still valid. Efficient animals on average are also the ones
that are more efficient in mid-lactation because this is the longest
stable period of the lactation without dramatic changes. But be-
cause energy sinks and their relative importance are changing
over the lactation, efficient animals in mid lactation are not neces-
sarily the ones which are the most efficient in early or late lactation.
It clearly appears that during the lactation, the first 5 to 7 weeks are
the most challenging period for the cows, during which they have to
face a huge increase of milk production associated with a pro-
nounced mobilization of their body reserves. This pronounced neg-
ative energy balance increases the risk of health issues and reduces
fertility (Esposito et al., 2014), and it is particularly important not to
increase the occurrence of health and fertility issues. Therefore, it is
critical to make sure that animals we are selecting for do not have a
deep body reserve mobilization in early lactation. Because BCS
(and/or its changes) was included in the RFI model, the two traits
are phenotypically independent but this may not be the case genet-
ically. In addition, the influence of some health events may have
been discarded in the analyses with the data filtering or because an-
imals with dramatic health or reproductive issues were culled.
Therefore, the importance of BCS in the overall efficiency (including
from an economic point of view) may have been underestimated in
the model. The relationship between efficiency and resilience or ro-
bustness needs to be further investigated in order to make better
choices in selection.

Furthermore, we now get individual trajectories of RFI which are
highly variable from one animal to another. An obvious next step
would be to perform a cluster analysis to see if we could identify specific
types of animals. This could help scientists and breeding companies to
determine what is the most suitable type of animals overall.

Other issues and future improvement

If this approach and the associated results are very promising,
numerous questions are still to be explored. First, the adaptability
of the model in more general situations needs to be tested. For in-
stance, the model should be able to perform well when mixing
data from different farms, or when animals are subject to diet

changes during the lactation. The diet is indeed a major component
of interest when studying feed efficiency because individual diges-
tive efficiency varies with diet composition. This influence of the
diet was observed, for example, by Tempelman et al. (2015) and
Lu et al. (2017), who determined substantial variability in partial
regression coefficients between different rations. Several studies
have already addressed this question. Durunna et al. (2011a) and
Cassady et al. (2016) found moderate phenotypic correlations of
0.33 and 0.40, respectively, between RFI of young beef cattle deter-
mined under both growing and finishing diets. Manafiazar et al.
(2015) reported a correlation of 0.30 for RFI between dry-lot condi-
tions and pasture. Estimate of genetic correlations between RFI of
the same animals in different situations is rare. They generally
give higher values than phenotypic correlations (0.50 in Durunna
et al., 2011b; 0.83 in Martin et al., 2019). These results suggest the
existence of a genotype-by-environment interaction, both in the
case of diet changes and variable systems of production, as
discussed by Berry and Crowley (2013) in their review. A test of
the ability of the present model to take different diets, or situations,
into account is therefore of primary interest.

Then, the mathematical modeling itself could be improved and
adapted to situations and data sets. For instance, one could use splines
instead of Legendre polynomials for the random regression or a hetero-
geneous variance for the residual effects.

A main point is also how animals with records for more than one
lactation should be taken into account. In its present form, the
model considers two lactations from the same animal as being inde-
pendent, which may have slightly influenced the figures, but with-
out changing any conclusion. As the first aim of this work was to
present the method, the focus was not on the RFI of individuals.
When the interest is on ranking animals, one can incorporate a
between-lactation (permanent environment) animal effect. More-
over, to be more accurate, all relationships among animals should
be taken into account when pedigree or genomic information is
available, which was not the case here. In addition, pedigree data
or genotypes are particularly interesting as it would allow to split
the animal effect into a genetic part and a permanent environmen-
tal part, as showed by Manzanilla Pech et al. (2014) in a model close
to ours. It is important to note that this new way of computing RFI
could indeed be used in genetic or genomic evaluations.

Conclusion

We have shown that it is possible to derive RFI from a multitrait ran-
dom regression model applied to DMI and its predictors, allowing for a
dynamic estimation of the traits, free from all time-related issues inher-
ent to the traditional RFI methodology. The model allows a better un-
derstanding of the correlations between the predictors during the
lactation, and it can be adapted and used in a genetic or genomic selec-
tion context.

**Supplementary materials**

Supplementary data to this article can be found online at https://doi.
org/10.1016/j.animal.2020.100101.

Ethics approval

During the experiment, all animals were kept indoors and handled
with care, in line with Aarhus University's ethics policy in compliance
with the European Union legislation for the protection of animals used
for scientific purposes.

Data and model availability statement

None of the data were deposited in an official repository.
## Author ORCIDs

Pauline Martin: https://orcid.org/0000-0003-2950-2728.

## Author contributions

P. Martin: Conceptualization, Formal analysis, Writing – original draft. V. Ducrocq: Conceptualization, Methodology, Writing – review & editing.

D.G.M. Gordo: Resources, Writing – review & editing.

N.C. Friggens: Conceptualization, Writing – review & editing, Project administration, Funding acquisition.

## Declaration of interest

The authors declare no conflicts of interest.

## Acknowledgments

GenTORE (www.gentore.eu) is a Horizon 2020 project running from 1 June 2017 to 31 May 2022. This research received funding from the European Union's H2020 Research and Innovation Program under agreement No.727213. The authors also would like to thank Amélie Fischer, Philippe Faverdin, Isabel Casasus and Daniel Villalba for their valuable advice on the work.

## Financial support statement

The research leading to these results has received funding from European Union's Horizon 2020 research and innovation programme - GenTORE - under grant agreement N° 727213.

## References

Archer, J.A., Reverter, A., Herd, R.M., Johnson, D.J., Arthur, P.F., 2002. Genetic variation in feed intake and efficiency of mature beef cows and relationships with postweaning measurements. Proceeding of the 7th world congress on genetic applied to livestock production, 19–23 August 2002, Montpellier, France, pp. 221–225.

Berry, D.P., Crowley, J.J., 2013. Cell biology symposium: genetics of feed efficiency in dairy and beef cattle. Journal of Animal Science 91, 1594–1613.

Byskov, M.V., Fogh, A., Løvendahl, P., 2017. Genetic parameters of rumination time and feed efficiency traits in primiparous Holstein cows under research and commercial conditions. Journal of Dairy Science 100, 9635–9642.

Cassady, C.J., Felix, T.L., Beever, J.E., Shike, D.W., 2016. Effects of timing and duration of test period and diet type on intake and feed efficiency of Charolais-sired cattle. Journal of Animal Science 94, 4748–4758.

Connor, E.E., 2015. Invited review: improving feed efficiency in dairy production: challenges and possibilities. Animal 9, 395–408.

Connor, E.E., Hutchison, J.L., Olson, K.M., Norman, H.D., 2012. TRIENNIAL LACTATION SYMPOSIUM: opportunities for improving milk production efficiency in dairy cattle. Journal of Animal Science 90, 1687–1694.

Durunna, O.N., Mujibi, F.D.N., Goonewardene, L., Okine, E.K., Basarab, J.A., Wang, Z., Moore, S.S., 2011a. Feed efficiency differences and reranking in beef steers fed grower and finisher diets. Journal of Animal Science 89, 158–167.

Durunna, O.N., Plastow, G., Mujibi, F.D.N., Grant, J., Mah, J., Basarab, J.A., Okine, E.K., Moore, S.S., Wang, Z., 2011b. Genetic parameters and genotype × environment interaction for feed efficiency traits in steers fed grower and finisher diets. Journal of Animal Science 89, 3394–3400.

van Elzakker, P.J.M., van Arendonk, J.A.M., 1993. Feed intake, body weight and milk production: genetic analysis of different measurements in lactating dairy heifers. Livestock Production Science 37, 37–51.

Esposito, G., Irons, P.C., Webb, E.C., Chapwanya, A., 2014. Interactions between negative energy balance, metabolic diseases, uterine health and immune response in transition dairy cows. Animal Reproduction Science 144, 60–71.

European Commission, 2018. EU agricultural and farm economics briefs. N°17. EU milk margin index estimate up to 2018. Retrieved on 12 November 2019 from. https://ec.europa.eu/agriculture/sites/agriculture/files/rural-area-economics/briefs/pdf/milk-margin-estimate-2018-no17_en.pdf.

FAO (Food and Agriculture Organization of the United Nations), 2010. Greenhouse gas emissions from the dairy sector: a life cycle assessment. Retrieved on 30 October 2019 from. http://www.fao.org/3/k7930e/k7930e00.pdf.

Fischer, A., Friggens, N.C., Berry, D.P., Faverdin, P., 2018. Isolating the cow-specific part of residual energy intake in lactating dairy cows using random regressions. Animal 12, 1396–1404.

Grummer, R.R., 2007. Strategies to improve fertility of high yielding dairy farms: management of the dry period. Theriogenology 68, S281–S288.

Hüttmann, H., Stamer, E., Junge, W., Thaller, G., Kalm, E., 2009. Analysis of feed intake and energy balance of high-yielding first lactating Holstein cows with fixed and random regression models. Animal 3, 181–188.

Islam, M.S., Jensen, J., Løvendahl, P., Karlskov-Mortensen, P., Shirali, M., 2020. Bayesian estimation of genetic variance and response to selection on linear or ratio traits of feed efficiency in dairy cattle. Journal of Dairy Science 103, 9150–9166.

Karacaören, B., Jaffrézic, F., Kadarmideen, H.N., 2006. Genetic parameters for functional traits in dairy cattle from daily random regression models. Journal of Dairy Science 89, 791–798.

Kennedy, B.W., van der Werf, J.H., Meuwissen, T.H., 1993. Genetic and statistical properties of residual feed intake. Journal of Animal Science 71, 3239–3250.

Koch, R.M., Swiger, L.A., Chambers, D., Gregory, K.E., 1963. Efficiency of feed use in beef cattle. Journal of Animal Science 22, 486–494.

Koenen, E.P.C., Veerkamp, R.F., 1998. Genetic covariance functions for live weight, condition score, and dry-matter intake measured at different lactation stages of Holstein Friesian heifers. Livestock Production Science 57, 67–77.

Li, B., Fikse, W.F., Lassen, J., Lidauer, M.H., Løvendahl, P., Måntysaari, P., Berglund, B., 2016. Genetic parameters for dry matter intake in primiparous Holstein, Nordic red, and Jersey cows in the first half of lactation. Journal of Dairy Science 99, 7232–7239.

Li, B., Berglund, B., Fikse, W.F., Lassen, J., Lidauer, M.H., Måntysaari, P., Løvendahl, P., 2017. Neglect of lactation stage leads to naive assessment of residual feed intake in dairy cattle. Journal of Dairy Science 100, 9076–9084.

Linamo, A.-E., Måntysaari, P., Måntysaari, E.A., 2012. Short communication: genetic parameters for feed intake, production, and extent of negative energy balance in Nordic red dairy cattle. Journal of Dairy Science 95, 6788–6794.

Lu, Y., Vandehaar, M.J., Spurlock, D.M., Weigel, K.A., Armentano, L.E., Staples, C.R., Connor, E.E., Wang, Z., Bello, N.M., Tempelman, R.J., 2015. An alternative approach to modeling genetic merit of feed efficiency in dairy cattle. Journal of Dairy Science 98, 6535–6551.

Lu, Y., Vandehaar, M.J., Spurlock, D.M., Weigel, K.A., Armentano, L.E., Staples, C.R., Connor, E.E., Wang, Z., Coffey, M., Veerkamp, R.F., de Haas, Y., Tempelman, R.J., 2017. Modeling genetic and nongenetic variation of feed efficiency and its partial relationships between component traits as a function of management and environmental factors. Journal of Dairy Science 100, 412–427.

Manafiazar, G., Basarab, J.A., Baron, V.S., McKeown, L., Doce, R.R., Swift, M., Undi, M., Wittenberg, K., Ominski, K., 2015. Effect of post-weaning residual feed intake classification on grazed grass intake and performance in pregnant beef heifers. Canadian Journal of Animal Science 95, 369–381.

Manzanilla Pech, C.I.V., Veerkamp, R.F., Calus, M.P.L., Zom, R., van Knegsel, A., Pryce, J.E., De Haas, Y., 2014. Genetic parameters across lactation for feed intake, fat- and protein-corrected milk, and liveweight in first-parity Holstein cattle. Journal of Dairy Science 97, 5851–5862.

Martin, P., Taussat, S., Vinet, A., Krauss, D., Maupetit, D., Renand, G., 2019. Genetic parameters and genome-wide association study regarding feed efficiency and slaughter traits in Charolais cows. Journal of Animal Science 97, 3684–3698.

Meyer, K., 2007. WOMBAT—a tool for mixed model analyses in quantitative genetics by restricted maximum likelihood (REML). Journal of Zhejiang University - Science B 8, 815–821.

Prendiville, R., Pierce, K.M., Delaby, L., Buckley, L., 2011. Animal performance and production efficiencies of Holstein-Friesian, Jersey, and Holstein × Holstein-Friesian cows throughout lactation. Livestock Science 138, 25–33.

Spurlock, D.M., Dekkers, J.C.M., Fernando, R., Koltes, D.A., Wolc, A., 2012. Genetic parameters for energy balance, feed efficiency, and related traits in Holstein cattle. Journal of Dairy Science 95, 5393–5402.

Strathe, A.B., Mark, T., Nielsen, B., Do, D., Kadarmideen, H.N., Jensen, J., 2014. Deriving genetic breeding values for residual feed intake from covariance functions of random regression models. Proceedings of the 10th world congress on genetic applied to livestock production, 18–22 August 2014, Vancouver, Canada, Abstract 116.

Tamminga, S., Luteijn, P.A., Meijer, R.G.M., 1997. Changes in composition and energy content of liveweight loss in dairy cows with time after parturition. Livestock Production Science 52, 31–38.

Tempelman, R.J., Spurlock, D.M., Coffey, M., Veerkamp, R.F., Armentano, L.E., Weigel, K.A., de Haas, Y., Staples, C.R., Connor, E.E., Lu, Y., Vandehaar, M.J., 2015. Heterogeneity in genetic and nongenetic variation and energy sink relationships for residual feed intake across research stations and countries. Journal of Dairy Science 98, 2013–2026.

Tetens, J., Thaller, G., Krattenmacher, N., 2014. Genetic and genomic dissection of dry matter intake at different lactation stages in primiparous Holstein cows. Journal of Dairy Science 97, 520–531.

Thorup, V.M., Højsgaard, S., Weisbjerg, M.R., Friggens, N.C., 2013. Energy balance of individual cows can be estimated in real-time on farm using frequent liveweight measures even in the absence of body condition score. Animal 7, 1631–1639.

Vallimont, J.E., Dechow, C.D., Daubert, J.M., Dekleva, M.W., Blum, J.W., Barlieb, C.M., Liu, W., Varga, G.A., Heinrichs, A.J., Baumrucker, C.R., 2010. Genetic parameters of feed intake, production, body weight, body condition score, and selected type traits of Holstein cows in commercial tie-stall barns. Journal of Dairy Science 93, 4892–4901.

Veerkamp, R.F., Brotherstone, S., 1997. Genetic correlations between linear type traits, food intake, live weight and condition score in Holstein Friesian dairy cattle. Animal Science 64, 385–392.

Veerkamp, R.F., Thompson, R., 1999. A covariance function for feed intake, live weight, and Milk yield estimated using a random regression model. Journal of Dairy Science 82, 1565–1573.