Milk, Fertility and Udder Health Performance of Purebred Holstein and Three-Breed Rotational Crossbred Cows within French Farms: Insights on the Benefits of Functional Diversity

Simple Summary When implementing dairy crossbreeding in purebred Holstein (HO) herds, farmers expect to improve the overall herd performance. However, they lack knowledge about how to manage and benefit from the diversity of genetic classes generated by three-breed rotational crossbreeding, which firstly refers to the cohabitation of purebred HO and first- and second-generation crosses (F1 and G2, respectively) within the herd. This study aimed to compare milk production, reproduction and udder health performance of HO, F1 and G2 cows, and to estimate how their combination in different proportions in the herd affects its profitability. We found that HO, F1 and G2 had different and complementary performance profiles, with two main trends. First, HO had higher milk yield, while F1 and G2 crosses had better fertility performance. Second, F1 had win-win trade-offs between milk production, fertility and udder health compared to HO and G2. We showed that HO-F1 or HO-F1-G2 (below 30%) mixed herds could be more profitable than purebred HO or fully crossbred herds with a conventional milk price. These findings can be used for advising purebred HO farmers who wonder about the benefits and the ways of managing the diversity of animal entities generated by the use of dairy crossbreeding in their herds. Abstract Using three-breed rotational crossbreeding in a purebred Holstein (HO) herd raises two questions: Do the different genetic classes of cows generated by crossbreeding perform differently? Are there any economic benefits of combining them within a herd? This study aimed at comparing the performance between the different genetic classes resulting from the use of three-breed rotational crossbreeding, and simulating the effect of combining them on herd profitability. Based on a dataset of 14 French commercial dairy herds using three-bred rotational crossbreeding from a HO herd over a 10-year period, we defined three genetic classes according to the theoretical value of heterosis and the percentage of HO genes. We performed linear models and estimated least square means to compare HO cows and the first and second generation of crosses (F1 and G2, respectively) on eight performance characteristics related to milk yield and solids, udder health and fertility. We used these to simulate profitability of five herd compositions differing according to HO, F1 and G2 proportions. We showed that HO, F1 and G2 cows had different and complementary performance profiles. HO had a win-lost trade-off between milk yield and fertility, G2 had the opposite trade-off and F1 had a win-win trade-off. Differences regarding milk solids and udder health were less clear-cut. We highlighted that combining HO with F1 or with both F1 and G2 (below 30%) could be more profitable than using purebred HO or crossbred herds in a conventional milk price scenario. These findings provide evidence on the benefits of functional diversity generated from the use of dairy crossbreeding in dairy herds.


Introduction
Using three-breed rotational crossbreeding in dairy herds is an appealing and relatively fast option to improve the functional traits of cows within high-yielding purebred herds and to help dairy farmers transitioning towards low-input dairy cattle systems [1][2][3].
Here, functional traits refer to characters of an animal that increase efficiency by reduced costs of input rather than by higher outputs of products (i.e., milk and meat) [4]. Functional traits that receive particular attention in dairy crossbreeding research are fertility, udder health and feed efficiency [5][6][7]. Despite its potential benefits through both breed complementarity and heterosis, using three-breed rotational crossbreeding still remains quite uncommon in France [8] and in many Western countries [2,9,10]. This situation has often been attributed to a lack of references on crossbreeding's benefits, which, however, has been gradually filled by research. Indeed, some studies have been estimating additive and non-additive genetic parameters of crossbreeding (i.e., heterosis and recombination loss) for production and functional traits [9,10] and assessing phenotypical cows' performance for specific three-breed crossbreeding programs [11][12][13]. To a lesser extent, others have also assessed the economic performance of crossbred herds for specific three-breed crossbreeding programs [14,15]. All these studies share a common approach: they focused on one or several specific crossbreeding programs, i.e., on specific combinations of dairy breeds. This built over the years an increasingly precise comparative of benefits and disadvantages of each program compared to purebred parental populations, according to the production systems [2,16], thus providing references toward dairy farmers on the available "genetic offer" for dairy crossbreeding. Finally, the few studies assessing the effect of three breed rotational crossbreeding on herd performance proceeded by simulation modeling and considered that one single crossbreeding program was introduced and managed integrally in the herd [14,15]. However, many dairy cattle farmers introduce and manage crossbreeding in a customized way and following a trial-and-error proceeding, which lead them to manage-temporarily or lastingly-several crossbreeding programs simultaneously [17], while keeping part of their herd purebred [14]. Consequently, their herd can be composed of diverse purebred and crossbred genotypes: introducing three-breed crossbreeding in a Holstein herd led to manage genetically diverse animal classes, regardless of crossbreeding programs. Hence, critical questions for dairy cattle farmers who manage the transition of a herd from pure-breeding to crossbreeding are: in which way the different genetic classes of cows within the herd differ in their performance? Do crossbred cows perform better than purebred cows in the herd? Do crosses of n + 1 generations maintain an advantage over those of n generations? Is it beneficial to herd profitability to combine these genetic classes? Answering these questions would contribute to address the issue of "functional diversity", defined as "the variety of functions played by individual production entities involved in overall performance" [14] and that may play a more critical role in overall performance than genetic or inter-specific diversity as such [14][15][16].
Therefore, the aim of the study was two-fold: (1) first, to assess and compare the performance for milk traits and fertility and udder health traits (referred to below as functional traits) between the different genetic animal classes resulting from the use of three-breed rotational crossbreeding, from a sample of French commercial dairy cattle farms practicing it over a 10-year period; (2) second, to simulate how different combinations of these genetic animal classes affect herd profitability.

Farm Sampling
The study draws on an original sample of 26 commercial dairy farms in France sampled as they managed three-breed crossbreeding and had at least 33% crossbred cows (all generations combined) and at least second-generation lactating cows in 2018 [17]. From these farms, we selected 14 in which herds were initially composed of purebred Holstein cows (n = 10) or with several pure breeds and prevailing Holstein breed (n = 4). All  available included milk, calving, fertility and udder health performance, as well as the breed and genealogy of females. These 14 dairy farms (8 organic and 6 conventional) were located in four main regions of dairy production in mainland France: Brittany (n = 4), Auvergne-Rhône-Alpes (n = 4), Hauts-de-France (n = 3) and Normandy (n = 3). The farms varied in size, forage systems and milk yield in 2008 and 2018 (Table 1).

Definition of Genetic Classes
Breed combinations varied greatly among the 14 herds (Figure 1), and more importantly, within each herd ( Figure A1). Due to the sampling criteria, purebred Holstein was the most common genotype among the sampled herds. Therefore, we used two indicators to define the genetic classes that resulted from using three-breed crossbreeding. The first one was the percentage of Holstein genes (HOg) of each female, which we calculated from the breed and genealogical data [7]. The second indicator was the coefficient of heterosis (H), which we calculated as follows [18]: where s i and d i are the proportions of sire genes and dam genes, respectively, from breed i. cows (n = 10) or with several pure breeds and prevailing Holstein breed (n = 4). All 14 herds had been enrolled in the French Milk Record Organization from 2008-2018. Data available included milk, calving, fertility and udder health performance, as well as the breed and genealogy of females. These 14 dairy farms (8 organic and 6 conventional) were located in four main regions of dairy production in mainland France: Brittany (n = 4), Auvergne-Rhône-Alpes (n = 4), Hauts-de-France (n = 3) and Normandy (n = 3). The farms varied in size, forage systems and milk yield in 2008 and 2018 (Table 1).

Definition of Genetic Classes
Breed combinations varied greatly among the 14 herds (Figure 1), and more importantly, within each herd ( Figure A1). Due to the sampling criteria, purebred Holstein was the most common genotype among the sampled herds. Therefore, we used two indicators to define the genetic classes that resulted from using three-breed crossbreeding. The first one was the percentage of Holstein genes (HOg) of each female, which we calculated from the breed and genealogical data [7]. The second indicator was the coefficient of heterosis (H), which we calculated as follows [18]: where and are the proportions of sire genes and dam genes, respectively, from breed i.   We then removed lactations for genetic classes that were insufficiently distributed among the sampled herds (i.e., third-, fourth-and fifth-generation crosses), which left only three genetic classes: purebred Holstein (HO), first-generation crosses (F 1 ) and secondgeneration crosses (G 2 ), which have 100%, 50% and 25% of Holstein genes, respectively ( Table 2). In doing so, however, the coefficient of heterosis became non-discriminatory in our database, since it was 100% for both F 1 and G 2 .

Data Editing
The French Milk Record Organization provided data for lactations of females of sampled herds. Available data included total and 305-day yields, somatic cell count (SCC), parity, lactation length, birth date, calving date, age at calving and length of the dry period. We calculated for each lactation from each cow milk yield (MY, in kg/cow) by either extrapolating or correcting for duration values of total and 305-day milk yield [19]. We then calculated corrected values for fat (FY, in kg/cow) and protein (PY, in kg/cow) yields, as well as fat (FC, in g/kg/cow) and protein (PC, in g/kg/cow) contents. We used three categories to express values for duration of the dry period of (n − 1) lactation: <50 days, 50-70 days and >70 days [20]. First, we calculated mean SCC for lactations with at least six SCC values available, then somatic cell score (SCS) as log 2 mean SCC 100,000 + 3. We calculated calving to first service interval (CFS, in days) as the number of days from the start of lactation to the first service, and calving interval (CI, in days) as the number of days from the start of (n − 1) lactation to the start of lactation n.
Since SCS and reproduction data were missing for some lactations (41% and 63%, respectively), for which production data were available, we created three separate datasets (Table 3). Dataset 1 consisted of milk performance data (i.e., MY, FY, PY, FC and PC) for 6672 lactations from 2730 cows. Dataset 2 consisted of SCS data for 3947 lactations from 2074 cows. Dataset 3 consisted of fertility performance data (i.e., CFS, CI) for 2449 lactations from 1411 cows. The breeds' combinations that were the most represented in datasets 1, 2 and 3 among F 1 and G 2 classes were the following: HO × Montbéliarde (MO), HO × Viking Red (VR), HO × Simmental and HO × Brown Swiss (BS) for F 1 , which represented about 54%, 9%, 8% and 6% of total sampled lactations in dataset 1, respectively, and (HO × MO) × BS, (HO × MO) × SR and (HO × SR) × MO for G 2 , which represented about 33%, 29% and 5% of total sampled lactations, respectively. The contributions to the three datasets of the diverse breed combinations involved in F 1 and G 2 crosses are detailed in Table A1. Table 3. Description of the three datasets used for the statistical analyses of milk production (dataset 1), udder health performance (dataset 2) and reproduction performance (dataset 3). For each, the number of lactations and the number of cows for Holstein (HO) purebred, first and second generation of crosses (F 1 and G 2 , respectively) obtained from three-breed rotational crossbreeding into purebred HO herds are indicated.

Statistical Analyses
We compared the three genetic classes on eight performance traits: five related to milk production (MY, FY, PY, FC and PC), one related to udder health (SCS) and two related to fertility (CFS and CI). We used linear models and estimated least square means (LSM) [21] for genetic classes (i.e., HO, F 1 and G 2 ) within parity (i.e., primiparous and multiparous). We then tested for differences in performance between pairs of genetic classes within parity using Tukey's tests (p < 0.05). We performed all statistical analyses using the emmeans package [22] in RStudio software (v. 4.0.4, RStudio Inc., Boston, MA, USA).
To estimate LSM of performance traits in dataset 1 (milk performance data), we used an adjusted version of the French genetic evaluation model for production traits [23]: where Y ijklmn = observation for the dependent variable, µ = overall population mean, P i × G j = fixed effect of the j-th genetic class (j = HO, F 1 , G 2 ) within the i-th parity (i = primiparous, multiparous), H k × Y l = fixed effect of herd-year for the k-th herd (k = 1-14) and l-th year (l = 2009-2018), P i × M m = fixed effect of the m-th calving month (m = January to December) nested within parity, P i × D n = fixed effect of length of the n-th previous dry period (n =< 50 days, 50-70 days, >70 days) within parity, P i × A o = fixed effect of o-th value for age at calving within the i-th parity (i = primiparous, multiparous), and ε ijklmn = residual error term.
To estimate LSM of performance traits in dataset 2 (somatic cells score data), we used an adjusted version of the French genetic evaluation model for cell count [23]: where Y ijklmno = observation for the dependent variable, µ = overall population mean, and ε ijklmno = residual error term. Finally, to estimate LSM of performance traits in dataset 3 (reproduction performance data), we used an adjusted version of the French genetic evaluation model for reproduction traits [23]: where Y ijklmo = observation for the dependent variable, µ = overall population mean, C i = fixed effect of i-th calving rank prior to insemination (i = primiparous, multiparous), and ε ijklmo = residual error term. The response variable used for the linear models of the reproduction performance (i.e., CFS and CI) were log-transformed. LSM were then transformed back.

Design of Simulated Herd Compositions and Estimation of Profitability
We considered a herd of 100 cows with several compositions of the three different genetic classes (i.e., combinations of HO, F 1 and G 2 ), based on the demographic dynamics of a crossbreeding herd [14]. Hence, we simulated five herd compositions (Table 4): one purebred Holstein herd (100HO) and four mixed herds including one HO-F 1 herd at 50% each (50HO-50F 1 ) and three mixed HO-F 1 -G 2 herds (32HO-63F 1 -5G 2 , 18HO-67F 1 -15G 2 and 10HO-60F 1 -30G 2 ). Table 4. Herd compositions (combinations of the three genetic classes HO, F 1 and G 2 ) simulated to estimate the economic performance at herd level.

Percentage of Cows of Each Genetic Class in the Herd
We calculated the performance of each herd composition for milk yield, fat and protein contents and calving interval as the mean value of LSM's weighted by the proportions of HO, F 1 , and G 2 in each herd composition simulated, as follows: where x is the weighted mean value for a given performance of the herd and x i is the LSM value estimated for the given performance of the i-th genetic class (i = HO, F 1 , G 2 ) and is the proportion of cows of the i-th genetic class in the herd composition. Then, for each herd composition simulated, we estimated the income generated by the volume of milk from the cows [24], the economic gains and costs related to the protein and fat contents, the reproduction costs of the cows [25] and finally the economic profitability generated by the milk production taking into account these bonuses/penalties. We considered a conventional milk price scenario in France.

Performance of Primiparous Cows
Primiparous cows of the three genetic classes had different performance profiles, with "win-win" trade-offs between milk production and fertility traits trade-offs for F 1 cows compared to HO and G 2 cows, for which trade-offs were to the detriment of fertility and milk yield, respectively (Table 5). F 1 had significantly higher FY than both HO and G 2 (310 vs. 296 and 294 kg, respectively). They also had higher MY (7655 vs. 7050 kg) and PY (245 vs. 240 kg) than G 2 , while not significantly differing from HO. G 2 primiparous cows had the highest performance for fat and protein contents of the three genetic classes: FC and PC were significantly higher for G 2 than for F 1 (+1.4 g/kg for FC and +1.0 g/kg for PC) and significantly higher for F 1 than for HO (+2.2 g/kg for FC and +0.8 g/kg for PC).
Regarding functional traits, both primiparous F 1 and G 2 crosses were more fertile than HO. Both CFS and CI were shorter for F 1 (+12 and +23 days, respectively) and G 2 (+17 and Animals 2021, 11, 3414 7 of 16 +26 days, respectively) than HO, while they did not differ significantly between F 1 and G 2 . Conversely, there was no significant difference between the three genetic classes for SCS. Table 5. Least square means estimates (standard error in parentheses) of 305-d milk, fat and protein yields, fat and protein contents, somatic cell score (SCS), calving to first service interval (CFS) and calving interval (CI) for primiparous cows of the three genetic classes: Holstein (HO), first-generation crossbreds with 50% of Holstein genes (F 1 ) and three-breed crossbreds of second generation with 25% of Holstein genes (G 2 ).

Performance of Multiparous Cows
As for primiparous cows, multiparous cows of the three genetic classes had different performance profiles, with F 1 maintaining a "win-win" trade-off-yet not as favorable as for multiparous cows-between milk production and functional traits compared to HO and G 2 cows (Table 6). Multiparous F 1 produced significantly less MY than HO (7596 vs. 7790 kg), but significantly more than G 2 (7596 vs. 6953 kg). They had higher FY than HO and G 2 (+14 and +16 kg, respectively) and higher PY than G 2 . Conversely, PY did not differ significantly neither between HO and F 1 , nor between HO and G 2 . FC and PC were significantly higher for F 1 than for HO (+1.9 g/kg for FC and +0.9 g/kg for PC) and significantly higher for G 2 than for F 1 (+0.8 g/kg for FC and +0.4 g/kg for PC). Table 6. Least square means estimates (standard error in parentheses) of 305-d milk, fat and protein yields, fat and protein contents, somatic cell score (SCS), calving to first service interval (CFS) and calving interval (CI) for multiparous cows of the three genetic classes: Holstein (HO), first-generation crossbreds with 50% of Holstein genes (F 1 ) and three-breed second-generation crossbreds with 25% of Holstein genes (G 2 ).  As for functional traits, CFS did not differ significantly between F 1 and G 2 . Conversely, CFS was shorter for F 1 and G 2 than for HO (−12 days and −13 days, respectively). F 1 also had shorter calving interval than HO (−18 days), while there was no significant difference neither between F 1 and G 2 nor between HO and G 2 . Finally, multiparous F 1 Animals 2021, 11, 3414 8 of 16 had significantly lower SCS than HO and G 2 , while the latter did not significantly differ for SCS.

Economic Performance of Herd Compositions Simulated
Mixed herd compositions generated higher profit per dairy cow than the 100% HO herd, with the exception of the mixed herd with a 30% G 2 share. Economic profitability was +29 €/cow for both the 50HO-50F 1 and 32HO-63F 1 -5G 2 mixed herds compared to purebred HO herd (Figure 2), while this was lower, equal to +17 €/cow, for the 32HO-63F 1 -5G 2 mixed herd (yellow bar in Figure 2).
Conversely, the mixed herd with the largest share of G 2 (10HO-60F 1 -30G 2 ) had a total net economic loss compared to purebred HO herd (−9 €/cow in Figure 2). Although it had the highest economic bonuses related to the increase of fat and protein contents and the improvement of the fertility of the cows (+39, +36 and +41 €/cow, respectively), it did not compensate for the penalty related to the lower milk yield of the crossbred cows compared to the purebred HO ones (−125 €/cow) in a conventional milk price scenario.

Economic Performance of Herd Compositions Simulated
Mixed herd compositions generated higher profit per dairy cow than the 100% HO herd, with the exception of the mixed herd with a 30% G2 share. Economic profitability was +29€/cow for both the 50HO-50F1 and 32HO-63F1-5G2 mixed herds compared to purebred HO herd (Figure 2), while this was lower, equal to +17€/cow, for the 32HO-63F1-5G2 mixed herd (yellow bar in Figure 2).
Conversely, the mixed herd with the largest share of G2 (10HO-60F1-30G2) had a total net economic loss compared to purebred HO herd (−9 €/cow in Figure 2). Although it had the highest economic bonuses related to the increase of fat and protein contents and the improvement of the fertility of the cows (+39, +36 and +41 €/cow, respectively), it did not compensate for the penalty related to the lower milk yield of the crossbred cows compared to the purebred HO ones (−125 €/cow) in a conventional milk price scenario.

Performance Profiles of Purebred HO, F1 and G2 Crossbred Cows
To discuss the different performance profiles between purebred HO and F1 and G2 crossbred cows, we based on studies that had focused on breeds' combinations that were the most represented in datasets 1, 2 and 3 among F1 and G2 classes (Table A1)

Performance Profiles of Purebred HO, F 1 and G 2 Crossbred Cows
To discuss the different performance profiles between purebred HO and F 1 and G 2 crossbred cows, we based on studies that had focused on breeds' combinations that were the most represented in datasets 1, 2 and 3 among F 1 and G 2 classes (Table A1) Regarding the milk performance of F 1 crosses compared to purebred HO in their first parity, our results are in accordance with those of many studies that reported HO × MO [26,27], HO × BS [28] and HO × SR [2,29,30] had similar milk yields and higher fat and protein Animals 2021, 11, 3414 9 of 16 yields compared to purebred HO in their first parity. Moreover, our results match those of many studies and conclude that F 1 have higher fat and protein contents compared to purebred HO in their first parity [29,30]. Conversely, they contrast with studies that also reported lower milk yield [29,31], lower fat and lower protein contents [29,30] for primiparous F 1 compared to purebred HO ones. Regarding the functional performance, our results contrast with studies [29,31,32] who found lower SCS for F 1 compared to purebred HO primiparous. Finally, regarding the fertility performance, our results are in accordance with those of many studies [2,29,30,33] that reported shorter DO for F 1 primiparous compared to purebred HO.
As for multiparous cows, our results are in accordance with most studies that reported that F 1 crosses had lower milk yield [27][28][29]33], higher fat [28,33] and protein yields [27,28] and higher fat and protein content [29] compared to HO. However, they contrast with fewer studies who reported lower fat and protein yields for multiparous F 1 compared to purebred HO [2,29]. Regarding SCS, our results contrast with those of studies that reported higher SCS for F 1 crosses compared to purebred HO [29,33]. Conversely, our results on the fertility performance match those of many studies that reported shorter days open for multiparous F 1 compared to purebred HO [2,29,32]. However, CFS and CI values for both primiparous and multiparous purebred HO were shorter in our results than the mean CFS and CI were in metropolitan France for purebred HO [25]. In fact, purebred HO in our sample had shorter CFS and CI than 50% of all French dairy farms: 85 vs. 87 days long for CFS and 402 vs. 403 days long for CI.
To our knowledge, very few studies compared G 2 with F 1 for the performance characteristics we investigated [12,32,34,35], whether in their first parity or the next ones: they more commonly compare G 2 with purebred HO [11,27]. Our results report lower milk, fat and protein yields for G 2 crosses compared to purebred HO in their first parity, and shorter DO, which is in accordance with those studies. Similarly, our results on multiparous cows report lower milk, fat and protein yields for G 2 crosses compared to purebred HO and shorter DO, which is in accordance with the same studies available [11,27].
Therefore, as most commonly reported in studies on dairy crossbreeding, our findings highlight that HO cows have a "win-lost" trade-off between milk yield and fertility, while G 2 cows have the opposite trade-off, and F 1 cows a "win-win" trade-off. By contrast, the advantages of both F 1 and G 2 crosses for milk solids (FY and PY) and SCS compared to purebred HO cows were less clear in our results, as they are in the literature on dairy crossbreeding: although F 1 crosses-and especially multiparous ones-had higher milk solids yields than HO cows, G 2 crosses do not really maintain such an advantage. This may be because these performance characteristics depend more on feeding and rearing conditions, and therefore, on the systems (commercial farms and experiments) in which the performance characteristics are evaluated.
The performances gaps that we observed between genetic classes are generally smaller than those reported in research on dairy crossbreeding. Comparing multiparous HO and F 1 cows, we observed a milk yield gap of −196 kg/cow, while other studies reported differences ranging from −347 to −1487 kg/cow [30,31,33]. Similarly, milk yield gap between multiparous HO and G 2 for MY was −837 kg/cow in our study, while it ranged from −850 to −1466 kg/cow in other studies [27]. The smaller gaps in milk yield and fertility performance that we observed between HO cows and the first-and second-generation crossbreds may be explained by the change in replacement and reform management by dairy farmers during their introduction of dairy crossbreeding in their herd. As the transition process proceeds from a purebred HO herd towards a partially crossbred herd, dairy farmers may reform HO cows with high milk yield but a low fertility performance [17]; consequently, our sample may be constituted accordingly by HO cows with a relatively higher fertility performance and lower milk yield than reported by most of the studies.

What Are the Benefits of Combining HO, F 1 and G 2 at the Herd Level?
Research on the benefits of using animal diversity to enhance the performance of livestock systems defines functional diversity as "the variety of functions played by individual production entities involved in overall performance" [36]. Many studies, thus, highlight that functional diversity plays a more critical role in overall performance than genetic or interspecific diversity [36][37][38]. Here, we showed that the performance profiles of the three genetic classes are different. Consequently, when considering the transition of one herd initially composed of purebred HO to rotational crossbreeding, its performance will change at the expense of milk yield and in favor of functional traits. Therefore, the challenge is to identify if combining these genetic classes is beneficial for the overall performance of the herd and which herd composition (i.e., the combination of HO, F 1 and G 2 ) allows to take advantage of the potential complementarity of their performance profiles.
Research on the benefits of using animal diversity to enhance the performance of livestock systems defines functional diversity as "the variety of functions played by individual production entities involved in overall performance" [36]. Our results provide empirical evidence of how such functional diversity is shaped and expressed within a herd and how it interacts with intra-specific diversity, which is still rare in animal production science [36,39].
Here, we showed that the performance profiles of the three genetic classes are different, which means that they are functionally diverse animal entities within the herd.
We also highlighted in our simulations that it may be beneficial to combine such functionally diverse entities. Thus, mixed herds HO-F 1 or HO-F 1 -G 2 with 15% G 2 or less would allow to benefit from complementary performance profiles of HO, F 1 and G 2 genetic classes: they combine proportions of the genetic classes in such a way that it both maintains the total milk yield in the herd and increases fat and protein contents, while reducing the operational costs related to cows' reproduction. Such a search for complementarity between milk yield, content and functional characteristics of genotypes within a herd has already been empirically demonstrated in the case of multi-breed herds in France [37].
Although consistent with previous studies that reported on the economic benefits of dairy crossbreeding in both organic and conventional farming systems [15] or regardless of the intensity level of the farming system [3], our first estimates are more nuanced. Two main reasons can explain that. Firstly, we considered milk price paid in conventional dairy farming in our simulations. However, a higher base milk price, as is the case in organic farming, could offset lower milk yield of G 2 crossbred cows in terms of the economic value of the total milk production of the herd. Moreover, this may explain why many dairy farmers who introduce dairy crossbreeding in their herds practice or convert to organic farming [1,7,19]. Secondly, we did not integrate the reduction in operational costs generated by the higher concentrate efficiency of crossbred cows [11] in our simulations, nor the added value generated by the valuation of meat products for specific crossbreeding programs [40].
Nevertheless, we provide initial insights about the compositions of crossbred or mixed herds that may be best to target while using rotational crossbreeding. They highlight that mixed HO-F 1 or HO-F 1 -G 2 herds with less than 30% of G 2 cows may be economically more beneficial than purebred HO herds or full-crossbred herds, which are yet often considered as the necessary target for dairy farmers using dairy crossbreeding, and consequently, that are the only herd compositions considered in studies assessing the economic benefits of dairy crossbreeding at the herd level [14,15]. To manage such mixed HO-F 1 or HO-F 1 -G 2 herds, dairy farmers can cross dairy breed females with both dairy and beef breed sires, in order to both increase the slaughter value of crossbred male calves [40] and regulate herd replacement for maintaining its mixed composition. Sexed semen can also be used in such mixed herds to have more flexibility in the replacement of the different genotypes [15].
In any case, mixing different genotypes within a dairy herd depends on what the farmer can accept losing and wants to gain. Indeed, based on our findings, three-breed rotational crossbreeding leads to a permanent loss of milk yield and a gain in the reproduction performance at herd level, as reported by most studies [9,14,29]. However, these gain and loss values may vary depending on the breeds used in crossbreeding programs [9,16,17], which cannot be taken into account with our method. Moreover, depending on the strategy dairy farmers adopt to increase the economic value of milk production, they can manage the herd composition differently: they may keep a small percentage of HO cows in the herd to complement loss of milk of F 1 crosses while benefiting from F 1 higher fat and protein contents. Alternatively, if economic bonuses for protein and fat contents or more generally base milk price are high enough, dairy farmers may mix F 1 and G 2 crosses and move towards a fully crossbred herd. However, this is not the only relevant strategy: other marketing strategies, such as Protected Designations of Origin (PDO) label for milk products, may encourage farmers to substitute HO cows with local pure breeds specified by PDO rules, and therefore, move towards a purebred herd rather than a crossbred one [41]. These results emphasize that an adequate herd composition (i.e., a combination of genetic classes with different trade-off profiles that meet dairy farmers' objectives), thus, depends on farmers' personal and business views of the performance of the herd [37,39].

Methodological Comments
In research that addresses the benefits of intra-specific animal diversity in dairy herds, one critical question is to develop a method to classify the different animal entities [36][37][38].
Here, we built genetic classes of cows based on theoretical values for heterosis and percentage of the initial pure breed to assess the effect of crossbreeding on cows' performance, both crossbred ones and purebred Holsteins that remain in the herd during the transition process to a crossbred herd. Our method of classification of cows' genotypes within crossbred herds is, therefore, similar to that developed by other studies [7,42]. However, we extended it to second-generation crosses of three-breed rotational crossbreeding, while previous studies only considered first-generation crosses and reciprocal crosses, i.e., twobreed crossbreeding. Our method is promising, as it allows to integrate dairy farmers' real ways of managing dairy crossbreeding, which are often less linear than introducing and managing one single crossbreeding scheme [17]. Such a method is still rare among studies that investigate dairy crossbreeding and its performance on both experimental [11,27] and commercial farms [2,29]. This raises questions about the applicability of their results under real farm conditions, where the strict application of a crossbreeding program may be limited by on-farm constraints but also the sociotechnical environment of the farms [1].
However, while our method enabled to provide results that are concordant with most of the studies focused on specific crossbreeding programs, such an approach of genetic categorization needs to be consolidated by testing it on a larger farm sample in France and beyond France to assess how robust our results are. We were not able to test the isolated effect of theoretical value for heterosis, as we did not have enough crossbred cows beyond second generation of three-breed rotational crossbreeding in our sample (e.g., G 3 and further). The heterosis value is yet another deciding factor of crossbreeding performance [18]. Moreover, the imbalance in lactation numbers between genetic classes (i.e., number of HO lactations > F 1 > G 2 ) makes it difficult to interpret the observed differences in standard error values between genetic classes for the estimated performance. It is, thus, difficult to know whether the higher values for G 2 are a statistical artefact or reflect a greater variability of performance for this genetic class. Most optimal would have been to follow an approach using samples of equivalent size for all crosses; unfortunately, this was practically unreachable and it would not have enabled us to integrate the real management practices of dairy farmers and the performance of cows while transitioning toward rotational crossbreeding. A larger sample should allow integrating more genetic classes, especially crosses beyond the second generation of rotational crossbreeding to test how relevant such genetic categorization is, and to provide a more comprehensive view of the performance of three-breed rotational crossbreeding in the long run. Moreover, a little more homogeneity in crossbreeding programs across sampled farms should allow to integrate the effects of complementarity between breeds into the method of categorizing genetic classes, which is currently lacking. For example, in F 1 crosses, 50% of genes comes from another breed than HO, but this other breed is not specified in our method. However, depending on whether it is Swiss Brown or Jersey, the characteristics provided by the breed may differ [37], and therefore, the performance of resulting cross with HO may differ too [28,43]. This may explain why we did not observe significant differences between F 1 and G 2 or between F 1 and HO for some traits. This may have been even more difficult if we had investigated performance characteristics, such as slaughter value for calves and cows, for which there is a particularly high variability between breeds [1]. A better consideration of the complementarity effect between breeds in the categorization method is, therefore, a challenge for future investigations.
In our study, we only investigated performance traits related to milk production, fertility and udder health, which is limited to accurately estimate the benefits of dairy crossbreeding on farms [15]. Therefore, there is a need to consider other production performance characteristics, such as those regarding meat production from dairy crossbred females [40] (since crossbred animals have "no real genetic value" [1]) and functional performance characteristics, such as feed conversion efficiency [6], general health [44], longevity or survival [44][45][46]. Finally, our economic estimates were rudimentary and need to be strengthened by integrating more items of production costs (e.g., feeding cost, veterinary and health costs) and by better modeling herd dynamics (e.g., by integrating stochastic herd management events that induce individual variability beyond genetic class [14,15]).

Conclusions
In this study, we developed a new method for classifying the different animal genetic entities generated by the use of three-breed rotational crossbreeding into purebred HO herds, which enables us to go beyond of the simple crossbreeding program. Based on the theoretical value of heterosis and the percentage of the initial HO purebred, this method was performed on a dataset of 14 French commercial dairy herds practicing three-breed rotational crossbreeding over a 10-year period. We showed that purebred HO and the first and second generation of crossbred cows had different and complementary performance profiles regarding milk production, reproduction and udder health. Thus, HO cows had a win-lost trade-off between milk yield and fertility, while G 2 cows had the opposite trade-off and F 1 cows had a win-win trade-off. Regarding milk solids and SCS, the advantages of F 1 and G 2 crosses compared to purebred HO were less clear-cut: however, multiparous F 1 crosses had a better performance than G 2 crosses. Based on these findings, we provided, for the first time, a few insights on the benefits of combining these functionally diverse animal genotypes at herd level. In a conventional milk price scenario, we highlighted that mixed HO-F 1 and HO-F 1 -G 2 herds with lesser that 30% of G 2 may be economically more beneficial than purebred HO herds than fully three-breed rotational crossbred herds, which are yet often considered as the necessary target for dairy farmers using dairy crossbreeding. However, future studies with greater size sample should verify these results and consolidate the method of classification of animal genetic entities developed here.