Genetic parameters for litter size in Black Slavonian pigs

The objective of this study was to estimate genetic parameters for litter size of Black Slavonian pigs using the repeatability, multiple trait, and random regression models, and to consider the possibility to increase litter size in Black Slavonian pigs by selection. A total of 4733 litter records from the first to the sixth parity from sows that farrowed between January 1998 and December 2010 were included in the analysis. Individual record consisted of the following variables: breeding organisation (eight regions), parity (1-6), service boar, and farrowing season (month-year interaction). Estimation of all the covariance components with three different models was based on the residual maximum likelihood method. Estimate of additive genetic variance and heritability for number of piglets born alive with repeatability model was 0.23 and 0.10, respectively. Estimates of additive genetic variance with multiple trait and random regression model were in a wider range from 0.05 to 0.65 across parities, and heritabilities were estimated in the range between 0.03 and 0.26. Estimates of phenotypic and additive genetic correlations were much smoother with random regression model in comparison with multiple trait model. Due to unexpected changes of variances along trajectory obtained with multiple trait and random regression model, the best option for genetic evaluation of litter size for now could be the use of repeatability model. With increasing number of data with proper data structure alternative modelling of litter size of Black Slavonian pig using multiple trait and random regression model could be taken into consideration.


Introduction
The Black Slavonian pig is a Croatian autochthonous breed developed at the end of the 19 th century. Until the middle of 20 th century this breed was the most widespread pig breed in the country. However, due to the introduction of modern breeds the Black Slavonian breed became less and less used up to the point of facing extinction in the middle of 1990s when the number of animals was reduced to only 60 sows and five boars (Uremovic et al., 2000). Since the year 1996 population of Black Slavonian pigs started to increase, primarily due to state subsidies and the efforts of individual pig breeders. Pigs of this breed are traditionally kept in outdoor system and its meat is used for the production of local dry meat products of high quality like "Slavonian kulen" sausage (Karolyi et al., 2008). During this period animal performance databa-se for reproduction traits of Black Slavonian pigs was established to provide information for management and selection. Namely, litter size of Black Slavonian pigs is lower compared to the modern pig breeds, on average six to seven piglets are born alive (Sencic et al., 2001). The increase of litter size of Black Slavonian breed is crucial for the market-oriented production of sufficient amounts of quality meat products. In the last few years, population of Black Slavonian pigs increased to more than 800 sows and 100 boars (Croatian Agricultural Agency, 2012), which opens a possibility to renew breeding programme for this local breed and to determine genetic parameters for litter size as one of the most important traits in pig breeding.
Knowledge of genetic parameters for litter size in pigs is necessary to estimate accurate breeding values, to optimize breeding schemes, and to predict response to selection (Roehe & Kennedy, 1995). There are several approaches regarding the method for genetic evaluation of litter size. They mainly depend on genetic correlations between litter size in subsequent parities, size of data set, capacity of computing etc. Due to simplicity many pig breeding organizations use repeatability model for the estimation of breeding values for litter size (e.g., Alfonso et al., 1997;Noguera et al., 1998). A repeatability model assumes that genetic correlations between litter size in subsequent parities are equal to one and that the variance is constant along this trajectory. On the other hand, multiple trait analysis supposes subsequent observations to be different traits. Estimates of genetic correlations between litter sizes in different parities are often lower than one (Alfonso et al., 1994;Irgang et al., 1994), especially between the first and later parities (Serenius et al., 2003). The low correlations between parities indicate that partly different genes are responsible for litter size in different parities. Therefore, multiple trait analysis is preferred in such situation in order to increase the efficiency of selection. Beside a simple repeatability model and a full multiple trait model, a possible approach of handling repeated observations would be to use covariance functions and random regression (RR) model. Main advantages of the RR approach in comparison to the multiple trait approach are: smaller number of parameters to describe the longitudinal process, smoother covariance estimates, as well as a possibility to estimate covariance components and predict breeding values at any point along the trajectory. Litter size in pigs is a repeated measure and could also be considered as a longitudinal trait (Schaeffer & Dekkers, 1994). Although litter size differs from typical longitudinal traits such as growth, primarily because parity order is not a continuous variable, RR could be also applied for the genetic evaluations (Lukovic et al., 2004;Schaeffer, 2004).
The objective of this study was to estimate genetic parameters for litter size of Black Slavonian pigs using different approaches, such as repeatability, multiple trait, and random regression models, and to consider possibility to increase litter size in Black Slavonian pigs by selection.

Material and methods
Data for this study were obtained from the Croatian Agricultural Agency. Litter records from the first to the sixth parity from sows that farrowed between Ja-nuary 1998 and December 2010 were included in the analysis. Litter size was presented as the number of piglets born alive (NBA). In preliminary analysis of pedigree around 10% of records were deleted due to inconsistencies (bisexual animals, wrong identification, etc.). Beside the identification number of a sow and NBA, individual record consisted of the following variables: breeding organisation (eight regions), parity (1-6), service boar, and farrowing season (monthyear interaction). All months with less than five records were joined to next closest month, while boars with less than five litters were joined in one group separately by region. After data editing, there were 132 seasons, and 164 service boars, and a total of 4,733 litter records (Table 1).
Average NBA was 6.12 piglets with minimum in the first parity of 5.77 piglets and maximum in the fourth parity of 6.42 piglets with slightly higher variation in parities with higher means. The distribution of NBA by litters showed that the majority of sows (> 60%) farrowed between five and seven piglets born alive. Further, more than nine piglets born alive was recorded in less than 5% of all litters.
Litter records were obtained from 1,649 animals with 2.9 litters per sow ( Table 2). The number of sows with complete records of all six parities was 272 which comprised 16% of the observed population. The pedigree file was prepared for three generations and inclu-  ded 1862 triplets of animal, sire and dam with 179 sires and 581 dams. Proportion of inbred animals in data set was 3%.
The choice of the fixed part of the model for NBA was based on significance, proportion of variation explained by effects, as well as the simplicity and interpretation of a model. Significance was tested at the level of p < 0.05 using the GLM procedure (SAS Inst. Inc., 2001). Preliminary analysis using the repeatability model showed that following fixed effects should be included in the statistical model for NBA: farrowing season as year-month interaction, service boar, and parity. Three classes were considered for parity effect: the first for parity one, the second for parity two, and the third for the later parities. The following repeatability model was used: where y ijkl is NBA, S i farrowing season, B j serviceboar, p k parity, p il permanent environmental effect, and a il additive genetic effect, and e ijkl residual term. In the matrix notation repeatability model can be written as: where y is a vector of observations, X is incidence matrix for fixed effects, β is a vector of unknown parameters for fixed effects, Z p and Z a are incidence matrices for permanent environmental and additive genetic effect, p and a are the corresponding vectors of parameters for random effects, and vector e presents residual. The following covariance structure was assumed: , where A is the numerator relationship matrix, while I p and I e are identity matrices. Covariances between random effects were assumed to be zero. The multiple trait (MT) analysis was performed in two ways. First, litter records by parity were treated as different traits resulting in the six-trait model [4], which implies that the model consisted only of additive genetic effect as a random effect. The appropriate MT model was: with the following covariance structure: where G 0 is a matrix of additive genetic covariances for NBA between parities and R 0i is a diagonal matrix of parity class specific residual variances for animal i. Residuals within and between animals were assumed independent. Second, litter records were grouped in three classes: the first parity, the second parity, and the third and later parities. The model for the third group of parities had permanent environmental effect included as in the repeatability model [1] but the whole model was treated as MT model [4]. The following RR model with different order of Legendre polynomials was fitted: The set of fixed effects used in the RR model analysis included the same effects as for the repeatability model [1]. Permanent environmental effect of a sow (s = 1) and direct additive genetic effect (s = 2) were fitted as random regressions on parity using Legendre polynomials (LG). Legendre polynomials from linear (LG1) with two terms to cubic (LG3) with four terms were fitted. Covariate in LG was standardized parity with a range from -1 to +1: where p min is the first parity and the p max sixth parity. The RR model in matrix notation is the same as [2], however the list of unknown parameters and corresponding assumptions for random effects differs: where K 0p and K 0a are covariance matrices between random regression coefficients for the permanent environmental effects and additive genetic effects, respectively. From the estimates of K 0a additive genetic covariances for NBA between parities C a were calculated as: where Φ is a matrix of covariates for Legendre polynomials on standardized parities. Estimation of all the covariance components with three different models was based on the residual maximum likelihood method using the VCE-5 software package (Kovac et al., 2002).

Results
Estimates of variance components and their ratios of phenotypic variance using a repeatability model are shown in Table 3. Total phenotypic variance as a sum of all estimated variance components was 2.39, whereas permanent environmental variance was 0.04 and additive genetic variance was 0.23. Ratio of permanent environmental effect of phenotypic variance was only 0.02, while heritability was 0.10.
Additive genetic and phenotypic correlations for NBA among the first six parities estimated using the MT analysis are presented in Table 4. Phenotypic correlations between different parities were generally low (< 0.40), in particular between the first and other parities (< 0.13) and between any parity and the last two (the fifth and sixth) parities. The highest phenotypic correlations were between the second, third, and fourth parity. Genetic correlations between the adjacent parities were higher and ranged between 0.01 and 0.96.
The highest values were observed among the pairs of the first four subsequent parities and remarkably lower for last two parities. Pattern of decrease of genetic correlation along the distance between parities was not smooth, due to oscillations for some pairs of parities.
Estimates of variance components and ratios of phenotypic variance by parities are presented in Table 5. Total phenotypic variance as a sum of all variance components by parities ranged between 1.59 and 2.35, with the increasing values up to the third parity and decreasing values afterwards. Direct additive genetic variances ranged between 0.06 in the first parity and 0.48 in the third parity with heritabilities between 0.03 in the first parity and 0.21 in the third parity. Similar trends were observed for the estimates of variance components and ratios of phenotypic variance by parities with the third and later parities grouped (Table 6). However, additive genetic variances and heritabilities tended to be lower with this analysis. Similar trends were observed also with the analysis of the first parity and the second and later parities grouped (results not shown). Additionally, using the second approach estimation of permanent environmental effect was obtained as very small proportions in the phenotypic variance.   Random regression models were fitted with Legendre polynomials up to the third order (from LG1 to LG3). Modelling the variance of the permanent effect along the trajectory of parities did not show considerable partitioning of variability into different components of individual production curves as described by the analysis of eigenvalues (Table 7). Namely, the zero-th term explained 99% of variability with the LG1 model and 96% of variability with the LG2 model. Only with the LG3 model 27% of variance of permanent environmental effects along the parity trajectory was explained by the first term describing the variation in the slopes of individual production curves of sows. Further two terms did not account for a significant amount of variation. For the additive genetic effect fitting the LG1 model resulted in 65% of variability explained by the zero-th term describing the level of individual production curves, while the first term explained 35% of variability describing the slope of individual production curves (Table 7). With further partitioning using the model LG2 and LG3 the amount of variation due to the first term remained practically constant, while the second and the third term describing the shape of individual production curves increased on the account of the decrease in the explained variability by the zero-th term.
Estimates of phenotypic correlations with the RR model (Table 8) were in agreement with the results obtained with the MT model. They ranged between 0.12 and 0.16 for all pairs of subsequent parities. Decrea-Genetic parameters for litter size in Black Slavonian pigs 93   sing tendency of phenotypic correlations with the increasing distance between parities was also similar as in the MT model. For the additive genetic effect correlations between parities were also similar to the MT model, although the genetic correlations changed more smoothly with the increasing distance when using the RR model. Between the subsequent parities genetic correlations ranged from 0.74 to 0.99. The smallest genetic correlations were between the first and the higher parities.
Estimates of variances as well as ratios of the phenotypic variance for random effects obtained with the RR model (Table 9) were generally in agreement with the estimates obtained using the MT model (Table 5). Estimates of permanent environmental variance ranged between 0.04 and 0.24 and were higher than with the repeatability model, except in the first and the sixth (last) parity. Estimates of additive genetic variance with the RR model differed the most from estimates from other models, with considerably higher values in the first and the sixth parity and lower values for the intermediate parities in comparison to the MT model. Estimate of residual variance was lower with the RR model for the first parity than with the MT model, while they were larger for later parities, with higher differences for the second, third, and fourth parity.

Discussion
Estimate of phenotypic variance for NBA using the repeatability model was lower (2.39) than many estimates obtained in the analysis of modern pig genotypes (e.g., Ferraz & Johnson, 1993;Crump et al., 1997;Hanenberg et al., 2001;Chen et al., 2003). Generally smaller litters characteristic for local pig breeds, as well as piglet losses during the first few days after farrowing reduce the phenotypic variance for NBA (Kaplon et al., 1991). Recording procedure where piglets are counted few days after farrowing at pig family farms additionally reduce phenotypic variance for NBA. Although, late recording reduces phenotypic variation, it can also be interesting alternative to improve litter size at weaning and piglet survival (Su et al., 2007). The first step in achieving higher litter size for this breed should be to ensure adequate production conditions, such as good health status of herds, and to overcome other problems related with extensive production systems. Moreover, small populations are often characterised by high level of inbreeding. However, increase of population of Black Slavonian pig in recent years caused decline in proportion of inbreed animals (Uremovic et al., 2000), which is confirmed by low  proportion (3%) of inbreed animals in our data set. Such a low value of inbreeding should be considered with caution due to partial knowledge of the pedigree. On the other hand, considerable increase of inbreeding might be expected if selection using animal model without restrictions is performed. Additional efforts through creation of mating schemes are necessary to maintain and increase the size of the population in order to avoid increase in the level of inbreeding.
Permanent environmental effect is a characteristic of each individual sow and the accuracy of estimating this effect depends on the number of litters per sows obtained. Estimate of permanent environmental variance was generally small in this study when using the repeatability model and in the first and the last parity when using the RR model. Chen et al. (2003) estimated the ratio of permanent environmental effect between 0.06 and 0.08, while even higher ratios (between 0.16 and 0.17) were obtained in the study of Ferraz & Johnson (1993).
Estimate of additive genetic effect using the repeatability model in this study (0.23) was twice lower than those estimated in study by Chen et al. (2003), but on other hand heritability estimate of 0.10 is in agreement with estimates obtained on modern pig genotypes (Haley et al., 1988;Rothschild & Bidanel, 1998). Smaller estimates of heritability than those found in this study were found in some other local European pig breeds, as Iberian or Cinta Senese pig (Pérez-Enciso & Gianola, 1992;Crovetti et al., 2005).
Low genetic correlations between litter size in different parities indicates that phenotype records from different parities should be treated as different traits (Irgang et al., 1994;Roehe & Kennedy, 1995). Direct additive genetic correlations obtained with MT analysis for six traits are lower than in Hanenberg et al. (2001), who reported an increase from 0.79 between the first and the second parity to 0.96 between the fifth and the sixth parity. Considerably lower values of additive genetic correlations between the pairs of the last three parities in this study suggest the potential issue with the data structure to apply the MT model. Possible reasons for such a low values for additive genetic correlations between pairs of last parities might be smaller amount of data in later parities and imprecise recording of parity order. Moreover, standard errors were several times larger than estimated correlations, indicating that additional efforts should be done in improving recording of parity order and increasing number of records in later parities.
Heritability estimates from the second to the fourth parity using the MT model were substantially higher than the estimate obtained with the repeatability model. In some previous studies heritabilities obtained with the MT model were generally smaller than with the repeatability model (Alfonso et al., 1997;Hanenberg et al., 2001). Small additive genetic variances and consequently low heritabilities were found in this study in the first and the last parity. Lower heritability estimated in the first parity with the MT model is in line with the work of Roehe & Kennedy (1995). Exceptionally high estimates of heritabilities in higher parities were found by Kovac & Sadek-Pucnik (1997), and they explained them as a consequence of numerical problems.
Random regression model provides a possibility of modelling the change in variances along the trajectory. To model only the variance in the level of NBA the repeatability model is sufficient. However, to model the variance in the slope and the general shape of individual production curves and their components RR models can be used. In the case of substantial genetic variation for the shape of production curve of litter size over parities selection could be enforced not only on the level but also on the shape of production curve (Lukovic et al., 2007). Fernández et al. (2008) applied RR model and found substantial genetic variability for the shape of production curves for litter size in Iberian pigs making such selection possible.
The proportion of variance explained by the constant (zero-th) term for permanent environmental effect in the RR model accounted almost all variability (96-99%) using Legendre polynomials LG1 and LG2. Only with using LG3 polynomial larger part of variation was captured by the first component which determines the slope of production curve for litter size. Results obtained for permanent environmental effect together with the fact that very small variation of additive genetic effect was explained with the LG3 model supports previous results from Lukovic et al. (2004) that LG2 with three regression coefficients might be sufficient to explain majority of variation in individual production curves. For this reason only results from the LG2 RR model were presented in this study.
Estimates of phenotypic and additive genetic correlations using the RR model with LG2 were smoother in relation to MT, which is a general property of RR models (Schaeffer, 2004). This is especially useful in the case where MT analysis give doubtful results due to worse data structure. However, with the RR model the estimates of variance components (especially for the additive genetic variance) in this study showed opposite trends as with the MT model -much higher values in the first and the sixth parity than in the intermediate parities. This suggests that RR model is not adequate for application in the Black Slavonian pig population based on the results using currently available data. In order to improve data quality and enable the use of more advanced statistical models, it is necessary to overcome problems which are linked to extensive production systems, such as imprecise recording of data (e.g. parity order). This could be achieved by increasing a level of control of data collection and emphasizing importance on data quality for results of the selection. Nevertheless, with the heritability of 0.10 obtained by repeatability model and high phenotypic coefficient of variation, there is opportunity to carry out successful selection for the litter size.
According to the practicality and computing requirements of genetic evaluation of litter size in pigs the repeatability model is the method of choice. With several measurements per animal, it requires less computational efforts and smaller number of parameters to estimate than the MT or RR models. However, if genetic correlations between litter sizes in different parities are lower, MT analysis is preferred (Alfonso et al., 1994). In such case RR model can be more robust compared as MT models often require several runs with different starting values to achieve convergence (Lukovic et al., 2004). RR model needs substantially less computing time to obtain comparable results to MT model in the case of proper data structure. This condition was not fulfilled in the analysis of Black Slavonian breed and although the estimates of genetic and phenotypic correlations were smoother the trend of change of additive genetic variance deviates from expectations making the use of RR model in this dataset questionable.
For estimation of genetic parameters of litter size in pigs different methods can be used. With the assumption of unit genetic correlations between parities and constant variance along parities the repeatability model is the model of choice. In the situation when genetic correlations between parities are lower, the multiple trait model is suggested. Random regression model can be used instead of the multiple trait model in such cases to lower the number of parameters to estimate. However, sufficient data with proper structure is needed to apply random regression model and if this requirement is not met as is often the case with less numerous local pig breeds, the best option is to use the repeatability model. However, additional efforts should be made in order to improve data quality. Improved data quality with emphasis on precise recording would result with more eff icient breeding programme for this breed. In this way, genetic improvement and economical effect might stimulate breeders to choose this particulate breed. Future steps might include obtaining information about connectedness between herds and breeding organisations included in breeding programme for Black Slavonian pig.