domain: sdv sdu sdv.mhep.csc

text:
MD, PhD João Pedro Ferreira
email: j.ferreira@chru-nancy.fr

MD John G F Cleland

MD, PhD Carolyn S P Lam

MD, PhD Dirk Van Veldhuisen

MD William M Byra

David La Police

MD, PhD Stefan D Anker

MD Mandeep R Mehra

Celine Leroy

PhD Valerie Eschwege

BS David A La Police

MSc Céline Leroy

PharmD, PhD Marie Toussaint-Hacquard

MD, PhD Patrick Rossignol

MD Barry Greenberg

MD, PhD Faiez Zannad

Impact of Geographic Region on the COMMANDER-HF Trial

Keywords: regional differences, rivaroxaban, heart failure HF, heart failure NPs, natriuretic peptides LVEF, left ventricular ejection fraction DAPT, dual anti-platelet therapy MI, myocardial infarction EE, Eastern Europe WE & SA, Western Europe and South Africa NA, North America AP, Asia-Pacific LA, Latin America HR, hazard ratio

à 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.

Introduction

Globalization of heart failure (HF) trials provides access to many research sites that enables a large number of patients to be enrolled in a timely fashion and potentially makes the results generalizable to a much broader population. However, there may be important regional differences in patient characteristics, outcomes, background treatment, treatment adherence practice patterns, health-care organization, social-support, level of income and inequality in income [START_REF] Tromp | Heart failure around the world[END_REF][START_REF] Ferreira | Geographic differences in heart failure trials[END_REF][START_REF] Pfeffer | Regional Variation in Patients and Outcomes in the Treatment of Preserved Cardiac Function Heart Failure With an Aldosterone Antagonist (TOPCAT) Trial[END_REF] . These differences may affect event-rates, the proportions of events contributing to a composite outcome (e.g., hospitalizations versus mortality), the adherence to treatment, and ultimately the observed treatment effect [START_REF] Pfeffer | Regional Variation in Patients and Outcomes in the Treatment of Preserved Cardiac Function Heart Failure With an Aldosterone Antagonist (TOPCAT) Trial[END_REF][START_REF] Dewan | Income Inequality and Outcomes in Heart Failure: A Global Between-Country Analysis[END_REF][START_REF] Ferreira | Income level and inequality as complement to geographical differences in cardiovascular trials[END_REF] .

Regional heterogeneity may exist despite the application of uniform trial entry criteria [START_REF] Mcmurray | Dapagliflozin in Patients with Heart Failure and Reduced Ejection Fraction[END_REF][START_REF] Packer | Evaluation of the effect of sodium-glucose co-transporter 2 inhibition with empagliflozin on morbidity and mortality of patients with chronic heart failure and a reduced ejection fraction: rationale for and design of the EMPEROR-Reduced trial[END_REF][START_REF] Armstrong | Vericiguat in Patients with Heart Failure and Reduced Ejection Fraction[END_REF] . Recent HF trials have required plasma natriuretic peptides (NP) as an entry criterion in order to exclude patients at low risk of events who may not even have HF [START_REF] Pfeffer | Regional Variation in Patients and Outcomes in the Treatment of Preserved Cardiac Function Heart Failure With an Aldosterone Antagonist (TOPCAT) Trial[END_REF][START_REF] Mcmurray | Dapagliflozin in Patients with Heart Failure and Reduced Ejection Fraction[END_REF][START_REF] Packer | Evaluation of the effect of sodium-glucose co-transporter 2 inhibition with empagliflozin on morbidity and mortality of patients with chronic heart failure and a reduced ejection fraction: rationale for and design of the EMPEROR-Reduced trial[END_REF][START_REF] Armstrong | Vericiguat in Patients with Heart Failure and Reduced Ejection Fraction[END_REF] .

Although this may increase event rates, it may fail to address other regional factors.

Accordingly, we compared patient characteristics, outcomes and treatment effect amongst prespecified regions in the COMMANDER-HF trial (Rivaroxaban in Patients with Heart Failure, Sinus Rhythm, and Coronary Disease), a contemporary, international, double-blind, randomized trial comparing the factor-Xa inhibitor rivaroxaban (2.5 mg twice daily) vs. placebo in addition to background antiplatelet therapy in patients with heart failure, coronary artery disease, reduced left ventricular ejection fraction and sinus rhythm following recent decompensation.

Methods

Study population

The study design of the COMMANDER-HF trial has been previously described [START_REF] Zannad | Rationale and design of a randomized, double-blind, event-driven, multicentre study comparing the efficacy and safety of oral rivaroxaban with placebo for reducing the risk of death, myocardial infarction or stroke in subjects with heart failure and significant coronary artery disease following an exacerbation of heart failure: the COMMANDER HF trial[END_REF][START_REF] Zannad | Rivaroxaban in Patients with Heart Failure, Sinus Rhythm, and Coronary Disease[END_REF] . Key inclusion criteria included history of chronic HF for three or more months, treatment for decompensated HF in the previous 30 days, left ventricular ejection fraction (LVEF) of 40% or less, history of coronary artery disease, and absence of atrial fibrillation or other indication for chronic anticoagulation. Decompensated HF was defined by symptoms of worsening dyspnea or fatigue, objective signs of congestion, and/or adjustment of HF medications requiring hospital admission. Rivaroxaban or placebo was given in addition to background single or dual anti-platelet therapy (DAPT).

Study outcomes

The primary efficacy outcome was the composite of myocardial infarction, stroke or all-cause mortality. Secondary efficacy outcomes included death from cardiovascular causes, rehospitalization for worsening heart failure, rehospitalization for cardiovascular events, and the composite of rehospitalization for worsening heart failure or death from cardiovascular causes. Rivaroxaban did not reduce the incidence of the primary efficacy outcome nor the rate of HF re-hospitalization; however, the incidence of ischemic stroke was reduced [START_REF] Mehra | A comprehensive analysis of the effects of rivaroxaban on stroke or transient ischaemic attack in patients with heart failure, coronary artery disease, and sinus rhythm: the COMMANDER HF trial[END_REF] , as well as the incidence of thrombo-embolic events [START_REF] Greenberg | Association of Rivaroxaban With Thromboembolic Events in Patients With Heart Failure, Coronary Disease, and Sinus Rhythm: A Post Hoc Analysis of the COMMANDER HF Trial[END_REF] . The principal safety outcome was the composite of fatal bleeding or bleeding into a critical space with a potential for causing permanent disability. Secondary safety outcomes included bleeding events requiring hospitalization and clinically overt major bleeding events as defined by the International Society on Thrombosis and Haemostasis (ISTH) (i.e., associated with a decrease in hemoglobin level of ≥2 g/dL, transfusion of 2 or more units of packed red cells or whole blood, a critical site [intracranial, intraspinal, intraocular, pericardial, intraarticular, intramuscular with compartment syndrome, or retroperitoneal], or a fatal outcome).

Investigators reported outcomes on detailed case report forms, which were verified by the sponsor's clinical operations team using source data. All participants provided written informed consent. The protocol was approved by the appropriate institutional review board or ethics committee at each site.

Based on a review of patient characteristics and event rates, blind to assigned treatment for the first 1,155 patients, the steering committee amended the enrolment criteria to require a plasma NT-pro BNP level ≥800 ng/L or BNP level ≥200 ng/L measured at any time between the index admission for decompensated HF and randomization in the local labs.

Simultaneous with the amendment, enrolment in the Asia-Pacific region and several additional countries began [START_REF] Zannad | Rivaroxaban in Patients with Heart Failure, Sinus Rhythm, and Coronary Disease[END_REF] . Overall, 5,022 patients were enrolled, 3,582 in countries that enrolled both before and after the amendment and 1,440 in countries that enrolled only after the protocol amendment [START_REF] Cunningham | Natriuretic Peptide-Based Inclusion Criteria in a Heart Failure Clinical Trial: Insights From COMMANDER HF[END_REF] . Most patients (included either before or after the amendment) were enrolled in Eastern Europe (64.2%; n = 3224) [START_REF] Cunningham | Natriuretic Peptide-Based Inclusion Criteria in a Heart Failure Clinical Trial: Insights From COMMANDER HF[END_REF] . Median follow-up was 21.1 months (percentile25-75, 12.9-32.8).

Plasma Concentrations of D-dimer and rivaroxaban

Treatment with rivaroxaban decreases plasma concentrations of D-dimer which can be used as a surrogate for treatment adherence [START_REF] Baglin | Effects on routine coagulation screens and assessment of anticoagulant intensity in patients taking oral dabigatran or rivaroxaban: guidance from the British Committee for Standards in Haematology[END_REF] . Samples were taken at 4 weeks after randomization in a random subset of patients in the COMMANDER-HF whom had been randomized to rivaroxaban (n =394), stored as plasma at -20°C and plasma concentrations of D-dimer (Siemens Healthcare Diagnostics INNOVANCE®) measured in a core laboratory and expressed in ng/mL after the trial had ended.

Samples were also taken in a smaller random subset of patients, in whom D-dimers at 4-weeks had also been measured, for the measurement of plasma concentrations of rivaroxaban (n =170) after the trial had ended. Rivaroxaban has dose-proportional bioavailability and a predictable, dose-dependent pharmacokinetics. In a fasted state, rivaroxaban pharmacokinetics is linear up to about 15 mg/day and is not affected by food [START_REF] Kubitza | Safety, pharmacodynamics, and pharmacokinetics of single doses of BAY 59-7939, an oral, direct factor Xa inhibitor[END_REF] .

The administered oral dose of rivaroxaban is absorbed rapidly, with the maximum plasma concentration (Cmax) occurring 2-4 hours after tablet intake and its elimination may take up to 13 hours [START_REF] Kubitza | Safety, pharmacodynamics, and pharmacokinetics of single doses of BAY 59-7939, an oral, direct factor Xa inhibitor[END_REF] . At total daily oral dose of 5 mg of rivaroxaban should provide a Cmax of 30 to 90 μg/l and a Cmin of 6 to 50 μg/l, depending on renal function, age and lean body mass [START_REF] Gosselin | International Council for Standardization in Haematology (ICSH) Recommendations for Laboratory Measurement of Direct Oral Anticoagulants[END_REF] . In our study, rivaroxaban was measured using an anti-Factor Xa chromogenic assay (Method DiXaI [Hyphen BioMed®] on STARMax [Diagnostica Stago®] automat).

Statistical analysis

Baseline characteristics were described using means ± SD for normally distributed continuous variables, median (percentile25-75) for skewed continuous variables, and number (proportion) for categorical variables, and compared across geographical regions. As prespecified, the studied regions were Eastern Europe (EE), Western Europe and South Africa (WE & SA), North America (NA), Asia-Pacific (AP), and Latin America (LA).

Number and proportion of events, event-rates (per 100 person-years), and hazard ratios (HR) were compared across regions. The proportional hazards assumption was assessed by testing for significant trend in the change of HR over time and was not violated. Treatment effect was studied by the intention-to-treat principle and heterogeneity by region was assessed through an interaction term between region and treatment within the Cox model. Cumulative incidence rates were compared using the Nelson-Aalen method. A 2-sided p value of <0.05 was considered significant. Stata® version 16 (StataCorp. 2019. College Station, TX:

StataCorp LLC) was used for the statistical analysis.

Results

Baseline patients` characteristics by region

In COMMANDER-HF, 3224 (64.2%) patients were from EE, 458 (9.1%) from WE & SA, 149 (3.0%) from NA, 733 (14.6%) from AP, and 458 (9.1%) from LA. Table 1. Compared with EE, patients from WE & SA, NA, and AP were older and those from LA had a similar age. Most patients were treated with angiotensin converting enzyme inhibitors/angiotensin receptor blockers (ACEi/ARBs), mineralocorticoid receptor antagonists (MRAs), and beta-blockers regardless of region. Digoxin was used more often in AP and LA, whereas angiotensin receptors neprilysin inhibitor (ARNi), percutaneous coronary intervention (PCI), coronary artery bypass grafting (CABG), and cardiac devices were more frequent in WE & SA and NA. The baseline use of dual antiplatelet therapy was highest in AP (46%). Patients from EE were more likely to have had a prior myocardial infarction (MI), whereas diabetes was more common in other regions. The proportion of patients with anaemia was higher in WE & SA and NA. Patients from EE had better renal function (higher estimated glomerular filtration rate [eGFR]), lower plasma concentrations of BNP/NT-pro BNP and higher left ventricular ejection fraction (LVEF). Randomization was well balanced amongst regions (as per stratification procedure).

Event rates by region

Patients from EE had lower event rates compared with other regions. For the composite of MI, stroke or death the event rates per 100 person-years in EE were 11.6 (10.8-12.5), in WE & SA 19.5 (16.5-23.0), in NA 14.2 (10.5-19.2), in AP 17.7 (15.4-20.3), and in LA 18.6 (15.6-22.1). The same pattern was observed for the composite of MI, stroke and cardiovascular death, its individual components, HF hospitalizations and "pump failure" death. Table 2. The rate of bleeding events was lower in EE. Bleeding requiring hospitalization was much less

common in EE: 0.8 (0.6-1.0) in EE, 2.0 (1.2-3.4) in WE & SA, 1.4 (0.5-3.6) in NA, 2.0 (1.3- 3.0) in AP, and 1.8 (1.1-3.2) in LA. Table 2.

Treatment effect by region

For the primary outcome of MI, stroke or death, the overall effect of rivaroxaban was neutral with no significant interaction by region (p for interaction =0.14), but point estimates for the HR ranged widely from 0.65 in NA to 1.0 in EE and 1.11 in AP. Differences were accentuated for the outcome of MI, stroke or cardiovascular death, where there was evidence of treatment heterogeneity (p for interaction =0.017) with potential benefit in WE & SA and NA. Similar patterns were observed for cardiovascular death alone (p for interaction =0.062) but not the composite of HF hospitalization and cardiovascular death (p for interaction =0.20). In EE the efficacy and safety of rivaroxaban was overall neutral. Figures 1,2 & Table 3.

Plasma concentrations of D-dimer and rivaroxaban

Plasma D-dimers were measured in 394 patients

(EE n =223, WE & SA n =49, NA n =14, AP n =56, LA n=52

Supplemental Table 2.

The totality of the findings is summarized in the Central Illustration.

Discussion

This analysis suggests that regional differences might have contributed to the overall neutral treatment effect in the COMMANDER-HF trial. Geographic region was the only factor used to stratify randomisation, and while the p-value for interaction did not reach statistical significance for the primary outcome, it was significant for the secondary composite, which substituted all-cause mortality with cardiovascular deaths. The rate of both efficacy and bleeding events was lower in EE and the overall treatment effect neutral. Plasma concentrations of rivaroxaban were lower in EE and 21% of patients had no detectable levels.

The fall in plasma D-dimer concentrations was also less pronounced with rivaroxaban in EE.

These observations suggest lower rates of adherence in EE, possibly due to investigators being less selective about enrolling patients who were likely to comply with investigational treatments and thereby enhancing enrolment rates. Alternatively, investigators may have had less time to spend on educating their patients about the protocol because of high recruitment rates.

Important differences were also found in other regions. Patients in AP might have experienced an increased risk of HF-related "pump failure" death and an excessive risk of major bleeding, findings that may help explain the neutral treatment effect in this region for very different reasons than in EE. The more frequent use of DAPT in AP may explain the higher observed bleeding rates in this region. LA patients had a highest rate of cardiovascular death but relatively few HF hospitalizations (the second lowest after EE), and the lowest MI rate, which might reflect important differences in health-care systems and HF treatment patterns. Regional differences in mortality rates and a "mismatch" between mortality and hospitalizations have also been documented in other trials. For example, in the PARADIGM-HF (angiotensin-neprilysin inhibition versus enalapril in heart failure) trial [START_REF] Mcmurray | Angiotensinneprilysin inhibition versus enalapril in heart failure[END_REF] , patients from LA had higher mortality rates than WE and NA but HF hospitalization rates were lowest in LA [START_REF] Kristensen | Geographic variations in the PARADIGM-HF heart failure trial[END_REF] . These findings may reflect the differences in the health-care systems and access to care and thus endpoints such as HF hospitalizations and other major vascular endpoints including MI/revascularization, may not reflect the same disease severity across regions.

The primary region of patient enrolment for contemporary HF trials lies outside WE and NA. A meta-analysis of 300 trials showed that the proportion with their primary region of enrolment in NA or WE decreased from almost 70% to just over 50% within a decade [START_REF] Vaduganathan | Globalization of heart failure clinical trials: a systematic review of 305 trials conducted over 16 years[END_REF] . An important driver of this trend is the greater speed of trial enrolment and lower costs in EE and LA [START_REF] Tromp | Heart failure around the world[END_REF] . This strategy may have drawbacks if the effect of treatment depends on background therapy and the organisation of health services and these are inhomogeneous. An intervention that is effective in one healthcare setting may not be effective in another. A recent and striking example was the TOPCAT (spironolactone for heart failure with preserved ejection fraction) trial, where patients from EE (Russia and Georgia) had low event rates (compared with patients from NA and LA), possibly due to different diagnostic standard being applied, and high rates of non-adherence leading to a neutral overall treatment effect [START_REF] Pfeffer | Regional Variation in Patients and Outcomes in the Treatment of Preserved Cardiac Function Heart Failure With an Aldosterone Antagonist (TOPCAT) Trial[END_REF][START_REF] De Denus | Spironolactone Metabolites in TOPCAT -New Insights into Regional Variation[END_REF] . COMMANDER-HF and TOPCAT indicate that greater scrutiny of patient enrolment and patient education about the protocol is required worldwide, but perhaps particularly in EE.

In lower income regions, patients may agree to participate in clinical trials in order to access better health-care. In some countries, payment for enrolling patients into trials are made directly to doctors rather than their hospital or university and employment and promotion of junior staff may depend on their ability to enrol patients into trials. These can act as perverse incentives to enrolling patients that are less likely to comply with the protocol. This accelerates recruitment to the detriment of good trial conduct. Once the patient is in the trial, they are there until the end (regardless of treatment adherence). The health-care professional is rewarded by inclusion and not by conduct. Cultural differences and beliefs may also have a major impact on therapeutic adherence and long-term compliance, with negative beliefs being more common among non-adherent patients [START_REF] Gagnon | Patient Beliefs Have a Greater Impact Than Barriers on Medication Adherence in a Community Health Center[END_REF][START_REF] Sjölander | The association between patients' beliefs about medicines and adherence to drug treatment after stroke: a cross-sectional questionnaire survey[END_REF] . Although we cannot ascertain if patients enrolled in EE have more negative beliefs, many of EE countries lack public funding and rely on private financing and co-payments, which may limit the reimbursement and availability of essential medicines [START_REF] Vogler | Inequalities in medicine use in Central Eastern Europe: an empirical investigation of socioeconomic determinants in eight countries[END_REF] . Hence, participating in the trial may serve as an incentive to access health-care regardless of therapeutic adherence. These issues need to be considered when designing future global trials [START_REF] Tromp | Heart failure around the world[END_REF] .

Global clinical trials should require objective inclusion criteria that are less prone to investigator interpretation, ensure that financial incentives to recruitment are carefully controlled and confirm that participants are aware of the consequences for others if they are non-adherent without good reason. In COMMANDER-HF a mid-trial protocol amendment was performed to include NP thresholds for inclusion in order to overcome the low eventrates observed in EE. Patients enrolled post-amendment in EE were "sicker" and had higher event rates. However, this protocol amendment did not modify the neutral effect of rivaroxaban and cannot overcome the low treatment adherence [START_REF] Cunningham | Natriuretic Peptide-Based Inclusion Criteria in a Heart Failure Clinical Trial: Insights From COMMANDER HF[END_REF] . For improving the scrutiny about treatment adherence, blood and urine samples should be collected to monitor treatment adherence in "real-time" during the trial (and not a posteriori), and the corresponding attitudes in the case of important non-adherence (e.g., >10%) in certain centres or regions should be prespecified in the protocol (e.g., stopping recruitment in those centres or regions, improve monitoring and education, and readapt sample size). This is particularly important in trials with a median follow-up duration superior to one year, such as COMMANDER-HF, as after the first year of follow-up nearly 40% of trial participants might have stopped taking their trial medication [START_REF] Blaschke | Adherence to medications: insights arising from studies on the unreliable link between prescribed and actual drug dosing histories[END_REF] .

Strengths and limitations

The inclusion of geographic region as a stratification variable in COMMANDER-HF, ensured a good balance between the rivaroxaban and placebo arms within regions. However, the trial was not powered to study the treatment effect within and between subgroups [START_REF] Pocock | Eastern Europe; WE & SA, Western Europe and South Africa; NA, North America; AP, Asia-Pacific; LA, Latin America; ACEi/ARB, angiotensin converting enzyme inhibitor/angiotensin receptor blocker; MRA, mineralocorticoid receptor antagonist; ARNi, angiotensin-receptor neprilysin inhibitor; DAPT, dual anti-platelet therapy; PCI, percutaneous coronary intervention; CABG, coronary artery bypass grafting; BMI, body mass index; SBP, systolic blood pressure; DBP, diastolic blood pressure; eGFR, estimated glomerular filtration rate; BNP, brain natriuretic peptide; NT-pro BNP, N-terminal pro brain natriuretic peptide; LVEF, left ventricular ejection fraction[END_REF] .

Notwithstanding, these findings along with those from TOPCAT, support that patients enrolled in EE countries may be less compliant with the trial good practices and future trials should take this information into account and implement some of the measures discussed in this manuscript (e.g., strict and objective entry criteria and real-time treatment adherence monitoring).

Conclusions

In the COMMANDER-HF trial important differences amongst geographic regions were observed in patient characteristics, the rate and pattern of events and treatment adherence that may have contributed to an overall neutral outcome. These differences were most marked in Eastern European countries that presented lower event rates and treatment adherence. These observations highlight the need for greater scrutiny of trial conduct and treatment adherence.

Clinical perspectives

-The rate of both efficacy and bleeding events was lower in Eastern Europe and the overall treatment effect neutral.

-Plasma concentrations of rivaroxaban were lower and 21% of patients had no detectable levels in Eastern Europe.

-These regional differences might have contributed to the overall neutral treatment effect in the COMMANDER-HF trial.

Translational outlook

-Global clinical trials should require objective inclusion criteria that are less prone to investigator interpretation, ensure that financial incentives to recruitment are carefully controlled and confirm that participants are aware of the consequences for others if they are non-adherent without good reason.

-Blood levels of treatments should be routinely determined and action should be taken during the trial if adherence rates are low in certain regions (e.g., stop patient inclusion in areas with lower treatment adherence).

Figure Legends: Figure 1. Forest plot representing the treatment effect on the time-to-first myocardial infarction, stroke or all-cause death and the same composite with cardiovascular death instead of all-cause death by the region of stratification

Rivaroxaban (vs. placebo) might have been beneficial in regions other than Eastern Europe and Asia-Pacific, with significant effect heterogeneity for the composite outcome of myocardial infarction, stroke or cardiovascular death. EE, Eastern Europe; WE & SA, Western Europe and South Africa; NA, North America; AP, Asia-Pacific; LA, Latin America; MI, myocardial infarction. Patients from Eastern Europe had lower event rates and no signal of efficacy or safety from rivaroxaban. EE, Eastern Europe; WE & SA, Western Europe and South Africa; NA, North America; AP, Asia-Pacific; LA, Latin America; MI, myocardial infarction; CV, cardiovascular; HF, heart failure.

Central Illustration. Main findings of the COMMANDER-HF trial according to the prespecified world region.

Patients from Eastern Europe presented lower event rates, lower treatment adherence and a neutral effect of rivaroxaban. TTx, treatment

Figure 2 .

 2 Figure 2. Nelson-Aalen cumulative incidence curves for the time-to-first myocardial infarction, stroke or cardiovascular death, heart failure hospitalization or cardiovascular death, and any hospitalization due to bleeding by the region of stratification.Patients from Eastern Europe had lower event rates and no signal of efficacy or safety from rivaroxaban. EE, Eastern Europe; WE & SA, Western Europe and South Africa; NA, North America; AP, Asia-Pacific; LA, Latin America; MI, myocardial infarction; CV, cardiovascular; HF, heart failure.

Table 1 .

 1 Numbers enrolled in each country is shown inSupplemental

Table 1 .

 1 Bulgaria, Poland, Russia, Ukraine and Romania contributed 83.2% of EE patients.

  ). D-dimer levels decreased over 4 weeks in all regions, with the smallest change (95%CI) found in EE -89 (-241, 62) vs. -202 (-558, 154) in WE & SA, -594 (-1417, 229) in NA, -372 (-824, 80) in AP, and -801 (-1877, 275) in LA. Furthermore, in a subset of 170 patients randomized to rivaroxaban (EE n =128, WE & SA n =14, NA n =2, AP n =14, LA n =12), there was no detectable blood concentration in 27 of 128 patients in EE (21.1%) at 4 weeks, compared to 2 of 42 patients from other regions (4.8%). The median concentrations were 37% lower in EE patients vs. rest of the world (28 vs. 45 μg/l, p =0.021).

Table 1 . Characteristics of the population by the region of stratification

 1

	Factor	EE	WE & SA	NA	AP	LA	p-value
	N.	3224	458	149	733	458	-
	Age, years	65.8 ± 9.7	70.5 ± 10.6	67.0 ± 11.1	67.3 ± 11.1	64.9 ± 10.2 <0.001
	Age ≥65	1773 (55.0%)	327 (71.4%)	82 (55.0%)	443 (60.4%)	238 (52.0%) <0.001
	Female	764 (23.7%)	97 (21.2%)	41 (27.5%)	146 (19.9%)	102 (22.3%) 0.11
	Race
	White	3216 (99.8%)	439 (95.9%)	101 (67.8%)	13 (1.8%)	359 (78.4%) <0.001
	Black or African
	American	0 (0.0%)	3 (0.7%)	40 (26.8%)	1 (0.1%)	21 (4.6%)
	Asian	0 (0.0%)	4 (0.9%)	5 (3.4%)	716 (97.7%)	2 (0.4%)
	Other	8 (0.2%)	12 (2.6%)	3 (2.0%)	3 (0.4%)	76 (16.6%)
	ACEi/ARB	3094 (96.0%)	378 (82.5%)	130 (87.2%)	639 (87.2%)	419 (91.5%) <0.001
	MRA	2555 (79.2%)	297 (64.8%)	68 (45.6%)	545 (74.4%)	375 (81.9%) <0.001
	Beta-blocker	3049 (94.6%)	426 (93.0%)	143 (96.0%)	617 (84.2%)	407 (88.9%) <0.001
	Digoxin	154 (4.8%)	8 (1.7%)	21 (14.1%)	160 (21.8%)	90 (19.7%) <0.001
	Loop diuretic	3223 (100.0%) 452 (98.7%)	148 (99.3%)	723 (98.6%)	453 (98.9%) <0.001
	ARNi	7 (0.2%)	24 (5.2%)	7 (4.7%)	1 (0.1%)	2 (0.4%)	<0.001
	Aspirin	3035 (94.1%)	414 (90.4%)	140 (94.0%)	639 (87.2%)	447 (97.6%) <0.001
	DAPT	1064 (33.0%)	138 (30.1%)	37 (24.8%)	334 (45.6%)	173 (37.8%) <0.001
	Diabetes	1187 (36.8%)	232 (50.7%)	86 (57.7%)	340 (46.4%)	207 (45.2%) <0.001
	Hypertension	2578 (80.0%)	317 (69.2%)	122 (81.9%)	442 (60.3%)	324 (70.7%) <0.001
	Myocardial
	infarction	2651 (82.2%)	317 (69.2%)	101 (67.8%)	416 (56.8%)	318 (69.4%) <0.001
	Stroke	289 (9.0%)	41 (9.0%)	21 (14.1%)	78 (10.6%)	24 (5.2%)	0.005
	Cardiac device	52 (1.6%)	19 (4.1%)	7 (4.7%)	10 (1.4%)	6 (1.3%)	<0.001
	PCI or CABG	1958 (60.7%)	331 (72.3%)	109 (73.2%)	493 (67.3%)	259 (56.6%) <0.001
	NYHA
	I	26 (0.8%)	20 (4.4%)	3 (2.0%)	74 (10.1%)	26 (5.7%)	<0.001
	II	1242 (38.5%)	267 (58.3%)	69 (46.6%)	309 (42.2%)	331 (72.3%)
	III	1895 (58.8%)	167 (36.5%)	66 (44.6%)	244 (33.3%)	90 (19.7%)
	IV	61 (1.9%)	4 (0.9%)	10 (6.8%)	106 (14.5%)	11 (2.4%)
	BMI, kg/m2	28.5 ± 4.9	27.8 ± 5.5	29.4 ± 6.7	23.6 ± 3.8	28.3 ± 5.4	<0.001
	BMI categories
	<25	777 (24.1%)	154 (33.8%)	43 (28.9%)	505 (68.9%)	123 (26.9%) <0.001
	25-29	1369 (42.5%)	172 (37.8%)	51 (34.2%)	192 (26.2%)	190 (41.5%)
	≥30	1077 (33.4%)	129 (28.4%)	55 (36.9%)	36 (4.9%)	145 (31.7%)
						120.2 ±
	SBP, mmHg	124.6 ± 14.0	120.3 ± 18.8	116.0 ± 17.4 119.1 ± 16.1	16.7	<0.001
	DBP, mmHg	75.0 ± 8.6	68.6 ± 10.7	69.3 ± 10.4	70.3 ± 11.0	72.7 ± 9.9	<0.001
						415.0
		330.0 (200.0,	475.0 (302.5,	430.0 (285.0,	385.0 (225.0,	(245.0,
	D-dimer, ug/L	590.0)	817.5)	795.0)	760.0)	710.0)	<0.001

Table 2 . Events, event rates and hazard ratios by the region of stratification

 2

	Region	EE	WE & SA	NA	AP	LA
	Outcome
	MI, Stroke, Death
	Events, n. (%)	777 (24.1)	139 (30.4)	42 (28.2)	198 (27.0)	128 (28.0)
	Event rates (per 100py)	11.6 (10.8 -12.5)	19.5 (16.5-23.0)	14.2 (10.5-19.2) 17.7 (15.4-20.3)	18.6 (15.6-22.1)
	Hazard ratio (95%CI)	Ref.	1.63 (1.36-1.95)	1.21 (1.18-1.64) 1.46 (1.25-1.71)	1.53 (1.27-1.85)
	MI, Stroke, CV death
	Events, n. (%)	694 (21.5)	114 (24.9)	34 (22.8)	168 (22.9)	111 (24.2)
	Event rates (per 100py)	10.4 (9.6-11.2)	16.0 (13.3-19.2)	11.5 (8.2-16.1)	15.0 (12.9-17.4)	16.1 (13.4-19.4)
	Hazard ratio (95%CI)	Ref.	1.48 (1.21-1.80)	1.09 (0.77-1.54) 1.37 (1.15-1.62)	1.47 (1.2-1.80)
	MI
	Events, n. (%)	132 (4.1)	32 (7.0)	9 (6.0)	33 (4.5)	10 (2.2)
	Event rates (per 100py)	2.0 (1.6-2.3)	4.4 (3.1-6.2)	3.0 (1.6-5.8)	2.9 (2.0-4.0)	1.4 (0.8-2.6)
	Hazard ratio (95%CI)	Ref.	2.13 (1.45-3.14)	1.54 (0.78-3.02) 1.39 (0.94-2.04)	0.68 (0.36-1.3)
	Stroke
	Events, n. (%)	69 (2.1)	14 (3.1)	4 (2.7)	25 (3.4)	15 (3.3)
	Event rates (per 100py)	1.0 (0.8-1.3)	1.9 (1.1-3.2)	1.3 (0.5-3.6)	2.2 (1.5-3.2)	2.1 (1.3-3.6)
	Hazard ratio (95%CI)	Ref.	1.77 (1.0-3.15)	1.29 (0.47-3.54) 2.01 (1.26-3.18)	1.99 (1.13-3.48)
	All-cause death
	Events, n. (%)	673 (20.9)	114 (24.9)	38 (25.5)	161 (22.0)	116 (25.3)
	Event rates (per 100py)	9.8 (9.1-10.6)	15.0 (12.5-18.0)	12.6 (9.2-17.3)	13.7 (11.8-16.0)	16.3 (13.5-19.5)
	Hazard ratio (95%CI)	Ref.	1.50 (1.23-1.83)	1.27 (0.91-1.75) 1.36 (1.14-1.61)	1.61 (1.32-1.96)
	CV death
	Events, n. (%)	587 (18.2)	87 (19.0)	29 (19.5)	128 (17.5)	98 (21.4)
	Event rates (per 100py)	8.6 (7.9-9.3)	11.4 (9.3-14.1)	9.6 (6.7-13.8)	10.9 (9.2-13.0)	13.7 (11.3-16.7)
	Hazard ratio (95%CI)	Ref.	1.30 (1.03-1.63)	1.10 (0.76-1.6)	1.22 (1.0-1.47)	1.53 (1.24-1.9)

Table 3 . Treatment effect (rivaroxaban vs. placebo) by the region of stratification

 3 .73-1.17) 0.61 (0.29-1.29) 0.72 (0.16-3.21) 0.80 (0.48-1.33) 1.09 (0.48-2.47) 0.82 "Pump failure" death 1.11 (0.84-1.47) 0.60 (0.31-1.13) 0.38 (0.12-1.21) 1.74 (0.97-3.11) 0.92 (0.51-1.64) 0.056 CV death or HF hosp. 1.04 (0.93-1.17) 0.86 (0.65-1.12) 0.62 (0.38-0.99) 1.04 (0.84-1.Eastern Europe; WE & SA, Western Europe and South Africa; NA, North America; AP, Asia-Pacific; LA, Latin America; MI, myocardial infarction; CV, cardiovascular; HF, heart failure; ISTH, International Society on Thrombosis and Haemostasis; HR, hazard ratio; CI,

	Region	EE	WE & SA	NA	AP	LA	Interaction P
	Outcome, HR
	(95%CI)
	MI, Stroke, Death	1.0 (0.87-1.15)	0.76 (0.54-1.06) 0.65 (0.35-1.21) 1.11 (0.84-1.46) 0.73 (0.52-1.04) 0.14
	MI, Stroke, CV death 1.0 (0.86-1.16)	0.61 (0.42-0.89) 0.46 (0.22-0.94) 1.07 (0.79-1.45) 0.72 (0.49-1.04) 0.017
	MI	0.94 (0.67-1.33) 0.99 (0.5-1.98)	0.28 (0.06-1.33) 0.74 (0.37-1.47) 0.25 (0.05-1.16) 0.27
	Stroke	0.81 (0.50-1.3)	0.16 (0.04-0.72) 0.32 (0.03-3.09) 1.09 (0.50-2.4)	0.36 (0.12-1.14) 0.14
	All-cause death	1.02 (0.88-1.19) 0.87 (0.6-1.25)	0.62 (0.32-1.19) 1.16 (0.85-1.58) 0.81 (0.56-1.17) 0.31
	CV death	1.02 (0.87-1.2)	0.67 (0.44-1.03) 0.43 (0.2-0.94)	1.15 (0.81-1.62) 0.81 (0.55-1.21) 0.062
	Sudden death	0.93 (03)	0.94 (0.69-1.27) 0.20
	HF hosp.	1.06 (0.92-1.22) 0.80 (0.59-1.08) 0.65 (0.39-1.09) 0.94 (0.74-1.2)	1.11 (0.77-1.6)	0.21
	Bleeding main safety	0.9 (0.37-2.21)	0.33 (0.03-3.16) 1.0 (0.12-8.11)	2.37 (0.61-9.16) 0.14 (0.02-1.16) 0.13
	ISTH-major bleeding 1.47 (0.89-2.42) 1.81 (0.61-5.39) 0*	2.73 (1.27-5.87) 0.70 (0.27-1.85) 0.19
	Bleeding req. hosp.	1.71 (0.94-3.12) 2.33 (0.60-9.01) 0*	2.22 (0.77-6.38) 1.0 (0.29-3.45)	0.55
	Legend: EE, confidence interval.

*Unprecise estimates.

The COMMANDER-HF trial was supported by Janssen Research and Development. Acknowledgements JF, PR and FZ are supported by the RHU Fight-HF, a public grant overseen by the French National Research Agency (ANR) as part of the second "Investissements d'Avenir" program (reference: ANR-15-RHUS-0004), and by the French PIA project "Lorraine Université d'Excellence" (reference: ANR-15-IDEX-04-LUE). The biobanking for the rivaroxaban dosages is managed by Biological Resource Center Lorrain BB-0033-00035.

Major bleeding is defined by the ISTH as overt bleeding associated with a decrease in hemoglobin level of at least 2 g/dL, transfusion of two or more units of packed red cells or whole blood, a critical site (intracranial, intraspinal, intraocular, pericardial, intraarticular, intramuscular with compartment syndrome, or retroperitoneal), or a fatal outcome.


Luciano Barreto Mendes
email: luciano.mendes@inrae.fr.

M Coppa

J Rouel

B Martin

B Dumont

A Ferlay

C Espinasse

F Blanc

Profiles of dairy cows with different productive lifespan emerge together from multiple traits assessed at first lactation: the case of a grassland-based dairy system

Keywords: clustering of dairy cows, robustness, trade-offs between functional and productive traits, discriminant analysis, survival, longevity

Profiles of dairy cows with different productive lifespan

emerge together from multiple traits assessed at first lactation: the case of a grassland-based dairy system 2

HIGHLIGHTS

• A trade-off was observed among productive lifespan of cows and milk yield, fertility and body weight at first lactation;

• Contrasting productive lifespan profiles defined as long, average and short were characterized by differing reproductive performance and milk production traits in first lactation;

• Individual placement of cows into groups of long, average or short productive longevity can be correctly predicted;

• Studying traits of cows at first lactation enables early identification of cows that have a greater probability of living longer.

INTRODUCTION

Recent studies have emphasized that robustness is a relevant characteristic of dairy cows that should be addressed in the future [START_REF] Roche | New considerations to refine breeding objectives of dairy cows for increasing robustness and sustainability of grass-based milk production systems[END_REF]. This is particularly true in extensive grassland-based production systems, where cows are exposed to limiting conditions and fluctuating environments [START_REF] Delaby | Robust animals for grass-based production systems[END_REF] and are able to express their genetic potential only partially. When facing nutritional constraints or challenges, modern dairy cows are expected to show a trade-off between life functions (growth, survival, health, reproduction, and lactation) [START_REF] Blanc | Adaptive abilities of the females and sustainability of ruminant livestock systems[END_REF][START_REF] Friggens | Towards a biological basis for predicting nutrient partitioning: the dairy cow as an example[END_REF]. Such trade-offs result from breeding programs mainly focused on lactation and, as a consequence, milk yield per dairy cow has doubled in the last 40 years [START_REF] Oltenacu | The impact of genetic selection for increased milk yield on the welfare of dairy cows[END_REF]. Such a priority given to milk production traits has been shown to be deleterious to other traits, such as fertility [START_REF] Albarrán-Portillo | The relationship between fertility and lactation characteristics in Holstein cows on United Kingdom commercial dairy farms[END_REF][START_REF] Berglund | Genetic Improvement of Dairy Cow Reproductive Performance[END_REF][START_REF] Pryce | Fertility in the high-producing dairy cow[END_REF], overall health [START_REF] Berry | Genetics of animal health and disease in cattle[END_REF], welfare [START_REF] Oltenacu | The impact of genetic selection for increased milk yield on the welfare of dairy cows[END_REF] and longevity [START_REF] Essl | Longevity in dairy cattle breeding: a review[END_REF], and it has yielded genotypes that show better dairy performance and efficiency in optimal environments (a rich environment where individuals are allowed to fully express their genetic potential [START_REF] Vandehaar | Harnessing the genetics of the modern dairy cow to continue improvements in feed efficiency[END_REF]), but are more susceptible to challenging environments (where animals are exposed to scarcities or restrictions [START_REF] Blanc | Adaptive abilities of the females and sustainability of ruminant livestock systems[END_REF]).

Some authors have suggested that profiles of cows defined from expression of trade-offs between life functions that play an important role in challenging situations are linked to robust properties [START_REF] Ollion | Assessing the diversity of trade-offs between life functions in early lactation dairy cows[END_REF]. Recently, [START_REF] Friggens | Review: Deciphering animal robustness. A synthesis to facilitate its use in livestock breeding and management[END_REF] defined animal robustness as "the ability, in the face of environmental constraints, to carry on doing the various things that the animal has to do to favor its future ability to reproduce". From this definition it clearly appears that robustness is linked to the ability of the animal to remain in the herd (not being culled or sold and staying alive) and to match the farmer's expectations [START_REF] Ollion | How livestock farmers define dairy cows[END_REF].

Thus, an integrative way to estimate animal robustness is by looking at longevity [START_REF] Friggens | Review: Deciphering animal robustness. A synthesis to facilitate its use in livestock breeding and management[END_REF]. Our study aimed to test the link between these two previous approaches to assess animal robustness, i.e. to connect longevity to the ways cows express trade-offs between life functions.

As is well known, the longevity of a dairy cow is ultimately determined by the culling rules applied by the farmer. Culling can be driven by economic reasons such as herd renewal needs, feed costs, market demand for milk production, or slaughter [START_REF] Edwards-Callaway | Culling Decisions and Dairy Cattle Welfare During Transport to Slaughter in the United States[END_REF], but most often it is first determined by biological reasons when cows are not able to meet the farmer's expectations anymore or fail to match the requirements of the production system, such as production amount and quality, for instance. Biological reasons for culling may be quite diverse like poor udder health or conformation, low fertility, low production, poor milking, leg problems, diseases (including metabolic diseases), poor temperament, injuries or old age. For example, poor udder health and poor fertility were reported to be the main reason for culling in organic and conventional Swedish dairy herds, respectively [START_REF] Ahlman | Culling reasons in organic and conventional dairy herds and genotype by environment interaction for longevity[END_REF].

Because culling decisions involve risks related to discarding a cow that would otherwise be productive for a few more lactations, or retaining in the herd an unproductive individual for another cycle, we considered that the earlier we could detect cows likely to have high versus short longevity the more helpful it would be for the farmer. As [START_REF] Adriaens | Productive life span and resilience rank can be predicted from on-farm first-parity sensor time series but not using a common equation across farms[END_REF] have recently shown that production and activity traits measured at first lactation have the potential to predict the cows' lifetime resilience ranking within a herd, we hypothesized that, within a herd, groups of cows that would express different trade-offs between life functions during first lactation will have different productive lifetime performances (i.e. different robustness).The first objective of the present study was to establish whether cows differing in productive lifespan have different ways to prioritize life functions at first lactation. This is an a posteriori hypothesis testing approach, by assessing the variability of all traits and allowing intrinsic trade-offs to emerge. We then tested a priori if it would be possible to predict whether or not individual cows would belong to groups of dairy cows that are statistically similar, according to their expression of health, growth, body reserves, milk production and reproductive traits at first lactation.

MATERIALS AND METHODS

Description of farming facilities and management

The data used in this study were collected between during the indoors season (from November to mid-April, between 1.0 and 1.5 tDMI•animal - 1 •year -1 of concentrate) and grazed for the rest of the year.

Housing facilities consist of a naturally ventilated barn subdivided into cubicles, with a maximum capacity of 110 cows. Storage of raw data [milk production and composition, body weight (BW), body condition score (BCS), fertility and health events] was performed using the data collection and management software Datalogic Aladdin TM (Datalogic S.p.A, Bologna, Italy). Milking was done year-round inside the barn with an automatic milking parlor twice a day (at 06:00 and 16:00). During the indoor season, forage was offered individually twice a day, and records of forage intake were obtained on an individual basis (forage offered minus refusal weighed three/four days a week). The amount of concentrate that was offered to the cows on a given week was determined to cover their nutritional requirements as quantified by the INRA feeding system (INRA, 1988) based on their milk production in the previous week. Concentrate intake was monitored on an individual basis using an automatic feeder (Feeding Station FSC-40 DeLaval, Elancourt, France).

Calving strategy was mostly seasonal, from cows artificially inseminated during the winter months. After a voluntary waiting period of at least two months, natural occurrence of heat (or estrus) was detected visually by trained staff in charge of the cows. Cases of anestrus were only dealt with for heifers, which were treated with hormones. Then, an artificial insemination campaign was conducted on the farm from November to February (4 months). Conception date was estimated from the size of the fetus obtained via ultrasound tests, and was further confirmed by comparing this date with the closest date of an observed heat event or an artificial insemination event. Cows that did not show signs of heat or were not confirmed as pregnant at the beginning of spring had a chance to be naturally impregnated by a bull, which was present in the herd during the entire grazing period. If a given cow was not confirmed pregnant at the end of the grazing/reproductive period, that cow was culled.

Animal data

Data on dairy performance, BW and BCS, fertility and health events along with feeding were collected at first lactation, and culling dates and reasons were recorded at the end of the productive career. Individual milk yield (MY, kg•d -1 ) was recorded daily, with the use of recording jars (DeLaval®, Ghent, Belgium) which rested on loading cells. Milk protein (MP, g•kg -1 ) and fat (MF, g•kg -1 ) contents, were analyzed weekly during 4 consecutive milkings (two consecutive days) using an infrared method (Milkoscan 4000, Foss System, Hillerod, Denmark). All cows were weighed with an electronic scale (Mettler-Toledo SAS, model Mettler Spider 2, Viroflay, France) after the morning milking for two consecutive days at least once a month. BCS was estimated through palpation using a scale from 0 to 5 points, with 0.25 interim steps [START_REF] Bazin | Grille de notation de l'etat d'engraissement des vaches pie noires[END_REF], by two trained staff members at least once a month. For each cow, reproductive status at first lactation was assessed by accounting for the total number of heat events (heatcount), the interval between calving and first heat (intcalving-heat, d), the total number of artificial inseminations or natural inseminations by a bull ((AI + NI)count), the number of successful inseminations (conceptioncount) and the interval from calving to conception also commonly known as 'days open' (DO, d). In this study, the variable conceptioncount also included possible abortions that the cow had during the lactation.

The number of natural inseminations was defined by adding up the total number of times when a heat event was observed at the field by a staff member (under the assumption that the cow mated with the bull every time a heat event happened) from turnout (start of natural 1 along with their definitions, formulas and units.

Herd composition and culling strategy

Data used in this study included individuals for which both first calving and end of the last lactation occurred during the study period. Only cows culled at the second lactation or later were kept. Outliers were defined as cows with productive performance levels and calving-tocalving intervals out of the interquartile range ± 1.5 and then were removed. This yielded a total of 185 cows (114 Montbéliarde (MO) and 71 Prim'Holstein (HO)) that were used for further analysis.

Reasons for culling were: (1) low fertility (including empty, abortion, late calving and calving troubles); (2) poor health (including mastitis, lameness, acidosis, other infections and overall health status); (3) poor production performance, (4) herd renewal (including old age), and (5) accidents on farm.

Statistical analyses

All analyses in this study were performed with the statistical software R © (version 3.5.2, R Development Core Team, 2020).

With the aim of reducing the number of variables to perform a principal component analysis (PCA), a first variable selection was done to avoid variance inflation due to collinearity, i.e., a multicollinearity test was performed among selected variables. We used the variance inflation factor (VIF) as criterion to judge the presence of multicollinearity. The variable with the highest VIF was excluded, and new VIFs were calculated again. These steps were repeated until only variables that had a VIF ≤ 10.0 remained [START_REF] Hair | Multivariate Data Analysis, 7th Editio[END_REF]. All remaining variables were then normalized and scaled (between -1 and 1) prior to the analysis steps described below in order to avoid the variance inflated effect due to heterogeneous scaling of variables. Using the R © package FactoMineR v2.3 [START_REF] Husson | Package "FactoMineR": Multivariate Exploratory Data Analysis and Data Mining[END_REF], we then performed PCA to determine which variables accounted for most of the variability found among individuals. With this, variables were further omitted by using a stepwise suppression criterion of significant correlation between a variable with any of the first three principal components. Only variables showing a significant r ≥ 0.50 (P ≤ 0.05) with at least one given component were kept.

Once identified through PCA, cluster analysis of individuals was performed using the Agglomerative Hierarchical Clustering (AHC) technique. The decision regarding the optimum number of clusters was based on an analysis of within-cluster variance, by detecting the minimum number of clusters after which an abrupt change in variance occurred.

Then, each variable (at original scale) was submitted to comparison tests in two steps: (1) the difference between the mean of each cluster and the respective mean obtained from the entire experimental herd was tested. For this process, the one-sample two-sided t-test was used for normally distributed variables, and the one-sample Wilcoxon signed rank test was applied to other variables; (2) between-cluster comparisons were performed using one-way ANOVA (considering cluster as a factor) for normally distributed variables and the Kruskal-Wallis test for the rest, followed by their respective post-hoc tests (Tukey and Wilcoxon signed rank tests, respectively).

We performed a trade-off analysis of groups of cows in each cluster by biologically interpreting the significant trait deviations (according to the cluster mean compared with the experimental herd mean) in each cluster. A trade-off was identified in a cluster when at least two traits presented statistically significant deviations with an opposite sign. When a trade-off occurred, we quantitatively profiled the group of individuals of that cluster according to the nature of the observed trade-off. Finally, we interpreted deviations in productive lifespan as a consequence of the combination of other traits and their trade-offs. The multi-trait profiling method used in this study was inspired by the approach proposed by [START_REF] Ollion | Assessing the diversity of trade-offs between life functions in early lactation dairy cows[END_REF]. For illustrative purposes, one paragon (representative cow) per cluster was selected based on two criteria: (1) the shortest geometrical distance (calculated with individual coordinates obtained from PCA) of the individual to the cluster barycenter, and (2) the individual representativeness of cluster characteristics, i.e. whether or not paragons resemble average individuals per cluster.

The dependency of contingents of cows per cluster were compared to the contingents per reason for culling (poor fertility, poor health and other causes including low milk yield, old age, herd renewal and accidents), using Fisher's exact test.

The ability to discriminate cows according to their lifespan profile obtained in the previous step was tested through discriminant analysis. We used the package 'psych' v2.0.9 R © [START_REF] Revelle | Package "psych": Procedured for Psychological, Psychometric and Personality Research[END_REF] to run our linear discriminant (LD) regression, using the cluster classification obtained from PCA as our dependent variable.

We first calibrated the LD model with the entire available dataset, and then applied the k-fold cross-validation method (k=4), also with the aid of the R © package 'psych' v2.0.9. The quality of both calibrated and cross-validated models was further assessed through sensitivity (proportion of correctly identified individuals per class, or true positives), specificity (proportion of individuals that are correctly identified as not belonging to a specific class, or true negatives), precision (sensitivity divided by the number of true negatives) and accuracy (overall ability of the model to predict class correctly), as defined by [START_REF] Fawcett | An introduction to ROC analysis[END_REF].

In this study, results from all statistical analyses were considered significant when P≤0.05.

RESULTS

Table 2 provides results for variables describing milk production, BW, BCS, health status, fertility and feeding at first lactation, and longevity traits. Although we initially kept the "breed" variable in our statistical analyses, it did not turn out to be significant. For that reason, data per MO and HO were pooled together, and breed effect was not further considered in this study.

Descriptive summary statistics on milk production, BW, BCS, health status, fertility and feeding at first lactation are provided in Table 3.

The main reason that we generated as many as 32 variables (Table 1) to compose our dataset was to increase the chances that, after a variable selection process, at least one variable representative of the most important life functions would remain. Indeed, 15 variables eventually remained after PCA (Table 3), representing the main life functions, as defined by [START_REF] Friggens | Towards a biological basis for predicting nutrient partitioning: the dairy cow as an example[END_REF]: survival (productive lifespan), growth (BWcalving, BWloss), body reserves (BCSmin), reproduction (heatcount, (AI+NI)count, intcalving-1st heat, DO), and lactation (through productive traits MY , cumulativeMYpeak, cumulativeMY4wk, Persistency, MP ), plus feeding (concentrate ).

Health variables were eliminated in the selection process for PCA to avoid inflation problems and do not appear in Table 3. However, they were considered again when performing profile analysis and LD regression, where inflation is not an issue.

PCA and clustering

Through the variable reduction process involved in PCA, we obtained 15 variables from the initial 32 variables considered in this study (Table 1). A total of four principal components (PC) with eigenvalues ≥ 1 accounted for 68% of the total variability of the dataset ( PC3 separated cows according to high protein content in milk (MP ) and their high ability to maintain body reserves (BWcalving, BWmin and BCSmin). Finally, PC4 mainly distinguished cows that presented longer times to restart estrus (int ! " ).

From the AHC results, individual cows were differentiated into three clusters. The partitioning of cows per cluster was 80, 28 and 77 individuals for clusters 1, 2 and 3, respectively.

Cluster 1 was characterized by cows that had a productive lifespan that was 11% longer than the experimental herd mean. They also produced 10 to 24% (for MY and cumulativeMYpeak, respectively) less milk at first lactation than the average. The average lengths of lactation and the dry period were significantly 5 and 12% shorter for cows in cluster 1, as compared to the herd mean, respectively. Cows from cluster 1 had a BWcalving 2% lower than the herd mean, but a BCSmin 8% higher. Cluster 1 was characterized by cows that had %TMF/MP>1.4 and lamenesscount 30 and 44% lower than the herd mean, respectively. These cows were more fertile than the herd mean (heatcount and (AI+NI)count decreased by 21 and 28%, respectively) and became pregnant earlier than the herd mean, as indicated by shorter int ! " , int 2!3 " # , int # !$" # , int ! " # !$" # and CI (6 to 23% shorter than the herd mean). Concentrate consumption was significantly lower in this cluster compared to the herd mean (-16%).

Cluster 2 was characterized by cows that had a productive lifespan 27% shorter than the herd mean. Average milk production of these primiparous cows was, in general, similar to that of the herd mean, except for variable MY , which was reduced by 13% compared to the herd mean. Moreover, cows from cluster 2 had significantly longer DIM (+16%) which was followed by a longer dry period (intdry, +29%) than the herd mean. They took a longer time to express first heat after calving (intcalving -1st heat, +27%), had lower fertility than the herd mean, presenting more heatcount and (AI + NI)count (76% and 91%, respectively) and finally showed longer DO than the herd mean (+59%). Health variables observed in this cluster did not differ from the herd mean.

Cows from cluster 3, which was characterized by a productive lifespan that did not differ from the herd mean, presented better milk production performance in relation to the herd mean, by presenting higher MY , MY , cumulativeMY305d, cumulativeMYpeak, cumulativeMY4wk and cumulativeMY24wk (increased by between 14 and 27%). They showed lower persistency of milk production than the herd mean. These cows also presented BWcalving and BWmin that were significantly higher than the herd mean by 3 and 2%, respectively.

Concerning health variables, cows in cluster 3 presented 11% lower %TMF/MP<1.0, but presented increased mastitiscount and lamenesscount (30 and 35%, respectively) compared to the herd mean. They presented an intcalving-1st heat 14% higher than the herd mean, but an int1st heat- conception 18% shorter, resulting in an int # !$" # similar to the herd mean. Cows in this cluster ate more concentrate than the herd mean (+18%).

Reasons for culling per cluster

The majority of cows in the herd were culled due to poor fertility (59%) or poor health (26%) and only 15% were culled due to other reasons like poor milk production, herd renewal, old age or accidents. This is further confirmed by the contingents of cows per cluster and reasons for culling that are presented in Table 5. Poor fertility was the leading reason for culling cows in clusters 1, 2 and 3, but was more frequent in cluster 2 (49%, 75%, 64% respectively). The second reason for culling was poor health with 33%, 25% and 22%, respectively. Culling for reasons other than reproduction or health was observed in clusters 1 and 3 (18% and 15%, respectively), but not in cluster 2.

Discrimination of cows according to clustered profiles

Results from LD analysis revealed that between-class variance explained by linear discriminants 1 and 2 was 56 and 44%, respectively. Figure 1 illustrates how individual cows per cluster are distributed in a Cartesian plane which has LD1 and LD2 as orthogonal axes.

Figure 2 depicts the variables' contributions to LD1 and LD2. The main positive contributions to LD1 came from MY , DO, concentrate and (AI+NI)count, while the main negative contributions to this LD came from persistency and cumulativeMY4wk.

The main positive variable contribution to LD2 came from MY , while the main negative contributions to this LD were from variables (AI+NI)count, intdry, heatcount, persistency and DO.

Figure 1 suggests that LD1 discriminates clusters 1 and 3 better, while LD2 discriminates cluster 2 from both clusters 1 and 3 better.

Between-cluster comparison tests (Table 4) confirm the effect of LD1 and LD2 in discriminating clusters. In other words, means for clusters 1 and 3 were significantly different for the main variables that describe LD1 (MY , concentrate , Persistency, CumulativeMY4wk and (AI+NI)count). Means from cluster 2 were significantly different from those of clusters 1 and 3 for the variables that mostly contribute to LD2 (mainly heatcount, (AI+NI)count, intdry and DO).

Results from cross-validation of the LD model indicate an overall accuracy of 86%. Other model quality measures are presented in Table 6.

DISCUSSION

Characterization of trade-offs and profiling clusters of dairy cows

Cluster analysis revealed different strategies of individual cows regarding priority given to lactation, reproduction and health performance at first lactation that ultimately were associated with significant differences in productive lifespan. A trade-off between milk production and functional traits (reproduction and health) at first lactation existed in cluster 1 in parallel to a trend to higher productive lifespan. Thus, productive lifespan was observed to be longer for cows in cluster 1 in which reproduction, health and BCS safeguarding were prioritized to the detriment of milk production at first lactation. If we suppose that such priorities remain in successive lactations, the higher longevity observed in this cluster is not surprising as reproduction and health have been reported to be the main reasons for culling in dairy herds [START_REF] Ahlman | Culling reasons in organic and conventional dairy herds and genotype by environment interaction for longevity[END_REF][START_REF] De Vries | Review: Overview of factors affecting productive lifespan of dairy cows[END_REF]. Other studies have reported positive correlations between longevity and reproduction traits (De Vries and Marcondes, 2020; O'Sullivan et al., 2020), good health [START_REF] Alvåsen | Farm characteristics and management routines related to cow longevity: a survey among Swedish dairy farmers[END_REF][START_REF] Shabalina | Influence of common health disorders on the length of productive life and stayability in German Holstein cows[END_REF] and ability to preserve body reserves [START_REF] Zavadilová | Effect of type traits on functional longevity of Czech Holstein cows estimated from a Cox proportional hazards model[END_REF]. Jimenez-Krassel et al.

(2015) found that young adult dairy heifers with higher concentrations of anti-Müllerian hormone (which is positively associated with follicle count, ovarian function, and fertility) had a longer productive lifespan. O'Sullivan et al. ( 2020) reported that survival of Holstein-Friesian cows in pasture-based systems has been significantly improved in the genetic divergent group of cows presenting the highest Irish Economic Breeding Index (EBI) compared to the average. These "top" EBI cows achieved greater pregnancy rates, presented fewer services per cow and shorter intervals between mating start date and pregnancy compared to the average, and, like cows in cluster 1, they maintained higher mean BCS during lactation. Our results are also in agreement with those of [START_REF] Shabalina | Influence of common health disorders on the length of productive life and stayability in German Holstein cows[END_REF], who estimated negative genetic correlations between longevity and occurrence of common health disorders in German herds of Holstein cows. Furthermore, [START_REF] Alvåsen | Farm characteristics and management routines related to cow longevity: a survey among Swedish dairy farmers[END_REF] found that Swedish farmers who received more visits from the veterinarian had dairy herds with shorter longevity.

Cows in cluster 2 succeeded in conceiving at first lactation (as all cows considered in this study), but with a large delay compared to the herd mean (+48 days). Consequently, these cows re-calved with a delay of more than 2 months (+72 days for CI) compared to the herd mean and started their subsequent breeding season later. In a pasture-based system, the reproductive season is often time-bound in order to make lactation match grass availability [START_REF] Delaby | Robust animals for grass-based production systems[END_REF]. In other words, seasonal or grouped calving derive from the strong constraints that characterize pasture-based dairy farms, with a reproductive season of only 3 months in order to synchronize the beginning of lactations with the period of high availability of forage [START_REF] Roche | New considerations to refine breeding objectives of dairy cows for increasing robustness and sustainability of grass-based milk production systems[END_REF]. Consequently, cows of cluster 2 had less chance of becoming pregnant at the second and subsequent lactations within the period defined by the farmer, first because their breeding season is shortened compared to the herd mean and second because they need more time to conceive. This is in agreement with the reasons for culling observed in this cluster with 75% of cows culled because of reproduction issues. [START_REF] Reimus | On-farm mortality and related risk factors in Estonian dairy cows[END_REF] also

showed that Estonian dairy cows with lower conception rate presented higher mortality, i.e.

unassisted death and euthanasia. No trade-off between reproduction and any other life function was observed in our study for cluster 2, except that they showed slightly higher milk solid contents (MF and MP) than the herd mean. In fact, Bedere et al. (2017a) found that, when underfed, high-yielding cows that presented higher milk solid content (MF + MP) were more prone to fertility failure, which might explain why they presented conception delays.

When reproduction is not a prioritized life function, the risk for the cow is lower longevity and for the farmer it is inability to use voluntary culling basing on economic reasons. No cows in cluster 2 were culled for reasons other than biological ones (poor fertility or poor health).

Cows in cluster 3 were able to manage high levels of milk yield regarding no trade-off with reproduction as they succeeded in conceiving with no delay compared to the herd mean. Such an ability to maintain reproduction was observed despite a higher calving to first heat interval (+10 days) and has also been reported in high-producing cows by Bedere et al. (2017b), who considered this to be due to delays in commencement of luteal activity and first service.

However, the productive lifespan of these cows did not differ from the herd mean, probably because of high milk yield traded-off with health. Higher occurrence of mastitis and lameness observed in this cluster at first lactation may have increased their culling risk due to poor health, as also described by [START_REF] Shabalina | Influence of common health disorders on the length of productive life and stayability in German Holstein cows[END_REF] for German Holstein cows. However, the proportion of cows culled for poor health reasons was not higher for cows in cluster 3 (22%) than for those in clusters 1 (33%) and 2 (25%), which probably means that greater health issues observed at first lactation for cows of cluster 3 would have not been amplified for subsequent lactations.

The difference in concentrate intake ( concentrate ) that was observed between clusters is explained by the fact that profiles were characterized by different average levels of milk production. In our study, the feeding system used to define amounts of concentrate to be offered per cow and per day was adjusted to the cow's milk yield the week before. The fact that we did not identify any difference in forage intake between clusters is probably because our analysis considered data from first lactation only. [START_REF] Grandl | Biological implications of longevity in dairy cows: 1. Changes in feed intake, feeding behavior, and digestion with age[END_REF] observed that the average organic dry matter intake of dairy cows increased with age.

Profiling reveals that high productive lifespan is mainly associated with reproductive success.

Poor fertility and poor health were the main reasons why cows were ultimately culled from the experimental herd. In this study, culling due to old age or economic reasons was present in clusters 1 and 3, but not in cluster 2. Cows culled for old age or herd renewal are probably the most robust cows because they went on "doing what they have to do" throughout their lives, without showing reproductive or health problems exceeding the tolerance threshold defined by the farmer. The dependence we found between contingents of cows per reason for culling (reproduction, health and economic) and per productive longevity profile (low, medium and high) indicates that the individuals that presented unfavorable milk production and reproductive or health traits at first lactation were ultimately culled for that same reason.

For instance, cows from the low productive longevity profile (cluster 2) were mostly culled due to reproductive issues, while they also presented the poorest fertility at first lactation.

These outcomes sustain the approach tested in this study and already suggested by [START_REF] Adriaens | Productive life span and resilience rank can be predicted from on-farm first-parity sensor time series but not using a common equation across farms[END_REF], that productive longevity ranking within a herd can be predicted by phenotyping several traits at first lactation. What our study adds to this idea is that clustering cows at first lactation from the way they express life functions could also help to predict the reason for their culling.

Study results are bound to grassland-based dairy farming systems with seasonal calving, where reproduction represents is a major cause of culling. This is particularly the case in France, Ireland and New Zealand, which is non-negligible. In other parts of the world it seems that health is the major cause for culling, while reproduction comes in second place [START_REF] De Vries | Review: Overview of factors affecting productive lifespan of dairy cows[END_REF].

Profiling of cows at first lactation with low, medium and high productive lifespan is predictable

Clustering and profiling of cows allowed us to observe the trade-offs between productive lifespan, productive and functional traits a posteriori, i.e., without previously hypothesizing that productive lifespan is explained by such trade-offs. Conversely, the LD model obtained in this study stems from an a priori testing of the hypothesis that productive longevity and functional traits measured at first lactation are linked.

Discrimination by LD1 and LD2 of groups of cows agrees well with the profile analysis performed above. Although the variable productive lifespan was not directly present in the discriminant analysis, it was indirectly represented in the classification of clusters. This was clearly seen by taking a look at the selected individual paragons (Figure 3), where all individuals had distinct productive lifespans.

The main contributions to LD1 and LD2 are shared between milk production and reproduction traits, which have opposing or polarizing effects. A distinction should be made between production variables that describe the entire lactation cycle, such as MY , and variables that describe the geometry or shape of the lactation curve, such as persistency and cumulativeMY4wk. Indeed, two cows with identical MY may have lactation curves of different shapes.

The polarizing effect seen in LD1 (high overall productivity, longer lactation and poorer fertility, as opposed to high productivity at the beginning of lactation, with a more accentuated drop in the lactation curve after the peak) is supported by [START_REF] Muir | Genetic Relationships Between Persistency and Reproductive Performance in First-Lactation Canadian Holsteins[END_REF], who showed that primiparous cows which had more calving or conception difficulties also presented more persistent first lactation. The same idea is further supported in this study, as illustrated by paragons for clusters 1 and 3 (Figure 3).

The inverse production vs. reproduction relationship comprised in LD2 is well known and has been mentioned before in several studies on dairy cows, such as [START_REF] Nebel | Interactions of High Milk Yield and Reproductive Performance in Dairy Cows[END_REF], [START_REF] Leroy | Reduced Fertility in High-yielding Dairy Cows: Are the Oocyte and Embryo in Danger? Part II Mechanisms Linking Nutrition and Reduced Oocyte and Embryo Quality in High-yielding Dairy Cows*[END_REF] and [START_REF] Pires | Physiological adaptations and ovarian cyclicity of Holstein and Montbéliarde cows under two low-input production systems[END_REF]. Further evidence in this study that explains the effect of LD2 is the contrast between the paragons for cluster 2 with longer lactation and higher milk yield associated with greater (AI+NI) counts and longer dry period, and the paragons from clusters 1 and 3 (Figures 3A and 3C,respectively).

Concerning the predictive quality of the LD model, its overall ability to distinguish individuals from clusters 1 and 2 (the clusters with high and low productive lifespan, respectively) with 97 and 86% correct classification is considered good for calibrated and cross-validated models, respectively. This shows that despite the fact that the variable productive lifespan was not directly included in the model, its influence on clustering individuals was significant.

The discriminant analysis with cross-validation suggests that our approach is relevant for moving from the average level per cluster/profile to the individual level in order to classify the animals. Our good results in cross-validation reinforce the relevance of the grouping approach provided by LD analysis and the effectiveness of the selected variables in predicting productive lifespan.

CONCLUSIONS

In this study, we revealed a clear association between productive longevity and trade-offs between life functions (milk production, reproduction, BW or body reserves lability and health) assessed at first lactation. Such an association was observed irrespective of the two breeds used in this study (MO and HO) and for cows managed in a grassland-based dairy farm. Cows with good reproductive performance but moderate milk production during first lactation had a longer productive lifespan, whereas cows with moderate milk production but poor reproductive performance had a shorter one. Productive lifespan can be predicted from a trade-off strategy at an individual level through discriminant analysis, by correctly classifying the group membership of individuals. Key variables to predict lifespan group membership were the overall milk production at first lactation, traits characterizing the shape of the lactation curve and fertility traits.

Our approach allows the prediction of longevity ranking of cows from traits measured at first lactation and provides a means for the early identification of those cows that have a greater probability of living longer and thus for the selection of robust cows within the herd. Cumulative MY (kg) during a reference lactation cycle of 305 d ( 7)

cumulativeMY $! = ∑ MY 6G @ABC 6
Cumulative MY (kg) up to the lactation peak at Tpeak ( 8)

cumulatuveMY ) = ∑ MY 6) 6
Cumulative MY (kg) up to the 4 th week of 1 st lactation ( 9)

cumulativeMY H) = ∑ MY 6H) 6
Cumulative MY (kg) up to the 24 th week of 1 st lactation ( 10) Average milk protein (MP) content for the entire lactation cycle (g•kg -1 ) (13)

Persistency = KLMY $! -MY M" N LT $! -36N Q K Persistency of
MF = ∑ MF 6 6
Average milk fat (MF) content for the entire lactation cycle (g•kg -1 ) ( 14)

%T >T >U ⁄ W . = Y if MF MP ⁄ < 1.0, ,100 n ⁄ 1 ∑ i 6 6 if MF MP ⁄ ≥ 1.0, 0
Percentage of time during the total lactation period when the ratio MF/MP was lower than 1.0 (higher risk of acidosis) (15)

%T >T >U ⁄ _ .) = Y if MF MP ⁄ > 1.4, ,100 n ⁄ 1 ∑ i 6 6 if MF MP ⁄ ≤ 1.4, 0
Percentage of time during total lactation period when the ratio MF/MP was higher than 1. MY is the average daily milk yield for the entire lactation cycle; MY is the average daily milk yield for up to a standard lactation cycle of 305 days; DIM stands for days in milk; MY $! is the milk yield corresponding to the time when the lactation cycle reaches its peak; T $! is the time in the lactation cycle when the lactation peak occurs; Int DE is the time interval during a lactation cycle of calving j-1 and before calving j when the cow was no longer being milked; cumulativeMY is the cumulative milk yield during a reference lactation cycle of 305 days ; cumulativeMY $! is the cumulative milk yield up to the lactation peak; cumulativeMY ) is the cumulative milk yield up to the 4 th week of lactation; cumulativeMY H) is the cumulative milk yield up to the 24 th week of lactation; Persistency is the persistency of the lactation curve, considering as a reference the 36 th week of the lactation cycle; MY rise rate is rate of rise in daily milk yield between the 1 st week of lactation and the lactation peak; MP is the average milk protein content over the entire lactation cycle; MF is the average milk fat content over the entire lactation cycle; %T >T >U ⁄ W . and %T >T >U ⁄ _ .) represent the percentages of time when the ratio MF/MP (milk fat content divided by milk protein content) is less than 1.0 and higher than 1.4, respectively, as indicators of acidosis and ketosis; BW is the average weekly body weight at calving; BW 3 is the minimum average weekly body weight over the lactation; BCS 3 is the minimum average weekly body condition score achieved during the lactation period; mastitis #+ " and lameness #+ " are the total counts of mastitis and lameness, respectively, that occurred during the lactation cycle; med. int #+ " is the total count of medical interventions that occurred during the lactation cycle; heat #+ " is the total count of observed heat events; ,AI + NI1 #+ " is the total count of artificial and/or natural inseminations per lactation; int ! " is the interval between calving and first heat; int CI stands for calving interval; forage is the average forage intake up to 10 weeks after calving; concentrate is the average concentrate intake up to 10 weeks after calving; Differences followed by "'", "*", "**" and "***" have P-values of <0.05, <0.01, <0.001 and <0.0001; "ns" means that cluster mean is not significantly different from the entire herd mean; Letters stem from results of post-hoc tests:

Tukey and Wilcoxon signed rank tests for normally and non-normally distributed variables, respectively (statistical significance at P≤0.05).

Variables presenting same letters per cluster denote that cluster means are not significantly different. 86) 97 ( 96) 96 (89) *Values in bold represent the percentage of individuals from a given cluster which were misclassified in other clusters by the linear discriminant method; **Accuracy value obtained from cross-validation is significantly different from zero at a confidence of 95%.

Figure 3 .

 3 FIGURE CAPTIONS

  the lactation curve (kg•d -1 ), considering as a reference the 36 th week of the lactation cycle (11)MY rise rate = LMY $! -MY 2" N T $! Q MY rise rate (kg•d -1) is the rate of rise in MY between the 1 st week of lactation and the lactation peak

  first insemination; DO stands for days open; int1stheat-conception stands for time interval between first heat and conception;

Table 3

 3

	),

  Raw data were presented on a daily basis, but to cope with missing data for the variables that were not measured daily, average values per week after calving were calculated. For the variables presented in this table, we still present data on a daily basis by dividing weekly based data by 7.

	int	2!3	" # = T		-T	2!3	" #	Time interval (d) between calving and first insemination	(26)
	DO = T	-T # !$" #				Days Open: time interval (d) between calving and conception	(27)
	int ! " # !$" # = T ! " -T # !$" #	Time interval (d) between first heat and conception	(28)
	CI = T	o -T	o					Calving Interval: interval (d) between two subsequent calvings (j-1 and j)	(29)
									Average concentrate intake (kgDM•d -1 ) up to 10
	concentrate	= ,1 10 ⁄ 1 ∑	6 6	,concentrate 1	weeks after calving; concentratei is the average concentrate intake (kgDM•d -1 ) at week i after	(30)
									calving
	forage	= ,1 10 ⁄ 1 ∑	6 6	,forage 1	Average forage intake (kgDM•d -1 ) up to 10 weeks (kgDM•d -1 ) at week i after calving after calving; foragei is the average forage intake	(31)
	productive lifespan = ageculling -age1st calving	Time period between first calving and culling event (d)	(32)
									4 (risk	(16)
	BW BW 3 = min,BW 1| 6 = BW 6G eBfghij 6						of ketosis) Average weekly body weight (kg) at calving Minimum average weekly body weight (kg) achieved during the 1 st lactation	(17) (18)
	BCS 3 = min,BCS 1| 6 6						between 0 and 5, with 0.25 interim steps. Minimum average weekly body condition score achieved during the lactation period. Ranging	(19)
	mastitis #+ " = ∑ mastitis #+ ", 6 6	Total count of mastitis that occurred during a lactation cycle of n weeks	(20)
	lameness #+ " = ∑ lameness #+ ", 6 6	Total count of lameness that occurred during n weeks after calving	(21)
	med. int #+ " = ∑ med. int #+ ", 6 6 heat #+ " = ∑ heat #+ ", 6 6 ,AI + NI1 #+ " = ∑ AI 6 6 + ∑ NI 6 6 int ! " = T -T ! "	Total count of medical interventions that occurred during n weeks after calving Total count of observed heat events Total count of artificial and/or natural inseminations per lactation Time interval (d) between calving and first heat	(22) (23) (24) (25)

*

Table 2 .

 2 Descriptive statistics for age at first calving, milk production, body weight, body condition score, fertility and health traits and longevity of primiparous cows. N=185 represent the percentages of time when the ratio MF/MP (milk fat content divided by milk protein content) is less than 1.0 and higher than 1.4, respectively, as indicators of acidosis and ketosis. AI and NI stand for artificial and natural inseminations, respectively.

	Mean ± SE Maximum Minimum

Table 3 .

 3 Summary of results from PCA for the variables that remained in the analysis. Bold numbers represent the highest significant correlations between variables and components. MY is the average daily milk yield for the entire lactation cycle; Int DE is the time interval during a lactation cycle of calving j-1 and before calving j when the cow was no longer being milked; cumulativeMY $! is the cumulative milk yield up to the lactation peak; cumulativeMY ) is the cumulative milk yield up to the 4 th week of lactation; Persistency is the persistency of the lactation curve, considering as a reference the 36 th week of the lactation cycle ;

	Eigenvectors of each variable, and correlations between variables with the first four principal
	components (PC)
		Significant correlations (95% probability)
	Variables	between variables and components	VIF
		PC1	PC2	PC3	PC4
	MY Int DE (d) (kg•d -1 ) cumulativeMY $! (kg) cumulativeMY ) (kg) Persistency (kg•d -1 ) MP (g•kg -1 ) BW (kg) BW 3 (kg) BCS 3 heat #+ " ,AI + NI1 #+ " int ! " (d) DO (d) (kgDM•d -1 ) concentrate productive lifespan (d)	0.69 0.54 0.77 0.86 -0.61 -0.31 0.32 0.21 -0.32 0.36 0.44 0.53 0.54 0.75 -0.30	-0.53 0.47 -0.23 -0.24 0.32 0.26 -0.40 -0.42 -0.21 0.69 0.75 0.18 0.70 -0.29 -0.34	-0.16 0.17 ---0.54 0.73 0.77 0.59 ---0.17 --	--0.24 --0.38 ------0.68 0.23 --	8.2 2.3 4.0 8.0 3.6 6.2 4.1 3.6 1.6 2.4 6.1 2.3 8.2 3.7 1.3
	Eigenvalue	4.4	2.9	1.9	1.0	-
	Cumulative variability (%)	29	48	61	68	-

MP is the average milk protein content over the entire lactation cycle; BW is the average weekly body weight at calving; BW 3 is the minimum average weekly body weight over the lactation; BCS 3 is the minimum average weekly body condition score achieved during the lactation period; heat #+ " is the total count of observed heat events; ,AI + NI1 #+ " is the total count of artificial and/or natural inseminations per lactation; int ! " is the interval between calving and first heat; DO stands for days open; concentrate is the average concentrate intake up to 10 weeks after calving.

Table 4 .

 4 Deviations of cluster means from entire experimental herd mean and post-hoc results for between cluster comparisons

	Variable	Herd average	Significant differences from population mean Cluster 1 Cluster 2 Cluster 3
	Sample size MY (kg•d -1 ) MY (kg•d -1 ) DIM (d) MY $! (kg•d -1 ) T $! (d) Int DE (d) cumulativeMY cumulativeMY $! (kg) (kg) cumulativeMY ) (kg) cumulativeMY H) (kg) Persistency (kg•d -1 ) MY rise rate (kg•d -1 ) MP (g•kg -1 ) MF (g•kg -1 ) %T >T >U ⁄ W . (%) %T >T >U ⁄ _ .) (%) BW (kg) BW 3 (kg) BCS 3 mastitis #+ " lameness #+ " med. int #+ " heat #+ " ,AI + NI1 #+ " int ! " (d) int 2!3 " # (d) DO (d) int ! " # !$" # (d) CI (d) forage (kgDM•d -1 ) (kgDM•d -1 ) concentrate productive lifespan (d)	185 18 18 320 26 8 78 5755 942 580 3640 -0.0573 11 38 33 8 8 628 564 1.2 0.4 0.5 1.5 1.8 1.8 72 88 113 41 398 9.6 5.7 1546	80 -3*** a -2*** b -17*** c -4*** c -1* b -9*** c -962*** c -228*** c -79*** c -601*** c 0.0120*** a ns a ns a 1' a ns a -2*** b -12' b ns b 0.1* a ns a -0.2' b ns a -0.4*** c -0.5*** c -17*** b -8** b -26*** c -10* b -26*** c ns a -0.9*** c 167* a	28 ns b -2*** b 50*** a ns b ns a 22*** a ns b ns b ns b ns b 0.0078' a ns a 1' a 1* a ns a ns a ns b ns b ns a ns a ns a ns a 1.4*** a 1.7*** a 20* a ns a 67*** a 48*** a 72*** a ns a ns b -424** c	77 3*** c 3*** a ns b 5*** a ns a ns b 995*** a 234*** a 84*** a 659*** a -0.0153*** b ns a ns a -1** a -1** b ns a 16* a 11' a ns b 0.1' a 0.2' a ns a ns b ns b 10* a ns a ns b -7' b ns b ns a 1.0*** a ns b

Table 6 .

 6 Quality assessment of calibrated and cross-validated (between brackets) linear

	discriminant models
	Cluster	1	Predicted classification (count)* 2	3	Sensitivity (%)	Specificity (%)	Precision (%)	Accuracy (%)
	1	77(70)	0 (0)	3 (10)	96 (94)	98 (89)	98 (96)
	2	1 (3)	27 (21)	0 (4)	96 (86)	100 (97)	100 (97)	97 (89**)
	3	1 (6)	0 (4)	76 (67)	96 (

ACKNOWLEDGEMENTS 23 This research received partial funding from GenTORE Project, a European Union H2020 Research and Innovation Program, under agreement No. 727213, and falls within the thematic area of the French government IDEX-ISITE initiative 16-IDEX-0001 (CAP 20-25). The authors would like to acknowledge Philippe Pradel and Matthieu Bouchon for their aid retrieving some of the data used in this study. Final thanks go to the anonymous reviewers for their helpful comments on earlier versions of the manuscript.

LIST OF TABLES

 Table 1. Definition of variables related to milk production, body weight, body condition score, health status, fertility and feeding at first lactation and longevity traits, their equations and units

Variable and equation

Definition and units

Average daily milk yield (MY) for the entire lactation cycle (kg•d -1 *); n is the number of weeks in a single lactation

(1)

Average daily milk yield (MY) for up to a standard lactation cycle of 305d (kg•d -1 *); 44 is the considered number of weeks comprised in a standard lactation

(2)

Mendes2

MY is the average daily milk yield for up to a standard lactation cycle of 305 days; Int is the time interval during a lactation cycle of calving j-1 and before calving j when the cow was no longer being milked; cumulativeMY is the cumulative milk yield up to the lactation peak; cumulativeMY is the cumulative milk yield up to the 4 th week of lactation; Persistency is the persistency of the lactation curve, considering as a reference the 36 th week of the lactation cycle; MP is the average milk protein content over the entire lactation cycle; BW ! "#$%& is the average weekly body weight at calving; BW '$% is the minimum average weekly body weight over the lactation; BCS '$% is the minimum average weekly body condition score achieved during the lactation period; heat !+,%-is the total count of observed heat events; .AI 0 NI2 !+,%-is the total count of artificial and/or natural inseminations per lactation; int ! "#$%&45 67 8 -is the interval between calving and first heat; DO stands for days open concentrate 5 is the average concentrate intake up to 10 weeks after calving.

Mendes3

A paragon is a representative cow selected per cluster based on two criteria: (1) the shortest geometrical distance (calculated with individual coordinates obtained from principal component analysis) of the individual to the cluster barycenter, and (2) the individual representativeness of cluster characteristics, i.e. whether or not paragons resemble average individuals per cluster.