Multivariate restricted maximum likelihood estimation of genetic parameters for growth, carcass and meat quality traits in French Large White and French Landrace pigs

ADG2, DP, ECLC and MQI were respectively 0.30, 0.64, 0.22, 0.52, 0.39, 0.60, 0.33 in the LW and 0.34, 0.56, 0.25, 0.46, 0.31, 0.68, 0.23 in the LR breed. Common litter effects ranged from 5% (ABT in LW breed) to 16% (ADG2 in LR breed) of phenotypic variance. Growth traits and FCR exhibited favourable genetic correlations, but were unfavourably correlated to DP and carcass lean content. MQI also showed unfavourable though generally low genetic correlations with all the other traits. These antagonisms were apparent in both breeds, but tended to be larger in the LW than in the LR breed.


INTRODUCTION
Best linear unbiased prediction (BLUP) applied to an individual animal model (IAM; Henderson, 1988) is becoming increasingly used to predict breeding values in most species of farm animals (Carabano and Alenda, 1990). In pigs, numerous studies have been conducted to evaluate the advantage of using BLUP instead of standard selection indexes in a variety of situations (Belonsky and Kennedy, 1988;Keele et al, 1988;Sorensen, 1988;Wray, 1989;Long et al, 1990;De Vries et al, 1990;Roehe et al, 1990). National genetic evaluation programmes based on the BLUP-IAM technique have been implemented in several countries since 1985 (Hudson and Kennedy, 1985;Van Hofstraeten and Vandepitte, 1988;Harris et al, 1989; Sorensen and Vernessen, 1991). In France, a national genetic evaluation programme based on BLUP-IAM methodology is being implemented for both production and reproduction traits. The first step of this project concerns the genetic evaluation of station tested animals, ie x5 6 000 animals/yr. The use of BLUP procedures requires knowledge of variance components in the unselected base population. In practice, these components have to be estimated from available data. The method of restricted maximum likelihood (REML; Patterson and Thompson, 1971) has been shown to be the method of choice for estimating variance components in selected populations, mainly because of its ability to account for selection bias (Gianola et al, 1986). The aim of this study is to estimate genetic parameters of French Large White and Landrace breeds for traits measured in central test stations using a multiple trait IAM-REML.

MATERIAL
Animals and data recording Genetic evaluation of centrally tested pigs in France is currently based on a combined selection index involving performance of the candidate for selection and of 1 full-sib. Full-sibs are females or barrows castrated before entering the station, ie shortly after weaning. Breeders are asked to choose centrally tested animals at random in litters of at least 8 piglets. The data analysed in the present study concerned both candidates for selection and relatives slaughtered at the end of test. In order to keep computing costs reasonable, 2 computationally manageable sets of data (ie 3 671 and 3 630 candidates, 3 039 and 2 695 relatives in, respectively, LW and LR breeds) were created by considering all the data collected in the 3 stations in which both candidates for selection and relatives were tested. The period of time considered was [1985][1986][1987][1988][1989][1990] for the LW and 1980-1990 for the LR. For computational reasons, only 2 generations of ancestors, ie the parents and grand-parents of tested animals, were considered. The structure of the 2 data sets is shown in table I.
Testing of both candidates for selection and their relatives was performed in discontinuous batches. A batch was defined by the year of test, the testing station and the 2-wk period of entering station (! 4 levels for each year x station combination) and will consequently be referred to as the year x station x batch (Y x S x B) effect hereafter.
Young boars were tested between 35 and 90 kg live weight. Until 1988, they were penned in groups of 4, but individually fed on a liberal feeding diet based on the voluntary feed intake of the animal during 2 daily meals of 20 min each. From 1988 they were allotted to pens of 10-12 animals, with ad libitum feeding. Animals were weighed twice at the beginning and at the end of the test. Dates of measurement were chosen so that the 2 initial and the 2 final weights flanked 35 and 90 kg, respectively. This allowed us to adjust the different traits to a constant initial and/or final weight. Feed intake was recorded individually during the whole test until 1988. Backfat thickness was measured twice at the same time as final weights. The ultrasonic measurements were taken on each side of the spine, 4 cm from the mid-dorsal line at the levels of the shoulder, the last rib and the hip joint, respectively.
Animals from the second group were tested between 35 and 100 kg live weight. They were allotted to pens of 2 animals until 1988, and to pens of 10-12 animals thereafter. Pigs were fed ad libitum, but feed intake was not individually recorded. Animals were weighed once ! 35 kg and twice ! 100 kg. They were slaughtered during the week following the last weight measurement. Standardized cutting of one half-carcass was performed as described by Ollivier (1970Ollivier ( ) until 1988, and since 1989 as described in Anonymous (1990).
Three measurements of meat quality were taken on the ham on the day after slaughter, ie: 1) ultimate pH (pHu) of Adductor femoris muscle; 2) water-holding capacity (WHC) as assessed by the time (in tens of s) for a piece of pH paper to become wet when put on the freshly cut surface of Biceps femoris (until 1988) or Gluteus superficialis (since 1989) muscle; 3) reflectance (REF) of Gluteus superficialis muscle at 630 nm, using a Manuflex reflectometer (scale 0 to 1000).

Traits analyzed
Seven different traits were defined from the above-mentioned measurements: average daily gain (ADG1) and feed conversion ratio (FCR) from 35 to 90 kg and backfat thickness at 90 kg (ABT) of young boars candidates for selection. Adjustments to a constant initial and/or final weight were made by interpolation between the 2 weights flanking 35 and/or 90 kg, respectively; -average daily gain from 35 to 100 kg (ADG2), dressing percentage (DP), estimated carcass lean content (ECLC) of the carcass with head (EC reference) and meat quality index (MQI) of candidates' relatives. DP was computed as the ratio of carcass weight with head and feet to live weight. ECLC was estimated from the relative weight of 6 joints expressed as percentage of half-carcass weight, according to the following prediction equations: [1] ECLCI = -3.539 + 0.751 (percentage of ham) + 1.216 (percentage of loin) -0.610 (percentage of backfat) -0.453 (percentage of leaf fat) + 0.328 (percentage of belly) ; [2] ECLC2 = -42.035 + 1.282 (percentage of ham) + 1.818 (percentage of loin) -0.678 (percentage of backfat) + 0.040 (percentage of leaf fat) + 0.701 (percentage of belly) + 0.616 (percentage of shoulder). Equation [1] was used until 1988 and was replaced by equation [2] simultaneously with the change in cutting procedure. Both equations have been shown to be highly correlated with the true carcass lean content (R Z = 0.911 for the first equation and 0.930 for the second) so that ECLCI and ECLC2 were considered as the same trait. The meat quality index (MQI), established as a predictor of the technological yield of Parisian ham processing, was computed as a linear function of the 3 meat quality measurements defined above (Guéblez et al, 1990): MQI = -35 + 8.329 pHu + 0.127 WHC -0.00744 REF.
Elementary statistics for the 7 traits studied are shown in table II.

Model
The model varied according to the trait, but had the following basic structure (in matrix notation): where y is the vector of observations for the 7 traits, b is the vector of fixed effects, p is the vector of litter effects, a is the vector of additive genetic values of animals, e is the vector of residuals and X, W, Z are incidence matrices relating observations to the effects in the model. Location and dispersion parameters for the random effects were as follows: where: The exact model used for each of the 7 traits in shown in table III. ADG2 was not pre-adjusted for initial weight, which was consequently included as a covariable in the model for that trait. Year x station x slaughter date has been shown to be the most important environmental factor affecting meat quality traits (Monin, 1983) and was therefore included as a fixed effect in the model for MQI instead of Y x S x B.

Computing strategies
The (co)variance components were estimated using the derivative free multiple trait restricted maximum likelihood procedure described by Groeneveld (1991).
It was not possible computationally to analyse all traits simultaneously and only 2-trait analyses reached convergence at a reasonable computing cost. Hence, 21 2trait analyses were performed for each breed. A derivative-free Quasi-Newton (DF--QN) algorithm was used to maximise the likelihood when possible because of its good convergence rate. The subroutine E04JAF from the NAG library (Numerical Algorithms Group, 1990) was used for this purpose. The convergence criterion was defined as the standardized norm of the changes in variance components between 2 consecutive iterations. The stopping criterion was set at 10-7 . However, convergence could not be reached in all cases. The DF-QN procedure approximates the matrix of the second derivatives (Hessian matrix) by finite difference and works relatively well for smooth portions of the likelihood function with continuous derivatives. But when one of the parameters is located at the border of the parameter space, adding a finite difference may generate numerical values outside the parameter space, in which case DF-QN does not work. When this situation occurred, the DF-QN algorithm was replaced by a Downhill Simplex based algorithm (Press and Flannery, 1986). Being a sampling procedure, Downhill Simplex can handle values outside the parameter space much better, thus finding the optimum where gradient based optimisers fail (Kovac, 1992). The convergence criterion was determined as the fractional range from the highest to the lowest likelihood value of the vertex (Groeneveld, 1991). The stopping criterion was set at 10-'!.
Approximate asymptotic standard errors of variance components and genetic parameters were obtained from the approximate Hessian matrix computed by the Quasi-Newton E04JBF subroutine from the NAG library (Numerical Algorithms Group, 1990). These approximate standard errors are lower bounds of the standard errors for the parameter estimates (Gianola, 1989).

RESULTS
The Quasi-Newton algorithm was used for all pairs of traits except the (ADG1-ADG2) pair, which had convergence problems. Computing time on an IBM 3090-17T varied from 31 to 165 min of CPU time in the LW breed and from 32 to 172 min of CPU time in the LR breed. The Downhill Simplex algorithm was used for the (ADG1-ADG2) pair of traits. Computing times were 366 and 591 min of CPU time for LW and LR breeds, respectively.
Six estimates of variance components were available for each trait. Variation among estimates was small (variation coefficient of heritabilities was lower than 5%, in both breeds), so that only average values are presented (table IV). ABT and ECLC had the largest heritabilities (h 2 ) in both breeds. The LW had a higher h 2 for ABT (0.64 vs 0.56), but a lower one for ECLC (0.60 vs 0.68) than the LR. ADG2 was much more heritable in both breeds than ADG1 (0.52 and 0.46 vs 0.30 and 0.34 in LW and LR, respectively). Heritability values were very similar in the 2 breeds for FCR, but larger in the LW breed for DP and MQI (0.39 vs 0.31 and 0.33 vs 0.23). Common environmental effects were not negligible for all the traits, with estimates ranging from 5 to 16% of the phenotypic variance. The largest values were obtained for growth rate and, in the LR breed, for MQI. Differences between breeds were small, except for ABT (0.05 vs 0.12) and MQI (0.06 vs 0.15).
Estimates of genetic and phenotypic correlations for the LW and LR breeds are shown in tables V and VI, respectively. Genetic correlations were well estimated for traits measured on the same animals but had larger standard errors for traits measured on different animals. Genetic correlations between ADG1 and ADG2 were close to unity in both breeds (0.97 and 0.99 in LW and LR breeds, respectively).
Similarly, traits predicting carcass lean content, ie ABT in candidates and ECLC in relatives were highly correlated (-0.86 and -0.90 in LW and LR breeds, respectively). Growth traits and FCR were negatively, ie favourably, correlated (-0.61 and -0.63 on average for ADG1 and ADG2, respectively) but showed unfavourable genetic correlations with carcass traits. Genetic correlations between MQI and growth or carcass traits were generally low, but were unfavourable in both breeds. Some correlations differed between breeds. This was particularly the case for the correlations between MQI and ECLC (-0.44 in LW V8 ! 0.02 in LR.), ADG2 and DP (0.08 in LW v.s -0.53 in LR) or ADG1 and ABT (0.48 in LW v.s 0.25 in LR). In general, the antagonisms between carcass lean content on one hand, growth rate and meat quality on the other hand, tended to be larger in the LW breed than in the LR breed.

Methodology
There is general agreement that BLUP and RENIL methodologies using an animal model are the methods of choice for estimating location and dispersion parameters for traits that can be described by linear models, because of their desirable statistical and genetic properties (Harville, 1977;Kennedy et al, 1988 ;Robinson, 1991). In particular, this method accounts for the effects of selection if all the information related to selection is included in the analysis.
Practical applications of BLUP-IAM in pig breeding are steadily increasing. Its use has been greatly facilitated over tlre last few years by the increasing power of computers and the appearance of general purpose software such as &dquo;PEST&dquo; (Groeneveld et al, 1990), or &dquo;PIGBLUP&dquo; (Brandt, 1990). On the other hand, the use of multivariate REML-IAM is still infrequent, mainly because of its huge computational requirements. Limited applications are possible but, in most cases, at the expense of some departures from the ideal situation. In the present case, a total of 72 (co)variances had to be estimated. A problem of comparable size (ie 60 covariance components) was treated by Groeneveld (1991), but with larger computing facilities and a more favourable data structure. The present data set had several drawbacks, such as different traits measured on different individuals, a low number of animals per parent and per litter, and a low number of performancetested parents. Each of them resulted in convergence problems, which were solved by: 1) limiting the number of generations of ancestors; 2) 1>utting the litter covariance component to zero ; 3) running 2-trait analyses. Limiting the number of ancestors and decreasing the number of (co)variance components resulted in a reduction of the number of likelihood functions to be computed and of the CPU time per likelihood. The impact of 1) was investigated with univariate models. Adding a third generation of ancestors considerably increased computing time, but did not change variance components at all. The effect of 2) was investigated for the (ADG1, ABT) and (ECLC, ABT) pairs of traits. Non-negligible correlation for the litter effect was obtained in botlr cases (0.21 and 0.75, respectively). However, this was mainly due to the low magnitude of litter variances, litter covariance being much smaller (< -1/10) than genetic or residual covariances. As a consequence, the impact on the other covariance components and on predicted breeding values was very limited. The impact of 3) was theoretically more critical. First, all traits selected upon should be included in the analysis to take into account properly the effects of selection. The effect of this simplification could not be tested, but the stability of estimates of variance components obtained from different 2-trait analyses tends to indicate that it should be rather limited, at least for variances. Then, tlre positive definiteness of the variance-covariance matrices is no longer guaranteed. Indeed, inconsistencies were obtained in both breeds. Yet they were mainly due to the (ADG1, ADG2) pair of traits, and positive definiterress was obtained when 1 of the 2 traits was removed from the matrix.
Another underlying assumption of the present REML analysis was the homogeneity of within-station variances. In fact, as first mentioned by Ollivier et al (1981), noticeable differences in variability may exist between stations. Methods have recently been developed to estimate variance components in situations of heteroskedasticity (Foulley et al, 1990(Foulley et al, , 1992. These methods are computationally much more demanding than standard RE1!IL and are currently intractable in applications such as the present one.

Genetic parameters
Genetic parameters in French LW and LR breeds using data from testing stations had previously been estimated by Ollivier (1970), Ollivier et at (1980Ollivier et at ( , 1981, Tibau i Font and Ollivier (1984), Sellier et at (1985) and Cole et at (1988). All these studies used classical estimation methods such as parent-offspring regression or Henderson's methods (Henderson, 1953 (1980,1981) and Tibau i Font and Ollivier (1984). This result suggests that growth rate is genetically more variable and more heritable under an ad libitum than under a restricted feeding diet. Variation in genetic parameter estimates due to feeding regime is well established in pigs (Wyllie et al, 1979;Cameron et al, 1988). As suggested by Cameron (1990), heritability estimates of ADG and ABT might be lower under a restricted feeding diet because competition effects would increase non-additive variance components. In spite of these heritabilities, feeding regime, sex and period of test differences, the 2 traits appeared as genetically very similar, the genetic correlation being very close to unity in both breeds. This result somewhat disagrees with the value of 0.55 reported by Tibau i Font and Ollivier (1984).
The signs of the genetic correlations are generally in good agreement with literature means (table VII). In general, the magnitude of these correlations is closer to previous French results and literature means in the LR than in the LW breed. Yet, several exceptions have to be mentioned. FCR and ECLC are less correlated than usually reported in the literature (-0.22 vs values ranging from -0.30 to -0.57: Lundeheim et al, 1980;Ollivier et al, 1981;Tibau i Font and Ollivier, 1984;Costa et al, 1986;Johansson et al, 1987;Van Hofstraeten and Vandepitte, 1988). The lack of antagonism between ECLC and MQI disagrees with most literature results in the LR breed (Malmfors and Nilsson, 1979;Lundeheim et al, 1980;Ollivier, 1983;Andersen and Vestergaard, 1984;Busse and Groeneveld, 1986;Merks, 1987;Schworer et at, i987;Cole et al, 1988). Differences between estimates may be related to variations in the frequency of the halothane gene in the populations studied (Sellier, 1988). ABT and ECLC are more closely correlated in both breeds than in previous French studies (-0.38;Tibau i Font and Ollivier, 1984) or than the literature average (-0.65). The large negative correlation between ADG2 and DP is also in disagreement with the only available literature estimates (0.93 ! 0.05; Cameron, 1990;0.36;Johansson et at 1987). There is no clear explanation for these large differences.
Discrepancies with literature means and previous French results are greater in the LW breed. The antagonism between growth rate and carcass composition, particularly between ADG1 and ABT, appears to be stronger than previously reported by Tibau i Font and Ollivier (1984), Kwon et at (1986), Van Steenbergen et at (1989), and Kaplon et at (1991). However, literature estimates are quite variable and large unfavourable correlations were also reported in LW or Yorkshire by Merks (1987), Savoie and Minvielle (1988), and Johansson et at (1987). The correlations between growth rate and FCR also differ from usual patterns; the relationship between ADG and FCR is generally closer under a restricted feeding diet than with an ad libitum diet (Ollivier et al, 1980;Tibau i Font and Ollivier, 1984;Cameron et al, 1988). The lack of antagonism between growth rate and meat quality is in disagreement with most studies in the LW breed (Ollivier, 1983;Johansson et al, 1987;Cole et al, 1988). In contrast, the unfavourable correlation between MQI and ECLC is larger than those reported by Lundeheim et al (1980), Johansson et at (1987), or Cole et al (1988).
CONCLUSION Accurate estimates of genetic parameters are essential to evaluate and compare alternative breeding plans as well as to predict breeding values. The genetic parameters estimated in the present study are likely to be more adequate than previous estimates, because multivariate REML-IAM procedures allow the fit of more realistic models, which are similar to those used for the genetic evaluation and give parameters which tend to be less biased by selection. However, multivariate REML-IAM procedures are computationally very demanding, so that only rather small data sets and only few traits (two in the present case) can be treated simultaneously. As a consequence, selection is probably not fully taken into account and variance-covariance matrices are less likely to be positive-definite. Better optimization algorithms, more powerful computers and better data structures will all be necessary to be able to run full multiple trait REML-IAM.