Genetic parameters of the twisted legs syndrome in broiler chickens

Summary - Genetic parameters of two types of angulations, described in the twisted legs syndrome as ’valgus’ and bilateral or unilateral ’varus’, were investigated in two commercial broiler strains. In the first line, 14 264 chickens of both sexes born from 111 sires, 76 maternal grandsires and 768 dams were studied. In the second line, corresponding figures were 8 164 chickens of both sexes born from 94 sires, 71 maternal grandsires and 553 dams. Chickens were classified at 6 weeks according to the type of pathology. Since deformities under study were unordered categorical traits, a generalized linear model using a multinomial logistic transformation as link function was applied. Location parameters were estimated by ’maximum a posteriori’, and variance components by ’maximum marginal likelihood’ using a tilde-hat approximation. The model of analysis took into account the fixed effects of the hatch and the sex as well as the random effects of the sire, maternal grandsire and dam within maternal grandsire. ’Pseudoheritability’ of latent susceptibility to valgus was equal to 0.16 and 0.29 in lines A and B respectively, when estimated from the sire/maternal-grandsire component, and to 0.40 and 0.35 respectively, when estimated from the dam component; for varus, estimates of the pseudoheritability were equal to 0.21 and 0.24 in lines A and B when estimated from the sire/maternal-grandsire component and to 0.30 and 0.26 when estimated from the dam component. Higher values of the dam heritabilities could suggest the existence of maternal or dominance effects. The average estimated genetic correlations between valgus and varus obtained from sire/maternal-grandsire and dam components were small to


INTRODUCTION
Selection of meat-type chickens has been aimed until now mainly at increasing growth rate. Phenotypic change of growth rate during the past 40 years has been spectacular: according to L'Hospitalier et al (1986), who compared eight commercial crosses from four countries, mean daily gain between hatching day and 42nd day increased from 20 to 47 g/day between 1962 and 1985. Even if this increase can be partly explained by improvements in nutrition and management, it seems that a great part of the progress is due to selection. Annual genetic gain for body weight measured at 6 weeks estimated in two commercial meat-type strains on large data sets, was equal to 94.6 g for the sire strain and to 72.6 g for the dam strain (Jego et al, 1995). However, leg disorders have appeared at higher frequencies concomitantly to this selection on growth performance. Hartmann and Flock (1979) compared the incidence of twisted legs in commercial lines between 1963-1968 and 1977-1978. Between these two periods, the incidence measured on male offspring at slaughter had increased from 20 to 32% (70% when including slight deformities).
Leg disorders have important economic consequences, such as a decrease of body weight (Hartmann and Flock, 1979;Leenstra et al, 1984;Leterrier and Nys, 1992), and culling of the most affected birds. Furthermore, as discussed by Sorensen (1989), decreasing leg disorders should contribute to improving animal welfare. Twisted legs are one of the most frequent deformities among leg disorders (Stuart, 1989). The goal of this study was to estimate in two meat-type strains, genetic parameters of the two main deformities observed in this syndrome, ie, 'valgus' and 'varus' angulations. Different angulations were scored as disjoint categories, and a multinomial sampling model was assumed. Since usual linear methods are not optimal for such traits, the generalized linear model theory (McCullagh and Nelder, 1989) was employed. Because scoring was considered as an unordered polytomy, a multivariate logit transformation (Cox, 1970), previously applied in a mixed-model context by Gianola (1980), provided the relevant link function between continuous latent variates and expected occurrences.

Animals
The present study was conducted on 14264 chickens born from 111 sires, 76 maternal grandsires and 768 dams in line A and 8 164 chickens born from 94 sires, 71 maternal grandsires and 553 dams in line B. Both male and female animals were considered. Chickens in the A and B lines were kept on the floor in 14 and 11 hatches respectively. Animals were examined for twisted legs at 6 weeks of age and the gravity of the deformity was recorded as mild or severe. According to the suggestions of Leterrier and Nys (1992), 'valgus' and 'varus' angulations of the tibiotarsal articulation were distinguished. The former is often bilateral, and displaced tendons are observed only in the most severe cases. The latter is most often associated with a medial tendinous displacement. 'Unilateral' or 'bilateral' varus were further distinguished as suggested by Riddell (1983) and Leterrier and Nys (1992), but these data were pooled for the analysis. All animals were thus classified as healthy, valgus or varus.

Statistical model
Let 7ri = (7ril, !i2, ... , !ri!)' be the vector of the probabilities of the different discrete categories in the ith (i = 1, ... , s) stratum (ie, combination of levels of the factors involved in the model). Because an animal can only be assigned to one category at a time, the probability of observing n 2k animals in the kth (k = 1, ... , c) category was assumed to be given by the multinomial distribution: where n i . was the total number of observations in the ith stratum. McCullagh and Nelder (1989) distinguished between ordinal and nominal polytomous data. In the ordinal case, various categorical responses can be classified and considered as expressions of a single latent variate in reference to several thresholds. In the nominal case, such a classification is impossible. This is clearly our situation because physiological studies on leg disorders have led to the conclusion that valgus and varus could correspond to different defects (Leterrier and Nys, 1992;Riddell, 1992), ie, each of them would be related to one specific susceptibility. In this case, the appropriate link between the latent risk variates and the observations is less easy to set up than in the threshold model. The aforementioned logit multinomial model corresponds to one possible situation, in which discrete expression corresponds to the result of a competition between latent variates: Gianola (1982), Judge et al (1985), and Albert and Chib (1993) remind us that the discrete observed code corresponds to the largest value among the c underlying logistic latent susceptibility variates. Let us assume these c standardized logistic latent variates to be yi ! Y2' ... , Y! with means Jli, Jl2, ... , Jl!, so that y* = !,2 + e2 , with Var(ei ) _ 7[ 2 /3 and Cov(Ei, E j) _ !r2/6. The differences between the variates Yj and a given variate Y i (j ! i) are still logistic variates with the same standardized variance and covariance (Johnson and Kotz, 1970). Assuming that y* is the largest amongst all the y * values, the probability of observing the ith category corresponds to: When p) and pg are known, the probability of such an event is given by the cumulative distribution function (denoted by F) of the c &mdash; 1 variates E! -E! (i 34 i), 1 following a multivariate standardized logistic distribution. Indeed, y) < yi implies: p) + éj < f -Li + El and E* -6 * < pgf -Lj. Therefore (Johnson and Kotz, 1970): Considering category c as a reference, and after reparameterization, one can write: Similar developments can be found in Bock and Jones (1968) and Gianola (1982). Since only c &mdash; 1 categories are independent, only c &mdash; 1 logits corresponding to the differences between the expectations of the various logistic latent variates and the expectation of the logistic latent variate chosen as the reference can be estimated.
Provided that the baseline 'risk' associated with the healthy category (noted here as the cth category) was chosen as this reference, the probabilities of response for the ith stratum are: where !,2k = log( 7 ri k/7 r ic ) was the kth (k = 1,..., c &mdash; 1) logit for stratum i.
Hence, and, as Inferences on location and dispersion parameters pertain to latent susceptibility variates related to each deformity (ie, 'latent valgus' and 'latent varus'), bearing in mind that these parameters depend on the relevant deformity and on the reference category. For this reason, reference to valgus and varus hereafter should be considered as applying to the corresponding latent variables.

Genetic model
The genetic model used for logits assumed additive genetic effects and no maternal effects. In this context, the statistical model used for logits, assuming three categories, namely valgus, varus and healthy, can be represented by: where !! is an (s x 1) vector, b k a (p x 1) vector pertaining to the fixed effects of the hatch (numbers of levels were 14 and 11 in lines A and B respectively) and sex, Ulk a (q, x 1) vector pertaining to random effects of sires and maternal grandsires and u 2k a (q 2 x 1) vector pertaining to random effects of dams (within maternal grandsires) on the kth logit. X, Z l and Z 2 are, respectively, (s x p), (s x q l ) and (s x q 2 ) known incidence matrices. Sire elements of the vector u ik represent one half of their additive direct genetic value. The maternal grandsire effect represents a quarter of his additive genetic value, so that it is expressed as 0.5 times the corresponding sire effect. Therefore, Z l = Zs+0.5 Z MGS , where Zs and Z MGS were incidence matrices pertaining to sires and maternal grandsires respectively, with the appropriate number of zero columns to give them the same (s x q l ) dimensions.
Elements of the dam within maternal-grandsire vector U2k represented one half of her additive direct genetic value for trait k deviated from the contribution of her sire effect, which itself was equal to one quarter of his genetic value (eg, see Manfredi et al, 1991).
As usual, genetic effects of male and female ancestors were assumed to follow a multivariate normal distribution with E(u l ) _ !J, E(u2) _ ! and with u í = ( U í l' U í 2 ) and u 2 = ( u 2 1 , u2 2 ). G 1 and G 2 are (2 x 2) matrices of the genetic variance-covariance components for the male ancestors (either sires or maternal grandsires) and dam within maternal grandsire effects respectively. From previous considerations, it can be shown that element (i,j) of G 1 (denoted g l ij) and element (i, j) of G 2 (denoted g 2ij ) have expectations respectively equal to 1/4 and 3/16 of the genetic variance (or covariance) pertaining to logits i and j. A l is the relationship matrix among sires and maternal grandsires of recorded animals. This was computed by considering relationships created by common male ancestors occurring in the pedigree (Henderson, 1975); there were totals of 376 and 357 male breeding animals in lines A and B respectively. A 2 was the relationship matrix among the dams created by considering relationships due to common female ancestors available in the pedigree information, ie, totals of 1625 and 1466 females in the two lines respectively. The linear logistic model presents peculiarities in contrast to the probit model for ordered categorical data. In the latter, the residual variance is equal to 1, and the marginal distribution of the underlying variate is normal. In the logistic distribution, the residual variance is !r2/3, as noted earlier. Although the conditional (given the random effects) distribution of the latent variate is logistic, the unconditional distribution is not, because the random effects are normal. However, the total variance in the latent scale decomposes additively. Because distribution of the unobserved latent variate corresponds to the sum of a normal logit and a standardized logistic residual, we shall use the term 'pseudoheritability'. In this study, pseudoheritabilities based on the variance components were: hi i = 4gi2i/a!2, and h2 i = (16/3)g 2 ida; i where the phenotypic variance a!2 = gl2i (ie, variance between sire groups) +1/4g lii (ie, variance between maternal grandsire groups) +g2!i (ie, variance between dam groups within maternal grandsire) + 7[ 2 /3 (ie, residual variance), ie, 2 . = 1.25 g ,ji + g 2 iz + !2/3. Genetic correlations were calculated as g l ij / ( g lii gljj )°. 5 and g2ij / (g2ii g2jj )°.5 respectively, from sire/maternal-grandsire and dam components. As mentioned above, the residual correlation is forced to be 0.5 in the logit model. Consequently, phenotypic correlations (not given) should be considered as pseudocorrelations.
Estimation of location parameters by maximum a posteriori (MAP) Location parameters were defined as 0' = (b', a'), with b' = (b', b') and a' = (ul 1 1, u2 1 , U!2' U22 ) where b and uij are defined in equation [3]. Following a Bayesian approach, they were estimated by maximizing the log of the posterior density for known dispersion matrices G I and G 2 according to Bayes theorem: Such an estimation is therefore MAP. As prior information about the distribution of b was considered to be vague, the a priori density of b, p(b), was flat. From !4!, the log of the a priori density of a is: where E a was obtained from £ u after sorting by trait.
For given b and a, the probabilities of each category in each population can be obtained from expressions [2a] or !2b!; equation [1] then gives the conditional probability of the observed data. Hence, the logarithm of the posterior density is equal to: Because finding the mode of L(0]£a) led to non-linear equations, the iterative Newton-Raphson algorithm was used; first and second derivatives of L(9) with respect to fixed and random effects are described in AP!e!cdix A. After rearrangement, the system of equations providing solutions is: where Z = [Z l Z 2] and E&dquo;' contained elements of E d corresponding to traits 1 and 1'. W kk (k = 1, 2) and W kk' (k = 1, 2; k' = 1, 2; k' 3 4 k) were (s x s) diagonal matrices: Vf l was obtained by where Vk was an (s x 1) vector, Conditional probabilities !rik(k = 1, 2) were calculated from !2a!, using estimates of b and a obtained at the round t. As already described in analyses of discrete traits with a threshold model (Gianola and Foulley, 1983), the system of equations providing MAP estimates can also be written in a form similar to linear mixedmodels equations; indeed, the right hand side of [6] can be expressed as: where the y j 's are the following working variates: Foulley (1993)  where D! is the inverse of the block diagonal part of the coefficient matrix pertaining to the genetic effects of sex k, which is obtained from the coefficient matrix in [6] after absorbing all the other effects and considering a block per breeding animal (so that approximation for u is better than when using a diagonal matrix since it takes into account correlation between the traits); Rhs k is a vector corresponding to the right-hand-side terms in [6] after absorption of the fixed effects and the effects of the ancestors of the other sex.
Expectation of the quadratic form Qkz! = GkiA §!3k; (k = 1, 2) is analogous to the form obtained by Van Raden and Jung (1988): where D kjm is the submatrix of D k (k = 1, 2) pertaining to traits j and m, and M!,&dquo;,,L is the submatrix of the absorption matrix pertaining to traits m and t. This algorithm did not recover standard errors, and methods based on second derivatives should be considered. An even better description of uncertainty could stem from a Monte-Carlo Markov chain implementation but would imply heavy calculations.

Numerical aspects
As described by Manfredi et al (1991), the complete analysis required three levels of nested iterations: outermost iterations for estimation of the variance and covariance components, the Newton-Raphson iterations for estimation of the location parameters of the model and innermost iterations for solving the system of linear equations corresponding to one iteration of [6]. In our case, this system was solved using a Gauss-Seidel algorithm; these innermost iterations, as well as iterations, for the calculation of the MAP estimates, were continued until the following condition was reached: Outermost iterations were stopped when the following condition was satisfied: where gkt2! was the estimate of the covariance component relative to sire/maternalgrandsire (k = 1) and dam within maternal grandsire (k = 2) effects and pertaining to traits i and j at the round t.

frequencies of the deformities
Valgus and varus deformities were first diagnosed at 3 weeks; incidences at this age are reported in table I (severity was not recorded at this early age). Frequencies at the age under study, 6 weeks, are reported in table II. While valgus incidence was already high at 3 weeks in both lines and sexes, varus deformity rarely appeared at this age. At 6 weeks of age, both sexes differed in the incidence of valgus deformity, which was twice as high in males (respectively equal to 63.0 and 63.8 in lines A and B) as in the females (respectively equal to 33.8 and 35.1 in lines A and B); this resulted from different frequencies of severe cases, which were more than 30% in the males and close to 6% in the females. Although prevalence of varus increased with age, this defect was markedly less developed than valgus deformity at 6 weeks. Total frequency of varus defects was rather homogeneous between sexes; incidences varied, according to line and sex, from 7.3 to 12.4% (table II). Moreover, very few severe cases of varus deformities were diagnosed at 6 weeks.

Estimates of variance and covariance components
Estimates of variance and covariance components are presented in

DISCUSSION
Deformities measured in the present study could not be considered as ordered categorical traits; therefore, use of the threshold model developed by Gianola and Foulley (1983) was precluded here. In contrast, the multiple logistic model is appealing for unordered categorical traits. Whereas many applications of the threshold model in animal breeding with a mixed model are now reported in the literature, use of the logistic model in this area is rare; it has been restricted until now to cattle breeding for analysis of survival data (measured as an 'all-or-none' trait) by De Lorenzo and Everett (1986) and for analysis of unordered categorical responses with a mixed model by Gianola (1980). In the latter case, logits pertained to ratios of observed counts and not to ratios of probabilities. Previous studies of the genetics of twisted legs syndrome have not taken advantage of advances in methods of analysis of discrete traits. Moreover few estimates of genetic parameters are available in the literature. However, all studies conclude that twisted legs are heritable. Hartmann and Flock (1979) estimated the heritability of twisted legs by the analysis of variance proposed by Robertson and Lerner (1949). Estimates ranged between 0.10 and 0.30 for the first period under study, and between 0.40 and 0.51 for the second period, for male and female offspring respectively. Leenstra et al (1984) compared three lines selected for increased body weight at 6 weeks or decreased incidence of twisted legs or for both. After three generations of selection, they obtained a significant' decrease of twisted legs in both lines selected against this disease. Mercer and Hill (1984) provided estimates of heritability of 'splay' and 'bow' deformities, analogous to valgus and varus respectively, in three meat-type strains. Both full and half-sibs analyses were conducted. Using the proband method (Falconer, 1965), heritability of splay (or valgus) deformity ranged between 0.14 and 0.29 when estimated from half-sib analysis; the mean value for this angulation was 0.21. When computed from dam components, heritability estimates ranged between 0.20 and 0.26 with a mean value of 0.23. For bow leg (or varus), heritability estimates were more variable between lines and ranged between 0.05 and 0.26, with a mean value of 0.11 when obtained from paternal half-sib analysis. When full-sibs were considered, estimates were between 0.18 and 0.24, with a mean value of 0.20. With the exception of valgus in line A based on the sire/maternal-grandsire component, our heritability estimates were slightly higher than those reported by Mercer and Hill (1984).
In their study, Mercer and Hill (1984) indicated that leg problems (including leg and keel defects) were positively correlated but noticed the possible exception of bow and splay deformities. The means of the estimates of the genetic correlation obtained by an analysis of variance were -0.12 and -0.06 when based on halfsibs and full-sibs respectively. Our estimates of the genetic correlation between valgus and varus also appeared to be very small or slightly negative. These values cannot be attributed to the fact that valgus and varus are mutually exclusive: the logistic model takes this into account. Indeed, when simulating data (on two traits) whose underlying variates were positively correlated, the mean of the estimates of this parameter was positive (and close to the true value, when incidences of the two traits were higher than 0.10). When considering leg defects, a positive genetic correlation would correspond to a favourable situation in which selecting against one of the deformities would allow a decrease in incidence of all the various deformities (and thus an increase in the incidence of healthy animals as shown by the expression for the probability of this category given in [2b]). This situation was not encountered because valgus and varus seemed, according to our estimates, to be approximately independent. This was also supported by studies on clinical and anatomical differences between the two deformities (Leterrier and Nys, 1992) which suggested that they correspond to different aetio-pathogeneses. However, the results obtained here cannot be considered as clear evidence that valgus and varus should be recorded separately for selection purposes. This should be established from considerations of selection responses which are beyond the scope of the present paper. It should be noted that, for the purpose of this study, unilateral and bilateral varus were grouped. Indeed, previous analyses suggested substantial genetic correlations between both defects: estimates obtained with a sire model ranged between 0.69 and 0.75. This result supports the observations of Leterrier and Nys (1992) of the clinical similarities between these two deformities. However, more intensive research on larger data sets and considerations of responses to selection should make it clear whether or not our grouping procedure was fully optimal.
A sensible question to be answered is whether mild or severe expressions of the same deformity should be considered as different traits and the list of logits modified accordingly. Preliminary analyses according to a simplified model (ie, a sire model) were carried out in the males for each of the two selected lines. Estimates of the genetic correlation between susceptibilities to mild or severe defects appeared very high, between 0.88 and 0.91 in line A, and between 0.73 and 0.96 in line B. These results justify classification of mildly affected birds as ill. It made it possible to group together mild and severe symptoms, which avoided very extensive computing costs and the consequences of extreme levels of incidence observed for varus deformities.
In this analysis leg defects of both males and females were considered as the same trait. Since the incidence of deformities and the balance between mild and severe forms varied significantly according to sex (as shown specially for valgus in table II), suggesting genetic determinisms possibly different between sexes, a preliminary analysis was conducted to estimate genetic correlations between susceptibilities of males and females. A sire model was used because convergence on genetic parameters could not be reached with a sire, maternal grandsire and dam model. Whatever the deformity, genetic correlation between susceptibilities of males and females was very high: estimates were between 0.71 and 0.88 in line A, and between 0.71 and 0.91 in line B. This result suggests that genes common to both sexes are involved in susceptibility to twisted legs syndrome, although the moderate accuracy of our estimates (standard deviations were between 0.05 and 0.16 in line A, and between 0.04 and 0.24 in line B) does not exclude the possibility of specific genetic effects for each sex, and particularly of sex-linked effects. Sexes were pooled together when considering genetic effects but distinguished when considering fixed effects.
A discrepancy between sire/maternal-grandsire and dam heritabilities was observed; the latter were higher. In the same way, in the study of Mercer and Hill (1984), heritability obtained from full-sibs correlations was usually greater than heritability given by half-sibs analysis. Because permanent environmental maternal effects are unlikely in present systems of chick rearing, this trend could suggest the possible existence of genetic maternal additive effects and/or dominance effects.
Genetic maternal effects could be assessed from a sire, maternal grandsire and dam within maternal grandsire model in which sire and grandsire variances would be distinguished, as in Manfredi et al (1991). However, the estimation of variancecovariance components for direct and maternal effects requires the estimation of the covariance between sire and grandsire effects, which is available only if some male ancestors are sire and grandsire at the same time. Here, there were very few (less than five) such ancestors. Further investigations on more comprehensive data sets should allow testing of the hypothesis of the presence of maternal or dominance effects.
No explanation has been found yet for the origin of the twisted legs syndrome. It is highly probable that some of the genes coding for bone, tendon or cartilage growth and quality may be involved in variations of susceptibility to these disorders. Polygenic determinism was assumed in the present study; testing the hypothesis of the involvement of a major gene would be interesting although no evidence of such stems from our analysis. CONCLUSION These results indicate that selection against the various types of twisted legs can be effective. It is likely that a simplified selection scheme based on the presence or absence of twisted legs would reduce valgus deformity because of its larger incidence, while changes in incidence of varus would most probably be small or even unfavourable because of the negative genetic correlation between the two defects. However, further work is now needed to derive the expected response to selection on unordered categorical traits. Indeed, genetic parameters on the logit scale are not connected to response to selection in the same way as in the Gaussian situation, primarily because the observed response arises from a competition between latent variates. Moreover, before defining the best breeding strategy to implement in practice, further investigations should be carried out in order to estimate genetic relationships between these deformities and traits of economical importance, such as growth rate and meat conformation traits.