L Tusell

A Legarra

M García-Tomás

O Rafel

J Ramon

M Piles

Genetic basis of semen traits and their relationship with growth rate in rabbits 1

Keywords: ejaculate quality, genetic parameter, growth, semen production, rabbit

ABSTRACT: This work aims to estimate the genetic parameters of seminal and production traits in a paternal line of rabbits selected for ADG during the fattening period. The considered traits were male libido (Lib) defi ned as successful mounting of an artifi cial vagina; presence of urine (Ur) and calcium carbonate deposits (Ca) in the ejaculate; semen pH; individual sperm motility (IM); the suitability for AI of the ejaculate (Sui), which involves the subjective combination of several quality traits; the average ejaculate volume (Vol); sperm concentration (Conc); and the average sperm production per ejaculate (Prod = Vol × Conc). The genetic relationship between all of these traits with ADG is also provided. Male libido and seminal data came either from routine evaluations of the ejaculates in an AI center or from 2 experiments in which bucks from the same population were used. Two consecutive ejaculates per male and per week were collected, leaving 7 d within weekly collections. A linear tri-trait model was used to analyze Conc, Vol, and ADG, whereas linear and threshold-linear 2-trait models were used to analyze male libido and the remaining seminal traits with ADG. A Bayesian approach was adopted for inference.

Approximately 38% of ejaculates were rejected for AI primarily due to low IM scores. Variables related to the quality of the ejaculate (Ur, Ca, pH, IM, Sui) and Lib were found to be lowly heritable (h 2 ranged from 0.04 to 0.11), but repeatable. This indicates performance of bucks for seminal quality traits and libido in AI centers would be more strongly affected by management practices rather than genetic selection. Semen production traits exhibited moderate values of h 2 (0.22, 0.27, and 0.23 for Conc, Vol, and Prod, respectively), suggesting the possibility of effective selection for these traits. A moderate to high negative genetic correlation (r g ; posterior mean; highest posterior density at 95%, HPD 95% ) was estimated between Conc and Vol (-0.53, HPD 95% = -0.76, -0.27). The ADG was estimated to have an h 2 of 0.16, to have a low, positive r g with Conc (0.21, HPD 95% = -0.03, 0.48), to have a low, negative r g with Vol (-0.19, HPD 95% = -0.47, 0.08), and to be genetically uncorrelated with all remaining traits analyzed. Therefore, selection for increasing ADG in paternal lines is expected to have no detrimental effects on Ur, Ca, pH, IM, Sui, and Lib and little to no effect on Conc, Vol, and Prod.

INTRODUCTION

Terminal bucks used for AI come from lines exclusively selected for growth traits (Baselga, 2004). Although they have good genetic merit for growth and feed effi ciency, no attention is placed on their reproductive performance and seminal characteristics. However, bucks used to produce doses to inseminate in commercial farms should also have good semen production and quality to potentially produce a large number of fertile doses. Optimal reproductive results are achieved with commercial AI in rabbits (e.g., 80.5% of conception rate and 10.34 number of total kits born per litter were obtained on average in French rabbit farms; [START_REF] Itavi | Gestion technico-économique des éleveurs de lapins de chair (programmes RENACEB et RENALAP) Résultats[END_REF], largely because the procedures used for AI in rabbits aim to maximize the probability of fertilization of oocytes. However, effi ciency of the AI centers seems to be far from optimal: the production and characteristics of semen are highly variable between collections (García-Tomás et al., 2006c), the ejaculate rejection rate in AI centers is large (38 to 52%; Brun et al., 2002a;[START_REF] Theau-Clément | Comparaison de la production spermatique de trois souches de lapins: Moyennes et variabilités. 10èmes[END_REF], and lack of libido and infertility are the most important causes for culling males in rabbit farms [START_REF] Rosell | Culling and mortality in breeding rabbits[END_REF].

To examine the possibility of genetic improvement for increasing AI dose production per buck, this work aims to estimate the genetic parameters of the following seminal and production traits: male libido, defi ned as successful mounting to an artifi cial vagina; presence of urine and calcium carbonate deposits in the ejaculate; semen pH; individual sperm motility; the suitability for AI of the ejaculate, which involves the subjective combination of several quality traits; the average ejaculate volume; sperm concentration; and the average sperm production per ejaculate. Additionally, genetic correlations between each semen characteristic with ADG were estimated to determine if selection for growth would have a correlated effect on male semen production and quality.

MATERIALS AND METHODS

The research protocol was approved by the Animal Care and Use Committee of the Institut de Recerca i Tecnologia Agroalimentàries.

Management

Bucks belonged to the Caldes line, which is selected for growth rate during the fattening period. They were bred and reared in the nucleus of selection in Caldes de Montbui (Barcelona, Spain). This farm has insulated roof and walls and cooling equipment to avoid animal exposure to extreme temperatures [e.g., the average temperature (SD) of the whole period in which the seminal data was collected was 19.5°C (3.04)]. After weaning at 32 d, males were housed in cages of 8 individuals with a photoperiod of 16 h of light/d. Animals were fed a commercial diet ad libitum (15.5% CP, 2.3% fat, 17.2% fi ber) until d 60. Average daily gain was measured during the fattening period (from 32 to 60 d old). After the fattening period, animals were individually housed and their feed was restricted to 180 g/d of another commercial diet (16% CP, 4.3% fat, 17% fi ber). Fresh water was always available.

Good health status (i.e., animals were free from illness, injury, or pain) was the sole criterion used for selecting young bucks to be used both in the AI center and in the experimental studies where males were evaluated for seminal production traits.

Selected males were trained at 4.5 mo of age to use an artifi cial vagina for ejaculate collection. A homemade polyvinyl chloride artifi cial vagina containing water at a temperature of 50°C was used. For ejaculate collection, the artifi cial vagina was hand-held beneath a doe with the open end pointed in a caudal direction and together they were placed inside the cage of the buck to allow the male to mount the doe. As the buck began to mount, penetration of the male penis into the artifi cial vagina was allowed and the ejaculate was collected in a tube connected to the artifi cial vagina (for further details of the procedure, refer to [START_REF] Morrell | Artifi cial insemination in rabbits[END_REF]. To train young males to use an artifi cial vagina for ejaculate collection, bucks were solicited for ejaculate collection once a week for the fi rst 2 wk. At 5 mo of age, males were considered sexually mature, bucks were solicited for ejaculate collection 2 times each week with an interval of 30 min between solicitations within a week. The time period between male solicitations to ejaculate collection of 2 consecutive weeks was 7 d. Semen data analyzed in this study are from sexually mature males.

Male Libido and Seminal Trait Data

All data used were collected in the IRTA AI center from routine evaluations to prepare AI commercial doses as well as specifi c experimental periods in which extra seminal measurements were collected. Male libido (sexual desire) was recorded as a binary trait (Lib = 1 if male successfully mounted; i.e., if male showed signs of falling off the doe when mating; 0 otherwise). Visual detection was used to identify and discard for further evaluation ejaculates contaminated with urine (Ur, assessed by the yellowish color of the semen instead of the typically translucent with white, gray color), calcium carbonate deposits (Ca, assessed by the presence of sandy sediments in the ejaculate), and blood (assessed by the pink or reddish color of the semen). Both the Ur and Ca traits were defi ned as binary: 1 = presence, 0 = absence. All gel plugs were removed. Ejaculates were stored at 37°C and were evaluated within 15 min after collection. The pH of the ejaculate was determined using a 507 Crison pH-meter (Crison Instruments, S.A., Alella, Barcelona, Spain). Ejaculate volume was determined using either a graduate tube or a micro-pipette.

Ejaculates were diluted 1:4 (vol/vol) in a commercial extender (Galap, IMV Technologies, Saint Ouen sur Iton, France). Aliquots (25 μL) were evaluated under a microscope with a phase-contrast optic (Nikon, Lewisville, TX) at 400× magnifi cation to assess individual sperm motility (IM). The IM of each sample was recorded and analyzed on a subjective scale that ranged from 0 to 5, which cor-responded to a percentage of sperm showing progressive movement: 0 to 10, 11 to 25, 26 to 50, 51 to 70, 71 to 90, or 91 to 100%, respectively [START_REF] Roca | Viability and fertility of rabbit spermatozoa diluted in Tris-buffer extenders and stored at 15 degrees[END_REF].

The suitability for AI of the ejaculate (Sui) was defi ned as a binary trait (1 = suitable, 0 = unsuitable for AI) and was assessed subjectively by a technician. An ejaculate was considered unsuitable for AI if presented any of the following characteristics: presence of Ur, blood, and Ca in the ejaculate, an IM score ≤2, the presence of a large number (approximately >50%) of dead spermatozoa, and the presence of clustered spermatozoa (both assessed by visual detection under microscope). Because an ejaculate could present more than one characteristic of unsuitability for AI, it could be discarded for AI for more than one reason at the same time.

Seminal data of Lib, Ur, Ca, IM, and Sui came from routine evaluations performed in the IRTA AI center from December 2001 to December 2009. Within this period, data of semen pH, ejaculate volume, and sperm concentration were collected in the same AI center in 2 experimental periods: January 2002 to December 2002 and January 2006 to September 2007. In the fi rst experimental period, sperm concentration was measured in each ejaculate using a Thoma Zeiss counting cell chamber. The analyzed trait was the sperm concentration per male on the day of collection (Conc). In cases where there were 2 ejaculates suitable for AI per male on a given day, Conc was calculated as (Conc 1 × Vol 1 + Conc 2 × Vol 2 )/(Vol 1 + Vol 2 ), where subscripts 1 and 2 are the fi rst and second ejaculate of the male on the day of ejaculate collection. In the second experimental period, Conc was measured using a counter of sperm cells (NucleoCounter SP-100, ChemoMetec A/S, Allerod, Denmark) after pooling semen from both suitable ejaculates for AI obtained from the same buck on the same day. No differences in Conc measurements have previously been encountered using either a hemocytometer or the Nucleocounter SP100 [START_REF] Theau-Clément | Evaluation de la concentration de semence de lapins selon 2 méthodes: Hématimètre et nucleocounter sp100[END_REF].

The average ejaculate volume per male each day (Vol) was calculated as (Vol 1 + Vol 2 /n), where n is the total number of ejaculates obtained per male on the day of collection (n = 1, 2). Collections in which the male successfully mounted the doe but no ejaculate was collected in the collection tube were considered as Vol = 0 mL. Whenever a nonzero Vol was obtained, the average sperm production per ejaculate (Prod) was determined by multiplying Vol by Conc.

Seminal traits Lib, Ur, Ca, IM, and Sui involved 883 males from approximately 740 litters, whereas pH, Vol, Conc, and Prod involved 541 males from approximately 470 litters. Total numbers of records obtained for each seminal trait are presented in Table 1. There were a total of 118,306 records of ADG, which correspond to all the data from the foundation of the line. The pedigree included 118,574 animals.

Model and Statistical Analysis

The seminal traits Lib, Ur, Ca, pH, IM, Prod, and Sui were each analyzed with ADG in a bivariate analysis to avoid selection bias in the estimates of the (co) variance components and to determine the relationship between ADG and the seminal traits. The Vol and Conc were analyzed jointly with ADG in a trivariate analysis to estimate the correlations among the 3 traits. Models were Gaussian for the continuous traits pH, IM, Vol, Conc, Prod, and ADG, whereas threshold models were used to analyze the binary traits Lib, Ur, Ca, and Sui. In this study, a multivariate analysis including all seminal traits and ADG was not attempted because not enough data of the seminal traits were available to allow the estimation of the increased number of parameters of the resultant model.

The threshold model assumes the observed binary responses are indicators of an underlying continuous random variable (l; liability; [START_REF] Falconer | The inheritance of liability to certain diseases, estimated from the incidence among relatives[END_REF] and a fi xed threshold, which divides the continuous scale into 2 intervals and forms the 2 response categories [START_REF] Wright | An analysis of variability in number of digits in an inbred strain of guinea pigs[END_REF]. 1 ADG = ADG during the fattening period (g/d).

2 Lib = male libido (binary trait: 1 = success to mount artifi cial vagina, 0 = failure).

3 Ur = presence of urine in the ejaculate (binary trait: 1 = presence, 0 = absence).

4 Ca = presence of calcium carbonate deposits in the ejaculate (binary trait: 1 = presence, 0 = absence).

5 pH = pH of the ejaculate.

6 IM = individual sperm motility (subjective scale from 0 to 5 corresponding to a percentage of sperm showing progressive movement of 0 to 10, 11 to 25, 26 to 50, 51 to 70, 71 to 90, or 91 to 100%, respectively).

7 Sui = ejaculate suitability for AI (binary trait: 1 = suitable, 0 = unsuitable).

The general model assumed for the analysis of the seminal traits was

= + + + + y X Z u W p W c e β
where y s is a vector of data for the corresponding continuous seminal traits (s = pH, IM, Vol, Conc, Prod) or a vector of liabilities corresponding to the binary seminal traits (s = Lib, Ur, Ca, and Sui), β s is a vector of systematic effects, u s is a vector of additive genetic effects, p s is a vector of random permanent environmental effects, c s is a vector of common litter environmental effects, and e s is a vector of residuals. Terms X s , Z s , W 1,s , and W 2,s are incidence matrices relating data with the systematic, random genetic, random permanent environmental effects, and common litter environmental effects, respectively. The systematic effects included in the model were order of the ejaculate (2 levels, fi rst or second; this effect was not included in Vol, Conc, and Prod models), year-season of collection (2-mo intervals: 33 levels for Lib, Ur, Ca, IM, and Sui, and 13 levels for pH, Vol, Conc, and Prod), and male age (7 levels for Lib, Ur, Ca, IM, and Sui: 4 to 6, >6 to 8, >8 to 10, >10 to 12, >12 to 16, >16 to 24, >20 mo and 3 levels for pH, Vol, Conc, and Prod: 4 to 6, >6 to 8, >8 to 10 mo).

The following model was assumed for ADG:

= + + + + y X Z u W p W c e β
where y adg is a vector of ADG measurements, β adg is a vector of systematic effects of ADG, u adg and c adg correspond to the same effects as those defi ned for seminal traits in the previous model. Because seminal traits and ADG were collected at 2 different time periods in the lifetime of a given buck and different records of seminal traits and male libido corresponded to only 1 ADG record, the residual (e*) was decomposed into 2 residual terms: e* = p adg + e adg , where p adg relates to the male permanent environmental effect of the corresponding seminal trait. This decomposition of the residual was carried out to increase data connectivity and to allow for estimation of a possible environmental correlation between ADG and each seminal trait and male libido. Covariance between p adg and p s can be estimated because there are repeated measures of the seminal traits and male libido, and therefore variance of p s can be estimated. Incidence matrices X adg , Z adg , W 1,adg , and W 2,adg relate growth data with the systematic, random genetic, random permanent environmental, and common litter environmental effects, respectively. The systematic effects included in the ADG model were yearseason of birth (105 levels), parity order (6 levels: 1, …, 5, >5 parities), and number of kits born alive in the litter to which the individual was born (8 levels: ≤5, 6, …, 11, ≥12 kits born alive). A Bayesian framework was adopted for inference. Denote The terms G, P, C are the different (co) variance matrices of the corresponding random effects defi ned above, and R is the residual (co)variance matrix. The joint posterior distribution of all parameters for the joint analyses of 2 continuous traits was The assigned prior distributions for the parameters of the models were

( ) ( ) ( )
( ) ( ) ( ) ( ) ( ) ( ) ( ) , ~| ~, | ~, | ~, m p kp N p N p N ⊗ ⊗ ⊗ u G 0 G A p P 0 P I c C 0 C I β
where k is a constant and A is the numerator relationship matrix of all the individuals. Bounded uniform prior distributions were assumed for β and the components of G, P, and C. For the binary traits, the threshold and the residual variance were fi xed to 0 and 1, respectively. For the continuous traits, the prior distribution for the residuals was ( ) ( )

| ~, , p N ⊗ e R
0 R I and bounded uniform priors were assumed for the elements of R.

The off-diagonal elements of R were always set to zero except in the trivariate analyses, where the elements corresponding to the residual covariance between Vol and Conc were estimated.

The marginal posterior distributions of the parameters of interest were derived from the joint posterior density of all the unknowns. The Gibbs sampler algorithm was used to estimate the marginal posterior distributions of the systematic effects and the (co)variance components using the TM software developed by [START_REF] Legarra | TM Threshold Model[END_REF]. Conditional distributions of the model parameters, necessary for the implementation of this algorithm, can be found in [START_REF] Sorensen | Likelihood, Bayesian, and MCMC Methods in Quantitative Genetics[END_REF]. Single chains of 1,000,000 and 3,000,000 iterations were run discarding the fi rst 250,000 and 1,000,000 iterations of each chain in the bivariate and trivariate analyses, respectively. Samples of the parameters of interest were saved every 100 rounds. The number of discarded samples was, in all cases, much larger than the required burn-in determined by the procedures of [START_REF] Raftery | How many iterations in the Gibbs sampler? Bayesian statistics 4[END_REF] and [START_REF] Geweke | Evaluating the accuracy of sampling-based approaches to the calculation of posterior moments[END_REF]. The sampling variance of the chains was obtained by computing Monte Carlo SE [START_REF] Geyer | Practical Markov chain Monte Carlo[END_REF]. Summary statistics from the marginal posterior distributions were calculated directly from the samples saved.

RESULTS AND DISCUSSION

Summary Statistics of ADG and Seminal Traits

Summary statistics for the seminal traits and ADG expressed as average percentage per buck are shown in Table 1. Means of the different seminal traits were within the range of values obtained for rabbits in previous studies reviewed by Alvariño (2000). The mean of ADG was also in accordance with those obtained in other rabbit lines [START_REF] Vogt | Selection experiment with domestic rabbits[END_REF][START_REF] Moura | Two-way selection for daily gain and feed conversion in a composite rabbit population[END_REF][START_REF] Lavara | Genetic parameter estimates for semen production traits and growth rate of a paternal rabbit line[END_REF]. A large portion of the seminal and Lib data came from routine evaluations performed in the AI center to prepare AI commercial doses; therefore, it is possible that there is some bias associated with these raw means. This potential for bias stems from the selection strategy followed by this AI center in which bucks are replaced based upon a lack of libido and poor sperm production. After 3 mo in production, buck replacement is only performed in older bucks, for reasons of lack of libido, poor sperm production and quality, and poor health status, the last being the most common reason for buck culling in the AI center. Because bucks with poor libido and poor sperm production and quality are only replaced after 3 mo in production; presumably this selection strategy will not have a large effect on the estimates of male variance components because enough data are collected per individual.

The percentage of successful collection rate with the artifi cial vagina was very high and in agreement with the one obtained in previous studies. This result contributes to refl ect the high adaptability of this species to be used in AI centers for dose production (Brun et al., 2002a[START_REF] Brun | Semen production in two rabbit lines divergently selected for 63-d body weight[END_REF][START_REF] Theau-Clément | Comparaison de la production spermatique de trois souches de lapins: Moyennes et variabilités. 10èmes[END_REF].

However, the effi ciency in the use of bucks in AI centers still seems to be less than optimal. A low percentage (62%) of ejaculates were considered suitable for AI, and this estimate is in agreement with similar studies performed on rabbits (Brun et al., 2002a[START_REF] Brun | Semen production in two rabbit lines divergently selected for 63-d body weight[END_REF][START_REF] Theau-Clément | Comparaison de la production spermatique de trois souches de lapins: Moyennes et variabilités. 10èmes[END_REF]García-Tomás et al., 2006a). These results suggest criteria for ejaculate rejection could be too restrictive or that, in some cases, it could be necessary to improve the qualitative characteristics of the ejaculate through management procedures or genetic selection. Low IM was the main reason for ejaculate rejection in this line (representing a 66% of total rejections), followed by the presence of Ca and Ur in the ejaculates (Figure 1). Sperm concentration and sperm motility are considered important variables infl uencing fertility [START_REF] Castellini | Effect of motile sperms inseminated on reproductive performance of rabbit does[END_REF]Brun et al., 2002b;[START_REF] Lavara | Do parameters of seminal quality correlate with the results of on-farm inseminations in rabbits?[END_REF]García-Tomás et al., 2006b). For that reason, individual or mass motility score or both are commonly used as criteria for ejaculate rejection in most AI centers. Conversely, individual sperm concentration is not commonly evaluated in rabbits because pooling ejaculates from several bucks is a common practice in rabbit AI centers. Thus, sperm concentration is only measured in the pooled ejaculate to perform the proper dilution to reach the standardized concentration of a heterospermic commercial dose.

Presence of Ur in the ejaculate has been found to be one of the principal factors for ejaculate rejection in rabbits (Brun et al., 2002a[START_REF] Brun | Semen production in two rabbit lines divergently selected for 63-d body weight[END_REF][START_REF] Theau-Clément | Comparaison de la production spermatique de trois souches de lapins: Moyennes et variabilités. 10èmes[END_REF]. Calcium carbonate deposits, which in this line were present in 17% of ejaculates, are also encountered in the bladder. However, sources of this Ca are still unknown and, to our knowledge, their presence has only been noted in a previous study carried out on the same experimental farm (García-Tomás et al., 2006a). In that study, presence of Ca was observed in another line but at a decreased frequency, indicating there could be differences between lines in the expression of this trait. To optimize the management of the bucks in the AI center, further research should be carried out to establish the origin and causes of the presence of Ca in the ejaculates. The presence of Ca in the ejaculates complicates the evaluation of the ejaculate for seminal quality in the laboratory and may have a detrimental effect on reproductive performance. The pH of the semen was within the range of values obtained in the same line and in other rabbit breeds (Brun et al., 2002b[START_REF] Brun | Paramètres génétiques des caractéristiques de la semence de lapin. 13èmes[END_REF]García-Tomás et al., 2006a). The pH could be considered a general semen quality indicator: the greater the concentration and motility of the spermatozoids in the ejaculate, the greater the production of lactic acid due to the greater metabolic activity and the lower the pH [START_REF] Coffey | Androgen action and the sex accessory tissues[END_REF][START_REF] Brun | Paramètres génétiques des caractéristiques de la semence de lapin. 13èmes[END_REF]. As a consequence, the value of this trait lies in its relationship with fertility: the lower pH class corresponded to the greater kindling rate in Brun et al. (2002b); and a negative and almost linear effect of pH on fertility was denoted by [START_REF] Tusell | Different ways to model biological relationships between fertility and the pH of the semen in rabbits[END_REF] when they analyzed the relationship of these 2 traits using recursive or classical univariate or bivariate models. Because semen pH is not expensive to measure, it could be of interest to include semen pH in routine measurements of the AI centers to select ejaculates to be used for AI to improve both conception rates and the selection of the most fertile males to use in the AI centers.

(Co)Variance Component Estimates

Table 2 shows features of the estimated marginal posterior distributions (EMPD) corresponding to the ratios between variance components and phenotypic variances for ADG and the seminal traits. Parameter estimates of the effects included in the ADG model did not change among the analyses. Therefore, for simplicity, only ADG estimates obtained in one of the models are presented (bivariate model with ADG and Lib). The heritability (h 2 ) for ADG was in accordance with estimates from previous studies obtained in other rabbit lines [START_REF] Larzul | Genetics of growth traits and meat quality in the rabbit[END_REF][START_REF] Lavara | Genetic parameter estimates for semen production traits and growth rate of a paternal rabbit line[END_REF] and in the same line [START_REF] Piles | Elliptical selection experiment for the estimation of genetic parameters of the growth rate and feed conversion ratio in rabbits[END_REF].

Seminal traits showed low to moderate values of h 2 and repeatability (r; calculated as the sum of the h 2 and the ratios of variance of permanent and common litter environmental effects). The ratios of variance of the permanent environmental effects ranged from 0.07 to 0.18 for all the seminal traits. The proportion of variance due to the common litter effect was almost null for the seminal traits, whereas it represented a high proportion of the total phenotypic variance for ADG in agreement with results obtained by [START_REF] Lavara | Genetic parameter estimates for semen production traits and growth rate of a paternal rabbit line[END_REF] and [START_REF] Piles | Elliptical selection experiment for the estimation of genetic parameters of the growth rate and feed conversion ratio in rabbits[END_REF].

Male Lib and Presence of Ur and Ca in the Ejaculate

Traits Lib, Ur, and Ca were found to be lowly heritable, which could be attributed in part to the great variability inherent in these traits due to factors involved in semen collection (i.e., variation in the temperature of the artifi cial vagina that could lead to a greater presence of Ur and Ca in the ejaculate or unsuccessful mountings; [START_REF] Morrell | Artifi cial insemination in rabbits[END_REF]. Therefore, genetic selection for increasing semen production by improving Lib and reducing Table 2. Mean (highest posterior density interval at 95% in brackets, and Monte Carlo SE in parentheses) of the marginal distribution of heritability (h 2 ), phenotypic variance (σ 2 ), ratio of variance of male permanent (pm) and common litter (c) environmental effects to phenotypic variance for ADG and seminal traits. 2 Lib = male libido (binary trait: 1 = success to mount artifi cial vagina, 0 = failure).

3 Estimates of the binary traits Lib, Ur, Ca, and Sui are given in the liability scale.

4 Ur = presence of urine in the ejaculate (binary trait: 1 = presence, 0 = absence).

5 Ca = presence of calcium carbonate deposits in the ejaculate (binary trait: 1 = presence, 0 = absence).

6 pH = pH of the ejaculate.

7 IM = individual sperm motility (subjective scale from 0 to 5 corresponding to a percentage of sperm showing progressive movement of 0 to 10, 11 to 25, 26 to 50, 51 to 70, 71 to 90, or 91 to 100%, respectively).

8 Sui = ejaculate suitability for AI (binary trait: 1 = suitable, 0 = unsuitable).

9 Conc = average sperm concentration of the ejaculate per male on the day of collection (×10 6 spermatozoa/mL).

10 Vol = average volume of the ejaculate per male on the day of collection (mL).

the number of rejected ejaculates may not be effective.

To our knowledge, there is no information in the literature concerning h 2 for presence of Ur and Ca in the ejaculates. Classifying Lib into 3 categories (no mating, collection after 5 min, and intermediate collection), [START_REF] Panella | Heritability of some male reproductive traits in rabbit[END_REF] reported an h 2 of Lib of 0.30 in rabbits. However, results obtained in that work should be taken with caution because all the genetic parameter estimates for seminal traits were unusually high, likely due to the absence in their model of a permanent environmental effect related to the male. [START_REF] Khalil | Genetic evaluation for semen characteristics in a crossbreeding project involving Saudi and Spanish V-line rabbits[END_REF] estimated the h 2 of Lib divided into 5 classes (from 1 for low libido to 5 for strong libido) using data from bucks coming from 2 rabbit lines and their reciprocal crosses. [START_REF] Khalil | Genetic evaluation for semen characteristics in a crossbreeding project involving Saudi and Spanish V-line rabbits[END_REF] found a greater estimate for this parameter than the one obtained in our study, due in part to the heterogeneity of the genetic type of the bucks. [START_REF] Flowers | Genetic and phenotypic variation in reproductive traits of AI boars[END_REF] concluded that it would be diffi cult to improve Lib and mating behaviors in boars through genetic selection due to the low phenotypic variation they observed for these traits in several studies. The posterior mean (high posterior density interval at 95%; PM [HPD95%]) of the repeatability (r) for Ca was 0.21 [0.19, 0.24], indicating a certain stability of the values of this trait over collections of the same male. However, the magnitude of this parameter is not high enough to make decisions concerning buck replacement according to this trait at the beginning of the production period of the male.

Individual Sperm Motility and Suitability for AI of the Ejaculate

The IM was also found to be lowly heritable and repeatable. The PM [HPD95%] of r was 0.25 [0.22, 0.28]. Estimates reported by other authors in rabbits are in agreement with these results for IM or mass motility score [START_REF] Bencheikh | The effect of frequency of ejaculation on semen characteristics and sperm output in the rabbit[END_REF][START_REF] Brun | Paramètres génétiques des caractéristiques de la semence de lapin. 13èmes[END_REF] and with the percentage of progressive motility in different breeds of boars [START_REF] Wolft | Genetic correlations between production and semen traits in pig[END_REF] and bulls [START_REF] Karoui | Time trends, environmental factors and genetic basis of semen traits collected in Holstein bulls under commercial conditions[END_REF]. This low heritability and repeatability could be due in part to the great variability of this trait originating from semen manipulation and time elapsed until evaluation as well as to the subjective manner in which this trait is measured because it is determined by the technician.

The low h 2 for Sui is in accordance with the estimate obtained for the subjective semen score used in bulls [START_REF] Knights | Estimates of heritabilities and of genetic and phenotypic correlations among growth and reproductive traits in yearling Angus bulls[END_REF], which, in turn, had a high r g with Conc and Vol, perhaps because they were the most important traits taken into consideration for scoring ejaculates. The effect of the technician [START_REF] Theau-Clément | Etude de facteurs de variation de la production spermatique chez le lapin. 13èmes Journées de la Recherche Cunicole[END_REF] probably introduces some error to the criterion that determines Sui, decreasing the obtained h 2 .

Semen pH

The PM of the pH h 2 was similar to the value previously reported in rabbits by [START_REF] Brun | Paramètres génétiques des caractéristiques de la semence de lapin. 13èmes[END_REF] and less than the one obtained in a previous study using a subset of the data analyzed here (0.18; [START_REF] Tusell | Different ways to model biological relationships between fertility and the pH of the semen in rabbits[END_REF]. The difference between the h 2 values obtained with these 2 sets of data could be due to the fact that in the previous study the analyzed trait was the pH corresponding to the pooled semen obtained from each male on the day of collection. Higher h 2 estimates are obtained if traits consist of means (or weighted means) of observations of 2 consecutive measurements or means of several records than if the trait comes from the data of individual records. This fact has been previously denoted for seminal traits by several authors [START_REF] Ducrocq | Genetic characteristics and evolution of semen production of young Normande bulls[END_REF][START_REF] Wolft | Genetic correlations between production and semen traits in pig[END_REF]. The PM (HPD95%) of repeatability of pH was 0.33 [0.26, 0.39], similar to the value obtained in a previous study where a subset of the data analyzed here was jointly analyzed with data corresponding to bucks from another paternal line and their reciprocal crosses (0.38; García-Tomás et al., 2006c).

Sperm Production Traits

The Conc and Vol traits showed moderate to high values of repeatability (PM [HPD95%]: 0.48 [0.42, 0.55] for Conc and 0.46 [0.40,0.53] for Vol), indicating the existence of important individual variation for both traits. Similar values were obtained by García-Tomás et al. (2006c). [START_REF] Bencheikh | The effect of frequency of ejaculation on semen characteristics and sperm output in the rabbit[END_REF] estimated a repeatability of about 0.38 for Vol and 0.35 for Conc whereas More O'Ferrall and Meacham (1968) obtained a value of 0.29 for Vol, all of them in rabbit.

The h 2 for these traits (Conc and Vol) were greater than estimates previously reported in rabbits by [START_REF] Brun | Paramètres génétiques des caractéristiques de la semence de lapin. 13èmes[END_REF] and [START_REF] Lavara | Genetic parameter estimates for semen production traits and growth rate of a paternal rabbit line[END_REF] analyzing records from single ejaculates (0.10 and 0.08 for Conc, 0.13 and 0.10 for Vol, respectively), whereas the ratios of permanent environmental effects for males were similar. In general, there is a wide range of values of heritability and repeatability of Vol and Conc in different species, varying from extremely low to high values [START_REF] Robinson | Impact of genetic selection on management of boar replacement[END_REF], and it is clear that the variation in magnitude found in the literature depends highly upon the defi nition of the trait (i.e., individual ejaculates or average values).

The EMPD of the genetic, nonadditive genetic plus permanent environmental male effect and phenotypic correlations between Conc and Vol are shown in Figure 2. All the Monte Carlo SE for the estimates of these correlations were less than 0.02.

The genetic correlation (r g ) between Conc and Vol was moderate and negative (PM: -0.53; HPD95%: [-0.76, -0.27]). Several studies also reported similar values of this r g in cattle (-0.4 on average; [START_REF] Ducrocq | Genetic characteristics and evolution of semen production of young Normande bulls[END_REF][START_REF] Basso | Estimation de paramètres génétiques et détection de QTL liés à des caractères de fertilité mâle, de production de semence et de qualité de la semence chez le taureau laitier[END_REF], sheep (-0.36 on average in adult males; [START_REF] Rege | Reproductive characteristics of Ethiopian highland sheep: II. Genetic parameters of semen characteristics and their relationships with testicular measurements in ram lambs[END_REF]David et al., 2007), goats (-0.33 on average;[START_REF] Furstoss | Genetic and non-genetic parameters of several characteristics of production and semen quality in young bucks[END_REF] and pigs (-0.60 on average; [START_REF] Smital | Estimation of genetic parameters of semen characteristics and reproductive traits in AI boars[END_REF][START_REF] Wolf | Quantifi cation of factors affecting semen traits in artifi cial insemination boars from animal model analyses[END_REF]. Conversely, [START_REF] Brun | Paramètres génétiques des caractéristiques de la semence de lapin. 13èmes[END_REF] obtained in rabbits an r g between Conc and Vol that could not be considered different from 0 (0.38 ± 0.45). Their estimate of the correlation between Conc and Vol did not agree with that in the current study, their negative correlation between permanent environmental effects of the male did (-0.47 ± 0.14 and PM: -0.51, HPD95%: [-0.78, 0.18] in [START_REF] Brun | Paramètres génétiques des caractéristiques de la semence de lapin. 13èmes[END_REF] and in our work, respectively). In boars, [START_REF] Smital | Estimation of genetic parameters of semen characteristics and reproductive traits in AI boars[END_REF] also estimated the r g between Conc and Prod (0.21) and between Vol and Prod (0.63), respectively. [START_REF] Smital | Estimation of genetic parameters of semen characteristics and reproductive traits in AI boars[END_REF] obtained an r g between Conc and Vol similar to that in the current study (-0.61), indicating that a greater ejaculate volume does not seem to be associated with a greater number of sperm in the ejaculate that could remain constant or even decrease. Alternatively, the cause of greater concentration of sperm in bull ejaculates seems to be the greater number of sperm, leading to obtain a negative but smaller r g between Conc and Vol because the r g between Conc and Prod and r g between Vol and Prod were 0.60 and 0.66, respectively [START_REF] Karoui | Time trends, environmental factors and genetic basis of semen traits collected in Holstein bulls under commercial conditions[END_REF].

Production exhibited a moderate h 2 and a moderate to high repeatability (PM [HPD95%]: 0.42 [0.35, 0.49]), both values being greater than the corresponding values obtained by [START_REF] Lavara | Genetic parameter estimates for semen production traits and growth rate of a paternal rabbit line[END_REF] and [START_REF] Brun | Paramètres génétiques des caractéristiques de la semence de lapin. 13èmes[END_REF] for individual ejaculates. The repeatability for Prod was estimated to be 0.33 in previous research by García-Tomás et al. (2006c) with both purebred and crossbred bucks. [START_REF] Smital | Estimation of genetic parameters of semen characteristics and reproductive traits in AI boars[END_REF] proposed either total number of spermatozoa or number of insemination doses produced at a certain dosage for developing a composite trait for breeding purposes.

Correlations Between ADG and Seminal Traits

The EMPD of the genetic, the nonadditive genetic permanent environmental effects, and the phenotypic correlations between ADG and the seminal traits are shown in Figures 3,4, and 5, respectively. All the Monte Carlo SE for the estimates of these correlations were less than 0.02.

Selection for increasing ADG would have a favorably correlated response on sperm concentration because the r g (PM [HPD95%]) between ADG and Conc, although low, was positive (0.21 [-0.03, 0.48]; Figure 3). This result is consistent with the low, but negative, r g between ADG and Vol (-0.19 [-0.47, 0.08]); Figure 3). The antagonism between Conc-ADG and Vol-ADG seems to infl uence the near-zero genetic relationship between 0.38]). In another paternal line of rabbits also selected for growth rate, the r g between Conc, Vol, and Prod with ADG was estimated to be -0.09, 0.36, and 0.17, respectively [START_REF] Lavara | Genetic parameter estimates for semen production traits and growth rate of a paternal rabbit line[END_REF]. However, in that work, estimated r g had wide HPD95%, making it diffi cult to draw unambiguous conclusions. Our estimates are in accordance with the phenotypic results obtained in an experiment of divergent selection for BW in rabbits in which the high line had greater sperm concentration and less ejaculate volume than the low line [START_REF] Brun | Semen production in two rabbit lines divergently selected for 63-d body weight[END_REF]. Other estimates, reported in pigs, confi rm that genetic correlations between growth and Vol, Conc, and Prod always had low magnitude [START_REF] Oh | Estimates of genetic correlations between production and semen traits in boar. Asian-australas[END_REF][START_REF] Wolft | Genetic correlations between production and semen traits in pig[END_REF].

None of the genetic correlations between ADG and the other seminal traits analyzed in this work can be considered different from zero. Because of the null or very low magnitude of the genetic correlation between ADG and seminal traits, we would not expect selection for growth rate in this paternal line of rabbits to dam- age semen production. In boars, [START_REF] Wolft | Genetic correlations between production and semen traits in pig[END_REF] reported low values of the genetic correlations between ADG and sperm motility, percentage of abnormal spermatozoa, and the composite trait, number of functional spermatozoa. No effect in Lib was encountered when 2 rabbit lines divergently selected for BW were compared in terms of rate of successful solicitations and in the average time until ejaculation [START_REF] Brun | Semen production in two rabbit lines divergently selected for 63-d body weight[END_REF].

All the correlations between permanent environmental effects for ADG and the seminal traits were of low magnitude (Figure 4) and a similar pattern was observed for the phenotypic correlations (Figure 5).

The correlation (r co ) between common birth litter effects for ADG and the seminal traits was very low, including zero in the HPD 95% except for Ca. For the last trait, the PM [HPD95%] of r co was 0.37 [0.01, 0.81], suggesting that favorable litter effects for growth could be detrimental for the production of ejaculates free of Ca. However, the accuracy of this correlation estimate is very low because of the low magnitude of this vari- (.)] between ADG during the fattening period and the following seminal traits: male libido (Lib), presence of urine (Ur) and calcium carbonate deposits (Ca) in the ejaculate, semen pH, individual motility score (IM), and suitability for AI (Sui) of the ejaculate, average sperm concentration (Conc), average volume (Vol), and average sperm production per ejaculate (Prod).

ance component for Ca. It is necessary to know the origin and causes of this trait to better understand the nature of this relationship.

Male libido and seminal traits related to the quality of the ejaculate (Ca, Ur, pH, and IM) seemed to be lowly heritable but repeatable. This indicates performance of bucks for seminal quality traits and libido in AI centers would be more strongly affected by management practices rather than genetic selection.

Moderate h 2 were estimated for Conc, Vol, and Prod, suggesting the possibility of effectively selecting for increased sperm production in paternal lines. Because of the existence of a genetic antagonism between Conc and Vol, selecting for Prod would be of most interest because it is a composite trait that encompasses both traits. (.)] between ADG during the fattening period and the following seminal traits: male libido (Lib), presence of urine (Ur) and calcium carbonate deposits (Ca) in the ejaculate, semen pH, individual motility score (IM), and suitability for AI (Sui) of the ejaculate, average sperm concentration (Conc), average volume (Vol), and average sperm production per ejaculate (Prod).

Selection for increasing ADG in paternal lines is not expected to have detrimental correlated effects on seminal traits involved in sperm production.

To improve the effi ciency of bucks for AI dose production, it would be necessary to evaluate the expected response to selection applying different strategies and under different scenarios. It would be also required to establish the genetic relationship between semen production and quality traits (in order not to impair the fer-tilizing ability of AI doses), and to know the economic weights for these traits.

LITERATURE CITED

Alvariño, J. M. R. 2000. Reproductive performance of male rabbits. Proc. 7th World Rabbit Congr., Valencia, Spain A:13-35. Baselga, M. 2004. Genetic improvement of meat rabbits. Programmes and diffusion. Proc. 8th World Rabbit Congr., Puebla, Mexico 1:1-13. .)] between ADG during the fattening period and the following seminal traits: male libido (Lib), presence of urine (Ur) and calcium carbonate deposits (Ca) in the ejaculate, semen pH, individual motility score (IM), and suitability for AI (Sui) of the ejaculate, average sperm concentration (Conc), average volume (Vol), and average sperm production per ejaculate (Prod). Phenotypic correlation between ADG and binary traits (Lib, Ur, Ca, and Sui) were based on the liability scale.

Figure 1 .

 1 Figure 1. Relative importance of several causes of ejaculate rejection for AI: presence of urine (Ur), presence of calcium carbonate deposits (Ca), and individual motility score less than or equal to 2 (IM ≤ 2) on a subjective scale from 0 to 5 corresponding to a percentage of sperm showing progressive movement of 0 to 10, 11 to 25, 26 to 50, 51 to 70, 71 to 90, or 91 to 100%, respectively.

Figure 2 .

 2 Figure 2. Estimated marginal posterior distributions of the genetic [r g (.)], male nonadditive genetic permanent environmental [r pm (.)], and phenotypic [r p (.)] correlation between average sperm concentration (Conc) and average volume (Vol) of the ejaculate on the day of collection.

Figure 3 .

 3 Figure 3. Estimated marginal posterior distributions of the genetic correlation [r g(.)] between ADG during the fattening period and the following seminal traits: male libido (Lib), presence of urine (Ur) and calcium carbonate deposits (Ca) in the ejaculate, semen pH, individual motility score (IM), and suitability for AI (Sui) of the ejaculate, average sperm concentration (Conc), average volume (Vol), and average sperm production per ejaculate (Prod).

Figure 4 .

 4 Figure 4. Estimated marginal posterior distributions of the male nonadditive genetic permanent environmental correlation [r pm(.)] between ADG during the fattening period and the following seminal traits: male libido (Lib), presence of urine (Ur) and calcium carbonate deposits (Ca) in the ejaculate, semen pH, individual motility score (IM), and suitability for AI (Sui) of the ejaculate, average sperm concentration (Conc), average volume (Vol), and average sperm production per ejaculate (Prod).

Figure 5 .

 5 Figure 5. Estimated marginal posterior distributions of the phenotypic correlation [r p(.)] between ADG during the fattening period and the following seminal traits: male libido (Lib), presence of urine (Ur) and calcium carbonate deposits (Ca) in the ejaculate, semen pH, individual motility score (IM), and suitability for AI (Sui) of the ejaculate, average sperm concentration (Conc), average volume (Vol), and average sperm production per ejaculate (Prod). Phenotypic correlation between ADG and binary traits (Lib, Ur, Ca, and Sui) were based on the liability scale.

Table 1 .

 1 Summary statistics of ADG and seminal traits expressed as average percentage per buck.

	Trait	Minimum Mean Maximum	CV	No. of records
	ADG 1	10.78	45.98	88.95	0.18	118,306
	Lib 2	0	0.80	1	1.58	13,914
	Ur 3	0	0.08	1	2.15	14,383
	Ca 4	0	0.22	1	1.18	14,449
	pH 5	6.35	7.45	8.66	0.05	1,386
	IM 6	0	2.69	4.67	0.28	13,727
	Sui 7	0	0.61	1	0.50	14,240
	Conc 8	0	331.82 1,616.04	0.67	1,290
	Vol 9	0.07	0.77	2.1	0.53	1,225
	Prod 10	29.875 670.83 2,450.50	0.49	1,077

at INRA Institut National de la Recherche Agronomique on July 23, 2013 www.journalofanimalscience.org Downloaded from

Conc = average sperm concentration of the ejaculate per male on the day of collection (×10 6 spermatozoa/mL).

Vol = average volume of the ejaculate per male on the day of collection (mL).

Prod = average sperm production per ejaculate (×10 6 spermatozoa). at INRA Institut National de la Recherche Agronomique on July 23, 2013 www.journalofanimalscience.org Downloaded from

Prod = the average sperm production per ejaculate (×10 6 spermatozoa). at INRA Institut National de la Recherche Agronomique on July 23, 2013 www.journalofanimalscience.org Downloaded from

1 This research was supported by the Instituto Nacional de Investigación y Tecnología Agraria y Alimentaria (INIA, ) project RTA2005-00088-CO2. Part of it was carried out during a visit of the fi rst author to INRA UR 631 Station d'Amélioration Génétique des Animaux (Castanet-Tolosan, France). Project partly supported by Toulouse Midi-Pyrénées bioinformatic platform (Toulouse, France). Llibertat Tusell Palomero received a fellowship from the INIA. The authors are grateful to the staff of Unitat de Cunicultura, IRTA (Caldes de Montbui, Barcelona, Spain) for their invaluable cooperation and


C Colombani
email: carine.colombani@toulouse.inra.fr

P Croiseau

S Fritz

F Guillaume

A Legarra

V Ducrocq

C Robert-Granié

A comparison of partial least squares (PLS) and sparse PLS regressions in genomic selection in French dairy cattle

Keywords: partial least squares regression, sparse partial least squares, genomic selection, French dairy cattle

Genomic selection involves computing a prediction equation from the estimated effects of a large number of DNA markers based on a limited number of genotyped animals with phenotypes. The number of observations is much smaller than the number of independent variables, and the challenge is to find methods that perform well in this context. Partial least squares regression (PLS) and sparse PLS were used with a reference population of 3,940 genotyped and phenotyped French Holstein bulls and 39,738 polymorphic single nucleotide polymorphism markers. Partial least squares regression reduces the number of variables by projecting independent variables onto latent structures. Sparse PLS combines variable selection and modeling in a one-step procedure. Correlations between observed phenotypes and phenotypes predicted by PLS and sparse PLS were similar, but sparse PLS highlighted some genome regions more clearly. Both PLS and sparse PLS were more accurate than pedigree-based BLUP and generally provided lower correlations between observed and predicted phenotypes than did genomic BLUP. Furthermore, PLS and sparse PLS required similar computing time to genomic BLUP for the study of 6 traits.

INTRODUCTION

Genomic selection relies on computing genomic estimated breeding values (GEBV) using high-density SNP marker data. [START_REF] Meuwissen | Prediction of total genetic value using genome-wide dense marker maps[END_REF] suggested a 2-step approach to calculate GEBV. First, the effects of SNP are estimated to obtain a prediction equation using a reference population in which the animals are genotyped and phenotyped. Then, GEBV are predicted for the genotyped animals (without phenotypes) from this equation.

In the past few years, the accuracy of GEBV provided by genomic selection has been assessed using different methods in dairy cattle populations in the United States, New Zealand, Australia, the Netherlands, and France, among others. A simple BLUP, as described in [START_REF] Meuwissen | Prediction of total genetic value using genome-wide dense marker maps[END_REF] and known as genomic BLUP (GBLUP) in subsequent literature, was used as the reference method. The simple BLUP assumes that all SNP have an effect sampled from the same normal distribution. [START_REF] Hayes | Invited review: Genomic selection in dairy cattle: Progress and challenges[END_REF] treated Australian Holstein-Friesian bull data using a method derived from BayesA, which exploits the prior knowledge that many SNP have small individual effects on the trait and only a few have moderate to large effects. The Bayesian method was shown to be slightly more reliable (+0.02 to +0.07 compared with the reliability of BLUP) for most traits. Using New Zealand dairy cattle, [START_REF] Harris | Genomic selection in New Zealand and the implications for national genetic evaluation[END_REF] also compared the BLUP approach with Bayesian methods (BayesA and BayesB), in which some SNP may have zero effect [START_REF] Meuwissen | Prediction of total genetic value using genome-wide dense marker maps[END_REF]. Bayesian methods slightly improved reliability (2%), whereas the use of regression methods such as least angle regression [START_REF] Efron | Least angle regression[END_REF] did not lead to any improvement. [START_REF] Vanraden | Invited review: Reliability of genomic predictions for North American Holstein bulls[END_REF] compared the reliability of GEBV in US and Canadian young bulls, using a method similar to GBLUP that fits the allelic effects of each SNP as random effects with a normal distribution with known variance [START_REF] Vanraden | Efficient methods to compute genomic predictions[END_REF], and a similar method to BayesA with a heavier tail distribution for the SNP effects. As in the Australian and New Zealand results, the Bayesian approach slightly increased reliability (1% compared with the reliability of GBLUP). [START_REF] Moser | A comparison of five methods to predict genomic breeding values of dairy bulls from genome-wide SNP markers[END_REF] compared 5 methods on dairy bull data including regression methods (least squares regression), shrinkage methods [Bayes regression (Bayes-R) similar to BayesA, and random regression BLUP (RR-BLUP), comparable to GBLUP], sup-port vector machine learning methods (nonparametric support vector regression), and dimension reduction methods such as partial least squares (PLS) regression. The accuracy of Bayes-R, RR-BLUP, PLS, and support vector regression was very similar for the 2 traits studied by these authors. However, PLS and RR-BLUP required substantially less computation time than the Bayesian method.

Using simulated data, [START_REF] Coster | Sensitivity of methods for estimating breeding values using genetic markers to the number of QTL and distribution of QTL variance[END_REF] demonstrated the superiority of PLS over Bayesian methods with regard to the stability of results according to the number of QTL or the distribution of QTL variance. They also showed that the computation time for the PLS method required to fit, cross validate, and evaluate the models was less than that for the Bayesian method. However, the Bayesian method was more accurate. [START_REF] Solberg | Reducing dimensionality for prediction of genome-wide breeding values[END_REF] also used simulated data to compare PLS and principal component regression with BayesB. They obtained the same results: BayesB was more accurate than other methods but PLS and principal component regression were computationally faster and simpler.

The PLS regression [START_REF] Wold | PLS-regression: A basic tool of chemometrics[END_REF] appears to be an efficient method to deal with genomic selection data, both in its capacity to handle large data sets and its prediction ability. This approach is particularly suitable when the matrix of predictors has more variables than observations, and when multicollinearity exists among variables. The sparse PLS regression (sPLS, Lê Cao et al., 2008) is a recent approach that combines variable selection and modeling in a one-step procedure. Dimension reduction methods and variable selection approaches may be an attractive way to deal with the increasing number of markers used in genomic selection in dairy cattle by limiting computing time [START_REF] Coster | Sensitivity of methods for estimating breeding values using genetic markers to the number of QTL and distribution of QTL variance[END_REF]. Furthermore, even though PLS has already been studied in a genomic evaluation context, the authors used simulated data and did not compare PLS accuracy with current genomic selection methods such as BLUP and GBLUP [START_REF] Solberg | Reducing dimensionality for prediction of genome-wide breeding values[END_REF][START_REF] Coster | Sensitivity of methods for estimating breeding values using genetic markers to the number of QTL and distribution of QTL variance[END_REF]. [START_REF] Long | Dimension reduction and variable selection for genomic selection: Application to predicting milk yield in Holsteins[END_REF] introduced sparsity in PLS and tested the predictive ability of sPLS versus principal component regression and PLS methods but did not apply other current genomic selection methods on their real data. They showed that combining dimension reduction and variable selection for accurate prediction of genomic breeding values was promising.

The aim of the present study was to compare PLS and sPLS on a real data set with other methods currently used in the evaluation of dairy cattle such as pedigree-based BLUP and GBLUP. Both PLS and sPLS regressions were compared based on their predictive abilities and then with pedigree-based BLUP and GBLUP results to evaluate their accuracy.

MATERIALS AND METHODS

Data

A data set of genotyped French Holstein bulls was split into a training data set and a validation data set using a cut-off birth date defined so that the validation set included the youngest 25% genotyped bulls. First, the prediction equation was estimated with the training data set, which comprised 2,976 genotyped and phenotyped Holstein bulls born before June 2002. Then, phenotypes were predicted for the bulls in the validation data set, which comprised 964 bulls (born between June 2002 and2004).

Genotypes for 39,738 polymorphic SNP were used as independent variables. The selected SNP, provided by the Illumina Bovine SNP50 Beadchip (Illumina, San Diego, CA), had minor allele frequencies >3%. Mendelian segregation was checked. Missing genotypes were inferred from large family information with a low error rate using DualPHASE software (Druet and Georges, 2009).

Six traits with different heritability were used as dependent variables: milk yield, fat yield, and protein yield (h2 = 0.3), fat content and protein content (h 2 = 0.5), and conception rate (h 2 = 0.02; [START_REF] Boichard | Genetic analysis of conception rate in French Holstein cattle[END_REF]. The bulls' phenotypes used in this study were daughter yield deviations (DYD, VanRaden and Wiggans, 1991;[START_REF] Mrode | Calculating cow and daughter yield deviations and partitioning of genetic evaluations under a random regression model[END_REF]) from a French national evaluation October 2009; that is, the average performance of the daughters of a sire, adjusted for fixed and nongenetic random effects and for the additive genetic value of their dam. For each DYD, a weighting was added in the form of the effective daughter contribution (EDC;[START_REF] Vanraden | Derivation, calculation, and use of national animal-model information[END_REF][START_REF] Fikse | Weighting factors of sire daughter information in international genetic evaluations[END_REF]. To be included in the analysis, each observation required an EDC >20.

The general statistical model in BLUP and GBLUP is

y = 1μ + Zg + e,
where y is a vector of phenotypes (DYD), μ is the mean, Z is a design matrix allocating observations to breeding values, g is a random vector of additive genetic values, and e is a vector of random normal errors. In BLUP, Var g A ( ) = σ g 2 , where A is the pedigree-based relationship matrix, and σ g 2 is the additive genetic variance. In GBLUP, Var g

( ) = Gσ g G WW = ′ - ( ) = ∑ 2 1 1 q q j j j p ,
where p is the number of loci considered, q j is the frequency of an allele of the marker j, and W is a centered incidence matrix of SNP genotypes. The SNP marker effects are assumed to have a prior normal distribution and mixed model equations are used with the genomic relationship matrix [START_REF] Cole | Distribution and location of genetic effects for dairy traits[END_REF][START_REF] Vanraden | Invited review: Reliability of genomic predictions for North American Holstein bulls[END_REF]. [START_REF] Wold | Estimation of principal components and related models by iterative least squares[END_REF] is a data analysis method that generalizes and combines principal component analysis and multiple regression. The method was mainly developed for industrial applications (petroleum and food processing industries) and the social sciences. The PLS method was designed to deal with the "p n problem"; that is, when the number of independent variables (p) is much larger than the number of observations (n). Partial least squares regression is very useful to predict dependent variables from a very large number of predictors that might be highly correlated.

PLS. The PLS regression introduced by

In its general form, the PLS regression replaces the initial independent variable space (X) and the initial response variable space (Y) by smaller spaces that rely on a reduced number of variables named latent variables, which are included one by one in an iterative process. These factors will be the new variables of a usual linear regression. The main idea is to perform successive regressions by projections onto latent structures to reveal hidden or latent underlying biological effects [START_REF] Wold | The PLS method-Partial least squares projections to latent structures and its application in industrial RDP[END_REF][START_REF] Lê Cao | A sparse PLS for variable selection when integrating omics data[END_REF].

Using the same notation as in Lê [START_REF] Lê Cao | A sparse PLS for variable selection when integrating omics data[END_REF], the PLS regression looks for a decomposition of centered data matrices X and Y in terms of component scores, called latent variables: (ξ 1 ,…, ξ h ,…, ξ H ) and (ω 1 ,…, ω h ,…, ω H ), which are linear combinations of the columns of X and Y respectively, and associated loading vectors:

(u 1 ,…, u h ,…, u H ) and (v 1 ,…, v h ,…, v H )
, where H is the number of latent variables retained in the final model. However, the regression coefficients that define these components are not linear, as they are solved via successive local regressions on the latent variables. The loading vectors are estimated to solve the following optimization problem:

max cov , , , u v 1 h h h h X u Yv = = - ( ) 1 1 h
where X h-1 is the residual X matrix in the regression of Y on (ξ 1 ,…, ξ h-1 ) for each dimension h = 1, …, H, and the associated latent variables are denoted ξ h = X h-1 u h and ω h = Yv h . As in principal component analysis, the loading vectors and the latent variables are directly interpretable. The loading vectors u h and v h indicate how the x j and y i variables explain the relationship between X and Y. The latent variables contain information regarding similarities or dissimilarities between individuals [START_REF] Wold | The PLS method-Partial least squares projections to latent structures and its application in industrial RDP[END_REF].

sPLS. The sPLS regression [START_REF] Lê Cao | A sparse PLS for variable selection when integrating omics data[END_REF]) aims at combining variable selection and modeling in a one-step procedure. It was first proposed to handle transcriptomic data and was adapted for genomic data in this study. To understand the sPLS approach, it is helpful to describe first the principle of the PLSsingular value decomposition [START_REF] Lorber | A theoretical foundation for the PLS algorithm[END_REF] that solves PLS problems efficiently by decomposing the X Y matrix into singular values and vectors.

For a real matrix M (p × q) of rank r, the singular value decomposition of M can be obtained as follows:

M = ΓΔΘ
where Γ (p × r) and Θ (q × r) are orthonormal and Δ (r × r) is a diagonal matrix with singular values δ k (k = 1…r).

The loading vectors u 1 and v 1 of X and Y, respectively, correspond to the first singular vectors γ 1 and θ 1 if M = X Y. Then, for h = 2, …, H, M h is directly deflated by its rank-one approximation, as explained in Lê Cao et al. ( 2008

): M M u v h h h h 1 = - - δ h ′ .
Sparsity of the loading vectors is introduced iteratively by penalizing both u h and v h with a soft-thresholding penalization, as for sparse principal component analysis [START_REF] Shen | Sparse principal component analysis via regularized low rank matrix approximation[END_REF]. The optimization problem becomes min ,

u,v M uv u v - + ( )+ ( ) F g g ′ λ λ 2 1 2 where g x sign x x λ λ 1 1 ( ) = ( ) - ( ) | | + is the soft-threshold- ing penalty function.
When no sparsity is required, the same results are obtained as in classical PLS. The details of this algorithm are presented in Lê [START_REF] Lê Cao | A sparse PLS for variable selection when integrating omics data[END_REF]. Although PLS and sPLS can perform multi-trait analyses, each trait in this study was considered independently with X the matrix of SNP and y the vector of phenotypes of one trait. In this case, the sPLS used here is similar to the sPLS introduced by Chun and Keles (2010) and used by [START_REF] Long | Dimension reduction and variable selection for genomic selection: Application to predicting milk yield in Holsteins[END_REF], in a genomic selection con-text. When one trait is studied, PLS and sPLS are easy to implement and not time consuming because they are based on successive regressions and do not require matrix inversion.

Parameter Tuning

In both PLS and sPLS, the optimal number H of dimensions has to be determined. The parameter H can be tuned by cross-validation as in the original PLS and as proposed by [START_REF] Chun | Sparse partial least squares regression for simultaneous dimension reduction and variable selection[END_REF]. [START_REF] Coster | Sensitivity of methods for estimating breeding values using genetic markers to the number of QTL and distribution of QTL variance[END_REF] also proposed to use cross-validation to find the number of dimensions that minimized the prediction error. In this study, the root mean squared error of prediction (RMSEP) was minimized with 10-fold cross-validation in the training data set and for each given dimension h [START_REF] Mevik | Mean squared error of prediction (MSEP) estimates for principal component regression (PCR) and partial least squares regression (PLSR)[END_REF]:

RMSEP k k k = - ( ) = ∑ 1 10 1 10 ˆ, y y
where ŷk is the vector of predicted values for the sample k. [START_REF] Solberg | Reducing dimensionality for prediction of genome-wide breeding values[END_REF] suggested keeping the number of dimensions leading to the highest correlation between predicted values and observed values in the validation data set. This approach was also tested in this study and the results of the 2 ways used to fix H are discussed in the Results and Discussion section.

In sPLS regression, in addition to H (that was selected as above), the number of variables selected in each dimension of the model has to be fixed. [START_REF] Long | Dimension reduction and variable selection for genomic selection: Application to predicting milk yield in Holsteins[END_REF] tested different values of H and different values of the number of variables selected in each dimension, based on results of their PCR study, to maximize the cross-validation correlation. In our study, we chose to minimize the RMSEP with 10-fold cross-validation in the training data set, for the previously fixed number of dimensions H. In practice, for each trait, several sPLS were performed depending on the number of selected SNP in each latent variable or dimension (assumed to be constant), as a percentage of the number of SNP in the whole data set. By construction, the same SNP could be selected in several dimensions. Ten sPLS regressions (keeping 0.2 to 10% of all SNP for each dimension considered) and the PLS regression were tested using dimensions 1 to 100.

Importance of SNP Effects

To enable better interpretation of the models, coefficients that represent the power of x j to explain y has to be defined. The "variable importance in projection" (VIP) coefficients measure the contribution of x j to the construction of y through latent variables ξ h , (h = 1,…, H) and is defined by

VIP p cor cor Hj h H hj h H = ( ) ( ) = = ∑ ∑ 2 1 2 2 1 y y h h , , , ξ ξ ω with VIP p Hj j p 2 1 = = ∑ .
The sum of squares of the VIP coefficients of all the SNP in one dimension of the PLS models is equal to the number of independent variables. Thus, the VIP coefficient of a SNP is related to the number of SNP that have a nonzero effect in the model.

The contribution of x j to the construction of ξ h is measured by its weight ω hj , provided by PLS or sPLS. Although the weight ω hj of x j is interpretable, it does not account for the contribution of the latent variable ξ h . The VIP coefficients are able to classify the variables x j according to their weight in each latent variable and the weight of each latent variable in the construction of y. So they could be considered as an evaluation of the effects of SNP on the prediction of y. Both PLS and sPLS were performed using the R package named "mixOmics" [START_REF] Lê Cao | integrOmics: An R package to unravel relationships between two omics datasets[END_REF].

Comparison of Methods and EDC

Prediction Equation. In this study, we compared 2 currently used methods for the evaluation of dairy cattle, BLUP and GBLUP, with PLS and sPLS regressions. The application of the different methods followed the same pattern, regardless of the method. The prediction equation was estimated using the training data set. The y phenotypes were DYD. One EDC was associated with each DYD, reflecting its uncertainty. This generates heterogeneity of variances, so that the ith DYD has a (pseudo-residual) variance σ e i EDC 2 [START_REF] Vanraden | Derivation, calculation, and use of national animal-model information[END_REF]). An equivalent model was constructed in all cases (BLUP, GBLUP, PLS, and sPLS), multiplying y i and the ith row of the incidence matrix by the square root of the EDC to obtain homogeneous variances.

Accuracy of the Methods. Two criteria were used to test the accuracy of the different methods and to compare their predictive ability: the correlation (ρ) between observed and predicted values and the regression slope (b) of observed to predicted values (a value of 1 is expected; [START_REF] Henderson | Selection index and expected genetic advance. Statistical genetics and plant breeding[END_REF]. The bulls in the validation data set were predicted by the prediction equations provided by the different methods. Bulls in the validation set were progeny tested, so that observed DYD and associated EDC were also available for this population. These weights were taken into account in the calculation of the correlation, using and w i are the weights (EDC) for each observation. In the regression of observed DYD onto predicted DYD, EDC were introduced in the same way as for the production of prediction equations. The Hotelling-Williams procedure was used to test for differences between the correlations obtained from the different methods. It tests the null hypothesis of equality between 2 dependent correlations that share a variable [START_REF] Steiger | Tests for comparing elements of a correlation matrix[END_REF][START_REF] Vansickle | Analysing correlations between stream and watershed attributes[END_REF]. Under the null hypothesis, the statistical test is distributed as t with n -3 degrees of freedom. All the correlations discussed in this study were compared with one another using the Hotelling-Williams test with a 5% threshold.

ρ xy i i w i w i i i w i i w i i w x x

RESULTS AND DISCUSSION

Parameter Tuning

Figure 1 presents RMSEP obtained by cross-validation in the training data set by PLS, according to the number of dimensions. For milk yield, fat yield, protein yield, and protein content, the pattern of the different curves became very stable after only 10 dimensions. For fat content and conception rate, RMSEP stabilized after about 30 dimensions. The minimum RMSEP was reached after around 20 dimensions for milk yield, protein yield, fat content, and conception rate and after around 30 dimensions for fat yield and protein content. However, the differences in RMSEP between 2 values for the number of dimensions were very small (the minimum of the curves was not accentuated), so cross-validation did not appear to be the best criterion to choose the number of dimensions. The same conclusions were reached using the sparse PLS approach, so the number of dimensions that led to the highest correlation between phenotypes and predicted values in the validation data set was kept. This was also done for practical reasons. In fact, creating pseudo-training and pseudo-validation data sets within the training population to calibrate H was difficult, especially if a time structure (old vs. young) had to be used. Partial least squares regression was tested up to dimension 100 but the correlations obtained after more than 50 dimensions no longer increased for most traits. Figure 2 shows the correlations between observed phenotypes and predicted values in the validation data set obtained by PLS as a function of the number of latent variables built for each trait. The pattern of the different curves was the same whatever the trait: the correlation continued to increase with the number of dimensions until a plateau was reached at around dimension 30 for conception rate, 40 for milk yield, fat yield, and protein yield, and 80 for fat content and protein content. The number of dimensions in the final model was fixed at these minimum values to avoid overfitting the data [START_REF] Abdi | Partial least regression and projection on latent structure regression (PLS regression)[END_REF]. As can be seen in Figure 2, the number of dimensions is not critical, and therefore the choice of H using the validation data set did not overestimate the predictive ability of PLS.

Sparse PLS required 2 parameters: the number of latent variables (H) and the proportion of SNP selected in each dimension. The same criterion as in PLS was used to fix the number of dimensions kept in the sparse PLS models. The pattern was the same as in PLS (results not shown): a plateau was reached at around dimension 25 for conception rate, 40 for milk, fat, and protein yields, and 50 for fat and protein contents.

The proportion of SNP selected by sPLS was tuned by cross-validation within the training data set. Table 1 shows the RMSEP provided by PLS and the different sPLS (according to the proportion of SNP selected in each dimension) for each trait and for the previously fixed number of dimensions. The selected proportion of SNP for sPLS was the one that minimized RMSEP (Table 1, in bold, based on 3 decimal places). The minimum RMSEP obtained from sPLS was close to that obtained from PLS. For example, for fat yield, the RMSEP obtained with sPLS with 2% of the total number of SNP for each dimension was the smallest (0.18) and was the same as the error of prediction obtained with PLS. The heritability of the trait played a role in the magnitude of the RMSEP: traits with the same heritability led to similar prediction errors. Milk, fat, and protein yield (h 2 = 0.3) obtained an RMSEP of around 0.15. Fat and protein content (h 2 = 0.5) gave an RMSEP of 0.56. Conception rate had the highest error (0.79) but this result was the same as the error with the PLS model (0.78). As DYD were pseudo-phenotypes calculated from the performance of the bull's daughters, low heritability traits were difficult to predict.

Table 2 shows the accuracy of PLS and sPLS regressions for the different traits, with the number of SNP in sPLS that led to the minimum RMSEP, and the number of dimensions that led to the maximum correlation. The PLS regression gave significantly higher correlations whatever the trait with an average increase of 0.04 compared with sparse PLS. The best correlations were obtained for fat and protein content (0.70 and 0.71 in PLS, 0.66 and 0.65 in sPLS, respectively). Using the Hotelling-Williams test with a 5% threshold, production traits with a higher heritability, such as milk, fat, and protein yield (h 2 = 0.3), and fat and protein content (h 2 = 0.5) were predicted more accurately (from 0.53 in PLS and 0.48 in sPLS for milk yield to 0.71 in PLS for protein content and 0.66 in sPLS for fat content) than traits with lower heritability, such as conception rate (h 2 = 0.02) with an accuracy of 0.33 in PLS and 0.29 in sPLS. Thus, accuracy of prediction and heritability of the trait are closely related. [START_REF] Moser | Accuracy of direct genomic values in Holstein bulls and cows using subsets of SNP markers[END_REF] processed data from 2,144 Holstein-Friesian bulls and compared accuracy between traits of different heritability. Production traits such as protein content, fat content, and milk yield, which have high heritability (0.56, 0.52, and 0.28, respectively), achieved higher accuracy than survival (h 2 = 0.03), which showed similar heritability to conception rate in the present study. For conception rate, a larger reference population is required to achieve the same level of accuracy as for production traits [START_REF] Hayes | Invited review: Genomic selection in dairy cattle: Progress and challenges[END_REF]. However, unlike the data sets used by [START_REF] Hayes | Invited review: Genomic selection in dairy cattle: Progress and challenges[END_REF], which contained only 332 Australian Holstein bulls for fertility and 798 for the other traits, the data sets of this study for conception rate used a similar number of bulls with as many daughters to evaluate DYD as the number of bulls and daughters used for production traits. Therefore, for conception rate, the power of the analysis was not reduced by a smaller data set but by low heritability. Regarding the regression slope b, both PLS and sPLS gave values below 1, with PLS values closer to 1 (from 0.60 for conception rate to 0.83 for protein content) than those for sPLS (from 0.53 for milk yield to 0.76 for protein content). Furthermore, the relationship between the heritability of the trait and slope was less clear than between heritability and the correlation. For example, the sparse PLS slopes for milk yield (b = 0.53) and for conception rate (b = 0.54) were similar but milk yield is a trait with a moderate heritability (h 2 = 0.3), whereas conception rate is a low heritability trait (h 2 = 0.02). Two traits with 2 different heritabilities would lead to 2 different correlations but not necessarily to 2 different regression slopes.

Sparse PLS gave significantly less accurate predictions than PLS. However, sPLS performed a variable selection by allowing the number of variables in the final model to be reduced by 50% (Table 2). The number of selected SNP was reduced to 9,832 for fat content. One explanation could be the presence of DGAT1 [START_REF] Grisart | Genetic and functional confirmation of the causality of the DGAT1 K232A quantitative trait nucleotide in affecting milk yield and composition[END_REF], a gene on bovine chromosome 14 that leads to a mutation that has a major effect on fat content in milk in Holstein dairy cattle. Therefore, a small number of SNP, of which many were located around this QTL, was sufficient to obtain accurate predictions. However, the number of selected SNP could have been smaller, but the large number of latent components ( 50) used in the model led to a high number of selected SNP, irrespective of the presence of a large QTL.

Indeed, the number of SNP in the final model was directly related to the percentage of SNP kept for each latent variable and to the number of latent variables. The large number of SNP selected for protein yield was the consequence of the large proportion of SNP (10%) selected for each of the 38 latent variables. Fifty latent variables were used for fat content but the number of SNP was reduced because the sPLS that led to the The number of latent variables required by PLS and sPLS was very high (between 27 and 83 dimensions depending on the trait). [START_REF] Long | Dimension reduction and variable selection for genomic selection: Application to predicting milk yield in Holsteins[END_REF] evaluated PLS, sPLS developed by [START_REF] Chun | Sparse partial least squares regression for simultaneous dimension reduction and variable selection[END_REF], and principal component regression in Holstein bulls. They showed that to predict milk yield in Holsteins by PLS, only 15 latent components were sufficient to obtain the largest predictive correlation, suggesting a strong predictive power of the latent variables. The lack of predictive power of the additional latent components in the present study seems to be due to the presence of highly EDC-weighted bulls in the training data set, which had a major effect on the distribution of the phenotypes.

The distribution of the phenotypes in the training data set was normal, but applying EDC disturbed the normal distribution of the phenotypes (results not shown). Furthermore, the correlation between observed DYD and observed weighted DYD was surprisingly small: about 0.20 for milk, fat, and protein yield, about 0.45 for fat and protein content, and about 0.35 for conception rate. Consequently, we investigated the effect of the use of EDC on the training population.

Figure 3 shows the distribution of EDC in the training data set for the 6 traits studied. The computation of EDC relies on the number of daughters per bull and trait parameters (heritability and repeatability). The EDC were the same for milk yield, fat yield, and protein yield and fat content and protein content. In order for the information content in conception rate to be consistent with that of production traits, a large number of daughters per bull was assumed both for production traits and conception rate. The graphs show that some bulls differed from others in their very high EDC. These bulls were generally older than the other bulls in the training population and did not obtain stronger DYD. The significant difference in weights between bulls resulted in a bias and had a major effect on the number of latent components introduced in the PLS and sPLS models. To test this hypothesis, the same study was performed without considering EDC either in the weighting of DYD or in the calculation of accuracy and of the regression slope.

Table 3 shows the accuracy of PLS and sPLS regressions with respect to the different traits without EDC. For production traits, with PLS, the correlations were very similar to the results of PLS with EDC (Table 2). With sPLS, accuracy was significantly better in the study without EDC (fat yield ρ = 0.59 and protein content ρ = 0.72, for example) than in the study including EDC (fat yield ρ = 0.54 and protein content ρ = 0.65). For conception rate, using EDC reduced significantly the accuracy of both PLS and sPLS. Sparse PLS and PLS gave no significantly different correlations without EDC, whereas sPLS was shown to be significantly less accurate than PLS with EDC. The regression slopes were differently affected by the use or nonuse of EDC, with regression slopes lower than or close to those in the previous study with PLS, and regression slopes greater than or equal to those in the previous study with sPLS. Irrespective of the trait, the number of dimensions was considerably reduced with both methods, which led to a stronger variable selection in sPLS and a restricted number of SNP in the prediction equation. Only 10 or fewer latent variables were needed to obtain the best correlations for fat yield, protein yield, and conception rate with both PLS and sPLS. These results are in agreement with the number of dimensions obtained by [START_REF] Long | Dimension reduction and variable selection for genomic selection: Application to predicting milk yield in Holsteins[END_REF]. For protein yield and protein content, the number of SNP remained high because the sPLS, with the minimum RMSEP by cross-validation, kept 10% and 5% of the SNP in each dimension. Therefore, introducing EDC did not have a major effect on the predictive ability of PLS but did affect the number of latent variables of the model. With PLS methods, it is probably wiser to have a more homogeneous distribution of the weights to favor the dimension-reducing ability of the PLS variants and hence to reduce compu-tation time. In the remainder of the study, EDC were not included.

Both PLS and sPLS seemed to fit well in the context of genomic selection, but sPLS led to slightly smaller correlations than PLS but with no significant differences when EDC were not considered. However, sPLS favored a variable selection that can highlight the most important SNP, if required. A secondary aim of the present study was to underline the explanatory power of PLS regarding biological processes. However, to interpret the model in a biological context, coefficients that represent the explanatory power of the variables in the construction of the response variable are needed.

Figure 4 shows the VIP coefficient computed for each variable according to the position of the SNP on the genome. All the SNP variables are shown in the graph. A large number of VIP coefficients were set to zero in sPLS, whereas all the coefficients differed from zero in PLS. Therefore, sPLS was able to select variables based on VIP coefficients. The scale of the y-axis was the same for all the traits except fat content, which had the highest VIP coefficients.

Excluding the use of EDC, the variable selection performed by sPLS highlighted areas of interest. For fat yield, some SNP on chromosomes 2 and 5 were already highlighted by PLS but were clearly weighted up in sPLS (VIP coefficients of around 2.2 and 2.5 in PLS and of 12 and 16 in sparse PLS, respectively). For fat content, the SNP located at the beginning of chromosome 14 were highlighted by both methods. This location corresponds to a region of the genome that hosts the QTL DGAT1 [START_REF] Grisart | Genetic and functional confirmation of the causality of the DGAT1 K232A quantitative trait nucleotide in affecting milk yield and composition[END_REF]. Furthermore, chromosomes 5 and 20 were more clearly revealed by sPLS than by PLS. The same comments can be made for most traits. The differences between PLS and sPLS were obvious for all traits with higher VIP coefficients in sPLS due to the reduced number of SNP considered. Conception rate showed many regions of interest with both methods. The aim of this study was to highlight differences and similarities between PLS and sPLS. We are not yet able to affirm that the different genome areas localized by PLS or sPLS correspond to QTL locations. This is currently under study [START_REF] Colombani | Could genomic selection be efficient to detect QTL? 15th QTLMAS Workshop[END_REF].

The aim of this study was to compare PLS variants with currently used methods in the evaluation of dairy cattle. Table 4 shows the correlation between observed and estimated DYD for all the traits, with 3 methods of genomic selection: PLS, sPLS, and GBLUP, compared with pedigree-based BLUP. On average, the correlations obtained by genomic selection methods were significantly higher than with pedigree-based BLUP for the 5 production traits concerned (0.426 for pedigree-based BLUP and 0.614, 0.612, and 0.630 for PLS, sPLS, and GBLUP, respectively). The differences between pedigree-based BLUP and genomic selection methods were not as clear for the conception rate trait, with a correlation of 0.28 for BLUP and 0.21 and 0.35 for sPLS and GBLUP. For conception rate, the BLUP correlation and the PLS correlation were not significantly different. Genomic BLUP (from 0.35 to 0.73), PLS (from 0.25 to 0.71), and sPLS (from 0.21 to 0.72) gave similar results for all traits concerned except for milk yield and conception rate, for which GBLUP gave significantly better results.

As expected, the genomic selection methods tested in this study were more efficient than pedigree-based BLUP. Genomic BLUP was accurate for use with French Holstein data with significantly higher accuracy for some traits. However, PLS and sPLS methods were comparable to GBLUP if we considered one trait at a time. Regarding computing time, GBLUP requires one inversion of the genomic relationship matrix for all traits, which took about 1 h. Then, once the genomic relationship matrix was inverted, computation was a matter of seconds. For each trait, PLS took about 10 min and sPLS took less than 10 min, depending on the number of SNP selected. The disadvantage of PLS and sPLS with respect to GBLUP for some traits could be overcome by building a multi-trait model.

CONCLUSIONS

Sparse PLS regression was compared with PLS and with 2 currently used methods in the evaluation of dairy cattle: pedigree-based BLUP and GBLUP. The results demonstrated that PLS variants were more efficient than pedigree-based BLUP but less accurate than GBLUP for 2 out of 5 traits. Furthermore, GBLUP provided a clear biologic interpretation by the use of a genomic relationship matrix that PLS and sPLS may lack, and PLS and sPLS do not provide reliabilities of GEBV, in contrast to GBLUP. Sparse PLS enabled easier identification of relevant variables than PLS, which are possibly associated with QTL regions. Currently, more and more markers are being genotyped, forcing the handling of increasing quantities of data and consequently a critical need exists for methods that perform well in this context. Sparse PLS could be used as a preliminary step in genomic selection to reduce the number of SNP used in the prediction equations provided by other genomic selection methods, such as Bayesian methods. Most importantly, the sPLS algorithm is fast to compute even with a large reference population and a large number of explanatory variables in a one-trait evaluation.

Figure 1 .

 1 Figure 1. Root mean square error of prediction (RMSEP) for the 6 traits studied plotted against the number of dimensions for partial least squares regression.

Figure 2 .

 2 Figure 2. Correlation between observed and predicted daughter yield deviations for the 6 traits studied plotted against the number of dimensions for partial least squares regression.

Figure 3 .

 3 Figure 3. Distribution of effective daughter contribution (EDC) in the training data set.

Table 1 .

 1 Root mean square error of prediction (RMSEP) in partial least squares (PLS) and each sparse PLS tested as a function of the percentage of SNP selected in each latent variable [minimum RMSEP (3 decimal places) in bold]

					Sparse PLS (% of the SNP data set selected)				PLS
	Variable	0.2	0.4	0.6	0.8	1	2	3	4	5	10	100
	Milk yield	0.17	0.16	0.17	0.16	0.16	0.16	0.15	0.15	0.15	0.15	0.15
	Fat yield	0.18	0.20	0.19	0.18	0.21	0.18	0.18	0.19	0.18	0.18	0.18
	Protein yield	0.14	0.15	0.14	0.14	0.14	0.14	0.13	0.13	0.13	0.13	0.13
	Fat content	0.60	0.67	0.58	0.56	0.57	0.58	0.57	0.59	0.63	0.60	0.57
	Protein content	0.64	0.62	0.59	0.60	0.58	0.57	0.56	0.59	0.56	0.56	0.58
	Conception rate	0.94	0.83	0.81	0.83	0.86	0.81	0.87	0.79	0.83	0.82	0.78

Table 2 .

 2 Effective daughter contribution-weighted correlations (ρ) and regression slopes (b) provided by partial least squares regression (PLS) and sparse PLS (sPLS)

		Milk	Fat	Protein	Fat	Protein	Conception
	Item	yield	yield	yield	content	content	rate
	PLS
	ρ	0.53	0.58	0.55	0.70	0.71	0.33
	b Dim 1	0.65 42	0.83 37	0.67 36	0.80 83	0.83 75	0.60 29
	sPLS
	ρ	0.48	0.54	0.51	0.66	0.65	0.29
	b	0.53	0.70	0.60	0.69	0.76	0.54
	Dim No. of SNP 2	44 22,948	43 16,296	38 32,578	50 9,832	51 26,034	27 20,150
	1 Number of latent variables or dimensions included in the final model.

Table 4 .

 4 Correlations between observed daughter yield deviations (DYD) and predicted DYD provided by partial least squares regression (PLS), sparse PLS (sPLS), pedigree-based BLUP (BLUP), and genomic BLUP (GBLUP)

		Milk	Fat	Protein	Fat	Protein	Conception
	Item	yield	yield	yield	content	content	rate
	BLUP	0.38	0.40	0.44	0.44	0.47	0.28
	PLS	0.52	0.58	0.55	0.71	0.71	0.25
	sPLS	0.50	0.59	0.53	0.72	0.72	0.21
	GBLUP	0.56	0.59	0.55	0.72	0.73	0.35

Table 3 .

 3 Correlations (ρ) and regression slopes (b) provided by partial least squares regression (PLS) and sparse PLS (sPLS) without effective daughter contribution Number of latent variables or dimensions included in the final model. 2 Number of SNP selected by sPLS.

		Milk	Fat	Protein	Fat	Protein	Conception
	Item	yield	yield	yield	content	content	rate
	PLS
	ρ	0.52	0.58	0.55	0.71	0.71	0.25
	b Dim 1	0.60 2 0	0.77 9	0.64 1 0	0.80 2 3	0.82 2 6	0.62 3
	sPLS
	ρ	0.50	0.59	0.53	0.72	0.72	0.21
	b	0.56	0.76	0.59	0.81	0.77	0.49
	Dim No. of SNP 2	14 17,408	10 6,779	10 26,458	11 2,870	23 27,191	3 4,592

1

, and G is the genomic relationship matrix as defined by VanRaden (2008):

Number of SNP selected by sPLS.

Journal of Dairy Science Vol. 95 No. 4, 2012

Figure 4. Variable importance in projection (VIP) coefficients from partial least squares regression (PLS) and sparse PLS for the 6 traits without effective daughter contribution (EDC).

ACKNOWLEDGMENTS

This work was supported by the French project AMASGEN, financed by the French National Research Agency (ANR, Paris France) and ApisGene (Paris, France). Labogena (Jouy-en-Josas, France) is gratefully acknowledged for providing the genotypes. We thank the reviewers for their constructive suggestions that enabled major improvement of the manuscript